Reactive Molecular Dynamics Simulations of ... - ACS Publications

Subscriber access provided by Bethel University

C: Surfaces, Interfaces, Porous Materials, and Catalysis

Reactive Molecular Dynamics Simulations of Atomic Oxygen Impact on Epoxies with Different Chemistries Chowdhury Ashraf, Aniruddh Vashisth, Charles E. Bakis, and Adri C.T. van Duin J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.9b03739 • Publication Date (Web): 30 May 2019 Downloaded from http://pubs.acs.org on May 30, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Reactive Molecular Dynamics Simulations of Atomic Oxygen Impact on Epoxies with Different Chemistries Chowdhury Ashraf a, Aniruddh Vashisth b, Charles E. Bakis b and Adri C.T. van Duin a,b* a

Department of Mechanical Engineering, The Pennsylvania State University, 240 Research East Building, University Park, PA, 16802, USA b

Department of Engineering Science and Mechanics, The Pennsylvania State University, 212 EarthEngineering Sciences Building, University Park, PA, 16802, USA

ABSTRACT Atomic oxygen (AO) is one of the most abundant species present in the lower earth orbit (LEO) and is responsible for the aggressive degradation of polymers used in spacecraft structures. In this investigation, we use ReaxFF reactive force field molecular dynamics simulations to evaluate the disintegration of several different thermosetting epoxy polymers subjected to hypervelocity AO impact. Our simulations indicate that epoxy with aromatic curative displays higher resistance to AO impact due to its stable benzene functionality. Decreased cross-linking density and increased simulation temperature both lead to faster disintegration of the polymer. Our simulation results indicate that ReaxFF force field simulations can be a useful tool to evaluate the response of various thermosetting epoxies to AO impact and identify promising candidate materials for spacecraft applications.

*Corresponding Author: Adri C. T. van Duin ([email protected])

1 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 33

1. INTRODUCTION Low earth orbit (LEO) satellites are commonly used for long distance communications as the orbit closer to the ground and they do not typically require high signal strength. The high cost of launching a unit of mass to LEO1 has motivated the satellite launch industry to use lightweight fiber-reinforced polymer composites in satellite structures. However, hypervelocity impact of atomic oxygen (AO), which is abundantly present in LEO, can severely degrade polymeric materials and impose lifetime limitations on LEO satellites. The usual speed of AO is approximately 7-8 km/s, resulting in an impact of equivalent energy of 4.5-5 eV at a high flux rate of 1014 – 1015 Oxygen atoms cm-2s-1.2 These atomic oxygens have enough energy to initiate surface modification affecting the chemical, electrical, thermal, optical and mechanical properties of the surface.3 Experimental examinations of various neat polymers and composite such as polyethylene,2 polytetrafluoroethylene,4 polyimide (with and without nano-ZrO2 reinforcements),5 carbon/epoxy composites,6,7 poly(p-phenylene benzobisoxazole) fibers8 have shown that AO erodes the exposed material significantly. Experimentally, parameters that affect the reactivity of AO with a composite are composite density, exposure area, and AO fluence.3 For carbon fiber reinforced polymer composites, the matrix erodes much faster than the graphitic fibers. This happens because polymer matrices have a lower mass density as compared to graphitic fibers.9 Studies have been carried out to improve the resistance of the matrix to AO impact.5,10–12 As the chemical structure of the epoxy can significantly change the thermo-mechanical properties of a polymer,13–16 it is intuitive that the AO impact response should be dependent on polymer structure. However, a comprehensive understanding correlating the chemical structure of polymer matrix to the disintegration due to AO impact is experimentally challenging. To this end, molecular dynamics (MD) simulations can be used as an alternative to performing experimental investigations of AO degradation mechanisms in polymers. AO impact and subsequent disintegration of polymers has been examined using non-reactive molecular dynamics17 as well as reactive molecular dynamics.18–21 Non-reactive molecular dynamics does not capture the reaction kinetics enabled by AO impact. On the other hand, reactive force fields such as ReaxFF 22,23 2 ACS Paragon Plus Environment

Page 3 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


provide insight into the resulting chemical reactions. Using ReaxFF, Zeng et al.18 found that poly(vinylidene fluoride) subjected to AO impact has a higher resistance to erosion when embedded with a polyhedral oligomeric silsesquioxane compound. Srinivasan et al.20 studied degradation mechanisms of graphene-oxide (GO) nanoparticles subjected to high energy AO impact and reported that the Eley-Rideal type reaction mechanism is responsible for O2 removal from the surface. Rahnamoun and van Duin19 showed that Teflon is more resistant than Kapton to AO impact. Rahmani et al. 21 investigated the damage mitigation efficacy of pristine and polyimide (PI)-grafted polyoctahedral silsesquioxane, graphene, and carbon nanotubes (CNTs) in a polyimide matrix exposed to high energy atomic oxygens. Despite the good insights these studies have provided, they are significantly limited in their scopes. One major limitation of these studies is that they were only conducted on slabs made with monomers, thus neglecting the effect of cross-linking. For example, one study 21 assumed that cross-linking would not affect the degradation results as the material was not subjected to an external load. However, crosslinked monomers are favored for demanding structural applications because they have superior mechanical properties and thermal stability, compared to monomers. When AO impacts a material, it is subjected to a highly repulsive force. Therefore, to capture the appropriate kinetics in reactive molecular dynamics, the force field needs to be trained for short-range interactions. The force fields used in the previously mentioned studies19,21 were not adequately trained for short-range repulsions, even though these forces are essential for AO impact simulations. In the absence of the repulsive force, AO reacts more readily with the surface, thereby overpredicting the damage. ReaxFF allows the force field to be trained for such short-range interactions. Therefore, by using a well-trained ReaxFF parameter set, one can achieve detailed insights into chemical degradation caused by AO impact. Such an approach can be a powerful tool for designing more durable composite materials for spacecraft. The accelerated ReaxFF method24 was used to virtually manufacture several thermosetting epoxy polymers with varying percentages of cross-linking. Previously, cross-linking of polymers has been simulated using a method based on cut-off distances of reactive sites,13,25–27 where a bond is created when two reacting



Page 4 of 33

