University of Southern California Dissertations and Theses

High throughput computational framework for synthesis and accelerated discovery of dielectric polymer materials using polarizable reactive molecular dynamics and graph neural networks

(USC Thesis Other)

High throughput computational framework for synthesis and accelerated discovery of dielectric polymer materials using polarizable reactive molecular dynamics and graph neural networks

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content
High Throughput Computational Framework for Synthesis and
Accelerated Discovery of Dielectric Polymer Materials using
Polarizable Reactive Molecular Dynamics and Graph Neural
Networks

By

Ankit Mishra

A Dissertation Presented to the
FACULTY OF THE GRADUATE SCHOOL
UNIVERSITY OF SOUTHERN CALIFRONIA

In Partial Fulfillment of the
Requirements for the Degree
DOCTOR of PHILOSOPHY
(Chemical Engineering)

August 2021

Copyright 2021 Ankit Mishra

ii
Acknowledgments
I would like to thank my advisor, Professor Priya Vashishta, for constantly guiding me through
the research and providing me his valuable insights. Without his excellent supervision and
persistent mentoring, I would not have been able to complete this thesis. I would also like to thank
my Ph.D. committee: Professor Rajiv Kalia and Aiichiro Nakano for finding time from their hectic
schedule to provide supervision, support, and insightful comments for the thesis.

I am also immensely grateful to Professor Ken-ichi Nomura for teaching me molecular
dynamics simulation and helping me with every problem I had in my projects. Additionally, I am
thankful to Dr. Pankaj Rajak, Aravind Krishnamoorthy and Subodh Tiwari for their mentorship
and supervision throughout my research work at USC. I would also like to thank my current
colleagues at USC, Tom Linker, Nitish Baradwaj, Anikeya Aditya, Liqiu Yang, Ruru Ma for all
for their help and support. I am immensely grateful to Ms. Patricia Wong, too, for taking care of
all the logistics and for being very patient with me. She has been extremely helpful and
understanding. The USC HPC center and the Argonne National Lab, where I ran all my
simulations, I am also thankful to their support team as they had resolved numerous issues related
to my programs. My journey at USC would not have been a joyful ride without the help and support
of my friends, thus, I would like to take this opportunity to thank my friends, Jitin Singla, Palash
Goyal, Ayush Jaiswal, Nitin Karma, Amal Thomas, Akshay Gaur, Akarsh Goyal, Rohan Rout,
Shivani Singh and Pranali Khare for all the wonderful time we have spent together at USC.

Finally, I am extremely thankful to my family for their unconditional support, constant
encouragement, and enduring curiosity and interest in my scientific endeavors. I am especially
grateful to my wife, Ankita, who has been there from the start of this journey celebrating my
successes and supporting me during the tough times, without her support this would have not been
possible.

iii
Table of Contents
Acknowledgments ......................................................................................................................... ii
List of Figures ............................................................................................................................... vi
List of Tables ................................................................................................................................. ix
Abstract .......................................................................................................................................... x
1. Introduction ............................................................................................................................... 1
1.1. Computational framework for polymer synthesis to study dielectric properties using
polarizable reactive molecular dynamics .................................................................................... 2
1.2. Dielectric Polymer Property Prediction Using Recurrent Neural Networks with
Optimizations ............................................................................................................................... 3
1.3. Accelerated Discovery of Dielectric Polymer Materials using Graph Convolution Neural
Network ....................................................................................................................................... 4
2. Molecular Dynamics .............................................................................................................. 5
2.1 Molecular dynamics simulation ............................................................................................. 5
2.2 Interatomic force fields .......................................................................................................... 6
2.3. Reactive forcefield (ReaxFF) ................................................................................................ 8
2.3.1. Bond Order and Bond Energy ....................................................................................... 8
2.3.2. Lone pair energy .......................................................................................................... 10
2.3.3. Overcoordination ......................................................................................................... 11
2.3.4. Undercoordination ....................................................................................................... 12
iv
2.3.5. Valence Angle Terms .................................................................................................. 12
2.3.6. Torsion angle terms ..................................................................................................... 14
2.3.7. Hydrogen bond interactions ......................................................................................... 15
2.3.8. Nonbonded interactions ............................................................................................... 16
2.3.9. Coulomb Interactions .................................................................................................. 17
3. Computation Framework for Polymer Synthesis to Study Dielectric Properties Using
Polarizable Reactive Molecular Dynamics ................................................................................ 18
3.1. Background ......................................................................................................................... 18
3.2. Methodology ....................................................................................................................... 21
3.2.1. Charge State Aware ReaxPQ+ Model ......................................................................... 22
3.2.2. Amorphous Model Creation ........................................................................................ 26
3.2.3. Estimation of Dielectric Constant ................................................................................ 30
3.2.4. Role of morphology on the Dielectric Constant .......................................................... 31
3.2.5. System Preparation for Dielectric Constant Calculation ............................................. 32
3.3. Results and Discussion ....................................................................................................... 35
3.3.1. Validation of high frequency dielectric constant for model polymer systems ............ 35
3.3.2. Dielectric behavior of TDI-EDR148 and PDTC-HK511 ............................................ 36
3.3.3. Dielectric behavior of DVS-BCB and FPE ................................................................. 37
3.4. Conclusion .......................................................................................................................... 38
4. Polynorbornenes Defect Incorporation and Analysis .......................................................... 39
4.1. Background ......................................................................................................................... 39
4.2. Methodology ....................................................................................................................... 41
v
4.2.1. PNB Design Space ....................................................................................................... 41
4.2.2. High Throughput Workflow Design ............................................................................ 42
4.3. Dataset Generation ............................................................................................................. 44
4.4. Defect Analysis ................................................................................................................... 45
4.4.1. Pendant group analysis in various POFNBs ................................................................ 46
4.4.2. Free Volume Analysis in PFNBs ................................................................................. 49
4.5. Conclusion .......................................................................................................................... 54
5. Attention based Graph Neural Network Analysis of PNB dataset ..................................... 56
5.1 Background .......................................................................................................................... 56
5.2. Graph Attention Network (GAT) Computations ................................................................. 58
6.3. Graph Attention Network (GAT) Architecture ................................................................... 63
5.4. Results ................................................................................................................................. 66
5.5. Conclusion .......................................................................................................................... 69
Conclusion .................................................................................................................................... 71
REFERENCES ............................................................................................................................ 74

vi
List of Figures

Figure 2. 1. Dimensionless Lennard-Jones potential for inert gases. ............................................. 7
Figure 2. 2. Energy terms included in ReaxFF where Ebond is the bond energy, Elp is the lone pair
energy, Eover and Eunder are the over- and under- coordination energies, Eval is the valence angle
energy, Epen is the penalty energy, Ecoa is the three-body conjugation energy, Etors is the torsion
angle energy, Econj is the four-body conjugation energy, Ehbond is the hydrogen bond interaction
energy, EvdWaals, ECoulomb and Elg are the van der Waals, Coulomb energies and dispersion
energy, respectively. ........................................................................................................................ 8

Figure 3. 1. Computational framework of polymer dielectric-behavior estimate. (a) Synthesis of
amorphous polymer system. SMILES string and Open Babel are used to create a single polymer
chain. Multiple polymer chains are placed at sufficiently large distance from each chain in a
simulation system. The simulation system is subjected to a number of consolidation and
relaxation steps until the system reaches to a desired density. (b) ReaxPQ+ model development
that involves various QM-based calculations and validations such as atomic polarization and
charges for constituent atomic species, and dielectric constant estimate. .................................... 22
Figure 3. 2. Schematic of the response of core (black) and shell (red) charges to an external
electric field in ReaxPQ+ model. .................................................................................................. 23
Figure 3. 3.Amorphous-polymer synthesis of TDI-EDR148 involves accurate SMILES
representation of the polymer followed by its generation using Open Babel and final
consolidation to achieve desired density. ...................................................................................... 27
Figure 3. 4. The consolidated structures for (a) PDTC-HK511, (b) FPE and (c) DVS-BCB are
generated from their SMILES string representation using Open Babel and ReaxPQ+. .............. 28
Figure 3. 5. (a) Weak scaling parallel efficiency of ReaxPQ+ with iso-granular workload scaling
of 25,576 atoms on P cores (P = 64, ..., 131,072) of Theta. (b) Strong scaling parallel efficiency
(correct this) of ReaxPQ+ with a fixed problem size 50,331,648-atom on P cores (P = 2,048, …,
32,678) of Theta. The idea and measured efficiency are shown in black and blue color lines
respectively. ................................................................................................................................... 30
Figure 3. 6. Total dipole moments (in Debye) at (a) 100K, (b) 200K, (c) 300K and (d) 400K for
TDI-EDR148 show different variance in total dipole moment across frames. Dipole moment is
computed over a large number of RMD trajectories, the variations in x, y and z directions shown
in red, green and blue respectively show an increasing trend. ..................................................... 33
Figure 3. 7. Total dipole moments (in Debye) at (a) 100K, (b) 200K, (c) 300K and (d) 400K for
HK511 show different variance in total dipole moment across frames. Dipole moment is
computed over a large number of RMD trajectories, the variations in x, y and z directions shown
in red, green and blue respectively show an increasing trend. ..................................................... 34
vii
Figure 3. 8. Dielectric constants of (a) PDTC-HK511 and (b) TDI-EDR148 polymer systems
(blue dots) as a function of temperature. The increasing trend as well as predicted values show a
good agreement with experimental values (shown in red). ........................................................... 36
Figure 3. 9. Dielectric constants of (a) DVS-BCB and (b) FPE polymer systems as a function of
temperature. The predicted values show a good agreement with experiments[66, 67, 83-87]. .... 37

Figure 4. 1. PNB with heteroatom X and various functional groups which can be incorporated as
defects. ........................................................................................................................................... 42
Figure 4. 2. High throughput framework for generating dielectric constant values for various
PNB derivatives. ............................................................................................................................ 43
Figure 4. 3. The computed dielectric constants are shown here as, (a) the overall dielectric
constant values as computed based on ReaxPQv forcefield and our high throughput workflow,
and (b) there characterization based on the presence and absence of different kind of backbones
in this polymer. .............................................................................................................................. 45
Figure 4. 4. Various substitutions in PNB base polymer are shown in this figure. CF3 groups are
substituted from left to right in meta, para and ortho positions respectively. ............................... 46
Figure 4. 5. POFNB chains are repetitively consolidated and thermalized by shrinking the
simulation box to attain desired density of 1.3 gm/cc thus finally obtaining an amorphous
structure at the end. ....................................................................................................................... 47
Figure 4. 6. Pendant group rotation analysis over a period of 10ps shows highest rotation wrt to
para POFNB and lowest for ortho POFNB. The extent of rotations in these polymers directly
correlates to their dielectric property. .......................................................................................... 49
Figure 4. 7.The polymerization process for poly fluoro norbornene (PFNB) involves co-
polymerization with benzene containing PNB in different ratios to obtain various polymers with
different fluorine concentration. .................................................................................................... 50
Figure 4. 8.The experimental measurements for various fluorine concentrations as shown here
are corresponding to 0, 25 and 100% fluorine groups in the pendant groups. The measurements
show a decreasing trend with increasing fluorine concentration ................................................. 51
Figure 4. 9. Various fluorine containing structures are shown in this table. The structures chosen
for this study are made by changing the fluorine containing group in 8 chain polymer from 0 to
6. .................................................................................................................................................... 52
Figure 4. 10. The histogram represents Voronoi volume distribution around the pendant group
containing fluorine atoms, the trend clearly shows increasing free volume with increasing
fluorine concentration. .................................................................................................................. 53
Figure 4. 11. High frequency dielectric constant shows a decrease in net magnitude with
increasing fluorine concentration ................................................................................................. 54
viii

Figure 5. 1. Input graph representation is shown on the left and the 2-layer graph embedding
computation are shown on the left. The Graph Network performs two steps at each layer i.e. (1)
compute average information from neighbors and (2) apply the non-linearity through the
application of a neural network. .................................................................................................... 59
Figure 5. 2. The attention mechanism can be represented as shown in this computation graph. (a)
Attention coefficients 𝒆𝒗𝒖 are learned for every pair of neighbors and (b) the computed
attention coefficients are used to obtain the weighing factors 𝜶𝒗𝒖. The weighing factors are then
used to compute the final embedding 𝒉𝑨(𝒍) ................................................................................. 62
Figure 5. 3. A single layer of graph attention network many components as shown in this figure.
These layers can then be stacked against one another to obtain deep Graph Attention Networks
....................................................................................................................................................... 64
Figure 5. 4. The graph network training pipeline involves graph embedding representation
followed by a fully connected layer to make predictions across labels using L2 loss and RMSE
evaluation metric. .......................................................................................................................... 65
Figure 5. 5. The input polymer represented as a graph along with 20-dimensional feature vector
enables the GAT to learn a 64-dimensional graph embedding which is used to make predictions
using fully connected layer. ........................................................................................................... 66
Figure 5. 6. The fitting obtained from graph embeddings generated using GAT are shown in (a)
and (b) demonstrates the superior power of graph networks in learning the overall polymer
property. ........................................................................................................................................ 67
Figure 5. 7. The 64-dimensional graph embedding obtained from final layer of GAT is reduced
dimensionally and plotted using PCA and t-SNE based dimensionality reduction techniques. T-
SNE plot clearly shows the proximity of oPOFNB and PFNB. ..................................................... 69
Figure 5. 8. Local feature visualization for various polymer shows the importance of CF3 group
attached at ortho in (a) oPOFNB, meta position also effect the polymer property as shown in (b)
mPOFNB and (c) PFNB also shows important contribution wrt to fluorine containing pendant
group. ............................................................................................................................................. 69

ix
List of Tables
Table 3. 1. Computed atomic polarizability (Å𝟑) for various charge states of constituent atom
species. ........................................................................................................................................... 25
Table 3. 2. Optimized parameters (charge radius rc and spring constant ks) taking into account
the valence charge states for the constituent atoms in PDTC-HK511 and TDI-EDR148. ........... 26
Table 3. 3. Dielectric constants for crystal and amorphous PE. ................................................... 32
Table 3. 4. 𝝐∞ values for model polymer systems compared against QM values obtained using
Berry-Phase Method. ..................................................................................................................... 35

x
Abstract
Polymer based dielectric materials are used in a wide range of electrical storage devices. However,
the morphological complexity and diversity in polymer space poses a challenge in design and
identification of new polymer materials. Here, I would like to use classical molecular dynamics to
introduce a high throughput scalable computational framework based on charge state aware
reactive molecular dynamics to study dielectric response under external electric field and further
aid in identification of new polymers capable of showing desirable dielectric properties. The
overall work is divided into three parts as follows (1) Computational framework for polymer
property prediction using polarizable reactive molecular dynamics, (2) Recurrent Neural Network
based study of dielectric property generated using the computational framework, (3) Graph
Attention based analysis of polymer dielectric data to study underlying global as well local features
important for design of new polymers.
The increased energy and power density required in modern electronics poses a challenge
for designing new dielectric polymer materials with high energy density while maintaining low
loss at high applied electric fields. Recently, an advanced computational screening method coupled
with hierarchical modelling has accelerated the identification of promising high energy density
materials. It is well known that the dielectric response of polymeric materials is largely influenced
by their phases and local heterogeneous structures as well as operational temperature. Such inputs
are crucial to accelerate the design and discovery of potential polymer candidates. However, an
efficient computational framework to probe temperature dependence of the dielectric properties of
polymers, while incorporating effects controlled by their morphology is still lacking. In the first
work, I propose a scalable computational framework based on reactive molecular dynamics with
a valence-state aware polarizable charge model, which is capable of handling practically relevant
polymer morphologies and simultaneously provide near-quantum accuracy in estimating dielectric
xi
properties of various polymer systems. The predictive power of this framework is demonstrated
on high energy density polymer systems recently identified through rational experimental-
theoretical co-design. This scalable and automated framework may be used for high-throughput
theoretical screenings of combinatorial large design space to identify next-generation high energy
density polymer materials.
Despite the growing success of machine learning for predicting structure-property
relationships in molecules and materials, such as predicting the dielectric properties of polymers,
is still in its infancy. In the second work, we report on the effectiveness of solving structure-
property relationships for a computer-generated database of dielectric polymers using recurrent
neural network (RNN) models. The implementation of a series of optimization strategies was
crucial to achieving high learning speeds and sufficient accuracy: (1) binary and non-binary
representations of SMILES (Simplified Molecular Input Line System) fingerprints; and (2) back
propagation with affine transformation of the input sequence (ATransformedBP) and resilient
backpropagation with initial weight update parameter optimizations (iRPROP
-
optimized). For the
investigated database of polymers, the binary SMILES representation was found to be superior to
the decimal representation in respect of the training and prediction performance. All developed
and optimized Elman-type RNN algorithms outperformed non-optimized RNN models in the
efficient prediction of non-linear structure-activity relationships.
There is a great deal of interest in designing and identifying new dielectric polymer
materials that exhibit high power density as they have many applications in modern electronics
system. However, designing such polymers with desired properties is a challenging task due to
combinatorial large search space. Polynorbornene (PNB) is one such important amorphous
polymer system, which has potential applications as a high energy density polymer due to its high
xii
breakdown strength with low dielectric loss and high thermal stability. Moreover, electrical
properties of PNB can be significantly enhanced by incorporation of defects or synthesis with
controlled crystallinity by hydrogenation reaction, which involves experimental synthesis and
characterization of combinatorial large number of polymer systems to identify potential
candidates. In the third work, I propose an attention-based graph convolutional neural network
(GAT) model that can identify polymer systems capable of exhibiting increased energy and power
density. The GNN model is trained to predict dielectric constant for a polymer, where the training
data for the high frequency dielectric constant of the PNB polymers are computed via ab-initio
molecular dynamics simulation. This model can significantly aid experimental synthesis of
potentially new dielectric polymer materials which is otherwise difficult using simplistic statistical
procedures

1

1. Introduction

Polymers are an important class of materials owing to their chemical resistance, light weight, high
processability and ability to act as thermal and electrical insulators. Because of their useful
properties, polymer products are ubiquitous in our daily lives and are used in plastics, packaging,
military, aerospace and various other fields.[1] In insulation, isolation and electrostatic energy
storage technologies, polymers dielectric materials have become mainstays, since they are
relatively inexpensive and lighter than ceramic material alternatives. For example, biaxially-
oriented polypropylene (BOPP) is currently employed as state-of-the-art polymer dielectric films
in electrostatic energy storage capacitors. Despite having high dielectric strength and low dielectric
losses, BOPP exhibits relatively low dielectric constant and energy density.[2-4] Therefore,
concerted experimental and theoretical efforts have been made to develop new class of dielectric
polymer materials which possess low conduction losses under high electric fields and electrical
energy density and can be used in a wide range of applications in modern electronics and electric
power systems.[5-7]
Due to the rapid advances in machine-learning technologies, the closed-loop cycle
consisting of materials synthesis, characterization and computational modeling has been attracting
great attention to accelerate design and discovery of novel materials.[8-15]. This has aided design
of new dielectric polymers capable of showing low conduction losses and high energy density
under high electric fields. However, there is a need to enhance this design strategy by providing
accurate description of unexplored polymer space using cheap computational techniques which
allow search of large number of potential polymer candidates. Reactive molecular dynamics based
2
on ReaxFF force field and Polarizable Charge Equilibration allows us to achieve this by efficiently
incorporating first principle-based information.
1.1. Computational framework for polymer synthesis to study dielectric
properties using polarizable reactive molecular dynamics