atoms come within a predefined range. These methods ignore the correct reaction kinetics of the crosslinking mechanism and the thermodynamics, where steric hindrance and path of approach play significant roles. If the predicted cross-linked structure is incorrect, the predicted thermal and mechanical properties of the polymer will be incorrect as well.13 A recently developed accelerated method in ReaxFF – named Accelerated ReaxFF24 - resolves this problem by accurately reproducing the reaction barriers as calculated using quantum mechanics. This method also uses specific cutoff distance criteria. However, when the criteria are fulfilled, and a desired spatial orientation of the interacting atoms is achieved, additional energies are provided as restraint energy (discussed in detail in the computational details section and in ref24) to the reacting sites to help them to pass the energy barrier for cross-linking. Thus, unlike the cutoff distance method, a bond is not always created when the cutoff distance criterion is satisfied. Therefore, the Accelerated ReaxFF method of cross-linking polymers is much closer to the actual cross-linking scenario. In this investigation, we have examined the influence of epoxy polymer backbone chemistry on the mechanisms of degradation caused by high energy AO impact. While there are other causes of polymer degradation in LEO, such as X-rays, ultraviolet radiation, gamma radiation, and micrometeorite impact, we have only considered AO interactions in the current investigation. This study primarily focuses on the effects of AO impact on various combinations of epoxides and curing agents while addressing the significant limitations of the previous studies mentioned earlier. In this study, bisphenol A (bisA) epoxides are cured with three different amine curing agents—aromatic (diethyltoluenediamine, DETDA), cycloaliphatic (isophorone diamine, IPDA) and aliphatic (JEFFAMINE T-403 polyetheramine, abbreviated here as T403)—to investigate the effect of high energy AO impact on the cross-linked polymers with variable chemistries. Additionally, to investigate the effect of the epoxide molecule on the AO impact response of polymer slab, bisphenol F (bisF) epoxide is cross-linked with DETDA. Figure 1 shows the different epoxides and curing agents considered in this study. As mentioned previously, in all the cases, cross-linking is performed using Accelerated ReaxFF simulations and structures with high cross-linking percentages are considered for investigation. Additionally, to determine whether the cross-linking


Page 5 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


percentage affects mass-loss and other properties upon AO impact, three structures with varying percentage of cross-linking are considered for bisA epoxide and DETDA curing agent system. We also investigated the effect of the polymer’s temperature at the onset of AO impact because the temperature in LEO varies from -170˚C to 123˚C .28 The force field used in this study is similar to our previous studies,14,24 however, additional parameters are trained to reproduce the short-range repulsion energies between atoms. Therefore, this study addresses both the limitations of the previous studies and further provides insights into the effect of the chemical structure of the epoxy on the AO impact response by accurately simulating physical and chemical processes. This manuscript is organized as follows – Section 2 outlines the computational details, including a brief description of the ReaxFF method. Results are presented in Section 3, while conclusions are drawn in Section 4.

Figure 1: Molecular structures of the epoxides and the curing agents considered in this study (a) bisphenol A (bisA) epoxides are cured with (b) diethyltoluenediamine (DETDA), (c) isophorone diamine (IPDA) and (d) JEFFAMINE T-403 polyetheramine (T403) curing agents to virtually manufacture thermoset polymers. In addition (e) bisphenol F (bisF) is also considered as an epoxide. Cyan, red, blue and white spheres represent carbon, oxygen, nitrogen and hydrogen atoms respectively.



Page 6 of 33

2. COMPUTATIONAL DETAILS The ReaxFF Method All the MD simulations reported in this study were performed using ReaxFF reactive force field method. ReaxFF is a bond order29,30 based empirical potential, which is capable of bond formation and breaking during the simulation. ReaxFF uses the following expression to derive the forces on each atom:

𝐸system = 𝐸𝑏𝑜𝑛𝑑 + 𝐸𝑜𝑣𝑒𝑟 + 𝐸𝑢𝑛𝑑𝑒𝑟 + 𝐸𝑙𝑝 + 𝐸𝑣𝑎𝑙 + 𝐸𝑡𝑜𝑟 + 𝐸𝑣𝑑𝑊𝑎𝑎𝑙𝑠 + 𝐸𝑐𝑜𝑢𝑙𝑜𝑚𝑏 + 𝐸𝑅𝑒𝑠𝑡

(1)

where, Ebond = bond energy, Eover = over-coordination energy penalty, Eunder = under-coordination stability, Elp = lone pair energy, Eval = valence angle energy, Etor = torsion angle energy, EvdWaals = van der Waals energy, Ecoulomb = Coulomb energy, ERest = restrain energy All the connectivity dependent interactions such as angle and torsion are bond order dependent. These bond orders are calculated and updated at every iteration based on the interatomic distances. The bond order dependency ensures a smooth transition to zero energy contribution due to bond breaking. Unlike non-reactive force fields, the non-bonded interactions in the current simulation approach are calculated between every pair of atoms irrespective of their connectivity. A shielding term is used in the van der Waals and Coulomb energy terms to eliminate any excessive non-bonded interactions, while a seventh-degree taper function is employed31 to eliminate any discontinuity in these non-bonded interaction energies. Additionally, atom charges are determined using a geometry dependent charge calculation scheme known as the Electronegativity Equalization Method (EEM).32 A complete description of each energy term is given by Chenoweth et al.22 The force field parameters are embedded in these energy terms, which are trained against quantum mechanical calculations. Senftle et al.33 provide a general overview of the ReaxFF method and its applications.


Page 7 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Simulation Setup Some of the polymers investigated in this study have been examined in previous publications,14,24 which included the simulation of thermo-mechanical properties such as glass transition temperature (Tg), coefficient of thermal expansions, elastic modulus, and mass density. The procedure and the simulation protocol for Accelerated ReaxFF, including cutoff distances, and force parameters, can be found in our previous publications as well.14,24 Cross-linked structures from our previous studies have been modified by annealing them to generate dense 2D slab structures (details below). Table 1 summarizes the various systems considered in this study. We placed multiple polymer chains with various percentages of crosslinking (as listed in Table 1) for each system in a cubic box and generated the 2D slab-like structure. For the bisA/DETDA, bisF/DETDA, and bisA/IPDA with more than 50% crosslinking, each chain consisted of 16 molecules of epoxide and 8 molecules of amine curing agent crosslinked through Accelerated ReaxFF.23 The bisA/T403 chains had 9 epoxide molecules and 3 molecules of T403 amine curative molecules. The bisA/DETDA system with 0% crosslinking has 256 bisA epoxide and 128 DETDA molecules that are unreacted. Finally, the 50% crosslinked bisF/DETDA system had 128 small chains with 2 epoxides and 1 amine molecule in each chain. The ratios of epoxide and amine molecules were chosen to preserve a 100% stoichiometric ratio of glycidyl groups and amine hydrogens. From Table 1, it is clear that all the cases have a similar slab size (cross-sectional area). The slabs are bombarded with a specific number of Oxygen atoms leading to approximately similar fluence values. Table 1: Details of the polymer systems

System

Slab crosssectional area

Number of chains

% Crosslinking

Number of atoms

Fluence after the end of simulation (Oxygen atoms/cm2)

Simulated Average Density (g/cc)

Investigation 1: Study the effect of cross-linking chemistry bisA/DETDA

62.5  62.5 Å2

16

84

16512

5.12  1014

1.159

bisA/T403

61.3  61.3 Å2

16

72

15296

5.32  1014

1.155

bisA/IPDA

61.8  61.8 Å2

16

82

16768

5.23  1014

1.147 7

ACS Paragon Plus Environment