To further accelerate the computational-screening step, it is crucial to develop new descriptors that
capture the effects of complex chemistries and microstructures on the temperature-dependent
dielectric behavior while retaining quantum-mechanical (QM) accuracy. First-principles
molecular dynamics (MD) simulation[16-18] is capable of describing distinct polymer chain
structures that consist of designed repeat units, although the computational cost prohibits to
investigate complex local structures such as domains and interfaces between crystalline,
amorphous, and semi-crystalline phases. It has been shown that crystalline polymer materials
exhibit lower dielectric constant values due to restrictive dipolar relaxations caused by the
surrounding media. On the other hand, significantly higher dielectric constant values have been
observed in amorphous phases, in which the vast difference has been attributed to the increased
spatial distributions in the amorphous samples.[19] Critical role of morphology has further been
investigated for many polyolefins such as polyethylene and polypropylene.[19, 20] However,
system sizes handled by the first-principles approach is limited to a few hundred atoms due to the
high computational cost.
All-atom MD simulations using an empirical force field can be used to scale the simulation
system size up to practically relevant polymer lengths. Such a force field can accurately describe
structural motifs and features such as repeat units, single chain morphology, chain stackings, local
domains and their interface, as well as temperature dependence and thermo-mechanical properties.
However, due to the absence of the polarizability of the constituent materials as a function of their
3
valence charge state[21] this approach suffers from poor predictions of dielectric response
subjected to an electric field. Hence, there is a need for polarizable reactive molecular dynamics
framework to circumvent this issue.
1.2. Dielectric Polymer Property Prediction Using Recurrent Neural Networks
with Optimizations

In recent years, studies of structure-activity or structure-property relationships have been a focus
of life science[22-25], chemistry[26], and materials research[27, 28]. The computer-aided
molecular design represents a milestone in industrial evolution and forces a shift from the time and
resource-intensive combinatorial search to in-silico high throughput screening sometimes in a
matter of seconds[29]. Compared to the vast success of machine learning (ML) for predicting
structure-property relationships in molecular and material research,[28] predicting the dielectric
properties of polymers (i.e., dielectric polymer genome) is still in its infancy[30-32]. Recently,
computationally synthesized polymer databases have been created,[33] which involved highly
accurate computation of dielectric constants based on the new generation of first principles-
informed polarizable reactive-forcefield methods[34-36]. A variety of ML methods have been
used to predict the properties of compounds to avoid the need for complex molecular
characterization.[37] Artificial neural networks of various architectures are a prominent platform
for property prediction as well as materials design. In multiple reports as well as in our studies, it
was shown that multi-layer perceptron (MLP)[38], convolutional (CNN)[39], and graph
(GNN)[40] neural networks perform well in mapping the physicochemical properties of complex
polymeric systems. While deep neural networks (DNN) generally exhibit advanced performance
compared to more classic ML models, a careful choice of the input representation, as well as
optimization of model parameters, have shown that recurrent neural network (RNN) to be an
equally effective alternative for polymer genome studies[41]. Hence, the dielectric property data
4
concerning polymers need to be analyzed using superior sequence models such as RNN which can
analyze the SMILES encoding of these polymers to extract meaningful features relevant to design
of new dielectric materials.
1.3. Accelerated Discovery of Dielectric Polymer Materials using Graph
Convolution Neural Network

Development of high-performance dielectric polymers is urgently needed for the advancement of
a wide range of energy technologies, including energy-storage and pulsed-power technologies [30,
42, 43]. Recent advent of materials genome (i.e., applying informatics to design new materials
significantly faster than the conventional trial-and-error approach) [44] has enabled rational design
of dielectric polymers [30], thereby heralding an exciting era of polymer genome [45]. A major
performance indicator of dielectric polymer is dielectric constants. While dielectric genome has
previously been proposed for inorganic crystals [46], high computational costs of quantum-
mechanically evaluating dielectric constants [47-49] make it still a challenge. The situation is even
worse for polymers, where computing is far more expensive since complex chemical and
morphological features essentially dictate their dielectric constants. As a result, the highly-
anticipated “dielectric polymer genome” is still elusive.
This computational challenge may be addressed by recent developments in first principles-
informed reactive molecular dynamics (RMD) simulations based on polarizable reactive force
fields (ReaxPQ) [35, 36, 50]

5

2. Molecular Dynamics
2.1 Molecular dynamics simulation

Molecular dynamics[51-53] (MD) is an initial value problem. It follows the physical movements
of particles at any given time from their initial positions using Newton’s equation of motion
under a given interatomic force field. A N-body molecular dynamics simulation consist of N
particle with 3N-element positions, and the velocities,
. The trajectories of each atom are predicted by
numerically solving the Newton’s equations of motion, where interatomic force between
particles and the associated potential energy are defined by a parameterized force fields.

𝑎 ⃗
!
(𝑡)=
𝐹
⃗
!
(𝑡)
𝑚
!

(2.1)
where mi is the mass of ith particle.
If Dt is a small increment in time, the velocity of i
th
particle can be derived from the change of
position with respect to time,

𝑣 ⃗
!
(𝑡) =
𝑑𝑟 ⃗
!
𝑑𝑡
= lim
∆#→%
𝑟 ⃗
!
(𝑡+∆𝑡)−𝑟 ⃗
!
(𝑡)
∆𝑡

(2.2)
The acceleration of the i
th
particle can be inferred from the change in velocity with respect to
time.
The acceleration of that particle is the change of velocity with respect of time,

𝑎 ⃗
!
(𝑡)=
𝑑
&
𝑟 ⃗
!
𝑑𝑡
&
=
𝑑𝑣 ⃗
!
𝑑𝑡
= lim
∆#→%
𝑣 ⃗
!
(𝑡+∆𝑡)−𝑣 ⃗
!
(𝑡)
∆𝑡

(2.3)
!
r
N
= (x
1
,y
1
,z
1
,x
2
,y
2
,z
2
,...,x
N
,y
N
,z
N
)
!
v
N
= (v
1x
,v
1y
,v
1z
,v
2x
,v
2y
,v
2z
,...,v
Nx
,v
Ny
,v
Nz
)
6
We use velocity-Verlet numerical method to integrate the Newton’s equation of motion for each
particle. In the velocity-Verlet algorithm, expression for the velocity of each particle in the
acceleration form is given by half the time step,

𝑎 ⃗
!
(𝑡) = lim
∆#→%
𝑣 ⃗
!
>𝑡+
∆𝑡
2
@−𝑣 ⃗
!
>𝑡−
∆𝑡
2
@
∆𝑡

= lim
∆#→%
𝑟 ⃗
!
(𝑡+∆𝑡)−2𝑟 ⃗
!
(𝑡)+𝑟 ⃗
!
(𝑡−∆𝑡)
∆𝑡

(2.4)
Using equation (2.1) to Eq. (2.4), we can relate the positions and velocities of a system of N
particles with given interatomic force fields. Thus, as seen from the spatio-temporal evolution of
a physical system represented by that system could be predicted with given initial values.
2.2 Interatomic force fields