bisF/DETDA

59.8  59.8 Å2

16

82

14944

Page 8 of 33

5.56  1014

1.167

Investigation 2: Study the effect of cross-linking density bisA/DETDA

62.5  62.5 Å2

16

84

16512

5.12  1014

1.159

bisA/DETDA

62.1  62.1 Å2

N/A

0

16512

5.12  1014

1.115

bisA/DETDA

61.3  61.3 Å2

128

50

16512

5.2  1014

1.151

One of the challenges in the current AO impact study is to generate a 2D slab-like structure which is periodic in two directions, but non-periodic in the third direction along which the impact will take place. To achieve this geometry, we initially placed the polymer chains in a large simulation box and slowly reduced the dimensions of the box from two directions while keeping the other dimension constant. Due to the van der Waals interactions, the polymer chains attract each other, thereby increasing the density of the polymer slab. This builds our initial 2D slab-like structure which is non-periodic in the z direction. However, the overall density of the system was too low because the z dimension was too high. To get an adequately packed 2D slab structure, we filled up the vacuum space in the z direction with Ar atoms and performed annealing simulations. Argon atoms were used in these simulations because they are inert and do not interact with the polymer (stay phase separated). These Ar atoms also prevent any cross-over of the polymer slab in the z direction, preserving the non-periodicity of the slab in the z direction. Because the force field is capable of reproducing the correct densities of the systems under consideration,14 we ended up with dense systems that are periodic in the x and y directions and non-periodic in the z direction. The thickness of each slab was between 40 Å and 45 Å, which was sufficiently high to sustain 200 AO impacts. Here, “sufficient” means the impacted AOs were not able to penetrate through the entire thickness of the slab, which will be discussed later in more detail. Figure 2 shows the resulting 2D structures which were used to study AO impacts. The details of these structures are listed in Table 1. The Ar atoms are removed from the system after generating the polymer slab after annealing.


Page 9 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 2: Representative initial structures of different annealed systems: Top and side views of (a) bisA/DETDA, (b) bisA/IPDA, (c) bisA/T403 and (d) bisF/DETDA systems generated from annealing. Cyan, red, blue and white spheres represent carbon, oxygen, nitrogen and hydrogen atoms respectively.

Once the systems was prepared, an NVT simulation was performed at 300 K for 30 ps to stabilize the system temperature. Next, to simulate AO bombardment along the z direction, 200 atomic oxygens were placed simultaneously in z direction in a way that two consecutive AOs have a z-spacing of ~15Å so that with the given velocity each impact takes place at ~200 fs time interval. In order to include some uncertainty of impact time interval, a standard deviation of 1 Å is introduced. Because all the AOs are added simultaneously, the box size increased by ~3000 Å (15 Å  200 AOs) in the z direction. The x and y coordinates of these Oxygen atoms were generated randomly using a Gaussian distribution with a mean of 𝑥̅ and 𝑦̅ (where 𝑥̅ = length of x-dimension /2 and 𝑦̅= length of y-dimension/2). Three simulations were performed for each case, each having a different width of the distribution (close to 1/6th of the x and y dimensions of the box). All the AOs were launched at once from their initial positions with an initial velocity vector of (0, 0, -7.4 km/s). Next, with an initial slab temperature of 300 K, an NVE simulation is



Page 10 of 33

performed for 40 ps with a time step of 0.00005 ps. Since the AOs are initially separated by ~15 Å in the z direction, each AO impact was separated by approximately 0.2 fs, resulting in the desired flux rate. A similar simulation protocol was followed in previous studies19,21 too. Trajectory files are saved in every 0.05 ps to perform post-processing of the simulations. All the simulations are performed in LAMMPS software package.34 To improve statistical sampling, three separate simulations were run for each case, and relevant data sets were averaged. All the three simulations of a specific case have the same slab surface; however, the positions of the Oxygen atoms were varied by generating random x-y coordinates, as mentioned earlier. Our preliminary simulations suggested that correct boundary conditions are essential for simulating realistic slab behavior during AO impact. If the non-periodic direction is free floating, the momentum of the AOs imparts rigid body velocity to the slab in the direction of impact. Because part of the AO impact energy is transferred to bulk velocity of the slab, the damage observed due to the impact will be much less compared to the actual cases. On the other hand, fixing the atoms at the bottom of the slab to prevent this movement will make the slab rigid, which is not practical. Therefore, we used the reflective wall feature in LAMMPS at the bottom x-y plane of the simulation box and placed the slab closed to that plane. If an atom moves a distance  outside the wall on a given time step, then the reflective feature puts the atom back inside the wall with the same distance . At the same time, the sign of the corresponding velocity component is reversed, thereby preventing atoms from crossing the boundary. Another approach for capturing the correct boundary condition is to use the “fix wall” command in conjunction with different potentials, whereby the wall will reflect the atoms by generating a force perpendicular to the wall direction. Since both approaches serve the same purpose, we decided to use the “fix wall/reflect” command in LAMMPS.

Force Field Training In this study, we used a modified version of the force field developed by Zhang et al.35 for functionalized hydrocarbon/water weak interactions. The this modified force field was used in our previous studies14,24 which indicated that the force field is suitable for such systems. However, during bombardment simulations, 10 ACS Paragon Plus Environment

Page 11 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


the AOs come in close vicinity with other atoms on the polymer surface resulting in very high repulsive inter-atomic forces. The force field we used previously24 is not trained for such applications and therefore requires further training for close range interactions. However, quantum mechanical based DFT calculations cannot handle such close-range interactions properly, therefore they are not suitable for force field

training.

One

widely

used

potential

for

high

energy

collision

modeling

is

the

Ziegler−Biersack−Littmark (ZBL) potential,36 which only depends on the nuclear charges/atomic numbers of the two interacting atoms. The potential energy due to a pair of atoms at a distance rij is given by Equations 2-4,

𝐸𝑖𝑗𝑍𝐵𝐿 =

1 𝑍𝑖 𝑍𝑗 𝑒 2 𝜙(𝑟𝑖𝑗 /𝑎) 4𝜋𝜀𝑜 𝑟𝑖𝑗

(2)

0.4685 + 𝑍𝑗0.23

(3)

𝑎=

𝑍𝑖0.23

𝜙(𝑥) = 0.18175𝑒 −3.1998𝑥 + 0.50986𝑒 −0.94229𝑥 + 0.28022 𝑒 −0.4029𝑥 (4)

+ 0.02817𝑒 −0.20162𝑥

where e is the electron charge, o is the electrical permittivity of vacuum, and Zi, and Zj are the nuclear charges of the two atoms under consideration. Once the short-range repulsive energies of various combinations of atoms are calculated, the force field is trained to reproduce the numbers. The ReaxFF formulation includes an additional energy term along with the energy components shown in Equation (1) to describe the potential energy contribution of the shortrange interactions, namely, Ecore, which is given by Equation (5).37

𝐸𝑐𝑜𝑟𝑒 = 𝑒𝑐𝑜𝑟𝑒 exp[𝑎𝑐𝑜𝑟𝑒 (1 −

𝑟𝑖𝑗 )] 𝑟𝑐𝑜𝑟𝑒

(5)



Page 12 of 33

where rij is the distance between the atoms and ecore, acore, and rcore are force field parameters. The force field is trained in a manner that the energy contribution of Ecore is only substantial when the inter-atomic distance between atom pairs is much less than equilibrium bond distance— i.e. when the bond is significantly compressed. Therefore, the training does not affect the overall quality of the force field as the short-range interactions only occur during the high-impact simulations. The training results and the resultant force field are given in the Supporting Information. It should be noted that that, ecore, acore, and rcore are force field parameters for each atom type (30th, 31st, and 32nd parameters in the atom section). They are typically set to zero when the inner-wall correction is not used. Non-zero values of these parameters activate the inner-wall correction function in LAMMPS which we typically used for this study.

3. RESULTS AND DISCUSSION Analysis of material degradation upon bombardment After training the force field for close range interactions, we performed AO bombardment simulations for each case following the procedure described above using the ReaxFF reactive force field. Figure 3 shows snapshots of final systems for the four different slabs with a high percentage of cross-linking. In addition to that, we have provided a short movie of one of our representative simulations as Supporting Information. Since the AOs move with very high energy, in a MD simulation, they can transfer their kinetic energy onto the slab’s surface, which in turn increases the slab temperature. This increase in temperature induces strain among the local bonds and enables disintegration of the polymer through chemical reactions. The disintegration process can be both physical and chemical. The advantage of using ReaxFF is that it can simulate the bond formation and bond breaking on the fly, providing more insight into the reactions that lead to the chemical disintegration of the polymer. Apart from chemical processes, physical disintegration processes leading to the disintegration of the polymer are also simulated by ReaxFF reactive force field simulations, in contrast to non-reactive force field methods which can only simulate the physical disintegration. The chemical disintegration process creates smaller molecules which then move into the vacuum. Therefore, to quantify the system’s mass and temperature evolution we defined an imaginary box 12 ACS Paragon Plus Environment

Page 13 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


around the polymer which includes roughly 10 Å above the material’s surface in the z direction. Any atom that leaves this virtual box (the red box in Figure 3(a)) is considered lost and is not included in the mass and temperature calculations for the system.

Figure 3: Representative final structures of different systems after bombardment: Side views of (a) bisA/DETDA, (b) bisA/IPDA, (c) bisA/T403 and (d) bisF/DETDA systems generated from annealing. Cyan, red, blue and white spheres represent carbon, oxygen, nitrogen and hydrogen atoms respectively, while the yellow spheres represent the incident Oxygen atoms. Analysis of the AO impact simulations were performed by defining a box similar to the red one in (a), where the top edge of the box is roughly 10Å above the atoms of the slab surface.

The system mass and temperature were analyzed as functions of time, which is proportional to the number of AOs that have impacted the system (Figure 4). Since different polymer slabs had different initial masses, we divided the mass of the system at any given time by the initial mass of the system to obtain the normalized mass. Therefore, all the systems have a mass of unity at the beginning of the impact simulations. As seen in Figure 4(a), for the four different slabs considered no significant mass loss is observed up to the first 10 ps of simulation. However, after 10 ps, the mass of all four systems starts to decrease, with bisA/IPDA and bisA/T403 experiencing much faster mass loss than the other two systems. At the end of the simulation, bisA/T403 exhibits the largest mass loss while bisA/IPDA follows it closely. On the other hand, bisA/DETDA and bisF/DETDA are subjected to similar mass losses which are almost ~33% less than that of the other two systems. This is expected, because both bisA/DETDA and bisF/DETDA contain the



Page 14 of 33

aromatic ring structure in their curing agent. Aromatic rings are more stable than the cycloaliphatic IPDA rings and the aliphatic T403 structure because C

C bonds are more stable as compared to C-C bonds in

a cyclic molecule.38 The authors expect that the mass loss trends would extrapolate accordingly for longer AO exposure times for each resin systems.

(a)

(b)

Figure 4: (a) Average normalized mass, (b) average system temperature as a function of simulation time for different polymer slabs. The shaded region in the plots show the standard deviation of each curve.

The stability of aromatic curatives as compared to aliphatic and cyclic curatives can be explained in terms bond dissociation energies (BDEs) required for breaking various bonds present in the polymer backbone. The BDE for C – C in aliphatic polymers is 370 – 380 kJ/mol,38 whereas the BDE of CH3 attached to benzene (in a methylbenzene structure) in the DETDA molecule is 426.8 ± 5.4 kJ/mol 38 making the bond much stronger than C – C bonds in the aliphatic chain. The BDE for CH3 molecule attached to cyclohexane (methylcyclohexane) in the IPDA structure is 377 ± 7.5 kJ/mol38 which is comparable to C – C bonding in aliphatic chains. These trends can also be observed using ReaxFF where BDE for C – C in an aliphatic chain is 310 kJ/mol, the BDE of CH3 attached to a benzene ring is 358 kJ/mol, and the BDE of CH3 bonded to cyclohexane is 309.5 kJ/mol. The C – C bonds in the cyclohexane have a BDE of 313 kJ/mol38 whereas the C

C in benzene have 518 kJ/mol39 BDEs, making benzene more stable to external stimuli that could


Page 15 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


break the structure. Using the current ReaxFF force field, the BDEs were computed to be 428 kJ/mol and 562 kJ/mol for cyclohexane and benzene rings, respectively. These BDE results show close agreement and similar trends as described in literature. It should be noted that the BDE for the ReaxFF simulations were calculated by finding the difference between the minimized structures before and after bond dissociation. Figure 4(b) shows that the average temperature of the polymer slab increases nearly linearly for each system. The final temperatures of the systems range between 600 K and 650 K. From the standard deviations associated with each profile in Figure 4(b), it is observed that the temperature differences between the cases are not significant. Since the system temperature is an average ensemble value, it is not an ideal metric for discerning differences among the materials that may exist only near the impacted surface. Therefore, we examined temperatures at four different z cross-sections near the impacted surface as shown in Figure 5 for bisA/DETDA. Temperatures of these four cross-sections along z-direction are shown at two different times – 15ps (Figure 5(a)) and 35ps (Figure 5(b)). To generate these images, we first time averaged the temperatures for 200 frames (100 before and 100 after the reported times) of all the atoms within ±2.5 Å of the specified z-coordinates. Then, we divided x and y dimensions in 100 bins each and calculated the spatial average temperature for each bin (total 100 × 100 bins) using the time-averaged results. As shown in figure 5(a), at t = 15ps, some small regions on the slab surface have temperatures as high as roughly 8000 K and relatively more regions near 3000 K. However, at t = 15 ps most of the slab is at its initial temperature as seen for depths of 20 Å and 30 Å. Therefore, the simulation indicates that 15 ps of AO impact mainly affect local areas of the slab’s surface region. On the other hand, at t = 35 ps, the maximum temperature on the impacted surface reaches around 10000 K, while almost all the surface area attains a temperature close to 5000 K. Additionally, the cross-sectional area located 10 Å below the surface exhibits a significant temperature rise too. Some regions in the area have a temperature close to 6000 K at this time. Within the time frame of the current simulations, the top of the slab reaches significantly high temperature while the atoms deeper in the polymer slab are at a much lower temperature as slab does not equilibrates. This phenomenon is evident as the cross-sectional area at z = -20 Å in the polymer slab remains at nearly the