To describe the state of each particle in the system, it is essential to compute the force applied on
each particle. In molecular dynamics, we consider atomic forces that are derived from the total
potential energy, which govern the motion of the atoms. The potential energy can be written as a
function of a 3N-dimensional position vector, where N is the total number of atom present in the
system.
𝑉 = 𝑉(𝑟
'
,…..𝑟
(
) = 𝑉(𝑥
'
,𝑦
'
,𝑧
'
,…..𝑥
!
,𝑦
!
,𝑧
!
,…..,𝑥
(
,𝑦
(
,𝑧
(
) (2.5)
Force on the i
th
atom is represented as
𝐹
)
HH⃗
= −
𝜕𝑉(𝑟
'
,…..𝑟
(
)
𝜕𝑟 ⃗
!
= −J
𝜕𝑉
𝜕𝑥 ⃗
!
,
𝜕𝑉
𝜕𝑦 ⃗
!
,
𝜕𝑉
𝜕𝑧 ⃗
!
K (2.6)
In general, the total interatomic potential energy could contain multiple terms corresponding to a
specific physical effect. A very commonly used interatomic potential used for inert gases is
7
Lennard-Jones potential. Figure 2.1 shows a graph between the potential energy and distance in
dimensionless units.

𝑉(𝑟) = L𝑉M𝑟
!*
N =
!+*
L4𝜖QR
𝜎
𝑟
!*
T
'&
−R
𝜎
𝑟
!*
T
,
U
!+*

(2.7)

Figure 2. 1. Dimensionless Lennard-Jones potential for inert gases.

In this dissertation, we have used a first principle informed many-body potential, ReaxFF[54].
ReaxFF allows the bond formation and breaking, which is suitable for large scale high explosive
material simulation. We have also used density functional theory (DFT)[55] to describe the
interaction between particles. DFT computation cost is extremely high compared to ReaxFF, thus
small system with different approach such as multiscale shock theory has been used to model the
shock.
8
2.3. Reactive forcefield (ReaxFF)

Reactive forcefields can describe chemical reactions and dynamic properties for large-scale
reactive system, with both bonded and non-bonded interactions. For bonded terms, bond order is
applied to take account of the formation and dissociation of the chemical bond. For non-bonded
terms, ReaxFF considers shielded van der Waals (Morse Potential) and Coulomb interactions
with a finite range[56].

Figure 2. 2. Energy terms included in ReaxFF where Ebond is the bond energy, Elp is the lone
pair energy, Eover and Eunder are the over- and under- coordination energies, Eval is the valence
angle energy, Epen is the penalty energy, Ecoa is the three-body conjugation energy, Etors is the
torsion angle energy, Econj is the four-body conjugation energy, Ehbond is the hydrogen bond
interaction energy, EvdWaals, ECoulomb and Elg are the van der Waals, Coulomb energies and
dispersion energy, respectively.

The overall system energy in ReaxFF can be described using eq 2.8 and the breakdown for
various terms is also shown in figure 2.2
2.8

Below follows a description of the partial energies introduced in eq 2.8.

2.3.1. Bond Order and Bond Energy

A fundamental assumption of ReaxFF is that the bond order BO’ ij between a pair of atoms can be
obtained directly from the interatomic distance rij as given in eq (2.9). In calculating the bond orders,
ReaxFF distinguishes between contributions from sigma bonds, pi-bonds and double pi bonds.
ReaxFF energy: E = E
lp
+E
over
+E
under
+E
bond
+ E
val
+E
pen
+E
coa
Non-bonded
+E
tors
+E
conj
+ E
hbond
+ E
vdWaals
+E
Coulomb
Bonded
Coulomb vdWaals bond H conj tors
coa pen val under over lp bond system
E E E E E
E E E E E E E E
+ + + +
+ + + + + + + =
−
9
2.9
Based on the uncorrected bond orders BO’, derived from eq 2.9, an uncorrected over coordination
D’ can be defined for the atoms as the difference between the total bond order around the atom and
the number of its bonding electrons Val as shown in eq. 2.10a.
2.10a
ReaxFF then uses these uncorrected overcoordination definitions to correct the bond orders BO’ij
using the scheme described in eqs (2.11a-f). To soften the correction for atoms bearing lone
electron pairs a second overcoordination definition D’
boc
(eq 2.10b) is used in equations 2.11e and
2.11f. This allows atoms like nitrogen and oxygen, which bear lone electron pairs after filling their
valence, to break up these electron pairs and involve them in bonding without obtaining a full bond
order correction.
2.10b
2.11a

2.11b
€
BO
ij
'
= BO
ij
σ
+ BO
ij
π
+ BO
ij
ππ
=exp p
bo1
⋅
r
ij
r
o
σ
%
&
'
(
)
*
p
bo2
+
,
-
-
.
/
0
0
+exp p
bo3
⋅
r
ij
r
o
π
%
&
'
(
)
*
p
bo4
+
,
-
-
.
/
0
0
+
exp p
bo5
⋅
r
ij
r
o
ππ
%
&
'
(
)
*
p
bo6
+
,
-
-
.
/
0
0
€
Δ
i
'
=−Val
i
+ BO
ij
'
j=1
neighbours(i)
∑
€
Δ
i
'boc
=−Val
i
boc
+ BO
ij
'
j=1
neighbours(i)
∑
ππ π σ
ππ ππ
π π
σ σ
ij ij ij ij
ij j ij i j i j i ij ij
ij j ij i j i j i ij ij
ij j ij i j i ij ij
BO BO BO BO
BO f BO f f f BO BO
BO f BO f f f BO BO
BO f BO f f BO BO
+ + =
Δ ⋅ Δ ⋅ Δ Δ ⋅ Δ Δ ⋅ =
Δ ⋅ Δ ⋅ Δ Δ ⋅ Δ Δ ⋅ =
Δ ⋅ Δ ⋅ Δ Δ ⋅ =
) , ( ) , ( ) , ( ) , (
) , ( ) , ( ) , ( ) , (
) , ( ) , ( ) , (
' '
5
' '
4
' '
1
' '
1
'
' '
5
' '
4
' '
1
' '
1
'
' '
5
' '
4
' '
1
'
!
!
"
#
$
$
%
&
Δ Δ + Δ Δ +
Δ Δ +
+
Δ Δ + Δ Δ +
Δ Δ +
⋅ = Δ Δ
) , ( ) , (
) , (
) , ( ) , (
) , (
2
1
) , (
' '
3
' '
2
' '
2
' '
3
' '
2
' '
2
1
j i j i j
j i j
j i j i i
j i i
j i
f f Val
f Val
f f Val
f Val
f
10
2.11c

2.11d

2.11e

2.11f

A corrected overcoordination Di can be derived from the corrected bond orders using eq (2.12).
2.12
Eq. (2.13) is used to calculate the bond energies from the corrected bond orders BOij.
2.13

2.3.2. Lone pair energy

Eq. (2.15) is used to determine the number of lone pairs around an atom. D
i
e
is determined in Eq.
(2.14) and describes the difference between the total number of outer shell electrons (6 for oxygen,
4 for silicon, 1 for hydrogen) and the sum of bond orders around an atomic center.
2.14
€
f
2
(Δ
i
'
,Δ
j
'
) = exp(−p
boc1
⋅Δ
i
'
) + exp(−p
boc1
⋅Δ
j
'
)
€
f
3
(Δ
i
'
,Δ
j
'
) =−
1
p
boc2
⋅ln
1
2
⋅ exp −p
boc2
⋅Δ
i
'
( )
+ exp −p
boc2
⋅Δ
j
'
( ) [ ]
%
&
'
(
)
*
€
f
4
(Δ
i
'
,BO
ij
'
) =
1
1+exp(−p
boc3
⋅(p
boc4
⋅BO
ij
'
⋅BO
ij
'
−Δ
i
'boc
) + p
boc5
)
€
f
5
(Δ
j
'
,BO
ij
'
) =
1
1+exp(−p
boc3
⋅(p
boc4
⋅BO
ij
'
⋅BO
ij
'
−Δ
j
'boc
) + p
boc5
)
€
Δ
i
=−Val
i
+ BO
ij
j=1
neighbours(i)
∑
€
E
bond
=−D
e
σ
⋅BO
ij
σ
⋅exp p
be1
1− BO
ij
σ
( )
p
be2
( )
%
&
'
(
)
*
−D
e
π
⋅BO
ij
π
−D
e
ππ
⋅BO
ij
ππ
€
Δ
i
e
=−Val
i
e
+ BO
ij
j=1
neighbours(i)
∑
11
2.15
For oxygen with normal coordination (total bond order=2, D
i
e
=4), eq (2.15) leads to 2 lone pairs.
As the total bond order associated with a particular O starts to exceed 2, eq. 2.15 causes a lone pair
to gradually break up, causing a deviation D
i
lp
, defined in eq 2.16, from the optimal number of lone
pairs n
lp,opt
(e.g. 2 for oxygen, 0 for silicon and hydrogen).

2.16
This is accompanied by an energy penalty, as calculated by eq. 2.17.
2.17

2.3.3. Overcoordination

For an overcoordinated atom (Di>0), equations (2.18a-b) impose an energy penalty on the system.
The degree of overcoordination D is decreased if the atom contains a broken-up lone electron pair.
This is done by calculating a corrected overcoordination (eq. 2.18b), taking the deviation from the
optimal number of lone pairs, as calculated in equation (2.16), into account.
2.18a

2.18b
€
n
lp,i
=int
Δ
i
e
2
#
$
%
&
'
(
+exp −p
lp1
⋅ 2 +Δ
i
e
−2⋅int
Δ
i
e
2
+
,
-
.
/
0
#
$
%
&
'
(
2
1
2
3
3
4
5
6
6
i lp opt lp
lp
i
n n
, ,
− = Δ
€
E
lp
=
p
lp2
⋅Δ
i
lp
1+exp −75⋅Δ
i
lp
( )
€
E
over
=
p
ovun1
⋅D
e
σ
⋅BO
ij
j=1
nbond
∑
Δ
i
lpcorr
+Val
i
⋅Δ
i
lpcorr
⋅
1
1+exp p
ovun2
⋅Δ
i
lpcorr
( )
&
'
(
(
)
*
+
+
€
Δ
i
lpcorr
=Δ
i
−
Δ
i
lp
1+ p
ovun3
⋅exp p
ovun4
⋅ Δ
j
−Δ
j
lp
( )
⋅(BO
ij
π
j=1
neighbours(i)
∑
+ BO
ij
ππ
)
'
(
)
*
)
+
,
)
-
)
.
/
0
0
1
2
3
3
12

2.3.4. Undercoordination

For an undercoordinated atom (Di<0), we want to consider the energy contribution for the
resonance of the p-electron between attached under-coordinated atomic centers. This is done by
eq. 2.19 where Eunder is only important if the bonds between under-coordinated atom i and its
under-coordinated neighbors j partly have p-bond character.
2.19
2.3.5. Valence Angle Terms

2.3.5.1 Angle energy
Just as for bond terms, it is important that the energy contribution from valence angle terms goes
to zero as the bond orders in the valence angle goes to zero. Eqs (2.20a-g) are used to calculate the
valence angle energy contribution. The equilibrium angle Qo for Qijk depends on the sum of p-
bond orders (SBO) around the central atom j as described in eq (2.20d). Thus, the equilibrium
angle changes from around 109.47 for sp
3
hybridization (p-bond=0) to 120 for sp
2
(p-bond=1) to
180 for sp (p-bond=2) based on the geometry of the central atom j and its neighbors. In addition
to including the effects of p-bonds on the central atom j, eq (2.20d) also takes into account the
effects of over- and under-coordination in central atom j, as determined by eq (2.20e), on the
equilibrium valency angle, including the influence of a lone electron pair. Val
angle
is the valency
of the atom used in the valency and torsion angle evaluation. Val
angle
is the same as Val
boc
used in
equation (2.10b) for non-metals. The functional form of eq. (2.20f) is designed to avoid
singularities when SBO=0 and SBO=2. The angles in eqs (2.20a)-(2.20g) are in radians.

( )
( )
!
!
"
#
$
$
%
&
'
(
)
*
+
,
+ ⋅ Δ − Δ ⋅ ⋅ +
⋅
Δ ⋅ − +
Δ ⋅ −
⋅ − =
∑
=
) ( exp 1
1
) exp( 1
exp 1
) (
1
8 7
2
6
5
ππ π
ij
i neighbours
j
ij
lp
j j ovun ovun
lpcor
i ovun
lpcor
i ovun
ovun under
BO BO p p
p
p
p E
13
2.20a

2.20b

2.20c

2.20d

2.20e

2.20f

2.20g

2.3.5.2 Penalty energy
To reproduce the stability of systems with two double bonds sharing an atom in a valency angle,
like allene, an additional energy penalty, as described in eqs 2.21a and 2.21b, is imposed for such
( ) ( ) [ ] { }
2
2 1 1 8 7 7
exp ) ( ) ( ) (
ijk o val val val j jk ij val
BO p p p f BO f BO f E Θ − Θ − − ⋅ Δ ⋅ ⋅ =
€
f
7
(BO
ij
) =1−exp −p
val3
⋅BO
ij
p
val4
( )
( )
( )
( ) ( )
angle
j val
angle
j val
angle
j val
val val j
p p
p
p p f
Δ ⋅ − + Δ ⋅ +
Δ ⋅ +
⋅ − − = Δ
7 6
6
5 5 8
exp exp 1
exp 2
1 ) (
€
SBO = BO
jn
π
+ BO
jn
ππ
( )
n=1
neighbors( j)
∑
+ 1− exp −BO
jn
8
( )
n=1
neighbours( j)
∏
&
'
(
)
*
+ ⋅ −Δ
j
angle
− p
val8
⋅n
lp, j ( )
€
Δ
j
angle
=−Val
j
angle
+ BO
jn
n=1
neighbours( j)
∑
2 if 2 2
2 1 if ) 2 ( 2 2
1 0 if 2
0 if 0 2
9
9
> =
< < − − =
< < =
≤ =
SBO SBO
SBO SBO SBO
SBO SBO SBO
SBO SBO
val
val
p
p
€
Θ
0
BO ( ) =π −Θ
0,0
⋅ 1−exp −p
val10
⋅ 2−SBO2 ( ) [ ] { }
14
systems. Eq. (2.16b) deals with the effects of over/undercoordination in central atom j on the
penalty energy.
2.21a

2.21b

2.3.5.3. Three-body conjugation term
The hydrocarbon ReaxFF potential contained only a four-body conjugation term (see section
2.3.6.2), which was sufficient to describe most conjugated hydrocarbon systems. However, this
term failed to describe the stability obtained from conjugation by the –NO2-group. To describe
the stability of such groups a three-body conjugation term is included (eq. 2.16).
2.16

2.3.6. Torsion angle terms

2.3.6.1. Torsion rotation barriers
Just as with angle terms we need to ensure that dependence of the energy of torsion angle wijkl
accounts properly for BO ® 0 and for BO greater than 1. This is done by eqs (2.17a)-(2.17c).
2.17a
€
E
pen
= p
pen1
⋅ f
9
(Δ
j
)⋅exp −p
pen2
⋅ BO
ij
−2
( )
2
[ ]
⋅exp −p
pen2
⋅ BO
jk
−2
( )
2
[ ]
( )
( ) ( )
j pen j pen
j pen
j
p p
p
f
Δ ⋅ + Δ ⋅ − +
Δ ⋅ − +
= Δ
4 3
3
9
exp exp 1
exp 2
) (
( )
( ) [ ] ( ) [ ]
2
4
2
4
2
) (
1
3
2
) (
1
3
2
1
5 . 1 exp 5 . 1 exp exp
exp
exp 1
1
− ⋅ − ⋅ − ⋅ − ⋅
#
#
$
%
&
&
'
(
)
)
*
+
,
,
-
.
+ − ⋅ − ⋅
#
#
$
%
&
&
'
(
)
)
*
+
,
,
-
.
+ − ⋅ − ⋅
Δ ⋅ +
⋅ =
∑
∑
=
=
jk coa ij coa
i neighbours
n
kn jk coa
i neighbours
n
in ij coa val
j coa
coa coa
BO p BO p BO BO p
BO BO p
p
p E
€
E
tors
= f
10
(BO
ij
,BO
jk
,BO
kl
)⋅sinΘ
ijk
⋅sinΘ
jkl
⋅
1
2
V
2
⋅exp p
tor1
⋅ BO
jk
π
−1+ f
11
(Δ
j
,Δ
k
)
( )
2
{ }
⋅ 1−cos2ω
ijkl
( )
+
1
2
V
3
⋅ 1+cos3ω
ijkl
( )
(
)
*
+
,
-
15

2.17b
2.17c

2.3.6.2. Four body conjugation term
Equations (17a-b) describe the contribution of conjugation effects to the molecular energy. A
maximum contribution of conjugation energy is obtained when successive bonds have bond order
values of 1.5 as in benzene and other aromatics.

(2.18a)

2.18b

2.3.7. Hydrogen bond interactions
Eq. 2.19 described the bond-order dependent hydrogen bond term for a X-H—Z system as
incorporated in ReaxFF.

2.19
€
f
10
(BO
ij
,BO
jk
,BO
kl
) = 1−exp −p
tor2
⋅ BO
ij ( ) [ ]
⋅ 1−exp −p
tor2
⋅ BO
jk ( ) [ ]
⋅ 1−exp −p
tor2
⋅ BO
kl
( ) [ ]
€
f
11
(Δ
j
,Δ
k
) =
2 +exp −p
tor3
⋅ Δ
j
angle
+Δ
k
angle
( ) [ ]
1+exp −p
tor3
⋅ Δ
j
angle
+Δ
k
angle
( ) [ ]
+exp p
tor4
⋅ Δ
j
angle
+Δ
k
angle
( ) [ ]
€
E
conj
= f
12
(BO
ij
,BO
jk
,BO
kl
)⋅ p
cot1
⋅ 1+ cos
2
ω
ijkl
−1
( )
⋅sinΘ
ijk
⋅sinΘ
jkl [ ]
€
f
12
(BO
ij
,BO
jk
,BO
kl
) =exp −p
cot 2
⋅ BO
ij
−1
1
2
$
%
&
'
(
)
2
*
+
,
-
.
/ ⋅exp −p
cot 2
⋅ BO
jk
−1
1
2
$
%
&
'
(
)
2
*
+
,
-
.
/ ⋅exp −p
cot2
⋅ BO
kl
−1
1
2
$
%
&
'
(
)
2
*
+
,
-
.
/
€
E
Hbond
= p
hb1
⋅ 1−exp p
hb2
⋅BO
XH
( ) [ ]
⋅exp p
hb3
r
hb
o
r
HZ
+
r
HZ
r
hb
o
−2
$
%
&
'
(
)
*
+
,
-
.
/ ⋅sin
8
Θ
XHZ
2
$
%
&
'
(
)
16

2.3.8. Nonbonded interactions

In addition to valence interactions which depend on overlap, there are repulsive interactions at
short interatomic distances due to Pauli principle orthogonalization and attraction energies at long
distances due to dispersion. These interactions, comprised of van der Waals and Coulomb forces,
are included for all atom pairs, thus avoiding awkward alterations in the energy description during
bond dissociation.
2.3.8.1. Taper correction
To avoid energy discontinuities when charged species move in and out of the non-bonded cutoff
radius ReaxFF employs a Taper correction, as developed by de Vos Burchart (1995). Each
nonbonded energy and derivative is multiplied by a Taper-term, which is taken from a distance-
dependent 7
th
order polynomial (eq. 2.20).
2.20

The terms in this polynomal are chosen to ensure that all 1
st
, 2
nd
and 3
rd
derivatives of the non-
bonded interactions to the distance are continuous and go to zero at the cutoff boundary. To that
end, the terms Tap0 to Tap7 in eq. 2.20 are calculated by the scheme in eq. 2.21, where Rcut is the
non-bonded cutoff radius.
2.21
€
Tap =Tap
7
⋅ r
ij
7
+Tap
6
⋅ r
ij
6
+Tap
5
⋅ r
ij
5
+Tap
4
⋅ r
ij
4
+Tap
3
⋅ r
ij
3
+Tap
2
⋅ r
ij
2
+Tap
1
⋅ r
ij
+Tap
0
€
Tap
7
=20/R
cut
7
Tap
6
=−70/R
cut
6
Tap
5
=84/R
cut
5
Tap
4
=−35/R
cut
4
Tap
3
=0
Tap
2
=0
Tap
1
=0
Tap
0
=1
17
2.3.8.2 van der Waals interactions
To account for the van der Waals interactions we use a distance-corrected Morse-potential (eqs.
2.22a-b). By including a shielded interaction (eq. 2.22b) excessively high repulsions between
bonded atoms (1-2 interactions) and atoms sharing a valence angle (1-3 interactions) are avoided.
2.22a

2.22b

2.3.9. Coulomb Interactions
As with the van der Waals-interactions, Coulomb interactions are considered between all atom
pairs. To adjust for orbital overlap between atoms at close distances a shielded Coulomb-potential
is used.
2.23
Atomic charges are calculated using the Electron Equilibration Method (EEM)-approach. The
EEM charge derivation method is like the QEq-scheme; the only differences, apart from parameter
definitions, are that EEM does not use an iterative scheme for hydrogen charges (as in QEq) and
that QEq uses a more rigorous Slater orbital approach to account for charge overlap.

€
E
vdWaals
=Tap⋅D
ij
⋅ exp α
ij
⋅ 1−
f
13
(r
ij
)
r
vdW
%
&
'
(
)
*
+
,
-
.
/
0 −2⋅exp
1
2
⋅α
ij
⋅ 1−
f
13
(r
ij
)
r
vdW
%
&
'
(
)
*
+
,
-
.
/
0
1
2
3
4
3
5
6
3
7
3
€
f
13
(r
ij
) = r
ij
p
vdW 1
+
1
γ
w
#
$
%
&
'
(
p
vdW 1
)
*
+
+
,
-
.
.
1
p
vdW 1
€
E
coulomb
=Tap⋅C⋅
q
i
⋅q
j
r
ij
3
+ 1/γ
ij
( )
3
[ ]
1/3
18
3. Computation Framework for Polymer Synthesis to Study
Dielectric Properties Using Polarizable Reactive Molecular
Dynamics

3.1. Background

Polymers are an important class of materials owing to their chemical resistance, light weight, high
processability and ability to act as thermal and electrical insulators. Because of their useful
properties, polymer products are ubiquitous in our daily lives and are used in plastics, packaging,
military, aerospace and various other fields.[1] In insulation, isolation and electrostatic energy
storage technologies, polymers dielectric materials have become mainstays, since they are
relatively inexpensive and lighter than ceramic material alternatives. For example, biaxially-
oriented polypropylene (BOPP) is currently employed as state-of-the-art polymer dielectric films
in electrostatic energy storage capacitors. Despite having high dielectric strength and low dielectric
losses, BOPP exhibits relatively low dielectric constant and energy density.[2-4] Therefore,
concerted experimental and theoretical efforts have been made to develop new class of dielectric
polymer materials which possess low conduction losses under high electric fields and electrical
energy density and can be used in a wide range of applications in modern electronics and electric
power systems.[5-7]
Due to the rapid advances in machine-learning technologies, the closed-loop cycle
consisting of materials synthesis, characterization and computational modeling has been attracting
great attention to accelerate design and discovery of novel materials.[8-15] Advanced
computational screening procedures also have been applied and have significantly helped identify
important high energy density polymeric materials.[60, 61] Since polymers exist as repeating
monomer units, it is possible to express them as a series of strings and store in relational
19
databases.[62] The database can be utilized to create models using machine learning tools to
predict relevant polymer properties.[32] This idea has enabled fast computational screening,
coupled with hierarchical modelling, to accelerate the identification of promising repeat units to
synthesize.[60, 63] For instance, polyurea and polythiourea, with -NH-CO-NH-C6H4- and -NH-
CS-NH-C6H4- repeat units, have been identified as high dielectric constant organic polymers
utilizing the rational co-design strategy.[60, 64] The identified repeat units have been used as
starting points for further optimization of the chemistry so as to be suitable for small- and large-
scale thin film processing.[64] One example of such an optimized polymer is PDTC-HK511, a
polythiourea-based material by polymerizing p-phenylene isothiocyanate and an etheramine
jeffamine HK511.[65] Similar protocols have been used to create TDI-EDR148, a polyurea-based
material by addition polymerization of 1,4-toluene diisocyanate (TDI) and an etheramine
jeffamine EDR-148. These materials have been further characterized and studied for establishing
a relationship between their microscopic molecular structure and macroscopic dielectric
properties.[60, 65, 66] They have been found to have high breakdown strength, energy density,
dielectric constant (> 5) and maintain low losses at high electric fields. Moreover, this advanced
screening procedure has identified promising high energy density polymers from the existing
materials such as divinyl siloxane benzo-cyclobutene (DVS-BCB) and flurorene polyester
(FPE).[67, 68] DVS-BCB and FPE exhibit high glass transition temperatures (Tg > 300
°
C) and
high dielectric constants (> 2.5) at room temperature.
First-principles molecular dynamics (MD) simulation[16-18] is capable of describing
distinct polymer chain structures that consist of designed repeat units, although the computational
cost prohibits to investigate complex local structures such as domains and interfaces between
crystalline, amorphous, and semi-crystalline phases. It has been shown that crystalline polymer
20
materials exhibit lower dielectric constant values due to restrictive dipolar relaxations caused by
the surrounding media. On the other hand, significantly higher dielectric constant values have been
observed in amorphous phases, in which the vast difference has been attributed to the increased
spatial distributions in the amorphous samples.[19] Critical role of morphology has further been
investigated for many polyolefins such as polyethylene and polypropylene.[19, 20] However,
system sizes handled by the first-principles approach is limited to a few hundred atoms due to the
high computational cost.
All-atom MD simulations using an empirical force field can be used to scale the simulation
system size up to practically relevant polymer lengths. Such a force field can accurately describe
structural motifs and features such as repeat units, single chain morphology, chain stackings, local
domains and their interface, as well as temperature dependence and thermo-mechanical properties.
However, due to the absence of the polarizability of the constituent materials as a function of their
valence charge state[21] this approach suffers from poor predictions of dielectric response
subjected to an electric field.
To circumvent this deficiency, we have developed a scalable computational framework
based on a charge-state aware ReaxPQ+ model. The original ReaxPQ model combines all-atom
reactive molecular dynamics (RMD) employing reactive force field (ReaxFF)[54] and a recently-
proposed polarizable charge model (PQEq)[50]. To study dielectric response under electric field,
ReaxPQ+ model was introduced to incorporate the effects from both internal and external electric
fields. ReaxPQ+ has been successfully applied to a number of inorganic and organic materials[69],
however, the absence of charge-state in the ReaxPQ+ model resulted in time-consuming and
laborious parameter tunings, presenting a serious challenge to enable a fully automated high-
throughput computational screening for polymer dielectrics. To this goal, we here introduce a
21
charge-state aware ReaxPQ+ model to accurately describe the dielectric responses of amorphous
polymers based on the quantum-mechanically informed atomic polarizability. Once model
parameters have been developed, our scalable framework predicts dielectric properties of polymers
within a small fraction of QM calculation time[69] and handles industrially relevant polymer chain
lengths in a highly-automated fashion. We have applied our model to computationally synthesize
amorphous polymer systems, and demonstrated the role of morphological complexity in evaluating
the dielectric constant of polymers.[19, 20] In addition, we successfully apply the framework
validated by the dielectric responses of PDTC-HK511 and TDI-EDR148 to newly identified FPE
and DVS-BCB within an experimental accuracy. The subsequent sections present the key
components of the framework, obtained dielectric properties and discussions.

3.2. Methodology

Our framework can be divided into following components: a) synthesis of polymer system, and b)
ReaxPQ+ model development and dielectric-constant estimate (Fig. 3.1). Details of these
components are described in subsequent subsections. The synthesis of the polymer system involves
use of various open-source libraries which involve first creation of the SMILES string, followed
by the 3D coordinate generation using babel library which involves use of Universal Forcefield
(UFF) to generate and optimize the structure. The next sections would then detail further into
model development and prediction of dielectric constant using the generated amorphized structure
and the developed model parameters to estimate the dielectric constant values using computed
dipole moments by utilizing the variance in the computed dipole moments at different temperatures
to estimate the final dielectric constant at high frequency as well as dipole moments at various
finite temperatures.
22

Figure 3. 1. Computational framework of polymer dielectric-behavior estimate. (a) Synthesis of
amorphous polymer system. SMILES string and Open Babel are used to create a single polymer
chain. Multiple polymer chains are placed at sufficiently large distance from each chain in a
simulation system. The simulation system is subjected to a number of consolidation and
relaxation steps until the system reaches to a desired density. (b) ReaxPQ+ model development
that involves various QM-based calculations and validations such as atomic polarization and
charges for constituent atomic species, and dielectric constant estimate.

3.2.1. Charge State Aware ReaxPQ+ Model

To take into account the dielectric response to both internal and external electric fields, the
ReaxPQ+ model employs the polarizable charge-equilibration scheme called PQEq,[50] which in
turn introduce core and shell charges to capture the complex interplay between the electronic and
ionic dynamics. The potential energy 𝐸MW𝐫
ij
Y,W𝐫
ijk
Y,W𝐫
ijkl
Y,{𝑞
!
},WBO
ij
YN of the system in ReaxPQ+
is represented as a function of relative positions of atomic pairs rij, triplets rijk and quadruplets rijkl,
23
atomic charges qi and bond orders BOij between different atomic pairs. In ReaxPQ+, the total
charge on i-th atom is computed based on the contribution from a massless shell (ρis) connected to
the core (ρi) by an isotropic harmonic spring constant (ks) as shown in Fig. 3.2. The massless shell
(ρis) has fixed total charge -Zi while the core (ρi) charge consists of two parts: 1) a dynamically
variable atomic charge (qi), and 2) a fixed charge (Zi) compensating the massless shell counterpart.