Page 16 of 33

initial temperature. The high surface temperatures, as observed in these simulations contribute to the disintegration of the polymer into smaller stable or unstable species that are discussed further in later sections.

Figure 5: Surface plot of temperature for bisA/DETDA system at different x-y cross sections of the slab at (a) 15 ps and (b) 35 ps. Here, z = 0 Å indicates the slab surface while increasing negative z values indicate the deeper cross sections of the slab, and x and y both are measured in Å.

In order to better understand the material degradation characteristics, we plotted the AO penetration depth for all the cases (Figure 6). The AO penetration depth defines how far the AOs were able to infiltrate into the slab, and in Figure 6 we demonstrate the number density profiles of the AOs at the end of the simulation. The number density profiles of AO for all the systems follow Gaussian distributions. The peaks for bisA/DETDA and bisF/DETDA (~ -15 Å) appear farther from the material surface with higher peak values compared to the other two systems (~ -10 Å). The peak positions for bisA/T403 and bisA/IPDA are much closer to the surface, which gives an initial impression that these materials are more resistant to degradation. However, the width of the AO distribution for these two cases are much wider compared to the other cases;


Page 17 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


while the higher peak values for bisA/DETDA and bisF/DETDA makes their corresponding AO distribution curve narrower. For example, the number densities of AOs at z = -30 Å for bisA/T403 and bisA/IPDA are almost double than the bisF/DETDA system. This means that the impacted AOs can penetrate further into the slab for these two systems. Since the AO penetration depth is defined on how deeper AOs can penetrate the surface, but not by the position of peak relative to the surface - we can say that these materials are less resistant to degradation. This indicates that AOs penetrate deeper into the systems with cyclo-aliphatic and aliphatic curing agents as compared to systems with aromatic curing agents. These observations are also in line with the earlier findings that polymers with higher material density are more resistant to AO impact,3 because bisA/DETDA has a higher density than bisA/IPDA and bisA/T403 polymer.14,24

Figure 6: Number density of AO for the four systems under consideration after the end of the simulation (t = 40 ps). Here, z = 0 Å indicates the slab’s surface. Standard deviations of the calculated values are shown as error bars.



Page 18 of 33

Table 2: Average damage propagation depth (DPD) with standard deviations for various systems after the end of the simulation (t = 40 ps) bisA/DETDA

bisA/IPDA

bisA/T403

bisF/DETDA

DPD (Å)

17.03

18.57

22.18

16.27

Std. Dev. (Å)

0.65

1.59

0.69

0.56

Similar to Rahmani et al.,21 we also calculated the average damage penetration depth (DPD) for all the four systems. Figure 7 shows a representative DPD for bisA/DETDA system after completion of the AO impact simulation, and a similar procedure was followed for all other systems. The damage penetration depth (DPD) was determined by comparing the initial and final mass density profiles of the system along the z axis (which is the direction of bombardment). DPD is defined as the distance in the z direction that has significant mass density drop due to AO impact. Table 2 lists the average DPD and their standard deviations for various systems under consideration. According to the calculation, DPD is lowest for bisF/DETDA and bisA/DETDA, however, bisA/T403 has the highest DPD value, which is consistent with our previous observations. We observed from the DPD calculations that the polymers with lower densities had higher DPDs and therefore suffered more damage.


Page 19 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 7: Example of damage penetration depth (DPD) calculation considering bisA/DETDA as the model system. The mass density drops up to certain point within the material which is defined as DPD. Here, z = 0 Å indicates the material surface.

The bisA/DETDA system has a higher Tg in comparison to bisF/DETDA and therefore may be more applicable for satellite applications. The molecular structures of bisA and bisF are reasonably similar, except that bisA has additional two methyl groups attached to the carbon which connects the aromatic rings. This slight difference in the structure does not affect the material degradation characteristics, according to the simulation results. The curing agent DETDA with aromatic rings performs much better than the curing agents with cyclo-aliphatic rings (IPDA) and linear-aliphatic chains (T403). To investigate this further, we compare the degradation of the number of six-member rings in these three structures (Figure 8). For comparison, we plotted the normalized number of six-member rings for each structure. Six-member rings lost due to bombardment is higher for bisA/IPDA slab compared to bisA/DETDA slab as observed in Figure 8. It is important to mention here that the initial number of six-member rings were the same for both systems and the only difference was the chemical nature of the ring, because the bisA/IPDA system has a mixture of aromatic and alicyclic rings whereas bisA/DETDA has only aromatic rings. As mentioned previously, the six-member rings in DETDA are aromatic, while they are cyclo-aliphatic in IPDA. Since the aromatic



Page 20 of 33

rings are thermodynamically more stable than the cyclo-aliphatic rings, they can resist the AO impact and undergo less degradation. Since the bisA/T403 system does not have any six-member ring in the curing agent structure, the initial number of six-member rings are much smaller compared to the other two systems. Therefore, we cannot directly compare the normalized number of six-member rings. However, because bisA/T403 has the least number of six-member rings to start with, it is subjected to the highest material degradation. Another factor that has not been detailed in this investigation is the probability of AO impacting the six-member rings. From analyzing the six-member rings, ReaxFF simulations can correctly predict the material degradation properties of various cross-linked polymer system and explain the observed results based on the chemical structure of the molecules. This is the primary reason for choosing ReaxFF for MD simulations as the results obtained from the simulations indicate that it can be used as a useful tool for investigating various candidate materials subjected to high-velocity AO impact.

Figure 8: Normalized number of six-member rings for three systems with same epoxy resin but different curing agents


Page 21 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Small molecule production The study of small molecules produced (H2O, CO, OH and CH2OH) is essential as it provides a deeper understanding the degradation of the polymer slab due to AO flux. At such high temperatures, any other gases present around the polymer can potentially react with the simulated polymer slab. Such cases were not considered in this investigation for the sake of simplicity. From our simulation trajectories, it is possible to identify and plot the evolution of small molecule production due to bombardment. Such a simulation capability is essential for two reasons. Firstly, it is experimentally challenging to identify these small molecules and, secondly, the production of these molecules is directly related to material degradation. Figure 9 shows the average number of four important intermediate and final species as a function of simulation time. These four species are H2O, CO, OH, and CH2O. These species were chosen for analysis because they appear as products in all the polymer systems investigated.21 In Figure 9(a), no apparent trend is observed for H2O; however, in Figure 9(c), OH production is higher for epoxides cured with cyclic amines. On the other hand, CO production (Figure 9(b)) in bisA/T403 is the highest, which can be explained based on the material structure. In all four cases (different combinations of epoxide and amine molecules), the molecular structures of the epoxy resins (bisA and bisF) are similar, whereas the molecular structure of the curing agents are different. T403 is the only curing agent which has Oxygen atoms in it; this increases the number of Oxygen atoms dramatically in the bisA/T403 system. Also, the BDE of O – C in an aliphatic chain is 320 – 350 kJ/mol,38 which is comparable to the C – C BDEs in the same curatives. This makes both C – C and O – C, equally viable candidates for bond breaking and thereby forming smaller species. From Figures 9(a) and (b), it is evident that the water molecules are produced earlier than CO. Therefore, the systems without Oxygen atoms in the curing agent molecule suffer a scarcity of Oxygen atoms to generate CO. Two competing factors are at play for formation of H2O molecules: (a) abundance of dangling H atoms with an increased probability of reacting with the AO and (b) the high BDE of H from parent chain. The BDEs of H bonded to an aliphatic, cycloaliphatic and an aromatic chain are more than 410 kJ/mol 38 depending on the polymer backbone. For our simulations, we



Page 22 of 33

observe that even though the C – O and C – C BDEs are lower than H – R BDEs (where R is the polymer backbone), the probability of a high energy AO finding an H is much higher. This causes a low CO production for all the systems except for bisA/T403, which has a more open structure in comparison to cyclic structures. This trend is also visible in the high bisA/T403 production of CH2O. However, one interesting trend is observed for bisA/DETDA system. Although it produces the smallest number of CO molecules, the CH2O production is quite comparable with bisA/T403. This indicates that bisA/DETDA favors CH2O molecules and, because these molecules are highly volatile, they leave the slab without converting to CO.


Page 23 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 9: Small molecule production for different systems upon bombardment (a) H2O, (b) CO, (c) OH and (d) CH2O

Although bisA and bisF have a small difference in their molecular structure, they have a clear difference in CH2O production (Figure 9). ReaxFF not only gives the reaction products but also illustrates the reaction pathways. Using this capability, we investigated the pathway of CH2O production in bisA/DETDA. Figure 10 shows snapshots of a reaction pathway of CH2O formation in bisA/DETDA system observed during the simulation. In the first step (Figure 10(a)), an O atom shown as a red sphere approaches towards a CH3 radical of bisA which is attached to the C atom located between two aromatic rings. Next (Figure 10(b)), the O atom reacts with the CH3 group and, after an atomic rearrangement, a -CH3O group is created. This CH3O group loses a hydrogen atom (shown in the red circle in Figure 10(c)). Lastly, the C-C bond breaks and CH2O is released (Figure 10(d)). Because the BDE for C – C for CH3-C6H5 is 426.8 ± 4.2 kJ/mol, while 23 ACS Paragon Plus Environment


Page 24 of 33

C – H and C = O have BDEs of ~440 kJ/mol and ~530 kJ/mol, respectively, C – C bond failure and formation of CH2O is the preferred degradation route. The only difference between bisA and bisF is that the carbon connecting the aromatic rings has two CH3 groups in bisA and our observation shows that these CH3 groups favor CH2O production in the bisA/DETDA system. This structural difference explains the observation made in Figure 9(d). Therefore, the molecular structure of the epoxy resin is directly related to product distribution.

Figure 10: CH2O formation mechanism in bisA/DETDA system (a) an O atom approaches a CH3 group of bisA which is attached to the C connected with two aromatic rings, (b) the O-atom reacts with the CH3 group and rearranges to -CH3O group, (c) the -CH3O group loses a H atom and (d) C-C bond breaks which releases CH2O

Effect of percentage cross-linking However, experimental investigations3 have shown that material disintegration is a function of density, which is further a function of percentage crosslinking. Therefore, we next evaluated the effect of percentage crosslinking in polymer slabs subjected to AO impact. We performed AO impact simulations on bisA/DETDA system with three smaller percentages of cross-linking (see Table 1 for details)—namely, 0%, 50% and 84% cross-linking to evaluate the effect of cross-linking on the slab disintegration. Figure 11 shows the results obtained from the simulation for the normalized mass and AO penetration depth.


Page 25 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 11: Effect of AO impact of bisA/DETDA system with different percent cross-linking: (a) normalized mass (b) AO penetration depth after the end of the simulation (t = 40 ps). Here, z = 0 Å indicates the slab’s surface. Standard deviations are represented by shaded areas in (a) and by error bars in (b).

In Figure 11(a), increased percent cross-linking clearly reduces mass loss. In Figure11(b), the final distribution of AO density for 0% cross-linking is much wider, indicating that AO penetrates further into the slab made with no cross-linking. Also, the peak for the 84% cross-linking case is much higher with a much narrower distribution compared to the 0% cross-linking case. These results prove that the cross-linked structure of the polymer imparts more resistance to AO penetration. Additionally, the damage penetration depth (DPD) is calculated for these different percentage cross-linked systems (Table 3). DPD for the system with 0% cross-linking is much higher than the cross-linked systems, which also supports our findings in Figure 11.

Table 3: Average damage propagation depth (DPD) with standard deviations for bisA/DETDA system with various percentage of cross-linking after the end of the simulation (t = 40 ps)

DPD (Å)

No Cross-linking

50% Cross-linking

84% Cross-linking

25.42

18.53

17.03



Std. Dev. (Å)

0.97

0.73

Page 26 of 33

0.65

The bisA/DETDA system with 0% cross-linking has a lower density in comparison to bisA/DETDA with a higher percentage of cross-linking.15 This makes the monomer system highly susceptible to disintegration. To test this hypothesis, we plotted the small molecule productions for different cross-link density system in Figure 12. This figure demonstrates that the production of small molecules containing carbon atoms (CO and CH2O, Figure 12(b) and (d)) are much higher for the monomer system compared to the cross-linked systems. Since the carbon atoms are the major component of the polymer backbone, generation of a higher number of small molecules containing carbon indicates more damage in the structure. On the other hand, the cross-linked structures produce more OH and H2O (Figure 12(a) and (c)), thus protecting the backbone structure of the polymer slab. Therefore, a higher percentage of cross-linking effectively resists material disintegration and our simulations can demonstrate the phenomenon. This is a significant improvement upon previous studies, and according to our best knowledge, this is the first time the effect of cross-linking on AO impact from atomistic-level has been investigated.


Page 27 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 12: Small molecule production upon bombardment for different bisA/DETDA systems with various percentages of cross-linking (a) H2O, (b) CO, (c) OH and (d) CH2O