Figure 3. 2. Schematic of the response of core (black) and shell (red) charges to an external
electric field in ReaxPQ+ model.
The atomic charges qi are determined by minimizing the total Coulomb energy ECoulomb at every
step until the electrochemical potentials, 𝜕𝐸
-./0.12
/𝜕𝑞
!
, become equal among all atoms. The
Coulombic energy is described as
𝐸
-./0.12
({𝑟 ⃑
!3
,𝑟 ⃑
!4
,𝑞
!
})
=L_𝐸
!
5
+𝜒
!
5
𝑞
!
+
1
2
𝐽
!!
5
𝑞
!
&
+
1
2
𝐾
4
𝑟
!3,!4
&
d
(
!
+LW𝐶M𝑟 ⃑
!3,*3
N𝑞
!3
𝑞
*3
−𝐶M𝑟 ⃑
!3,*4
N𝑞
!3
𝑍
*
−𝐶M𝑟 ⃑
!4,*3
N𝑞
*3
𝑍
!
!7*
+𝐶M𝑟 ⃑
!4,*4
N𝑍
!
𝑍
*
Y, (3.1)
24
where 𝑟 ⃑
!3
, 𝑟 ⃑
!4
, 𝜒
!
5
and 𝐽
!!
5
are the positions of core and shell charges, the electronegativity and the
hardness of i-th atom, respectively. ria,jb (i, j = 1,…, N; a, b = core(c) or shell(s)) represents the
charge-charge distances. The electrostatic energy between two Gaussian charges is given in terms
of the error function 𝐶
ia,jb
M𝑟
!8,*2
N, and the Coulombic interaction is screened using a taper function
T(r).[50] The shell positions ris subjected to an external electric field is determined
𝐹
!9#:;
+𝐹
:<#:;980
= 𝐹
!9#;8
(3.2)
𝐅
inter
= −
𝜕
𝜕𝐫
!=
i L 𝑇M𝑟
!8,*2
N𝐶
ia,jb
M𝑟
!8,*2
N𝑞
ia
𝑞
jb
ia7jb
k,
𝐅
intra
= −
𝜕
𝜕𝐫
!=
J
1
2
𝐾
=
𝑟
!>,!=
&
K 𝑎𝑛𝑑 𝐹
:<#:;980
= L 𝑞
!8
𝜀
!8
.
ReaxPQ+ model describes the system polarization via the displacement of Gaussian-shaped shell-
charge relative to the core-charge position, which is controlled by two key parameters: the core-
charge radius Rc and the spring constant Ks between the core-shell charges. A critical component
of the ReaxPQ+ model is the choice of atomic polarizability dictated by Ks, which substantially
changes depending on their valence charge state. Therefore, the model must be aware of the atomic
state to realize accurate prediction of dielectric response as a system. For determining the
polarizability of constituent atomic species, we use Hartree-Fock calculations with def2-svpd basis
sets[70]. Table 3.1 summarizes the obtained atomic polarizability values with various charge
states. The polarizability of neutral atoms is presented as validation in the right-most column,
which shows good agreement with values in literature[71]. These polarizability values can be used
for initial estimation of 𝐾
4
:[50]
𝐾
4
= 𝐶
?
!
@
, (3.3)
25
where C is a constant conversion factor, 𝑍
!
is the charge on shell attached to core, and 𝛼 is the
atomic polarizability (Å
A
), as reported in Table 3.1.
In addition to the polarizability, an accurate description of atomic charges may improve
the estimation of dielectric constants. In this study, we use QM charges throughout our simulations
for prediction of dielectric properties. Table 3.2 presents optimized ReaxPQ+ parameters for the
constituent atoms, which are validated using the high-frequency dielectric constant 𝜖
%
obtained
by Berry-phase method[47, 72].
charge state +2 +1 0 -1 -2 literature[71]
C 0.70 0.97 1.70 3.07 NA 1.71
N 0.41 0.61 0.99 1.59 NA 1.10
O 0.27 0.40 0.69 1.81 4.12 0.80
S 1.09 1.60 2.77 5.77 11.40 2.90
H NA NA NA 2.45 2.651 0.670
Table 3. 1. Computed atomic polarizability (Å
𝟑
) for various charge states of constituent atom
species.
ReaxPQ+ has several parameters for each atom, out of which two key parameters are the core-
shell charge radius Rc and spring constant Ks. The initial value of these parameters is chosen based
on atomic polarizability corresponding to its valence state, then optimized to reproduce the high
frequency dielectric constant 𝜖
%
using the quantum-mechanical (QM) Berry-phase method. Table
3.2 shows the obtained parameters after parameter optimization.

26
Element charge radius (Å)
spring constant
(kcal/mol/Å
2
)
C 0.60 276.84
H 0.37 500.20
N 0.35 80.04
O 0.30 90.87
S 0.30 66.54
Table 3. 2. Optimized parameters (charge radius r
c
and spring constant k
s
) taking into account
the valence charge states for the constituent atoms in PDTC-HK511 and TDI-EDR148.

3.2.2. Amorphous Model Creation
The three major steps in the amorphous polymer model creation using TDI-EDR148 as an example
system is illustrated in Fig. 3.3. The initial step involves generating a Simplified Molecular Input
Line Entry System (SMILES) string for a polymer structure under consideration. SMILES is a
specification in the form of a line notation for describing the structure of chemical species using
short ASCII strings.[73] In the second step, Open Babel[74] is used to generate a low-energy
structure based on the input SMILES string. Open Babel is a computer software, a chemical expert
system mainly used to convert chemical file formats. The obtained polymer-chain structure is
relaxed using conjugate-gradient method, followed by RMD simulation at low temperature to relax
the bonds, angles, dihedrals, and monomer linkages connecting different repeat units together. The
final step involves stacking the polymer chains, a series of RMD simulations to consolidate and
relax the system to achieve a desired density at a given temperature. Fig. 3.4 shows amorphous
structures produced for PDTC-HK511, FPE and DVS-BCB.
27

Figure 3. 3.Amorphous-polymer synthesis of TDI-EDR148 involves accurate SMILES
representation of the polymer followed by its generation using Open Babel and final
consolidation to achieve desired density.

Snapshots of amorphous synthesis for PDTC-HK511, FPE and DVS-BCB are shown in Fig. 3.4.
We follow the same procedure described in the main text; see Fig. 3.3. Here, multiple SMILES[73]
strings are generated corresponding to various polymer systems followed by their relaxation and
consolidation RMD simulations to attain the final system at desired density.
28

Figure 3. 4. The consolidated structures for (a) PDTC-HK511, (b) FPE and (c) DVS-BCB are
generated from their SMILES string representation using Open Babel and ReaxPQ+.

Scalable polymer synthesis is realized using a parallel MD engine called RXMD[75]. RXMD
employs first principles-informed ReaxFF[54] and a divide-and-conquer strategy to achieve the
29
scalability from desktop computer to supercomputing platforms (Fig. 3.5). Message Passing
Interface (MPI) library is used in data exchanges between processes and inter-process
communications[76]. RXMD supports ReaxPQ+[50] model, thus allowing accurate study of the
polarization and dielectric response of materials subjected an electric field. Our all-atom RMD
based scheme enables a flexible design of polymer chains consisting of any repeat unit. To create
industrially relevant polymer morphologies,[77-80] it is crucial that the all-atom algorithm is
highly scalable with optimized time-to-solution. We performed strong and weak scaling tests on
the Theta supercomputer at Argonne Leadership Computing Facility (ALCF). We demonstrated
high scalability of ReaxPQ+ by achieving a weak-scaling parallel efficiency of 0.989 with iso-
granularity of 25,576 atoms per core (N/P), for a 3,221,225,472-atom system. This efficiency is
obtained by calculating the ratio between the speed of P cores to that of 64 cores reference system.
Strong scaling test yields a parallel efficiency of 0.766, obtained by a factor of 12.26 speed up
on 32,678 cores compared to 2048 cores as shown in Fig. S2(b). Comparison of weak and strong
scaling efficiencies in Fig. S2 reveals the difficulty to achieve high strong scaling parallel
efficiency relative to weak scaling parallel efficiency. This decrease arises due to surge in
communication to computation ratio as the workload per rank shrinks proportionally.
30

Figure 3. 5. (a) Weak scaling parallel efficiency of ReaxPQ+ with iso-granular workload
scaling of 25,576 atoms on P cores (P = 64, ..., 131,072) of Theta. (b) Strong scaling parallel
efficiency (correct this) of ReaxPQ+ with a fixed problem size 50,331,648-atom on P cores (P =
2,048, …, 32,678) of Theta. The idea and measured efficiency are shown in black and blue color
lines respectively.

3.2.3. Estimation of Dielectric Constant
Maybe you can change the space between subtitles and following contents, because this space is
larger than the space between the subtitles and above contents. But not sure whether it is the
required format. Just my opinion
The temperature-dependent dielectric constant 𝜖
5
(𝑇) is calculated using Eq. (4), which
consists of the high-frequency dielectric constant (𝜖
%
) and the temperature-dependent term,
𝜖
5
−𝜖
%
=
+CD
"
7
E
#$%
AFG
&
H
, (3.4)
where 𝛺 is the volume of cell, kB is the Boltzmann constant, T is the temperature, and 〈Δ𝑀
&
〉 ≡
〈𝑀
&
〉−〈𝑀〉
&
is the variance of the total dipole moment M(t) over time t. Each component of the
total dipole moment is computed as follows:
31
𝐌(𝑡)= >𝑀
<
(𝑡),𝑀
I
(𝑡),𝑀
J
(𝑡)@ = ∑ ∑ w𝐫
@*
(𝑡)−𝐫
5*
(𝑡)x𝑞
@*
(
$'()
@
(
)(*
*
(3.5)
𝑀
<
(𝑡) = ∑ ∑ w𝑥
@*
(𝑡)−𝑥
5*(#)
x𝑞
@*
(
$'()
@
(
)(*
*
(3.5a)
𝑀
I
(𝑡)= ∑ ∑ w𝑦
@*
(𝑡)−𝑦
5*(#)
x𝑞
@*
(
$'()
@
(
)(*
*
(3.5b)
𝑀
J
(𝑡) = ∑ ∑ w𝑧
@*
(𝑡)−𝑧
5*(#)
x𝑞
@*
(
$'()
@
(
)(*
*
(3.5c)
Here, 𝑁
MNOP
is the number of atoms in molecule j, 𝑞
@*
is the charge of atom 𝛼, 𝐫
@*
(𝑡) is
the position vector of atom 𝛼 with components 𝑥
@*
(𝑡), 𝑦
@*
(𝑡), 𝑧
@*
(𝑡) and 𝑟
5*
(𝑡) is the position of
a reference atom in molecule j with coordinates 𝑥
5*
(𝑡), 𝑦
5*
(𝑡), 𝑧
5*
(𝑡). The high-frequency
dielectric constant 𝜖
%
is obtained from the polarization due to the shell charge displacement at
fixed atom positions subjected an electric field. We obtain the best estimates for 𝜖
%
by further
tuning the initial parameter estimates of the ReaxPQ+ model to fit the first-principles Berry-phase
method.[47, 49, 81] Once the model has been developed, the framework described in this section
allows researchers to investigate a large sets of repeat unit structures and quickly evaluate the
morphology effect such as crystal versus amorphous and system temperature in their dielectric
properties. These insights are crucial to develop new descriptors to further accelerate the
theoretical screening process for novel high energy density polymer design.
3.2.4. Role of morphology on the Dielectric Constant

To demonstrate the role of morphology on the dielectric constant of polymers, we have created
crystalline and amorphous polyethylene (PE) containing 3,040 atoms. A small uniform electric
field was applied in x, y and z directions to estimate the dielectric-constant value in each direction.
Table 3.3 shows the compute dielectric constants averaged over x, y and z directions, which reveal
high dielectric constants for amorphous state as compared to the crystal. Hence, amorphous state
32
plays an important role in enhancing dielectric response of material under external applied electric
field.

Crystal Amorphous
density
(g/cc)
1.0 0.98
dielectric
constant
2.07
± 0.085
2.12 ± 0.06
Table 3. 3. Dielectric constants for crystal and amorphous PE.
3.2.5. System Preparation for Dielectric Constant Calculation

The dielectric constants of TDI-EDR148 and PDTC-HK511 were calculated following the
procedures described in sections 3.2.5.1 and 3.2.5.2.
3.2.5.1. TDI-EDR148 System

The system dimensions of this system are 39.89×39.89×39.89 (Å
3
), and the total number of atoms
is 7,096. The final amorphous structure at 1.3 g/cc was prepared as shown in Fig. 3.6 in the main
text. This system was prepared by generating a polymer chain using Open Babel[74] from their
corresponding SMILES[73] string. The chains are then orthogonally arranged followed by
thermalization and relaxation of the system at 300K. The final system size is chosen to have
individual box dimension greater than twice the monomer length in any direction to allow
unwrapping for the entire chain when periodic boundary condition is applied. The simulation steps
are as follows:
1) The system is brought to desired temperature by slowly heating and quenching from 300K.
33
2) Once the system is at a given temperature, it is given a velocity distribution corresponding
to that temperature and run under micro-canonical ensemble for 20,000 steps.
3) Multiple parallel trajectories (exactly 20) are spawned for the last 10,000 steps.
4) For every trajectory a new velocity distribution is provided and run under micro-canonical
ensemble for 50,000 steps.
5) For computation of 𝜖
5
−𝜖
%
discard the initial 10,000 steps are discarded. The remaining
40,000 steps for every trajectory are aggregated to compute < Δ𝑀
&
>.
6) 80,000 total frames for a given temperature are considered for estimating < Δ𝑀
&
> and
subsequently 𝜖
5
−𝜖
%
.
Fig. 3.6 shows the total dipole moment

Figure 3. 6. Total dipole moments (in Debye) at (a) 100K, (b) 200K, (c) 300K and (d) 400K for
TDI-EDR148 show different variance in total dipole moment across frames. Dipole moment is
computed over many RMD trajectories, the variations in x, y and z directions shown in red,
green, and blue respectively show an increasing trend.
34
3.2.5.2. PDTC-HK511 System
The dimensions of this system are 39.89×39.89×39.89 (Å
3
), and the total number of atoms is
6,480. The final amorphous system was prepared as discussed earlier in the section 3.2.4.1. The
schedule for this calculation was also identical to previous section. Fig. 3.7 shows the obtained
total dipole moments of this system.

Figure 3. 7. Total dipole moments (in Debye) at (a) 100K, (b) 200K, (c) 300K and (d) 400K for
HK511 show different variance in total dipole moment across frames. Dipole moment is
computed over many RMD trajectories, the variations in x, y and z directions shown in red,
green, and blue respectively show an increasing trend.
35

3.3. Results and Discussion

In this section, we demonstrate our framework on 𝜖
5
(𝑇) of known high density polymers (PDTC-
HK511 and TDI-EDR148) and newly identified polymers (DVS-BCB and FPE[67, 68]). All
amorphous systems are created at experimental density[60, 65] following the procedures discussed
in Section 3.2.2.

3.3.1. Validation of high frequency dielectric constant for model polymer systems

The high frequency dielectric constant 𝜖
%
for various model polymer systems (PDTC-HK511,
TDI-EDR148, DVS-BCB and FPE) were estimated by tuning the ReaxPQ+ parameters to fit the
values obtained using Berry Phase method[47, 49, 81]. The estimated values of 𝜖
%
for various
systems are shown in Table 3.4.

𝝐
%

PDTC-HK511 4.28
TDI-EDR148 4.1
DVS-BCB 2.6
FPE 3.1
Table 3. 4. 𝝐
%
values for model polymer systems compared against QM values obtained using
Berry-Phase Method.

36

3.3.2. Dielectric behavior of TDI-EDR148 and PDTC-HK511

Fig. 3.8 (a) and (b) show temperature dependence of the dielectric constants of PDTC-HK511 and
TDI-EDR148, respectively. The high-frequency dielectric constant 𝜖
%
is estimated by ReaxPQ+
scheme using an optimized core-shell radius and atomic polarizability of the constituent atoms.
Overall the dielectric constants for these polymer systems reasonably agree well with
experiments.[65] The high dielectric constant of these systems may be attributed to the interplay
between thermally activated local atomic motions that results in the enhanced fluctuation of
molecular dipole moments.[65] Both PDTC-HK511 and TDI-EDR148 are rigid-base polymer
containing a benzene ring and a flexible ether-amine in the backbone. This flexible segment within
the backbone gives an extra ability of local chain movement thus increasing the dielectric constant,
which may facilitate further polarization when subjected to an external electric field.[65]

Figure 3. 8. Dielectric constants of (a) PDTC-HK511 and (b) TDI-EDR148 polymer systems
(blue dots) as a function of temperature. The increasing trend as well as predicted values show a
good agreement with experimental values (shown in red).
37

3.3.3. Dielectric behavior of DVS-BCB and FPE

DVS-BCB and FPE are promising materials that exhibit high dielectric constant values and high
thermal stability.[67, 68] DVS-BCB is used in preparation of commercially available resins for
multichip modules and four-level GaAs interconnect circuits.[67] DVS-BCB and FPE both have
high glass transition temperature, thus allowing graceful failure and roll-to-roll production
capabilities.[82] Following the previously described procedures, we have computed the
temperature-dependence in dielectric constant 𝜖
5
(𝑇) based on samples taken from a large number
of RMD trajectories. Fig. 3.9 shows the temperature-dependent dielectric constant of amorphous
DVS-BCB and FPE systems. Overall, the obtained dielectric constants agree well with
experiments.[67, 68, 83-87]

Figure 3. 9. Dielectric constants of (a) DVS-BCB and (b) FPE polymer systems as a function of
temperature. The predicted values show a good agreement with experiments[66, 67, 83-87].

38
3.4. Conclusion

In conclusion, we have developed a scalable computational framework to evaluate the
temperature-dependent dielectric responses of amorphous polymer. The seamless creation of
amorphous polymer structure realized by SMILES strings, Open Babel, and RXMD will allow
researchers to explore the vast parameter space of high energy polymer design in a highly
automated fashion. Equipped with the scalable parallel all-atom RMD simulation engine, it is
possible to incorporate complex repeat units with industry-relevant chain lengths. The model
parameters, such as atomic polarizability, core-shell charges, and radii, are extensively optimized
and validated by experiments and first-principal calculations using the high-frequency as well as
temperature-dependent dielectric constants. Our computational framework has been successfully
applied to the recently identified four high energy density polymers, PDTC-HK511, TDI-EDR148,
DVS-BCB and FPE, and revealed that the thermally activated flexible segment motion in polymer
chains facilitates the fluctuation of the molecular dipole moments. The obtained dielectric
constants as well as its increasing trend quantitatively agree with the experiments. The
computational framework presented here for the first time realizes high throughput theoretical
screenings incorporating several key parameters in the high energy density polymer design.
Therefore, combined with efficient experimental synthesis and characterization, our framework is
expected to accelerate the discovery of next-generation high energy density polymer and device
development.

39
4. Polynorbornenes Defect Incorporation and Analysis

4.1. Background

With the increase in the growth of renewable energy, transportation, and space exploration there
is an increased demand of dielectric material with ultra-high energy and power density[6, 68, 88-
90]. Therefore, there is a need for greater understanding of materials pertinent to application under
high electric field environment while maintaining high energy storage and insulation[63].
Moreover, advances in wide bandgap semi-conductors can enable design of ultra-high power
density materials if novel dielectric materials with ability to operate at high temperature and
electric fields could be identified[68, 91, 92]. Inorganic materials such as ceramics have been used
for capacitor applications due to their extremely large dielectric constants, however suffer from
low breakdown strength and non-graceful failure[93]. Hence, ceramic capacitors fail in
applications involving medium to high energy density and therefore exhibit low energy
density[93]. Due to such problems in ceramic capacitors, polymer materials have been identified
as potential candidates for such applications due to their ability to operate under high electric fields
and graceful failure while operating at such conditions[93, 94]. However, design and identification
of polymer dielectrics at higher temperature beyond 100 °C is still challenging. For example,
Biaxially-oriented polypropylene (BOPP) is currently employed as state-of-the-art polymer
dielectric films in electrostatic energy storage capacitors. Despite having high dielectric strength
and low dielectric losses, BOPP exhibits relatively low dielectric constant and energy density[2-
4]. Advances in synthesis of polymer films via nanocomposites or coating modification for high
energy storage applications at high temperature have been encouraging however, their industrial
scale processing is quite expensive [68, 95, 96]. Hence, BOPP is still used with unwieldly cooling
process.
40

The rational design of dielectric polymers using a combination of organic synthesis
expertise along with high throughput density functional theory and machine learning has resulted
in several interesting classes of polymers that exhibit high energy densities, four times that of
biaxially oriented propylene (BOPP), at temperature as high as 170
0
C [63, 65, 93, 97, 98].
Specifically, the polymer to achieve this is polynorbornene (PNB)[99]. This polymer was prepared
by using Grubb’s generation II catalyst. In analogy to BOPP, this PNB is a polyolefin made using
the only the monomer and a very small amount of catalyst and, as such, only must go through
catalyst removal to obtain pristine purity. Dielectric and breakdown testing on this polymer in the
temperature range of 160
0
C to 170
0
C yielded spectacular results[99]. To briefly summarize, this
PNB maintained a dielectric constant of approximately 2.8 with losses at less than 0.1% across the
entire temperature range tested. Having a Tg of approximately 180
0
C, this polymer is amorphous,
non-semicrystalline in nature, due to structural irregularity of the cis/trans double bond isomers
being at a ratio of approximately 50:50, and this polymer retains its glassy state across the entire
temperature range tested.

However, this polymer can have many different derivatives and hence has an extremely
large combinatorial search space. This is a huge challenge for experimental synthesis and
characterization of important candidates among these many derivatives. Due to the rapid advances
in machine-learning technologies, the closed-loop cycle consisting of materials synthesis,
characterization and computational modeling has been attracting great attention to accelerate
design and discovery of novel materials.[8-15] Advanced computational screening procedures also
have been applied and have significantly helped identify important high energy density polymeric
41
materials.[60, 61]. To this goal there is a need to introduce a high throughput framework based on
ReaxPQ+ and Graph Neural Network model to accurately describe the dielectric responses of
amorphous polymers based on the quantum-mechanically informed atomic polarizability. Once
model parameters have been developed, our scalable framework predicts dielectric properties of
polymers within a small fraction of QM calculation time[69]. The subsequent section discusses
the polymer search space, description of framework, generation of dataset, defect analysis of
POFNBs and fluorine containing PNBs.
4.2. Methodology

This section describes the search space for PNB, and framework description for our high
throughput study based on ReaxPQ+.

4.2.1. PNB Design Space

The homopolymer shown in Fig. 4.1 has a vinyl group and benzene rings as endcaps. Starting with
this homopolymer structure, defects can be incorporated into the system using various means as
follows:
1. The X atom can be substituted with an imide, methylene, sulfur and various nitrogen containing
groups.
2.The alkene stereoisomers could be varied from 100% trans to 100% cis in which some cis
stereoisomer could be introduced as defects within a perfectly trans polymer and vice versa.
3. The double bonds could be hydrogenated acting as saturated defects within the backbone of the
PNB.
4. Defects can be incorporated by functionalizing the monomer units with a small amount of
specific R, R
’
groups. This groups can vary as shown in figure 1. Varying these functional groups
42
can result in different behaviors. For example, introduction of polar atoms such as fluorine,
trifluoro methylene and Cyano groups can be modelled and understood.
The incorporation of defects into the system as shown in Fig. 4.1 can generate many
possible configurations. In this study, we consider ~1276 such derivatives to screen potential
candidates.

Figure 4. 1. PNB with heteroatom X and various functional groups which can be incorporated
as defects.

4.2.2. High Throughput Workflow Design

We generate amorphous polymer structures using our computational synthesis framework [100].
The initial step as shown in Fig. 4.2 involves generating a Simplified Molecular Input Line Entry
System (SMILES) string for a polymer structure under consideration. SMILES is a specification
in the form of a line notation for describing the structure of chemical species using short ASCII
strings[73]. In the second step, Open Babel is used to generate a low-energy structure based on the
Parent Polymer Structure
CH
2
, O, S, N-R-R=H
X
R=
or
43
input SMILES string. Open Babel is a computer software, a chemical expert system mainly used
to convert chemical file formats[74]. The obtained polymer-chain structure is relaxed using
conjugate-gradient method, followed by RMD simulation at low temperature to relax the bonds,
angles, dihedrals and monomer linkages connecting different repeat units together. Multiple
polymer chains are then placed at sufficiently large distances from each other in a simulation box.
The simulation system is subjected to a number of consolidation and relaxation steps, until the
system reaches a desired density. Following this procedure, we thus constructed a synthetic dataset
of structurally diverse ~7500 polymers. Dielectric constants of theses polymers are then computed
with the ReaxPQ-v model.

Figure 4. 2. High throughput framework for generating dielectric constant values for various
PNB derivatives.

44

4.3. Dataset Generation

We represent every atom in the polymer chain with a 20-dimensional feature vector for
representing their atomic information along with adjacency list information for the entire polymer
to encode the connectivity information such as bond length and bond order. We use following
atomic properties to represent every polymer chain: Pauling electronegativity, Electron affinity,
group in periodic table, covalent radius, ionisation energy, s electrons, p electrons, atomic
polarizability, ionization affinity ratio, liquid range, liquid ratio, ratio of electron affinity to electro
negativity.
The training dataset was generated for all the possible structures as discussed in Fig. 4.1.
using the workflow as shown in Fig. 4.2. The high frequency dielectric constant reported in Figure
4.3 have been computed from ReaxPQ-v force field using our high throughput framework. The
dielectric constant values as shown in this figure vary between 1.5-3.5. These values have been
calculated for various PNB derivatives with different hetero-atoms and pendant groups
substitutions as discussed above. The characterization of entire dataset based on the type of their
backbone as shown in Fig. 4.3b, shows 3 distinct regions corresponding to imides, imides
copolymerized with non-imides and non-imides backbone. The dielectric constant values here,
agree with theoretical trend due to the presence of symmetric groups in imide which reduces the
net dipole moment of the overall system. Although, theoretically non-imide backbone will have
higher dielectric constant, but they can suffer from low glass transition temperatures and hence are
not always suitable for high energy electrical applications. Due to this complex interplay of various
factors in this polymer, a more complex non-linear analysis is needed to selectively screen out
essentially important polymer candidates suitable for high energy electrical applications. In the
45
next section, defect analysis is carried out on some polymers of the PNB family to demonstrate
the role that various kind of substitutions play in determining the overall dielectric property of the
polymer. The substitutions with respect to the pendant groups in Poly-oxa-fluoro-norbornene
(POFNB) are analyzed first by substituting the pendant group with various benzene ring
substitutions containing CF3 groups at ortho, meta and para positions followed by analysis of
various groups to that can lead to change in the overall free volume around the substitution and
hence finally affect the dielectric constant of the polymer.

Figure 4. 3. The computed dielectric constants are shown here as, (a) the overall dielectric
constant values as computed based on ReaxPQv forcefield and our high throughput workflow,
and (b) there characterization based on the presence and absence of different kind of backbones
in this polymer.

4.4. Defect Analysis

In this section, we would analyze the pendant groups in various poly oxo fluoro norbornenes
(POFNBs) to study the effect of introducing different pendant groups with various kind of CF3
group substitutions at ortho, meta and para positions. These groups are hypothesized to have
different effect on the overall rotational capabilities on the pendant group and hence affecting the
Structure ID
Imides-co-
Nonimides
Imides
Nonimides
Dielectric constant (!
!
)
Structure ID
(a)
(b)
46
overall dipole density thereby realizing different dielectric constants. We will also look into free
volume analysis in Poly fluoro norbornenes (PFNBs) which are special polymers derived out of
the PNB family by substituting the pendant group with various fluorine substituted benzene rings
at different concentrations to introduce net free volume around the pendant groups and thus
causing the difference in the overall dielectric property of the polymer.

4.4.1. Pendant group analysis in various POFNBs

Here we study, four different polymers which are derived out of the PNB base class by substituting
different kind of benzene rings with CF3 groups. The different polymers considered for this study
are shown in Fig. 4.4. These structures vary from each other primarily due to the difference in the
attached pendant group. The effect of this pendant group on the dielectric property is of great
interest and is the focus of our analysis in this section.

Figure 4. 4. Various substitutions in PNB base polymer are shown in this figure. CF
3
groups are
substituted from left to right in meta, para and ortho positions respectively.
These groups are further transformed into long chain polymers by successively repeating them end
to end and finally capping the end by vinyl group on one and benzene group on the other side. The
polymer chain as well as the successive amorphization schematics is shown in Fig. 4.5. The
polymer chains are first stacked orthogonally to each other in a simulation box and after repetitive
POFNB
M-POFNB
P-POFNB O-POFNB
(meta-POFNB) (para-POFNB) (ortho-POFNB)
47
consolidation and thermalization. These amorphized structures contain large number of pendant
groups in each chain and hence provide good samples for estimating dielectric constants which is
done using the fluctuation formula as detailed in previous chapter.

Figure 4. 5. POFNB chains are repetitively consolidated and thermalized by shrinking the
simulation box to attain desired density of 1.3 gm/cc thus finally obtaining an amorphous
structure at the end.

The obtained final structures corresponding to different POFNBs are subjected to electric field at
low temperatures to obtain the high frequency dielectric constant (𝝐
%
) as reported in table 4.1. The
data shown in the second column corresponds to the computed MD results while the last column
shows the experimentally obtained dielectric constant values at room temperatures. The values
shown in second and third column differ from each other since former corresponds to high
frequency part of the dielectric constant while the latter corresponds to the total dielectric constant.
However, the trend of the two agrees with each other and follows this trend,
POFNB < 𝒑−𝑷𝑶𝑭𝑵𝑩 < 𝒎−𝑷𝑶𝑭𝑵𝑩 ≈ 𝒐−𝑷𝑶𝑭𝑵𝑩
C H O N F
Consolidation &
thermalization
POFNB
Obtained amorphous
structure
48

Polymer 𝝐
!
(MD,T=0K) 𝝐 (Expt,T=300K)
POFNB 2.20 2.5
p-POFNB 2.30 3.0
o-POFNB 2.40 -
m-POFNB 2.44 3.5

Table 4. 1.Various POFNBs values are reported in this table, the computed values from our MD
framework follow similar trend as experimentally reported values at room temperature.

The primary reason owing to this trend is attributed to the difficulty in rotation
corresponding to the pendant group due to addition of CF3 groups at ortho, meta and para positions.
This hypothesis can be intuitively understood by visually looking at the structures of various
POFNBs. It can be clearly seen that in ortho-CF3 case the proximity of CF3 groups to the imide
groups in the backbone should be sufficient to sterically hinder its rotation and thus not allowing
sufficient dipole relaxation and thereby would lead to higher dipole moments per polymer chain
and overall higher dipole density which would lead to higher dielectric constant. Similarly,
movement of CF3 groups away from the imide groups in the backbone would lead to more freedom
for degree of rotation thus leading to relaxed dipoles and hence reduced dielectric constants. This
hypothesis can be clearly verified from the results in table 4.1. The results shown in Fig. 4.6.
further analyze this hypothesis by measuring the rotation of these pendant groups across an
extended time scale of 10 ps. It can be clearly seen from this curve that ortho-POFNB is most
49
restricted and para-POFNB is least restricted since the freedom of rotation is highest in para-
POFNB and its rotations are highest in the observed time period.

Figure 4. 6. Pendant group rotation analysis over a period of 10ps shows highest rotation wrt to
para POFNB and lowest for ortho POFNB. The extent of rotations in these polymers directly
correlates to their dielectric property.

4.4.2. Free Volume Analysis in PFNBs

In the last section we looked at defect incorporation in PNB by introducing various pendant groups
containing similar stoichiometry but different chirality. Here, we analyze other kind of possible
defects which are created by introducing groups that can have more free volume around them
owing to bulkier atoms in their composition. Fig. 4.7. shows the synthesis process involved in
creating poly fluoro norbornene with different fluorine concentrations. These polymers can be
created by co-polymerizing two different kind of polymers and varying their composition to obtain
different net fluorine concentration in the overall polymer.

!"#$#%"&(()
oPOFNB
mPOFNB
pPOFNB
Time (ps)
50

Figure 4. 7.The polymerization process for poly fluoro norbornene (PFNB) involves co-
polymerization with benzene containing PNB in different ratios to obtain various polymers with
different fluorine concentration.

These co-polymers thus obtained introduce different free volumes around the fluorine containing
pendant groups and thus are believed to change the overall dipole density and hence resulting in
varying dielectric properties at similar conditions. Fig. 4.8. shows three polymers with varying
fluorine concentration from 0-100% by co-polymerizing them with benzene containing groups as
shown in Fig. 4.7 to obtain these polymers. The results reported in Fig. 4.8. are obtained
experimentally at room temperature by varying the frequency from 1-10000 Hz. The polymers
show a relatively steady profile across the frequency spectrum as can be seen from this figure. The
dissipation factor for these polymers also maintains a relatively steady profile and hence shows
the importance of these polymers for high temperature and extreme applications. It can be clearly
seen here that increase in F concentration clearly decreases the overall dielectric constant.
1
ROMP
Copolymer
+
m
n
Homopolymers
Copolymer Synthesis-
m: n – 25 : 75
50 : 50
75 : 25
51

Figure 4. 8.The experimental measurements for various fluorine concentrations as shown here
are corresponding to 0, 25 and 100% fluorine groups in the pendant groups. The measurements
show a decreasing trend with increasing fluorine concentration

The results shown in Fig. 4.8. show a correlation between increasing F concentration and
decreasing dielectric constant in these polymers. We try to validate this correlation by studying
different systems with varying Fluorine concentration as shown in Fig. 4.9. For this study, we
consider 4 different co-polymers of fluorine containing pendant groups and benzene pendant
groups. These polymers are created at 25, 37.5, 50 and 75% fluorine concentration by varying the
number of fluorine containing pendant group in 8 chain polymer from 2 in 25% case to 6 in 75%
case. These structures will then be amorphized using our framework as discussed in Chapter 3 and
finally subjected to Voronoi volume analysis and dielectric property computation to validate the
correlation as presented in previous figure.
1
Dielectric Constant Measurement of ‘F’ containing Polymers
Polymer F100
Polymer F25
Polymer F0
52

Figure 4. 9. Various fluorine containing structures are shown in this table. The structures
chosen for this study are made by changing the fluorine containing group in 8 chain polymer
from 0 to 6.

The amorphized structures created corresponding to various polymers shown above have
different number Fluorine containing homopolymers in their backbone and thus should have
different free volumes around them. Here we label them as PNB1-4 in order of increasing Fluorine
concentration as discussed in previous figure. We then perform a Voronoi volume analysis to find
out the net free volume introduced by various Fluorine containing pendant groups in these
polymers. The Voronoi volume analysis is performed using the open source free software
Ovito[101]. This analysis computes the net free volume distribution corresponding to all the
fluorine containing pendant group by computing the total free volume around this group and then
1
Structures F %
PNB1 !"
PNB2 #$."
PNB3 "&
PNB4 $"
53
removing its own internal volume to yield the final distribution of free volume around these groups
in the overall system.

Figure 4. 10. The histogram represents Voronoi volume distribution around the pendant group
containing fluorine atoms, the trend clearly shows increasing free volume with increasing
fluorine concentration.

The results shown in Fig. 4.10 show that the free volume in various systems analyzed here
increases with increasing Fluorine concentration as is evident by the free volume calculation
performed using the Voronoi volume-based analysis. We further compute the high frequency
dielectric constant to validate the correlation shown in Fig. 4.8, where we could see a decreasing
trend of dielectric constant with increasing F in the system. Here, we compute only the high
frequency part of the dielectric constant as shown in the fluctuation formula in chapter 3 and not
the overall dielectric constant, as high frequency dielectric constant serves as a very good proxy
to the total dielectric constant and it’s also relatively inexpensive to compute. The results shown
below clearly show a decreasing trend with increasing F in the system thus validating the
hypothesis as presented by experimental findings in Fig. 4.8.
PNB1 PNB2
PNB3 PNB4
F % in various
PNB
Volume(Å
!
)
PNB1 (25%) "#±%#
PNB2 (37.5%) &'±%(
PNB3 (50%) &'±%"
PNB4 (75%) )*±%*
54

Figure 4. 11. High frequency dielectric constant shows a decrease in net magnitude with
increasing fluorine concentration
4.5. Conclusion

In this chapter we presented a polymer base class Polynorbornene (PNB) and the combinatorial
large derivative space around this polymer created by incorporation of various defect in this
polymer backbone. Various substitutions for PNB results in ~3800 polymer candidates that need
to be screened to identify promising polymers suitable for synthesizing in laboratories. These
polymers are created using the high throughput framework presented in chapter 3 and amorphized
using the same framework to finally compute the dielectric constants. The dielectric constant
values thus obtained follow certain theoretically established rules as imide containing groups
altogether have a lower mean dielectric value compared to co-polymerized imides with non-imides
followed by higher mean values for non-imide groups.
Fluorine (%)
!
!
55
We further analyze a few defects out of this large space to study the effect they can have on the
overall dielectric property of the polymer system. Here, we analyzed primarily two kinds of
defects, (i) defects introduced by introduction of various pendant groups in POFNBs and (ii) free
volume introduction in PFNBs. Both defects occur in various forms in the polymer space being
considered here for subsequent data driven studies. These analyses verify certain experimental
hypothesis and provide an overview of the interplay among various complex derivatives in the
PNB base class to alter dielectric property in so many ways. This complexity makes it very difficult
to screen desired polymers from this complex polymer space. Thus, we present neural network
based analysis of this dataset to learn the underlying complexity of this dataset and perform useful
predictions given a polymer or predict polymer similar to other given polymers.

56

5. Attention based Graph Neural Network Analysis of PNB dataset

5.1 Background

Machine learning has gained great popularity particularly when broad exploration of material
space is needed[102-104]. Conventionally, ML models in materials chemistry are descriptor-
based, where the key descriptors representing the system must first be defined prior to fitting a
suitable ML model for prediction. General examples of these descriptors include stoichiometry,
the elemental properties such as group, period, electronegativity and radius, and electronic
properties such as partial charges and s, p, d-band positions[105]. A number of structural
descriptors have also been proposed satisfying translation and rotational invariance, including but
not limited to the Coulomb matrix[105, 106], atom-centered symmetry functions (ACSFs)[105,
107, 108], and smooth overlap atomic potentials (SOAP)[105, 109]. However, finding effective
descriptors can prove challenging for problems with a large number of compositionally diverse
structures[105].

Most Machine learning methods convert the molecules to an appropriate representation
utilizing physiochemical properties from experimental or computational measurements[110, 111]
or by using molecular fingerprints[108, 112, 113]. Popular molecular conversion approaches
usually involve generating molecular fingerprints which encodes the molecule into a vector of
appropriate dimension suitable for carrying out further neural network related modelling[114].
Fixed fingerprint feature extraction related methods are good for accurately reflecting underlying
chemical substructures, but are not very general in their approach[115]. Static fingerprints are also
difficult to develop as they require expertise in the area and often coming up with physical
57
descriptors for task at hand becomes difficult. Simplified Input Line Entry System (SMILES)[116]
also have been popular in various deep learning implementations concerning various applications
as they can automatically learn molecular structure[117, 118]. SMILES representation can be
vague as same structure can be represented by different smiles representation depending upon the
choice of package. For example, the canonical smiles code for caffeine is
CN1C=NC2=C1C(=O)N(C)C(=O)N2C according to RDKit, Cn1cnc2c1c(=O)n(C)c(=O)n2C
according to oBabel, and CN1C=NC2=C1C(=O)N(C(=O)N2C)C according to PubChem[115].
Moreover, other popular deep learning frameworks such as Convolutional Neural Networks and
Sequence Models also do not generalize well to molecular data as these are designed for simple
datatypes such as images, texts etc. Molecules can be easily represented as graphs which are
difficult to model using CNN, RNN type architectures as they have arbitrary size, complex
topological structure, no spatial locality i.e., there is no reference point or fixed node ordering
compared to sequence of text or grid of pixels for an image. Hence there is a need for special
architectures to represent graph data which can account for underlying complex of this datatype as
mentioned above.

Recently, Graph based Neural Networks have shown promising applications by their
capability to learn the inherent features from the graph data directly rather than relying on static
descriptors[99, 105, 119]. Several GNN models have been proposed for chemistry related
problems mostly involving focus on molecular systems[105, 120-123]. Subsequently, GNN have
been used in material predictions involving periodic crystals as well as surfaces[39, 124-127].
These systems involve molecules and hence atoms and bonds can be represented as graphs where
atoms form the nodes while bonds represent edges naturally[105]. An appropriate graph network
58
can then take these node and edge information to learn the underlying node as well as edge
embeddings by taking into account the overall surroundings for every atom in the systems[128].
These kinds of embeddings are thus better than the static descriptors-based techniques since they
naturally adapt to graph structure and hence can be used to represent many complex graphs and
perform good predictions[129, 130].

In this work we utilize the complexity of the Polynorbornene dataset presented in the previous
chapter to build an attention based Graph Convolution Network[129]. This networks builds on top
of the vanilla graph convolutional networks[128] by aggregating the neighborhood information in
order of relative importance of node neighborhood one hop away from the node under
consideration. This method of aggregating node features based on relative importance of its
neighborhood gives a very accurate generalized representation on graph data and hence performs
well compared to vanilla graph convNets. In the next sections we would discuss about attention-
based graph neural network architecture, followed by model performance on the dataset and
subsequently the learning and prediction obtained from global and local feature representation for
these polymers.

5.2. Graph Attention Network (GAT) Computations

Here, we would first describe the general computational framework for a Graph Convolutional
Network (GCN) which involves two essential steps involving message passing and message
aggregation. We would then incorporate the attention-based message aggregation step into this
framework to take into account relative importance of node neighborhood and represent more
general framework capable to learning the underlying graph complexity better than GCN.

59

Given a general graph G with V nodes and E edges as shown in Fig. 5.1. the graph network
computes the embeddings for each node as shown in the Fig. 5.1. by aggregating the neighborhood
information from every node 1-hop away followed by optimization of the embedding thus obtained
using a neural network. This step is propagated through multiple layers and hence at layer k, any
node will be able to aggregate information k-hops away from it.

Figure 5. 1. Input graph representation is shown on the left and the 2-layer graph embedding
computation are shown on the left. The Graph Network performs two steps at each layer i.e. (1)
compute average information from neighbors and (2) apply the non-linearity through the
application of a neural network.

Let us assume each node in the graph has initial features 𝑋
Q
corresponding to node 𝑣 and all the
nodes in the neighborhood of 𝑣 are represented by 𝑁(𝑣). Then we can have the following
representation for every embedding ℎ
Q
(0)
, corresponding to node 𝑣 at layer 𝑙
A
D
B
C
E
F
Input graph
Layer 0
A
C
A
B
E
F
A
B
C
D
A
Layer 1 Layer 2
(1) Avg from
neighbor
(2) Apply NN
(1) Avg from
neighbor
(2) Apply NN
60
ℎ
Q
5
= 𝑋
Q
5.1
ℎ
Q
(0R')
= 𝜎𝑊
0
L
ℎ
/
(0)
|𝑁(𝑣)|
+ 𝐵
0
ℎ
Q
(0)
/E( (Q)
,∀ 𝑙 ∈ {0,1,2,….,𝑙−1} 5.2
𝑍
Q
= ℎ
Q
(S)
5.3
Here, ℎ
Q
5
is the initial feature for the node which is normalized by its in-degree given by 𝑁(𝑣) and
summed over all the neighbors followed by application to ReLU non-linearity to learn the
embedding vector at every layer 𝑙. The equation show above (5.2) encodes two very important
steps into a single equation as it involves first message computation followed by aggregation
computation. Each node first creates a message which needs to be sent to its neighboring nodes,
in this case the message involves linear transformation by performing a multiplication of 𝑊
0
with
ℎ
/
(0)
to yield
𝑚𝑒𝑠𝑠𝑎𝑔𝑒 𝑓𝑢𝑛𝑐𝑡𝑖𝑜𝑛,𝑚
/
(0)
= 𝑊
0
ℎ
/
(0T')
5.4
This step transforms the embeddings from previous layer to appropriate dimension using the
matrix of learnable parameters 𝑊
0
. The next step involves aggregation of these message functions
using an aggregation function 𝐴𝐺𝐺
(0)
( ∙ ), this can be 𝑆𝑢𝑚(∙),𝑀𝑒𝑎𝑛(∙),𝑀𝑎𝑥(∙) 𝑒𝑡𝑐 aggregator.
For graph convolution network this aggregation function is a simple summation over all the
neighbors of respective node as,
𝑎𝑔𝑔𝑟𝑒𝑔𝑎𝑡𝑖𝑜𝑛 𝑓𝑢𝑛𝑐𝑡𝑖𝑜𝑛,ℎ
Q
0
= 𝑆𝑢𝑚>𝑚
/
(0T')
, 𝑢 𝜖 𝑁(𝑣) @ 5.5
An issue that can happen at this aggregation steps is the message loss corresponding to the node
that is aggregating the information from its neighbor and hence another message term is added to
eq 5.5 as,
𝑎𝑔𝑔𝑟𝑒𝑔𝑎𝑡𝑖𝑜𝑛 𝑓𝑢𝑛𝑐𝑡𝑖𝑜𝑛, ℎ
Q
0
= 𝐴𝐺𝐺
(0)
>𝑚
/
(0T')
,𝑢 𝜖 𝑁(𝑣) ,𝑚
Q
(0T')
= 𝐵
(0)
ℎ
Q
(0T')
@ 5.6
61
Further an activation function adds more expressivity to the overall embedding. In the case of GCN
the non-linear activation applied to the net aggregated embedding. Further a normalization term is
added to the embedding which in GCN case is the in-degree of every node, this essential step keeps
the overall calculation stable and does not cause exploding or vanishing gradient problem.

Graph Attention Networks implement eq 5.2 by using a different weighing term for
messages received from every neighbor node, the general form for embedding obtained for a node
v at layer l can be written as shown in eq. 5.6
ℎ
Q
(0 )
= 𝜎R𝑊
0
L 𝛼
Q/
𝑊
(0)
ℎ
/
(0T')
/E( (Q)
T,∀ 𝑙 ∈ {0,1,2,….,𝑙−1} 5.7
The weighing factor 𝛼
Q/
is learned as a byproduct of an attention mechanism 𝑎, where
could be simple neural network, single layer deep. The attention mechanism involves computation
of attention coefficients as shown in eq 5.8-5.9 and figure 5.2. The learned attention coefficient
(𝑒
Q/
) obtained through the attention mechanism (𝑎) are used to compute the final weighing factors
𝛼
Q/
by applying a softmax function to the learned attention coefficients (𝑒
Q/
)
𝑒
Q/
= 𝑎>𝑊
(0)
ℎ
/
(0T')
,𝑊
(0)
ℎ
Q
(0T')
@ 5.8
𝛼
Q/
=
exp (𝑒
Q/
)
∑ exp (𝑒
QG
)
GE( (Q)
5.9

62

Figure 5. 2. The attention mechanism can be represented as shown in this computation graph.
(a) Attention coefficients (𝒆
𝒗𝒖
) are learned for every pair of neighbors and (b) the computed
attention coefficients are used to obtain the weighing factors (𝜶
𝒗𝒖
). The weighing factors are
then used to compute the final embedding 𝒉
𝑨
(𝒍)

The calculation of attention weights (𝛼
Q/
) involves choosing an appropriate attention
mechanism 𝑎. In our present study, we use a smile single layer neural network with trainable
parameters. These parameters are trained jointly with weight matrices and other parameters of
neural net 𝑊
(0)
,𝐵
(0)
in an end-to-end fashion. Hence the attention coefficients can be finally
represented as shown in eq 5.10
𝑒
Q/
= 𝑎>𝑊
(0)
ℎ
/
(0T')
,𝑊
(0)
ℎ
Q
(0T')
@ = 𝐿𝑖𝑛𝑒𝑎𝑟J𝐶𝑜𝑛𝑐𝑎𝑡>𝑊
(0)
ℎ
Z
(0T')
,𝑊
(0)
ℎ
[
(0T')
@K 5.10

A
D
B
C
E
F
(a)
!
!
( #$%)
!
'
( #$%)
"
!'
(b)
A
D
B
C
E
F
!
(
( #$%)
!
'
( #$%)
!
)
( #$%)
#
!(
#
!'
#
!)
!
!
( #)
!
!
#
= % #
!'
&
( #)
!
'
( #$%)
+#
!)
&
( #)
!
)
( #$%)
+#
!(
&
( #)
!
(
( #$%)
63
The attention weights thus computed can be made better by another scheme called multi
head attention computation. This scheme differs from the original scheme by creating multiple
attention scores where each score is a replica with a different set of parameters. The node
embedding thus obtained for different attention scores as shown in eq. 5.11-12 can be aggregated
by using an aggregate function of choice to obtain the final node embedding.
ℎ
Q
(0)
[𝑖] = 𝜎>∑ 𝛼
Q/
(')
𝑊
(0)
ℎ
/
(0T')
/E( (Q)
@ 5.11
ℎ
Q
(0)
= 𝐴𝐺𝐺>ℎ
Q
(0)
[𝑖]@ 5.12

The node embeddings thus obtained from the attention mechanism having many benefits
such as computational efficiency due to parallelization capabilities across all the edges. They are
storage efficient since these computations can be represented easily using sparse matrix
representations hence can be stored in no more than 𝑂(𝑉+𝐸) space. Also, the attention
mechanism is localized since it attends over local network neighborhood and hence has inductive
capability therefore is independent of the global graph structure.

6.3. Graph Attention Network (GAT) Architecture

The message and aggregation function discussion in the previous section clearly outlines the
computations of attention weights as well the various aggregation function computations involved
to obtain node embedding and therefore by extending it to graphs, we can obtain an overall graph
embedding by performing some aggregation over all the nodes of the graph. However, the
architecture of GAT layers in practice is also quite complex involving number of components. In
Fig. 5.3 we show different components involved in a single layer of GAT. The linear layer shown
64
in the figure corresponds to the computation of the message function. The weight matrix used in
linear layer is batch normalized to stabilize the neural network training. Batch Normalization is
followed by DropOut layer which drops certain nodes from the graphs with a probability p such
that every node aggregates its neighbor’s embeddings in a balanced way thus preventing
overfitting. The DropOut layer is succeeded by activation, attention and aggregation layers
respectively. Activation function can be any non-linear function of choice, however for the current
application we have used the ReLU activation. The attention and aggregation layers perform the
computations as detailed in section 5.2.

Figure 5. 3. A single layer of graph attention network many components as shown in this figure.
These layers can then be stacked against one another to obtain deep Graph Attention Networks

These layers can now be stacked against each other to obtain deep Graph Attention
Networks. However, we need to answer a question about how deep is too deep or what should be
an ideal depth for these networks. To address this issue, we need to understand that in a k-layer
deep graph network every node at k
th
layer has information about all neighbors within k hop from
this node. This is in a problem in practice as it would lead to oversmoothing since every node in
the graph will have similar embedding to every other node as information would be propagated
across almost all the nodes in a deep implementation and hence would render the embedding
meaningless. This problem is often overcome by either computing the diameter of graph and
setting the depth of the network to be less than the diameter or empirically choosing an appropriate
depth suitable for the problem. If the problem usually requires many layers implementation to
Single Layer of Graph Attention Network in Practice
Linear BatchNorm DropOut Activation Attention Aggregation
65
enhance the expressiveness of the overall network, it is usually recommended to add shortcuts in
the network by introducing skip connections. These connections create mixture models as N skip
connections can introduce 2
N
possible paths in the network and hence we automatically get a mix
of shall and deep networks as output. However, for our current implementation we restrict our
network to be 3 layers deep without adding any skip connections.

Finally, training is performed on this architecture to learn polymer graphs and predict
dielectric properties for various polymers in the training as well as test dataset. The overall training
pipeline is as shown in Fig. 5.4 and involves finding appropriate graph embeddings using graph
attention networks followed by using a fully connected regular neural network to predict the
dielectric property.

Figure 5. 4. The graph network training pipeline involves graph embedding representation
followed by a fully connected layer to make predictions across labels using L2 loss and RMSE
evaluation metric.

GNN Training Pipeline
Input
Graph
GAT
Graph
Embedding
Fully
connected
layer
Predictions Labels
Evaluation
metrics
(RMSE)
Loss function
(L2 loss)
66
5.4. Results

The computational details as well the layout description for various components in GAT can be
summarized as shown in Fig. 5.5. We represent every polymer chain with a 20-dimensional feature
vector for representing their atomic information along with adjacency list information for their
entire polymer to encode the connectivity information such as bond length and bon order. We use
following atomic properties to represent every polymer chain: Pauling electronegativity, Electron
affinity, group in periodic table , covalent radius, ionisation energy, s electrons, p electrons, atomic
polarizability, ionization affinity ratio , liquid range, liquid ratio, ratio of electron affinity to electro
negativity. This input feature is fed to 1
st
layer of Graph Neural Network (GNN) and subsequently
to 2
nd
and 3
rd
layers as shown in Fig. 5.5. At each layer convolution is performed as shown in eq.
5.7. The final embedding from 3
rd
layer is fed to a fully connected layer which learns the
predictions on dielectric constants by minimizing the L2 distance between the predicted and actual
values.

Figure 5. 5. The input polymer represented as a graph along with 20-dimensional feature vector
enables the GAT to learn a 64-dimensional graph embedding which is used to make predictions
using fully connected layer.

67
The model is effectively able to capture the global as well as local features dominating
responsible for enhancing or degrading the polymer performance. The learning efficacy of the
model is well represented by studying the evolution of error over successive epochs during the
training process as shown in Fig 5.6. The errors quickly drop down below 10% after 300 epochs
and subsequent model training does not lead significant performance increase as evident by the
error evolution over subsequent epochs. The prediction capability of the model is also clearly
shown in Fig. 5.6(b) where the predicted values for various structures and their ground truths value
agree well with each other.

Figure 5. 6. The fitting obtained from graph embeddings generated using GAT are shown in (a)
and (b) demonstrates the superior power of graph networks in learning the overall polymer
property.
The trained GNN model is effectively able to capture global as well as local features which
is evident from the excellent predictive power of this model. The local feature information further
provides information about complex interplay between various functional groups and atoms within
the polymer structure as shown in Fig. 5.7. By visualizing the global feature space as shown in
Fig. 5.7 we were able to locate areas of high and low performance in our dataset. We further
visualized all the local features of various polymers in this region and screened out a few polymer
Training Data: 3200
Test Data: 612
(a) (b)
68
candidates. Some of these have been validated as high performance polymer in a previous
study[99]. We identified Poly oxo fluoro norbornene (POFNB), Poly fluoro norbornene (PFNB),
Poly reduced fluoro norbornene, and Poly oxo chloro norbornene (POClNB) as high-performance
polymer suitable for applications in high energy dielectric systems. A recent study has also
corroborated our work by analyzing POFNB as a very high-performance polymer with high
dielectric constant (~2.5 at 10kHz, 300 K), high glass transition temperature (Tg) of 410 K and
high electric field energy storage.
We identified two other polymers PFNB and PRFNB with high dielectric constant as
reported in figure 5.8(c) we present here PFNB which has been synthesized in labs and has shown
superior dielectric property. POFNB and PFNB share almost exactly same structure except for
presence of oxygen hetero atom in the POFNB and absence of Oxygen hetero atom for Carbon in
PFNB. This small change in the polymer changes the properties of these two polymers as O in
POFNB is large atom with high atomic radius this contributing more to the net atomic volume and
hence resulting in reduced net dipole density as compared to PFNB (where O is replaced with C
in the backbone).

PCA
t-SNE
oPOFNB
PFNB
69
Figure 5. 7. The 64-dimensional graph embedding obtained from final layer of GAT is reduced
dimensionally and plotted using PCA and t-SNE based dimensionality reduction techniques. T-
SNE plot clearly shows the proximity of oPOFNB and PFNB.

Further, we identified another polymer POClNB as a high-performance candidate whose
dielectric. POClNB also exhibits high dielectric constant throughout the frequency range. This
polymer’s high performance in terms of dielectric constant can be explained by the loss of
symmetric bulky CF3 groups in favor of single Cl group which reduced the overall system volume
and adds more to the net dipole moment of the system thus resulting in higher performance as
compared to its CF3 based derivatives.

Figure 5. 8. Local feature visualization for various polymer shows the importance of CF3 group
attached at ortho in (a) oPOFNB, meta position also effect the polymer property as shown in (b)
mPOFNB and (c) PFNB also shows important contribution wrt to fluorine containing pendant
group.

5.5. Conclusion

In this work, we computed the dielectric constant of various PNB related derivatives using
ReaxPQ-v forcefield and our high throughput workflow. The screening of potential candidates
from the computed values which vary between 1.5-3.5 for various PNB derivatives is complex due
oPOFNB
(a) (b) (c)
mPOFNB PFNB
70
to the interplay of various factors such as backbones, heteroatom, pendant groups etc. Hence, we
employed a Graph Neural Network based analysis to screen out potentially useful candidates for
high energy dielectric applications. Our screened candidates (POFNB, PFNB, PRFNB, POClNB)
have been experimentally validated to be high performance of polymer with high dielectric
constant over the frequency range of 10-10KHz. Moreover, the dissipation factor over the same
frequency range for these polymers maintain an almost constant profile, thus these polymers can
store high electrical energy and exhibit low losses at increased electric fields. The GNN based
screening procedure can thus greatly help to identify promising polymer candidates from a
complex and combinatorically large space exhibiting high performance by studying the attention
score based local features made available by GNN.

71
Conclusion
In this dissertation, polarizable reactive molecular dynamics was used to study high energy
dielectric polymers. These polymer materials are essential for various applications including
electronics, lasers etc. However, estimation of their dielectric property is complex due to the long-
time scales associated with relaxation of dipole moments in long chain polymer systems. In this
thesis, a scalable parallel framework was developed for calculating the dielectric constant for any
polymer system. Another important aspect identified for these polymer systems is the difficulty to
design new polymers based on information available for the already existing polymers. A known
base polymer structure can have various derivatives based on several groups that can be attached
to it, hence producing various defects in the original system which can affect the dielectric
properties in various complex ways. To address this problem, we used this scalable framework for
estimating the dielectric properties to calculate properties for a combinatorically large polymer
space. Further, attention-based graph neural networks (GAT) were used to identify important
polymers based on their ease of synthesizability by assessing their proximity to already existing
polymers. This was done by obtaining the graph embedding for various polymers using the GAT
architecture. Once the polymers were identified, local feature analysis was done to obtain more
information about the role played by various groups attached to the polymer in affecting the final
property.

In the first part of this dissertation, we developed a scalable computational framework as
mentioned above to evaluate the temperature-dependent dielectric responses of amorphous
polymer. The seamless creation of amorphous polymer structure realized by SMILES strings,
Open Babel, and RXMD will allow researchers to explore the vast parameter space of high energy
72
polymer design in a highly automated fashion. Equipped with the scalable parallel all-atom RMD
simulation engine, it is possible to incorporate complex repeat units with industry-relevant chain
lengths. The model parameters, such as atomic polarizability, core-shell charges, and radii, are
extensively optimized and validated by experiments and first-principal calculations using the high-
frequency as well as temperature-dependent dielectric constants. Our computational framework
has been successfully applied to the recently identified four high energy density polymers, PDTC-
HK511, TDI-EDR148, DVS-BCB and FPE, and revealed that the thermally activated flexible
segment motion in polymer chains facilitates the fluctuation of the molecular dipole moments. The
obtained dielectric constants as well as its increasing trend quantitatively agree with the
experiments. The computational framework presented here for the first time realizes high
throughput theoretical screenings incorporating several key parameters in the high energy density
polymer design. Therefore, combined with efficient experimental synthesis and characterization,
our framework is expected to accelerate the discovery of next-generation high energy density
polymer and device development.
In the second part we presented a polymer base class Polynorbornene (PNB) and the
combinatorial large derivative space around this polymer created by incorporation of various
defect in this polymer backbone. Various substitutions for PNB results in ~3800 polymer
candidates that need to be screened to identify promising polymers suitable for synthesizing in
laboratories. These polymers are created using the high throughput framework and amorphized
using the same framework to finally compute the dielectric constants. The dielectric constant
values thus obtained follow certain theoretically established rules as imide containing groups
altogether have a lower mean dielectric value compared to co-polymerized imides with non-imides
followed by higher mean values for non-imide groups. We further analyze a few defects out of this
73
large space to study the effect they can have on the overall dielectric property of the polymer
system. Here, we analyzed primarily two kinds of defects, (i) defects introduced by introduction
of various pendant groups in POFNBs and (ii) free volume introduction in PFNBs. Both defects
occur in various forms in the polymer space being considered here for subsequent data driven
studies. These analyses verify certain experimental hypothesis and provide an overview of the
interplay among various complex derivatives in the PNB base class to alter dielectric property in
so many ways. This complexity makes it very difficult to screen desired polymers from this
complex polymer space.
In the last part, we performed screening of potential candidates from the computed values
which vary between 1.5-3.5 for various PNB derivatives. This step was complex due to the
interplay of various factors such as backbones, heteroatom, pendant groups etc. Hence, we
employed a Graph Neural Network based analysis to screen out potentially useful candidates for
high energy dielectric applications. Our screened candidates (POFNB, PFNB, PRFNB, POClNB)
have been experimentally validated to be high performance of polymer with high dielectric
constant over the frequency range of 10-10KHz. Moreover, the dissipation factor over the same
frequency range for these polymers maintain an almost constant profile, thus these polymers can
store high electrical energy and exhibit low losses at increased electric fields. The GNN based
screening procedure can thus greatly help to identify promising polymer candidates from a
complex and combinatorically large space exhibiting high performance by studying the attention
score based local features made available by GNN.

74

REFERENCES

[1] T. Wei et al., "Aromatic Polyamide Reverse-Osmosis Membrane: An Atomistic
Molecular Dynamics Simulation," The Journal of Physical Chemistry B, vol. 120, no. 39,
pp. 10311-10318, 2016/10/06 2016, doi: 10.1021/acs.jpcb.6b06560.
[2] H. W. Starkweather, P. Avakian, R. R. Matheson, J. J. Fontanella, and M. C. Wintersgill,
"Ultralow temperature dielectric relaxations in polyolefins," Macromolecules, vol. 25, no.
25, pp. 6871-6875, 1992/12/01 1992, doi: 10.1021/ma00051a023.
[3] E. J. Barshaw et al., "High Energy Density (HED) Biaxially-Oriented Poly-Propylene
(BOPP) Capacitors For Pulse Power Applications," IEEE Transactions on Magnetics,
vol. 43, no. 1, pp. 223-225, 2007, doi: 10.1109/TMAG.2006.887682.
[4] J. Ho, R. Ramprasad, and S. Boggs, "Effect of Alteration of Antioxidant by UV
Treatment on the Dielectric Strength of BOPP Capacitor Film," IEEE Transactions on
Dielectrics and Electrical Insulation, vol. 14, no. 5, pp. 1295-1301, 2007, doi:
10.1109/TDEI.2007.4339492.
[5] M. Rabuffi and G. Picci, "Status quo and future prospects for metallized polypropylene
energy storage capacitors," IEEE Transactions on Plasma Science, vol. 30, no. 5, pp.
1939-1942, 2002, doi: 10.1109/TPS.2002.805318.

[6] B. Chu et al., "A Dielectric Polymer with High Electric Energy Density and Fast
Discharge Speed," Science, vol. 313, no. 5785, p. 334, 2006, doi:
10.1126/science.1127798.

75
[7] P. Kim et al., "High Energy Density Nanocomposites Based on Surface-Modified
BaTiO3 and a Ferroelectric Polymer," ACS Nano, vol. 3, no. 9, pp. 2581-2592,
2009/09/22 2009, doi: 10.1021/nn9006412.

[8] F. Häse, L. M. Roch, and A. Aspuru-Guzik, "Next-Generation Experimentation with
Self-Driving Laboratories," Trends in Chemistry, vol. 1, pp. 282-291, 2019, doi:
10.1016/j.trechm.2019.02.007.
[9] D. M. Parry, "Closing the Loop: Developing an Integrated Design, Make, and Test
Platform for Discovery," ACS Medicinal Chemistry Letters, vol. 10, pp. 848-856, 2019,
doi: 10.1021/acsmedchemlett.9b00095.
[10] S. Lu, Q. Zhou, Y. Ouyang, Y. Guo, Q. Li, and J. Wang, "Accelerated discovery of stable
lead-free hybrid organic-inorganic perovskites via machine learning," Nature
Communications, vol. 9, p. 3405, 2018, doi: 10.1038/s41467-018-05761-w.
[11] S. K. Saikin, C. Kreisbeck, D. Sheberla, J. S. Becker, and A. Aspuru-Guzik, "Closed-loop
discovery platform integration is needed for artificial intelligence to make an impact in
drug discovery," Expert Opinion on Drug Discovery, vol. 14, pp. 1-4, 2019, doi:
10.1080/17460441.2019.1546690.
[12] L. Bassman et al., "Active learning for accelerated design of layered materials," npj
Computational Materials, vol. 4, no. 1, p. 74, 2018/12/10 2018, doi: 10.1038/s41524-
018-0129-0.
[13] G. X. Gu, C.-T. Chen, and M. J. Buehler, "De novo composite design based on machine
learning algorithm," Extreme Mechanics Letters, vol. 18, pp. 19-28, 2018/01/01/ 2018,
doi: https://doi.org/10.1016/j.eml.2017.10.001.
76
[14] J. Yeo et al., "Materials-by-design: computation, synthesis, and characterization from
atoms to structures," Physica Scripta, vol. 93, no. 5, p. 053003, 2018/04/17 2018, doi:
10.1088/1402-4896/aab4e2.
[15] Y. Li et al., "Machine Learning Force Field Parameters from Ab Initio Data," Journal of
Chemical Theory and Computation, vol. 13, no. 9, pp. 4492-4503, 2017/09/12 2017, doi:
10.1021/acs.jctc.7b00521.
[16] R. Car and M. Parrinello, "Unified Approach for Molecular Dynamics and Density-
Functional Theory," Physical Review Letters, vol. 55, no. 22, pp. 2471-2474, 11/25/
1985, doi: 10.1103/PhysRevLett.55.2471.
[17] M. C. Payne, M. P. Teter, D. C. Allan, T. A. Arias, and J. D. Joannopoulos, "Iterative
minimization techniques for ab initio total-energy calculations: molecular dynamics and
conjugate gradients," Reviews of Modern Physics, vol. 64, no. 4, pp. 1045-1097, 10/01/
1992, doi: 10.1103/RevModPhys.64.1045.
[18] H.-P. Chen, R. K. Kalia, A. Nakano, P. Vashishta, and I. Szlufarska, "Multimillion-atom
nanoindentation simulation of crystalline silicon carbide: Orientation dependence and
anisotropic pileup," Journal of Applied Physics, vol. 102, no. 6, p. 063514, 2007, doi:
10.1063/1.2781324.
[19] M. Misra, A. Mannodi-Kanakkithodi, T. C. Chung, R. Ramprasad, and S. K. Kumar,
"Critical role of morphology on the dielectric constant of semicrystalline polyolefins,"
The Journal of Chemical Physics, vol. 144, no. 23, p. 234905, 2016, doi:
10.1063/1.4953182.
77
[20] L. Chen, T. D. Huan, and R. Ramprasad, "Electronic Structure of Polyethylene: Role of
Chemical, Morphological and Interfacial Complexity," Scientific Reports, vol. 7, no. 1, p.
6128, 2017/07/21 2017, doi: 10.1038/s41598-017-06357-y.
[21] S. Kwon, S. Naserifar, H. M. Lee, and W. A. Goddard, "Polarizable Charge Equilibration
Model for Transition-Metal Elements," The Journal of Physical Chemistry A, vol. 122,
no. 48, pp. 9350-9358, 2018/12/06 2018, doi: 10.1021/acs.jpca.8b07290.
[22] K. Hasegawa, T. Deushi, O. Yaegashi, Y. Miyashita, and S. Sasaki, "Artificial neural
network studies in quantitative structure-activity relationships of antifungal azoxy
compounds," Eur. J. Med. Chem, vol. 30, no. 7-8, pp. 569-574, 1995, doi: 10.1016/0223-
5234(96)88271-7.
[23] T. Le, V. C. Epa, F. R. Burden, and D. A. Winkler, "Quantitative Structure–Property
Relationship Modeling of Diverse Materials Properties," Chem. Rev., vol. 112, no. 5, pp.
2889-2919, 2012, doi: 10.1021/cr200066h.
[24] G. Sliwoski, S. Kothiwale, J. Meiler, E. W. Lowe, and E. L. Barker, "Computational
Methods in Drug Discovery," Pharmacol. Rev., vol. 66, no. 1, pp. 334-395, 2013, doi:
10.1124/pr.112.007336.
[25] G. Idakwo, S. Thangapandian, J. Luttrell, Z. Zhou, C. Zhang, and P. Gong, "Deep
Learning-Based Structure-Activity Relationship Modeling for Multi-Category Toxicity
Classification: A Case Study of 10K Tox21 Chemicals With High-Throughput Cell-
Based Androgen Receptor Bioassay Data," Front. Physiol., vol. 10, pp. 1-13, 2019, doi:
10.3389/fphys.2019.01044.
78
[26] O. A. von Lilienfeld and K. Burke, "Retrospective on a decade of machine learning for
chemical discovery," Nat Commun, vol. 11, p. 4895, Sep 29 2020, doi:
https://doi.org/10.1038/s41467-020-18556-9.
[27] A. Jain, K. A. Persson, and G. Ceder, "Research update: the materials genome initiative:
Data sharing and the impact of collaborative ab initio databases," (in English), APL
Mater, vol. 4, no. 5, p. 053102, May 2016, doi: 10.1063/1.4944683.
[28] K. T. Butler, D. W. Davies, H. Cartwright, O. Isayev, and A. Walsh, "Machine learning
for molecular and materials science," (in English), Nature, vol. 559, no. 7715, pp. 547-
555, Jul 26 2018, doi: 10.1038/s41586-018-0337-2.
[29] J. H. Van Drie, "Computer-aided drug design: the next 20 years," J. Comput. Aided
Mater. Des., vol. 21, no. 10-11, pp. 591-601, 2007, doi: 10.1007/s10822-007-9142-y.
[30] V. Sharma et al., "Rational design of all organic polymer dielectrics," (in English), Nat
Commun, vol. 5, p. 4845, Sep 2014, doi: 10.1038/ncomms5845.
[31] C. Kim, G. Pilania, and R. Ramprasad, "From organized high-throughput data to
phenomenological theory using machine learning: the example of dielectric breakdown,"
(in English), Chem. Mater., vol. 28, no. 5, pp. 1304-1311, Mar 8 2016, doi:
10.1021/acs.chemmater.5b04109.
[32] C. Kim, A. Chandrasekaran, T. D. Huan, D. Das, and R. Ramprasad, "Polymer Genome:
A Data-Powered Polymer Informatics Platform for Property Predictions," J. Phys. Chem.
C, vol. 122, no. 31, pp. 17575-17585, 2018, doi: 10.1021/acs.jpcc.8b02913.
[33] A. Mishra et al., "Materials Genome Software Framework: Scalable Parallel Simulation,
Virtual Reality Visualization and Machine Learning," in Int. Conf. Sci. Comput., Las
Vegas, NV, USA, 2019 2019, pp. 125-131.
79
[34] S. Naserifar, D. J. Brooks, W. A. Goddard, and V. Cvicek, "Polarizable charge
equilibration model for predicting accurate electrostatic interactions in molecules and
solids," (in English), J. Chem. Phys, vol. 146, no. 12, p. 124117, Mar 28 2017, doi:
10.1063/1.4978891.
[35] K. Liu et al., "Shift-collapse acceleration of generalized polarizable reactive molecular
dynamics for machine learning-assisted computational synthesis of layered materials,"
Proc ScalA18, pp. 41-48, Nov 12 IEEE, 2018, doi: 10.1109/ScalA.2018.00009.
[36] Y. Li et al., "Scalable reactive molecular dynamics simulations for computational
synthesis," Comput Sci Eng, vol. 21, no. 5, pp. 64-75, Sep/Oct 2019, doi:
10.1109/MCSE.2018.110150043.
[37] T. Deng and G.-z. Jia, "Prediction of aqueous solubility of compounds based on neural
network," Mol. Phys., vol. 118, no. 2, p. e1600754, 2019, doi:
10.1080/00268976.2019.1600754.
[38] E. Heidari, M. A. Sobati, and S. Movahedirad, "Accurate prediction of nanofluid
viscosity using a multilayer perceptron artificial neural network (MLP-ANN)," Chemom.
Intell. Lab. Syst., vol. 155, pp. 73-85, 2016, doi: 10.1016/j.chemolab.2016.03.031.
[39] T. Xie and J. C. Grossman, "Crystal Graph Convolutional Neural Networks for an
Accurate and Interpretable Prediction of Material Properties," Physical Review Letters,
vol. 120, no. 14, p. 145301, 04/06/ 2018, doi: 10.1103/PhysRevLett.120.145301.
[40] K. C. Liu, K.-i. Nomura, R. K. Kalia, A. Nakano, P. Vashishta, and P. Rajak, "Graph
Neural Network Analysis Of Layered Material Phases," 2019.
80
[41] L. Bernazzani et al., "Predicting Physical−Chemical Properties of Compounds from
Molecular Structures by Recursive Neural Networks," J. Chem. Inf. Model., vol. 46, no.
5, pp. 2030-2042, 2006, doi: 10.1021/ci060104e.
[42] Q. M. Zhang, V. Bharti, and X. Zhao, "Giant electrostriction and relaxor ferroelectric
behavior in electron-irradiated poly(vinylidene fluoride-trifluoroethylene) copolymer,"
(in English), Science, vol. 280, no. 5372, pp. 2101-2104, Jun 26 1998, doi:
10.1126/science.280.5372.2101
[43] B. J. Chu et al., "A dielectric polymer with high electric energy density and fast discharge
speed," (in English), Science, vol. 313, no. 5785, pp. 334-336, Jul 21 2006, doi:
10.1126/science.1127798
[44] A. Jain, K. A. Persson, and G. Ceder, "Research update: the materials genome initiative:
Data sharing and the impact of collaborative ab initio databases," (in English), APL
Mater, vol. 4, no. 5, May 2016, doi: 10.1063/1.4944683.
[45] C. Kim, A. Chandrasekaran, T. D. Huan, D. Das, and R. Ramprasad, "Polymer genome: a
data-powered polymer informatics platform for property predictions," (in English), J
Phys Chem C, vol. 122, no. 31, pp. 17575-17585, Aug 9 2018, doi:
10.1021/acs.jpcc.8b02913.
[46] K. Andersen, S. Latini, and K. S. Thygesen, "Dielectric genome of van der Waals
heterostructures," (in English), Nano Lett, vol. 15, no. 7, pp. 4616-4621, Jul 2015, doi:
10.1021/acs.nanolett.5b01251.
[47] P. Umari and A. Pasquarello, "Ab initio molecular dynamics in a finite homogeneous
electric field," (in English), Phys Rev Lett, vol. 89, no. 15, p. 157602, Oct 7 2002, doi:
10.1103/PhysRevLett.89.157602.
81
[48] I. Souza, J. Iniguez, and D. Vanderbilt, "First-principles approach to insulators in finite
electric fields," (in English), Phys Rev Lett, vol. 89, no. 11, p. 117602, Sep 9 2002, doi:
10.1103/PhysRevLett.89.117602.
[49] S. Fukushima et al., "Effects of chemical defects on anisotropic dielectric response of
polyethylene," AIP Advances, vol. 9, no. 4, p. 045022, 2019, doi: 10.1063/1.5093566.
[50] S. Naserifar, D. J. Brooks, W. A. Goddard, and V. Cvicek, "Polarizable charge
equilibration model for predicting accurate electrostatic interactions in molecules and
solids," The Journal of Chemical Physics, vol. 146, no. 12, p. 124117, 2017/03/28 2017,
doi: 10.1063/1.4978891.
[51] S. Yip, Handbook of materials modeling. Springer, 2005.
[52] H. Stanley, S. Buldyrev, O. Mishima, M. Sadr-Lahijany, A. Scala, and F. Starr,
"Unsolved mysteries of water in its liquid and glassy phases," Journal of Physics:
Condensed Matter, vol. 12, no. 8A, p. A403, 2000.
[53] A. Rahman, "Correlations in the motion of atoms in liquid argon," Physical Review, vol.
136, no. 2A, p. A405, 1964.
[54] A. C. T. van Duin, S. Dasgupta, F. Lorant, and W. A. Goddard, "ReaxFF: A Reactive
Force Field for Hydrocarbons," The Journal of Physical Chemistry A, vol. 105, no. 41,
pp. 9396-9409, 2001/10/01 2001, doi: 10.1021/jp004368u.
[55] P. Hohenberg and W. Kohn, "Inhomogeneous Electron Gas," Physical Review, vol. 136,
no. 3B, pp. B864-B871, 11/09/ 1964.
[56] A. C. T. van Duin, S. Dasgupta, F. Lorant, and W. A. Goddard, "ReaxFF: A reactive
force field for hydrocarbons," (in English), J. Phys. Chem. A, vol. 105, no. 41, pp. 9396-
9409, Oct 18 2001.
82
[57] L. Liu, Y. Liu, S. V. Zybin, H. Sun, and W. A. Goddard, "ReaxFF-lg: Correction of the
ReaxFF Reactive Force Field for London Dispersion, with Applications to the Equations
of State for Energetic Materials," The Journal of Physical Chemistry A, vol. 115, no. 40,
pp. 11016-11022, 2011/10/13 2011, doi: 10.1021/jp201599t.
[58] B. M. Rice, J. P. Larentzos, E. F. C. Byrd, and N. S. Weingarten, "Parameterizing
Complex Reactive Force Fields Using Multiple Objective Evolutionary Strategies
(MOES): Part 2: Transferability of ReaxFF Models to C–H–N–O Energetic Materials,"
Journal of Chemical Theory and Computation, vol. 11, no. 2, pp. 392-405, 2015/02/10
2015, doi: 10.1021/ct5007899.
[59] T. Zhou, H. Song, Y. Liu, and F. Huang, "Shock initiated thermal and chemical responses
of HMX crystal from ReaxFF molecular dynamics simulation," Physical Chemistry
Chemical Physics, 10.1039/C4CP00890A vol. 16, no. 27, pp. 13914-13931, 2014, doi:
10.1039/C4CP00890A.
[60] V. Sharma et al., "Rational design of all organic polymer dielectrics," Nature
Communications, Article vol. 5, p. 4845, 09/17/online 2014, doi: 10.1038/ncomms5845
https://www.nature.com/articles/ncomms5845#supplementary-information.
[61] A. J. Liu et al., "Opportunities in theoretical and computational polymeric materials and
soft matter," Soft Matter, 10.1039/C4SM02344G vol. 11, no. 12, pp. 2326-2332, 2015,
doi: 10.1039/C4SM02344G.
[62] A. Mannodi-Kanakkithodi et al., "Scoping the polymer genome: A roadmap for rational
polymer dielectrics design and beyond," Materials Today, vol. 21, no. 7, pp. 785-796,
2018/09/01/ 2018, doi: https://doi.org/10.1016/j.mattod.2017.11.021.
83
[63] G. M. Treich et al., "A rational co-design approach to the creation of new dielectric
polymers with high energy density," IEEE Transactions on Dielectrics and Electrical
Insulation, vol. 24, no. 2, pp. 732-743, 2017, doi: 10.1109/TDEI.2017.006329.
[64] R. Ma et al., "Rational design and synthesis of polythioureas as capacitor dielectrics,"
Journal of Materials Chemistry A, 10.1039/C5TA01252J vol. 3, no. 28, pp. 14845-
14852, 2015, doi: 10.1039/C5TA01252J.
[65] Z. Li, G. M. Treich, S. K. Scheirey, A. Gregory, Sotzing, and Y. Cao, "A novel aromatic
polyurea for high energy density capacitors," in 2017 IEEE Conference on Electrical
Insulation and Dielectric Phenomenon (CEIDP), 22-25 Oct. 2017 2017, pp. 78-81, doi:
10.1109/CEIDP.2017.8257508.
[66] Z. Li, C. Xu, H. Uehara, S. Boggs, and Y. Cao, "Transient characterization of extreme
field conduction in dielectrics," AIP Advances, vol. 6, no. 11, p. 115025, 2016/11/01
2016, doi: 10.1063/1.4971158.
[67] M. E. Mills, P. Townsend, D. Castillo, S. Martin, and A. Achen, "Benzocyclobutene
(DVS-BCB) polymer as an interlayer dielectric (ILD) material," presented at the
Proceedings of the symposium J of the 1996 E-MRS Spring meeting conference on
Advanced materials for interconnections, Strasbourg, France, 1997.
[68] Q. Li et al., "Flexible high-temperature dielectric materials from polymer
nanocomposites," Nature, vol. 523, p. 576, 07/29/online 2015, doi: 10.1038/nature14647
https://www.nature.com/articles/nature14647#supplementary-information.
[69] K. Liu et al., "Shift-Collapse Acceleration of Generalized Polarizable Reactive Molecular
Dynamics for Machine Learning-Assisted Computational Synthesis of Layered
Materials," in 2018 IEEE/ACM 9th Workshop on Latest Advances in Scalable Algorithms
84
for Large-Scale Systems (scalA), 12-12 Nov. 2018 2018, pp. 41-48, doi:
10.1109/ScalA.2018.00009.
[70] A. Hellweg and D. Rappoport, "Development of new auxiliary basis functions of the
Karlsruhe segmented contracted basis sets including diffuse basis functions (def2-SVPD,
def2-TZVPPD, and def2-QVPPD) for RI-MP2 and RI-CC calculations," Physical
Chemistry Chemical Physics, 10.1039/C4CP04286G vol. 17, no. 2, pp. 1010-1017, 2015,
doi: 10.1039/C4CP04286G.
[71] T. M. Miller and B. Bederson, "Atomic and Molecular Polarizabilities-A Review of
Recent Advances," in Advances in Atomic and Molecular Physics, vol. 13, D. R. Bates
and B. Bederson Eds.: Academic Press, 1978, pp. 1-55.
[72] F. Shimojo et al., "A divide-conquer-recombine algorithmic paradigm for large
spatiotemporal quantum molecular dynamics simulations," The Journal of Chemical
Physics, vol. 140, no. 18, p. 18A529, 2014, doi: 10.1063/1.4869342.
[73] D. Weininger, "SMILES, a chemical language and information system. 1. Introduction to
methodology and encoding rules," pp. 31--36, 1988, doi: 10.1021/ci00057a005.
[74] N. M. O'Boyle, M. Banck, C. A. James, C. Morley, T. Vandermeersch, and G. R.
Hutchison, "Open Babel: An open chemical toolbox," Journal of Cheminformatics, vol.
3, no. 1, p. 33, 2011/10/07 2011, doi: 10.1186/1758-2946-3-33.
[75] K.-i. Nomura, R. K. Kalia, A. Nakano, P. Rajak, and P. Vashishta, "RXMD: A scalable
reactive molecular dynamics simulator for optimized time-to-solution," SoftwareX, vol.
11, p. 100389, 2020/01/01/ 2020, doi: https://doi.org/10.1016/j.softx.2019.100389.
[76] F. MPI, "MPI: A Message-Passing Interface," Oregon Graduate Institute School of
Science \& Engineering, 1994.
85
[77] M. H. Nafar Sefiddashti, B. J. Edwards, and B. Khomami, "Steady shearing flow of a
moderately entangled polyethylene liquid," Journal of Rheology, vol. 60, no. 6, pp. 1227-
1244, 2016/11/01 2016, doi: 10.1122/1.4963800.
[78] M. H. Nafar Sefiddashti, B. J. Edwards, and B. Khomami, "Evaluation of reptation-based
modeling of entangled polymeric fluids including chain rotation via nonequilibrium
molecular dynamics simulation," Physical Review Fluids, vol. 2, no. 8, p. 083301, 08/09/
2017, doi: 10.1103/PhysRevFluids.2.083301.
[79] M. H. Nafar Sefiddashti, B. J. Edwards, and B. Khomami, "Communication: A coil-
stretch transition in planar elongational flow of an entangled polymeric melt," The
Journal of Chemical Physics, vol. 148, no. 14, p. 141103, 2018/04/14 2018, doi:
10.1063/1.5026792.
[80] C. N. Edwards, M. H. Nafar Sefiddashti, B. J. Edwards, and B. Khomami, "In-plane and
out-of-plane rotational motion of individual chain molecules in steady shear flow of
polymer melts and solutions," Journal of Molecular Graphics and Modelling, vol. 81, pp.
184-196, 2018/05/01/ 2018, doi: https://doi.org/10.1016/j.jmgm.2018.03.003.
[81] I. Souza, J. Íñiguez, and D. Vanderbilt, "First-Principles Approach to Insulators in Finite
Electric Fields," Physical Review Letters, vol. 89, no. 11, p. 117602, 08/26/ 2002, doi:
10.1103/PhysRevLett.89.117602.
[82] N. Pfeiffenberger, F. Milandou, M. Niemeyer, T. Sugawara, M. Sanner, and J. Mahood,
"High temperature dielectric polyetherimide film development," IEEE Transactions on
Dielectrics and Electrical Insulation, vol. 25, no. 1, pp. 120-126, 2018, doi:
10.1109/TDEI.2018.006806.
86
[83] Q. Li, F.-Z. Yao, Y. Liu, G. Zhang, H. Wang, and Q. Wang, "High-Temperature
Dielectric Materials for Electrical Energy Storage," Annual Review of Materials
Research, vol. 48, no. 1, pp. 219-243, 2018, doi: 10.1146/annurev-matsci-070317-
124435.
[84] Q. Tan, P. Irwin, and Y. Cao, "Advanced Dielectrics for Capacitors," IEEJ Transactions
on Fundamentals and Materials, vol. 126, no. 11, pp. 1153-1159, 2006, doi:
10.1541/ieejfms.126.1153.
[85] D. Tan, L. Zhang, Q. Chen, and P. Irwin, "High-Temperature Capacitor Polymer Films,"
Journal of Electronic Materials, journal article vol. 43, no. 12, pp. 4569-4575, December
01 2014, doi: 10.1007/s11664-014-3440-7.
[86] J. Xu and C. P. Wong, "Effect of the polymer matrices on the dielectric behavior of a
percolative high-k polymer composite for embedded capacitor applications," Journal of
Electronic Materials, journal article vol. 35, no. 5, pp. 1087-1094, May 01 2006, doi:
10.1007/bf02692571.
[87] A. Modafe, N. Ghalichechian, B. Kleber, and R. Ghodssi, "Electrical characterization of
benzocyclobutene polymers for electric micromachines," IEEE Transactions on Device
and Materials Reliability, vol. 4, no. 3, pp. 495-508, 2004, doi:
10.1109/TDMR.2004.830289.
[88] W. J. Sarjeant, J. Zirnheld, and F. W. MacDougall, "Capacitors," IEEE Transactions on
Plasma Science, vol. 26, no. 5, pp. 1368-1392, 1998, doi: 10.1109/27.736020.
[89] J. Jiang et al., "Polymer Nanocomposites with Interpenetrating Gradient Structure
Exhibiting Ultrahigh Discharge Efficiency and Energy Density," Advanced Energy
Materials, vol. 9, no. 15, p. 1803411, 2019/04/01 2019, doi: 10.1002/aenm.201803411.
87
[90] H. Pan et al., "Ultrahigh–energy density lead-free dielectric films via polymorphic
nanodomain design," Science, vol. 365, no. 6453, p. 578, 2019, doi:
10.1126/science.aaw8109.
[91] J. Watson and G. Castro, "A review of high-temperature electronics technology and
applications," Journal of Materials Science: Materials in Electronics, vol. 26, no. 12, pp.
9226-9235, 2015/12/01 2015, doi: 10.1007/s10854-015-3459-4.
[92] C. Buttay et al., "State of the art of high temperature power electronics," Materials
Science and Engineering: B, vol. 176, no. 4, pp. 283-288, 2011/03/15/ 2011, doi:
https://doi.org/10.1016/j.mseb.2010.10.003.
[93] S. Nasreen et al., "Polymer Dielectrics for Capacitor Application," Kirk‐Othmer
Encyclopedia of Chemical Technology, pp. 1-29, 2020/05/02 2017.
[94] A. F. Baldwin et al., "Poly(dimethyltin glutarate) as a Prospective Material for High
Dielectric Applications," Advanced Materials, vol. 27, no. 2, pp. 346-351, 2015/01/01
2015, doi: 10.1002/adma.201404162.
[95] J. S. Ho and S. G. Greenbaum, "Polymer Capacitor Dielectrics for High Temperature
Applications," ACS Applied Materials & Interfaces, vol. 10, no. 35, pp. 29189-29218,
2018/09/05 2018, doi: 10.1021/acsami.8b07705.
[96] Y. Zhou et al., "Superexchange Effects on Oxygen Reduction Activity of Edge-Sharing
[CoxMn1−xO6] Octahedra in Spinel Oxide," Advanced Materials, vol. 30, no. 11, p.
1705407, 2018/03/01 2018, doi: 10.1002/adma.201705407.
[97] M. Tefferi, R. Ma, G. Treich, G. Sotzing, R. Ramprasad, and Y. Cao, "Novel dielectric
films with high energy density," in 2015 IEEE Conference on Electrical Insulation and
88
Dielectric Phenomena (CEIDP), 18-21 Oct. 2015 2015, pp. 451-454, doi:
10.1109/CEIDP.2015.7352142.
[98] A. Mannodi-Kanakkithodi et al., "Rational Co-Design of Polymer Dielectrics for Energy
Storage," Advanced Materials, vol. 28, no. 30, pp. 6277-6291, 2016/08/01 2016, doi:
10.1002/adma.201600377.
[99] C. Wu et al., "Flexible Temperature-Invariant Polymer Dielectrics with Large Bandgap,"
Advanced Materials, vol. n/a, no. n/a, p. 2000499, 2020/04/06 2020, doi:
10.1002/adma.202000499.
[100] A. Mishra et al., "Computational framework for polymer synthesis to study dielectric
properties using polarizable reactive molecular dynamics," ACS Central Sci, p.
submitted, 2020.
[101] S. Alexander, "Visualization and analysis of atomistic simulation data with OVITO–the
Open Visualization Tool," Modelling and Simulation in Materials Science and
Engineering, vol. 18, no. 1, p. 015012, 2010.
[102] C. Chen, Y. Zuo, W. Ye, X. Li, Z. Deng, and S. P. Ong, "A Critical Review of Machine
Learning of Energy Materials," Advanced Energy Materials, vol. 10, no. 8, p. 1903242,
2020, doi: https://doi.org/10.1002/aenm.201903242.
[103] R. Batra, L. Song, and R. Ramprasad, "Emerging materials intelligence ecosystems
propelled by machine learning," Nature Reviews Materials, 2020/11/09 2020, doi:
10.1038/s41578-020-00255-y.
[104] K. Alberi et al., "The 2019 materials by design roadmap," Journal of Physics D: Applied
Physics, vol. 52, no. 1, p. 013001, 2018/10/24 2018, doi: 10.1088/1361-6463/aad926.
89
[105] V. Fung, J. Zhang, E. Juarez, and B. G. Sumpter, "Benchmarking graph neural networks
for materials chemistry," npj Computational Materials, vol. 7, no. 1, p. 84, 2021/06/03
2021, doi: 10.1038/s41524-021-00554-0.
[106] M. Rupp, A. Tkatchenko, and K. R. a. Mller, "Fast and accurate modeling of molecular
atomization energies with machine learning," Physical Review Letters, vol. 108, no. 5,
pp. 1--5, 2012, doi: 10.1103/PhysRevLett.108.058301.
[107] J. Behler and M. Parrinello, "Generalized neural-network representation of high-
dimensional potential-energy surfaces," Physical Review Letters, vol. 98, no. 14, pp. 1--4,
2007, doi: 10.1103/PhysRevLett.98.146401.
[108] J. Behler, "First Principles Neural Network Potentials for Reactive Simulations of Large
Molecular and Condensed Systems," Angewandte Chemie - International Edition, vol.
56, no. 42, pp. 12828--12840, 2017, doi: 10.1002/anie.201703114.
[109] S. De, A. P. Bartók, G. Csányi, and M. Ceriotti, "Comparing molecules and solids across
structural and alchemical space," Physical Chemistry Chemical Physics,
10.1039/C6CP00415F vol. 18, no. 20, pp. 13754-13769, 2016, doi:
10.1039/C6CP00415F.
[110] Y. H. Tang, D. Zhang, and G. E. Karniadakis, "An atomistic fingerprint algorithm for
learning ab initio molecular force fields," Journal of Chemical Physics, vol. 148, no. 3,
pp. 1--30, 2018, doi: 10.1063/1.5008630.
[111] C. Hansch and T. Fujita, "$\rho$-$\sigma$-$\pi$ Analysis. A Method for the Correlation
of Biological Activity and Chemical Structure," Journal of the American Chemical
Society, vol. 86, no. 8, pp. 1616--1626, 1964, doi: 10.1021/ja01062a035.
90
[112] M. Rupp and R. a. Ramakrishnan, "Machine Learning for Quantum Mechanical
Properties of Atoms in Molecules," Journal of Physical Chemistry Letters, vol. 6, no. 16,
pp. 3309--3313, 2015, doi: 10.1021/acs.jpclett.5b01456.
[113] M. Rupp, A. Tkatchenko, K.-R. Müller, and O. A. von Lilienfeld, "Fast and Accurate
Modeling of Molecular Atomization Energies with Machine Learning," Physical Review
Letters, vol. 108, no. 5, p. 058301, 01/31/ 2012, doi: 10.1103/PhysRevLett.108.058301.
[114] S. Riniker and G. A. Landrum, "Open-source platform to benchmark fingerprints for
ligand-based virtual screening," Journal of Cheminformatics, vol. 5, no. 1, p. 26,
2013/05/30 2013, doi: 10.1186/1758-2946-5-26.
[115] B. Tang, S. T. Kramer, M. Fang, Y. Qiu, Z. Wu, and D. Xu, "A self-attention based
message passing neural network for predicting molecular lipophilicity and aqueous
solubility," Journal of Cheminformatics, vol. 12, no. 1, p. 15, 2020/02/21 2020, doi:
10.1186/s13321-020-0414-z.
[116] D. Weininger, "SMILES, a chemical language and information system. 1. Introduction to
methodology and encoding rules," Journal of Chemical Information and Computer
Sciences, vol. 28, no. 1, pp. 31-36, 1988/02/01 1988, doi: 10.1021/ci00057a005.
[117] M. Olivecrona, T. Blaschke, O. Engkvist, and H. Chen, "Molecular de-novo design
through deep reinforcement learning," Journal of Cheminformatics, vol. 9, no. 1, p. 48,
2017/09/04 2017, doi: 10.1186/s13321-017-0235-x.
[118] X. Li, X. Yan, Q. Gu, H. Zhou, D. Wu, and J. Xu, "DeepChemStable: Chemical Stability
Prediction with an Attention-Based Graph Convolution Network," Journal of Chemical
Information and Modeling, vol. 59, no. 3, pp. 1044-1049, 2019/03/25 2019, doi:
10.1021/acs.jcim.8b00672.
91
[119] J. Zhou et al., "Graph neural networks: A review of methods and applications," AI Open,
vol. 1, pp. 57-81, 2020/01/01/ 2020, doi: https://doi.org/10.1016/j.aiopen.2021.01.001.
[120] J. Gilmer, S. S. Schoenholz, P. F. Riley, O. Vinyals, and G. E. Dahl, "Neural Message
Passing for Quantum Chemistry," p. arXiv:1704.01212
[121] K. T. Schütt, P.-J. Kindermans, H. E. Sauceda, S. Chmiela, A. Tkatchenko, and K.-R.
Müller, "SchNet: A continuous-filter convolutional neural network for modeling quantum
interactions," p. arXiv:1706.08566
[122] D. Duvenaud et al., "Convolutional Networks on Graphs for Learning Molecular
Fingerprints," p. arXiv:1509.09292
[123] Z. Wu et al., "MoleculeNet: a benchmark for molecular machine learning," Chemical
Science, 10.1039/C7SC02664A vol. 9, no. 2, pp. 513-530, 2018, doi:
10.1039/C7SC02664A.
[124] A. Dunn, Q. Wang, A. Ganose, D. Dopp, and A. Jain, "Benchmarking materials property
prediction methods: the Matbench test set and Automatminer reference algorithm," npj
Computational Materials, vol. 6, no. 1, p. 138, 2020/09/15 2020, doi: 10.1038/s41524-
020-00406-3.
[125] S.-Y. Louis et al., "Graph convolutional neural networks with global attention for
improved materials property prediction," Physical Chemistry Chemical Physics,
10.1039/D0CP01474E vol. 22, no. 32, pp. 18141-18148, 2020, doi:
10.1039/D0CP01474E.
[126] C. W. Park and C. Wolverton, "Developing an improved crystal graph convolutional
neural network framework for accelerated materials discovery," Physical Review
92
Materials, vol. 4, no. 6, p. 063801, 06/01/ 2020, doi:
10.1103/PhysRevMaterials.4.063801.
[127] M. Karamad, R. Magar, Y. Shi, S. Siahrostami, I. D. Gates, and A. Barati Farimani,
"Orbital graph convolutional neural network for material property prediction," Physical
Review Materials, vol. 4, no. 9, p. 093801, 09/08/ 2020, doi:
10.1103/PhysRevMaterials.4.093801.
[128] T. N. Kipf and M. Welling, "Semi-Supervised Classification with Graph Convolutional
Networks," p. arXiv:1609.02907
[129] P. Veličković, G. Cucurull, A. Casanova, A. Romero, P. Liò, and Y. Bengio, "Graph
Attention Networks," p. arXiv:1710.10903
[130] J. A. Keith et al., "Combining Machine Learning and Computational Chemistry for
Predictive Insights Into Chemical Systems," p. arXiv:2102.06321

Abstract (if available)

Abstract Polymer based dielectric materials are used in a wide range of electrical storage devices. However, the morphological complexity and diversity in polymer space poses a challenge in design and identification of new polymer materials. Here, I would like to use classical molecular dynamics to introduce a high throughput scalable computational framework based on charge state aware reactive molecular dynamics to study dielectric response under external electric field and further aid in identification of new polymers capable of showing desirable dielectric properties. The overall work is divided into three parts as follows (1) Computational framework for polymer property prediction using polarizable reactive molecular dynamics, (2) Recurrent Neural Network based study of dielectric property generated using the computational framework, (3) Graph Attention based analysis of polymer dielectric data to study underlying global as well local features important for design of new polymers. ? The increased energy and power density required in modern electronics poses a challenge for designing new dielectric polymer materials with high energy density while maintaining low loss at high applied electric fields. Recently, an advanced computational screening method coupled with hierarchical modelling has accelerated the identification of promising high energy density materials. It is well known that the dielectric response of polymeric materials is largely influenced by their phases and local heterogeneous structures as well as operational temperature. Such inputs are crucial to accelerate the design and discovery of potential polymer candidates. However, an efficient computational framework to probe temperature dependence of the dielectric properties of polymers, while incorporating effects controlled by their morphology is still lacking. In the first work, I propose a scalable computational framework based on reactive molecular dynamics with a valence-state aware polarizable charge model, which is capable of handling practically relevant polymer morphologies and simultaneously provide near-quantum accuracy in estimating dielectric properties of various polymer systems. The predictive power of this framework is demonstrated on high energy density polymer systems recently identified through rational experimental-theoretical co-design. This scalable and automated framework may be used for high-throughput theoretical screenings of combinatorial large design space to identify next-generation high energy density polymer materials. ? Despite the growing success of machine learning for predicting structure-property relationships in molecules and materials, such as predicting the dielectric properties of polymers, is still in its infancy. In the second work, we report on the effectiveness of solving structure-property relationships for a computer-generated database of dielectric polymers using recurrent neural network (RNN) models. The implementation of a series of optimization strategies was crucial to achieving high learning speeds and sufficient accuracy: (1) binary and non-binary representations of SMILES (Simplified Molecular Input Line System) fingerprints; and (2) back propagation with affine transformation of the input sequence (ATransformedBP) and resilient backpropagation with initial weight update parameter optimizations (iRPROP-optimized). For the investigated database of polymers, the binary SMILES representation was found to be superior to the decimal representation in respect of the training and prediction performance. All developed and optimized Elman-type RNN algorithms outperformed non-optimized RNN models in the efficient prediction of non-linear structure-activity relationships. ? There is a great deal of interest in designing and identifying new dielectric polymer materials that exhibit high power density as they have many applications in modern electronics system. However, designing such polymers with desired properties is a challenging task due to combinatorial large search space. Polynorbornene (PNB) is one such important amorphous polymer system, which has potential applications as a high energy density polymer due to its high breakdown strength with low dielectric loss and high thermal stability. Moreover, electrical properties of PNB can be significantly enhanced by incorporation of defects or synthesis with controlled crystallinity by hydrogenation reaction, which involves experimental synthesis and characterization of combinatorial large number of polymer systems to identify potential candidates. In the third work, I propose an attention-based graph convolutional neural network (GAT) model that can identify polymer systems capable of exhibiting increased energy and power density. The GNN model is trained to predict dielectric constant for a polymer, where the training data for the high frequency dielectric constant of the PNB polymers are computed via ab-initio molecular dynamics simulation. This model can significantly aid experimental synthesis of potentially new dielectric polymer materials which is otherwise difficult using simplistic statistical procedures.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Simulation and machine learning at exascale

PDF

Neural network for molecular dynamics simulation and design of 2D materials

PDF

Large-scale molecular dynamics simulations of nano-structured materials

PDF

Physics informed neural networks and electrostrictive cavitation in water

PDF

Quantum molecular dynamics and machine learning for non-equilibrium processes in quantum materials

PDF

Nanomaterials under extreme environments: a study of structural and dynamic properties using reactive molecular dynamics simulations

PDF

Multi-scale quantum dynamics and machine learning for next generation energy applications

PDF

Controlling electronic properties of two-dimensional quantum materials: simulation at the nexus of the classical and quantum computing eras

PDF

Deep learning for subsurface characterization and forecasting

PDF

Exploring properties of silicon-carbide nanotubes and their composites with polymers

PDF

Molecular dynamics simulations of nanostructures

PDF

Physics-aware graph networks for spatiotemporal physical systems

PDF

Theoretical and computational foundations for cyber‐physical systems design

PDF

An FPGA-friendly, mixed-computation inference accelerator for deep neural networks

PDF

Molecular dynamics studies of protein aggregation in unbounded and confined media

PDF

Scaling up temporal graph learning: powerful models, efficient algorithms, and optimized systems

PDF

Dynamics of water in nanotubes: liquid below freezing point and ice-like near boiling point

PDF

Hardware-software codesign for accelerating graph neural networks on FPGA

PDF

Shock-induced poration, cholesterol flip-flop and small interfering RNA transfection in a phospholipid membrane: multimillion atom, microsecond molecular dynamics simulations

PDF

Modeling of bicontinuous nanoporous materials and investigation of their mechanical and hydraulic properties

Asset Metadata

Creator Mishra, Ankit (author)

Core Title High throughput computational framework for synthesis and accelerated discovery of dielectric polymer materials using polarizable reactive molecular dynamics and graph neural networks

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Chemical Engineering

Degree Conferral Date 2021-08

Publication Date 07/23/2021

Defense Date 06/14/2021

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag dielectric property,graph neural networks,machine learning,molecular dynamics,OAI-PMH Harvest,polarizable reactive molecular dynamics,polymer,polymer discovery,reactive molecular dynamics

Format application/pdf (imt)

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Vashishta, Priya (committee chair), Kalia, Rajiv (committee member), Nakano, Aiichiro (committee member)

Creator Email ankitmis@usc.edu,m.ankit011@gmail.com

Permanent Link (DOI) https://doi.org/10.25549/usctheses-oUC15619459

Unique identifier UC15619459

Legacy Identifier etd-MishraAnki-9846

Document Type Dissertation

Format application/pdf (imt)

Rights Mishra, Ankit

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the author, as the original true and official version of the work, but does not grant the reader permission to use the work if the desired use is covered by copyright. It is the author, as rights holder, who must provide use permission if such use is covered by copyright. The original signature page accompanying the original submission of the work to the USC Libraries is retained by the USC Libraries and a copy of it may be obtained by authorized requesters contacting the repository e-mail address given.

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA

Repository Email cisadmin@lib.usc.edu