Effect of slab temperature The temperature of the LEO vehicles varies between -170oC to 123oC (103K – 396K).28 Therefore, in addition to our regular simulations at 300 K, we performed two additional simulations at 103 K and 396 K, representing the two extreme temperature conditions, in order to investigate the effect of slab temperature when high energy AO impacts take place. We performed this investigation only on the bisA/DETDA system. The simulation setup remained the same as described in Table 1 (with the same % of cross-linking)



Page 28 of 33

except for the initial temperature of the slab. Similar to Figure 11, we plotted normalized mass as a function of time and AO density profile in Figure 13 for three different initial slab temperatures. It is clear from Figure 13(a) that when the slab temperature is lowest (103 K), the normalized mass is highest at the end of the simulation, i.e., the mass loss is lowest. This is expected because, at that low temperature (103K), the polymer molecules have less kinetic energy which makes the bonds, angles, and dihedrals much stiffer as compared to 300 K. Therefore, a significant amount of impact energy is needed to be transferred into the slab at a lower temperature to initiate the material degradation. This is reflected in Figure13(a) where a mass loss for the lowest temperature system is lowest. This can potentially be due to lower rate of reactions at lower temperatures. However, no significant difference can be observed for 300 K and 396 K cases. This is because the highest temperature case (396 K) is still below the glass transition temperature of the system as shown in our previous study14. A similar trend is also observed for the final AO density profile of the system upon bombardment where the AO distribution at 103 K has the sharpest peak with the lowest spread. That means more AOs get trapped close to the material surface and cannot penetrate further at 103 K. All these observations made in Figure 13 are consistent with the physical and chemical intuition.


Page 29 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 13: Effect of slab temperature upon bombardment on bisA/DETDA system (a) normalized mass (b) AO penetration depth after the end of the simulation (t = 40 ps), here z = 0 Å indicates the slab’s surface. Standard deviations are represented by shaded areas in (a) and by error bars in (b).

4. CONCLUSIONS In this study, we used ReaxFF to investigate the effect of cross-linking chemistry on material degradation due to Atomic oxygen (AO) impact. Three different curing agents – aromatic, cyclo-aliphatic and aliphatic curing agents were combined with two different epoxide monomers to generate polymer slabs with high cross-linking density. Our simulation results indicate that the predicted reduced mass loss for the aromatic curing agent—DETDA— is roughly 30% smaller (according to Figure 7(a)) compared to the other two curing agents simulated under similar condition, primarily due to the stable aromatic ring present in the structure. Additionally, we analyzed the small molecules produced during the simulations and demonstrated that ReaxFF can be used to investigate the possible formation pathways for those molecules. Therefore, through a ReaxFF simulation, one can identify the functional group which is responsible for chemistry that results in accelerated material degradation. Our results demonstrate the importance of higher cross-linking density of polymers benefits visibly to resists AO impact and reduce the material degradation. At the same time, we also investigated the effect of slab temperature and found that at higher simulation temperature the



Page 30 of 33

material degrades more due to the already available kinetic energy, which is according to expectations. Additionally, to validate our observation, we also calculated the BDEs using ReaxFF, which compare reasonably well with the literature. Similar studies can be performed for nano-particle enhanced polymer matrices to investigate the joint effect of nano-particles with high cross-linking density, which will lead to more robust material optimization procedure using ReaxFF MD simulations. Due to the recent advances in material sciences, more candidates for the aerospace application will emerge, so it is essential to carefully screen the properties of the materials and composites to ensure safe operation.

5. Supporting Information The polymerization process in ReaxFF, the force field training results for short-range interaction and the resulting force field are available as supporting material.

6. Acknowledgment CA and ACTvD acknowledge funding from DoE grant DE-EE008195. We would also like to thank Penn State Institute of Cyber Science (ICS) and ACI-cluster for computational resources where all the simulations were performed.

REFERENCES (1) (2)

(3) (4)

(5)

London III, J. R. LEO on the Cheap. Methods for Achieving Drastic Reductions in Space Launch Costs.; Air Univ Maxwell AFB AL Airpower Research Inst, 1994. Intrater, R.; Lempert, G.; Gouzman, I.; Grossman, E.; Cohen, Y.; Rein, D. M.; Khalfin, R. L.; Hoffman, A. Simulated Low Earth Orbit Environment Interaction with Different Types of Polyethylene. High Perform. Polym. 2004, 16, 249–266. Reddy, M. R. Effect of Low Earth Orbit Atomic Oxygen on Spacecraft Materials. J. Mat. Sci. 1995, 30, 281–307. Zhao, X.-H.; Shen, Z.-G.; Xing, Y.-S.; Ma, S.-L. An Experimental Study of Low Earth Orbit Atomic Oxygen and Ultraviolet Radiation Effects on a Spacecraft Material – Polytetrafluoroethylene. Polym. Degrad. Stabil. 2005, 88, 275–285. Lv, M.; Wang, Q.; Wang, T.; Liang, Y. Effects of Atomic Oxygen Exposure on the Tribological Performance of ZrO2-Reinforced Polyimide Nanocomposites for Low Earth Orbit Space Applications. Compos. Pt B Eng. 2015, 77, 215–222. 30 ACS Paragon Plus Environment

Page 31 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(6) (7)

(8)

(9) (10)

(11)

(12)

(13)

(14)

(15) (16)

(17) (18)

(19)

(20)

(21)

(22) (23)

Han, J.-H.; Kim, C.-G. Low Earth Orbit Space Environment Simulation and Its Effects on Graphite/Epoxy Composites. Compos. Struct. 2006, 72, 218–226. Zheng, N.; He, J.; Zhao, D.; Huang, Y.; Gao, J.; Mai, Y.-W. Improvement of Atomic Oxygen Erosion Resistance of Carbon Fiber and Carbon Fiber/Epoxy Composite Interface with a Silane Coupling Agent. Mater. Design. 2016, 109, 171–178. Chen, L.; Wang, C.; Wu, Z.; Wu, G.; Huang, Y. Atomic Oxygen Erosion Behaviors of PBO Fibers and Their Composite: Microstructure, Surface Chemistry and Physical Properties. Polym. Degrad. Stabil. 2016, 133, 275–282. Young, P. R.; Slemp, W. S. LDEF Materials Workshop 91. NASA Con& Pub 1992, 3162, 747– 748. Jiao, L.; Gu, Y.; Wang, S.; Yang, Z.; Wang, H.; Li, Q.; Zhang, Z. Atomic Oxygen Exposure Behaviors of CVD-Grown Carbon Nanotube Film and Its Polymer Composite Film. Compos. Pt. A Appl. Sci. Manuf. 2015, 71, 116–125. Minton, T. K.; Wright, M. E.; Tomczak, S. J.; Marquez, S. A.; Shen, L.; Brunsvold, A. L.; Cooper, R.; Zhang, J.; Vij, V.; Guenthner, A. J. Atomic Oxygen Effects on POSS Polyimides in Low Earth Orbit. ACS Appl. Mater. Inter. 2012, 4, 492–502. Vashisth, A.; Ashraf, C.; Bakis, C. E.; Van Duin, A. Reactive Molecular Dynamics Simulation of Accelerated Cross-Linking and Disintegration of Bisphenol F/DETDA Polymer Using ReaxFF. In American Society for Composites (ASC) 33rd Annual Technical Conference, the 18th US-Japan Conference on Composite Materials, and the ASTM D30; DEStech publications: Seattle, WA, USA, 2018. Bandyopadhyay, A.; Valavala, P. K.; Clancy, T. C.; Wise, K. E.; Odegard, G. M. Molecular Modeling of Crosslinked Epoxy Polymers: The Effect of Crosslink Density on Thermomechanical Properties. Polymer 2011, 52, 2445–2452. Vashisth, A.; Ashraf, C.; Bakis, C. E.; van Duin, A. C. T. Effect of Chemical Structure on Thermo-Mechanical Properties of Epoxy Polymers: Comparison of Accelerated ReaxFF Simulations and Experiments. Polymer 2018, 158, 354–363. Chanda, M. Advanced Polymer Chemistry: A Problem Solving Guide; Marcel Dekker Inc.: New York, 2000. Vashisth, A.; Bakis, C. E. Characterization of Nanosilica Filled Bis f Epoxide with Diamino Diphenyl Sulfone Curing Agents. In 31st Annual Technical Conference of the American Society for Composites, ASC 2016; DEStech publications: Williamsburg, VA, 2016. Chen, L.; Lee, C.-H. Interaction Potential between Atomic Oxygen and Polymer Surfaces in Low Earth Orbit. J. Spacecr. Rockets. 2006, 43, 487–496. Zeng, F.; Peng, C.; Liu, Y.; Qu, J. Reactive Molecular Dynamics Simulations on the Disintegration of PVDF, FP-POSS, and Their Composite during Atomic Oxygen Impact. J. Phys. Chem. A 2015, 119, 8359–8368. Rahnamoun, A.; van Duin, A. C. T. Reactive Molecular Dynamics Simulation on the Disintegration of Kapton, POSS Polyimide, Amorphous Silica, and Teflon during Atomic Oxygen Impact Using the Reaxff Reactive Force-Field Method. J. Phys. Chem. A 2014, 118, 2780–2787. Goverapet Srinivasan, S.; van Duin, A. C. T. Molecular-Dynamics-Based Study of the Collisions of Hyperthermal Atomic Oxygen with Graphene Using the ReaxFF Reactive Force Field. J. Phys. Chem. A 2011, 115, 13269–13280. Rahmani, F.; Nouranian, S.; Li, X.; Al-Ostaz, A. Reactive Molecular Simulation of the Damage Mitigation Efficacy of POSS-, Graphene-, and Carbon Nanotube-Loaded Polyimide Coatings Exposed to Atomic Oxygen Bombardment. ACS Appl. Mater. Inter. 2017, 9, 12802–12811. Chenoweth, K.; van Duin, A. C. T.; Goddard, W. A. ReaxFF Reactive Force Field for Molecular Dynamics Simulations of Hydrocarbon Oxidation. J. Phys. Chem. A 2008, 112, 1040–1053. van Duin, A. C. T.; Dasgupta, S.; Lorant, F.; Goddard, W. A. ReaxFF: A Reactive Force Field for Hydrocarbons. J. Phys. Chem. A 2001, 105, 9396–9409.



(24)

(25) (26) (27) (28) (29) (30) (31)

(32) (33)

(34) (35)

(36) (37)

(38) (39)

Page 32 of 33

Vashisth, A.; Ashraf, C.; Zhang, W.; Bakis, C. E.; van Duin, A. C. T. Accelerated ReaxFF Simulations for Describing the Reactive Cross-Linking of Polymers. J. Phys. Chem. A 2018, 122, 6633–6642. Li, C.; Strachan, A. Molecular Simulations of Crosslinking Process of Thermosetting Polymers. Polymer 2010, 51, 6058–6070. Lin, P.-H.; Khare, R. Molecular Simulation of Cross-Linked Epoxy and Epoxy- POSS Nanocomposite. Macromolecules 2009, 42, 4319–4327. Yarovsky, I.; Evans, E. Computer Simulation of Structure and Properties of Crosslinked Polymers: Application to Epoxy Resins. Polymer 2002, 43, 963–969. Pisacane, V. Fundamentals of Space Systems; Oxford University Press: New York, 1994. Brenner, D. W. Empirical Potential for Hydrocarbons for Use in Simulating the Chemical Vapor Deposition of Diamond Films. Phys. Rev. B 1990, 42, 9458–9471. Tersoff, J. Empirical Interatomic Potential for Carbon, with Applications to Amorphous Carbon. Phys. Rev. Lett. 1988, 61, 2879–2882. van Duin, A. C. T.; Strachan, A.; Stewman, S.; Zhang, Q.; Xu, X.; Goddard, W. A. ReaxFFSiO Reactive Force Field for Silicon and Silicon Oxide Systems. J. Phys. Chem. A 2003, 107, 3803– 3811. Mortier, W. J.; Ghosh, S. K.; Shankar, S. Electronegativity-Equalization Method for the Calculation of Atomic Charges in Molecules. J. Am. Chem. Soc. 1986, 108, 4315–4320. Senftle, T. P.; Hong, S.; Islam, M. M.; Kylasa, S. B.; Zheng, Y.; Shin, Y. K.; Junkermeier, C.; Engel-Herbert, R.; Janik, M. J.; Aktulga, H. M. The ReaxFF Reactive Force-Field: Development, Applications and Future Directions. npj Comput. Mater. 2016, 2, 15011. Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117, 1–19. Zhang, W.; van Duin, A. C. T. Improvement of the ReaxFF Description for Functionalized Hydrocarbon/Water Weak Interactions in the Condensed Phase. J. Phys. Chem. B 2018, 122, 4083–4092. Ziegler, J. F.; Biersack, J. P.; Littmark, U. The Stopping Power and Ranges of Ions in Matter. The Stopping and Range of Ions in Solids 1985, 1. Kański, M.; Maciążek, D.; Postawa, Z.; Ashraf, C. M.; van Duin, A. C. T.; Garrison, B. J. Development of a Charge-Implicit ReaxFF Potential for Hydrocarbon Systems. J. Phys. Chem. Lett. 2018, 9, 359–363. Luo, Y.-R. Handbook of Bond Dissociation Energies in Organic Compounds; CRC press, 2002. George, A.; Hambley, T. The Essentials of Organic Chemistry; Prentice Hall, 1996.


Page 33 of 33 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


TOC Graphic


Reactive Molecular Dynamics Simulations of ... - ACS Publications

Recommend Documents