Intermolecular Energy Transfer Dynamics at a Hot-Spot Interface in

Subscriber access provided by UNIV OF CALIFORNIA SAN DIEGO LIBRARIES

Article

Intermolecular Energy Transfer Dynamics at a Hot Spot Interface in RDX Crystals Kaushik L. Joshi, Martin Losada, and Santanu Chaudhuri J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.5b06359 • Publication Date (Web): 07 Jan 2016 Downloaded from http://pubs.acs.org on January 19, 2016

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 free 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 accessible to all readers and 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.

The Journal of Physical Chemistry A 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 35

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

Intermolecular Energy Transfer Dynamics at a Hot Spot Interface in RDX Crystals Kaushik Joshi1, Martin Losada2, and Santanu Chaudhuri1* 1

Applied Research Institute, University of Illinois at Urbana-Champaign, Champaign, IL 61820 2

Applied Sciences Laboratory, Washington State University, Washington, 99210-1495, USA

The phonon mediated vibrational up-pumping mechanisms assume an intact lattice and climbing of a vibrational ladder using strongly correlated multi-phonon dynamics under equilibrium or near-equilibrium conditions. Important dynamic processes far from-equilibrium in regions of large temperature gradient after the onset of decomposition reactions in energetic solids are relatively unknown. In this work, we present a classical molecular dynamics (MD) simulation based study of such processes using a non-reactive and a reactive potential to study a fully reacted and unreacted zone in RDX (1,3,5-trinitro-1,3,5-triazocyclohexane) crystal under nonequilibrium conditions. The energy transfer rate is evaluated as a function of temperature difference between the reacted and unreacted regions, and for different widths and cross-

*To whom correspondence should be addressed. E-mail: [email protected], Office: 217-300-9419

1

ACS Paragon Plus Environment


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

sectional area of unreacted RDX layers. Vibrational up-pumping processes probed using velocity autocorrelation functions indicate that the mechanisms at high-temperature interfaces are quite different from the standard phonon-based models proposed in current literature. In particular, the up-pumping of high-frequency vibrations are seen in presence of small molecule collisions at the hot spot interface with strong contributions from bending modes. It also explains some major difference in order of decomposition of C-N and N-N bonds as seen in recent literature on initiation chemistry.

1. Introduction Molecular solids which constitute different classes of high explosives (HE) are unique by virtue of their high density of chemical energy and rapid rates of energy release. The behavior of HE under extreme conditions (high pressure and temperature) is of interest for a wide range of both civilian and military applications. The ability to engineer safer HE which are insensitive in nature will benefit from a better understanding of the events before detonation.1 In the past, the focus of research concentrated on identifying the atomistic scale pathways for shock induced initiation of chemistry in HE crystals. Upon passage of the shock, low-frequency lattice vibrations, phonons, are primarily excited. In order to systematically analyze the reactive events during the initiation and growth of high-temperature regions within the energetic material, the first step is to understand the mechanisms of vibrational energy transfer from fully reacted hot molecular fragments to the unreacted solid. It is generally accepted that the initial energy in the low-frequency vibrations has to be deposited into higher-frequency vibrations, normally, bond stretching modes. Dlott, Fayer, and Tokmakoff2-3 have studied such multiphonon up-conversion

2


Page 2 of 35

Page 3 of 35

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


processes associated with shock-induced chemistry. They derived a master equation for phonon and vibron temperatures. The master equation depends on the phonon-vibron energy transfer rate, which was estimated from experimental results on naphthalene, arguing that the measured anharmonic coupling matrix elements in a variety of solid materials composed principally of C, N, H, and O atoms lie within about a factor of two of the naphthalene values. Kim and Dlott4 have further studied thermal conduction and vibrational cooling in naphthalene through molecular dynamics simulations.5-6 Their results show that the energy transfer from the hot molecules to its nearest neighbors is anisotropic. The results also show that the redistribution of vibrational energy occurs via an indirect mechanism where vibrational energy is transported by phonon intermediates. Besides identifying the gateway modes, few other attempts were made to investigate how energy is distributed into various vibrational modes. Fried et al. showed that the energy transfer rate should be measured as a function of frequency and temperature.7 Using normal mode analysis, Dawes et al.8 studied energy transfer rates between translational, rotational and vibrational modes in shocked solid nitromethane. Boyd and Boyd9 correlated center of mass motions with the occupation of internal degrees of freedom in αRDX. The authors found that the nitro-group rotations play a significant role in transferring the energy from lattice modes to intramolecular motions. On this basis, a sequence of energy transfer process is likely to begin with

the lowest frequency modes, then passing to mid-frequency;

normally bond bending modes, and finally exciting the high frequency bond stretching vibrations Apart from the decomposition of HE under phononic excitations, the excitation of intramolecular modes under the influence from a dense medium is also of particular interest for initiation chemistry that leads to the growth of so-called hot spots. Such excitations are generally

3



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

a precursor to triggering of the chemistry at the interface between a hot spot and unreacted HE crystal. Derivation of a probability of reaction due to vibrational excitation assumes the partitioning of energy along some commonly known dissociation pathways. A recent study using a classical potential for α-RDX analyzed the partitioning of energy among different modes and identified that there is a large contribution of N-N bond vibrations in lower temperatures but the contribution drops from 70% to nearly 20% in temperatures in 500-2000 K range.10 Pressure dependence of such modes have also been studied both for α- and γ-RDX phases using the same classical potential for infrared spectra in the THz region.11 We have recently demonstrated that the ignition of the high-pressure γ-RDX from nanoscale hot spots are possible at temperatures around 1700 K.12 Experimental and computational investigation of phonon-vibron coupling and multiphonon densities of states have been performed for different HE crystals showing difference in vibrational lifetime.13 The implication of lower vibrational lifetime of an HE crystal is in localization of energy which can cause chemical decomposition and lead to detonation. Therefore, up-pumping and energy transfer rates at an interface with large temperature gradient are important for determining the sensitivity of HE materials. In the present paper, we will focus on the intermolecular energy transfer from high temperature gas phase molecules to the unreacted solid α-RDX (1,3,5-trinitro-1,3,5triazocyclohexane) molecules at a hot spot interface using non-reactive and reactive molecular dynamics(MD). Under the very high-density and confined conditions that exist in these hot spot interfaces, the basic tenets of collisional energy transfer from strong colliders are not necessarily valid. Therefore, for the hot spot consisting of small molecules, the unreacted layers of 1,3,5trinitro-1,3,5-triazocyclohexane (RDX) molecules treated in this paper can be considered using 4


Page 4 of 35

Page 5 of 35

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 concepts from condensed phase theories of transfer of vibrational excitations. Here, a layer of polyatomic molecular oscillators in a dense hot medium is subjected to a chaotic phonon field. The theoretical background has been discussed in detail by Nitzan et al. in identifying different classical and quantum effects related to vibrational excitation process in condensed phase14-17. The temperature of the dense medium studied in the current work is in the range of 1000-3000 K. 2. Simulation Methods For studying the energy transfer process at a gas-solid interface, the simulated system is composed of two regions, the reacted region represented by the hot gas-phase molecules and the unreacted region that is represented by solid α-RDX (Figure 1). Since the molecules in the unreacted region can undergo chemical reactions during the energy transfer process, we have performed the MD simulations using both reactive and non-reactive interatomic potentials. In the non-reactive MD simulations, the high temperature gas-phase molecules are represented by the gas products that are likely to form in a fully reactive RDX system. For the reactive MD simulations, the high temperature gas-phase is represented by inert Argon atoms. The choice of inert gas ensures that the unreacted solid molecules at the interface do not undergo any chemical reactions with the high temperature gas-phase molecules. The details of non-reactive and reactive simulation methodologies are described in the next two sections. 2.1 Non-reactive MD simulations For the non-reactive MD simulations, we first took a 7×3×3 α-RDX cell (14 layers in xdirection) to represent the unreacted molecular system and added a second cell along the x direction, containing H2O, NO, CO2, N2, CO, and NH3 molecules to represent the high temperature gas products describing a fully reacted RDX system18. The amount of each gas 5



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

product in the cell was adjusted to match the mass balance with respect to the unreacted system. Thus a new cell of 14×3×3 size (184.548 Å ×34.722 Å ×32.127 Å) with 21168 atoms was obtained as shown in Figure 1. The configuration shown in Figure 1 is termed as basic supercell in rest of the article. For obtaining the estimates of the energy transfer rates at the gas-solid interface, we simulated three systems, the single layer interface, the double layer interface and the entire crystal as interface (14 layers) (Figure 1). The single and double layers contain 36 and 72 RDX molecules, respectively. During the simulations, the gas products and the RDX molecules at the interface were allowed to interact while the remaining RDX molecules (noninterface) representing the solid crystal were kept fixed.

6


Page 6 of 35

Page 7 of 35

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 1. Schematic diagrams showing the showing the different configurations used in nonreactive MD simulations to study energy transfer from hot gas-phase to RDX interface. Color scheme: red – hot gas phase, green – interface RDX, blue – unreacted RDX crystal. The intramolecular bonding in RDX was represented by the fully flexible, nonreactive molecular potential developed by Smith et al.

19

(SB). The non-bonded interactions were

described using both Lennard-Jones and Buckingham potentials, with parameters taken from Thompson et al. 20-21 and Bedrov et al. 22, respectively. The Coulomb interactions with the fixed partial charges on the nuclei were also included. The long-range electrostatic interactions were calculated using the PPPM method

23

with a cutoff of 10 Å. The k-space solver threshold for

computing the forces was set to 10-6. At elevated temperatures, the chemical bonds can get stretched towards the bond dissociation region of energy surface. In order to take into account this anharmonic effect, we also performed the non-reactive MD simulations in which the chemical bonds in RDX were described by the Morse potential. The anharmonic force field parameters for the stretching mode were taken from Wallis and Thompson24. For the single and double layer non-reactive simulations, the gas-phase molecules and the interface RDX molecules were equilibrated at 550K using isothermal-isochoric (NVT) ensemble for 150 ps duration. The equilibrium temperature of 550 K was used to represent the temperatures observed in bulk RDX crystals after passage of a strong shock wave using SB potential. The higher temperature was picked to ensure that the interfacial energy transfer process is not convoluted with melting of the crystal. The hot spot growth will require vibrational excitation of the molecules at the molten RDX interface and reactions at a very fast rate. The results discussed in this paper focus on potential pathways at the hot spot interface that can lead

7



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 8 of 35

to ultrafast excitation of vibrons such as C-N, C-H and N-N bonds which are assumed to be key for the initiation pathways. The interface between the hot spot and condensed-phase RDX is molten. Thus, the mechanisms discussed are valid for high temperature dynamics and activation of RDX molecules. To mimic these vibrational excitations, starting from the equilibrated spatial configuration, the gas phase molecules were provided the thermal pulse by assigning the initial temperatures of 1000K, 2000K and 3000K. The desired temperatures were achieved by reassigning the velocities of the gas-phase products according to the Maxwell-Boltzmann velocity distribution for the corresponding temperatures. Since the simulated interfaces are relatively small with only 36 (single layer) and 72 (double layer) molecules, it is not required to maintain the gas-phase molecules at elevated temperatures. Hence, after assigning instantaneous thermal pulse, both the interface and the gas-phase molecules evolved using an NVE ensemble for next 100 picoseconds to calculate the temperature and energy transfer rates from the gasphase to the RDX interface. The simulations were performed using an integration time steps of 0.1 fs. The non-reactive MD simulations described here were performed using the LAMMPS 26

25-

package. For the interface with 14 layers (entire crystal), the system was initially equilibrated at

550K. Then the gas-phase molecules were assigned the initial temperature of 1000K, 2000K and 3000K. Using NVT ensemble, the gas-phase molecules were maintained at those respective temperatures until the unreacted RDX crystal attains the thermal equilibrium with the gas-phase. This NVT-based temperature control is essential for gas-phase molecules because the entire portion of unreacted RDX crystal requires considerable amount of thermal energy to reach the

8


Page 9 of 35

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


elevated temperatures. During this energy transfer process, the entire unreacted RDX crystal was allowed to evolve using NVE ensemble. For studying the effect of gas-phase/interface temperature gradient, we initially equilibrated the system at 300K. Then the gas-phase molecules were given an initial thermal pulse of 350K, 375K, 400K and 450K respectively and were maintained at those respective temperatures using NVT. The rest of the unreacted RDX crystal was integrated using NVE ensemble. The entire unreacted RDX crystal was allowed to interact in these temperature gradient cases. These simulations were performed with harmonic, non-reactive SB potential with a time step of 0.25fs.

Figure 2: Schematic of larger supercell (84,672 atoms) with double layer RDX interface used in no-reactive MD simulations. Interface cross section is 4x bigger than the basic supercell (Figure 1). Color scheme: red – hot gas phase, green – interface RDX, blue – unreacted RDX crystal We also performed non-reactive simulations to investigate the effect of cell size on the energy transfer process. The double layer system with cross-sectional area that is four times 9



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 10 of 35

larger than the basic supercell was used for investigating the cell size dependence on the energy transfer process (Figure 2). The MD simulations were performed with harmonic, non-reactive SB potential. The rest of MD set-up used for this interface is similar to the previous double layer calculations on basic supercell.

2.2 Reactive MD simulations The reactive simulations were performed only for single and double layer configurations of basic supercell. The solid unreacted portion of the reactive simulations is exactly same as the non-reactive simulations. The gas products that were used in the non-reactive MD simulations have a finite probability of triggering chemical reactions in the interface RDX molecules. Since the main emphasis of this article is to study energy transfer rates, we replaced the gas-phase molecules with inert Argon atoms in high temperature gas-phase portion of the simulation cell. The number of Argon atoms added was adjusted to obtain the density of the entire system that is same as α-RDX. The Argon atoms were added using PACKMOL package.27 Similar to the nonreactive simulations, the single layer and double layer interface RDX molecules were allowed to interact with the high temperature gas-phase portion of the simulation cell. The non-interface RDX atoms were kept fixed during the reactive MD simulations. The interatomic interactions were described using the ReaxFF28 reactive force field method. The ReaxFF is a reactive potential that can model bond breaking and bond formation on the fly during molecular dynamics simulations. To model chemical reactions, it combines bond order/bond distance relationship with a polarizable charge description and bond-order dependent 10


Page 11 of 35

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


3- and 4-body interactions. The combination of covalent and columbic interactions has made ReaxFF applicable to a wide range of materials. A more detailed description of the force field method can be found in ref 27. The ReaxFF has been extensively used to obtain atomic-scale description of the decomposition of high energy materials.29-40 Recently, Wood, van Duin and Strachan re-fitted the existing force field parameters to merge it with the combustion force field for greater transferability.41 So we have used those parameters for performing the reactive MD simulations. The Argon parameters are from Kamat et al. study of soot formation in presence of various monatomic and polyatomic gas molecules like He, Ne, Ar, N2, CO2, and CH4.42 Before performing the actual reactive simulations, the system was energy minimized so that the atoms can relax into local energy minima. The remaining MD set up is same as nonreactive simulations. It should be noted that only the chemical reactions between the gas-phase Argon atoms and the interface RDX molecules were suppressed. Although Argon atoms cannot react with interface RDX, the thermal energy up-pumping from the gas-phase molecules can still trigger the chemical reactions in the interface RDX molecules due to thermal excitations. Such thermally excited chemical reactions within the interface RDX molecules were allowed to evolve during these reactive MD simulations. 3. Results and Discussion 3.1 Heating of RDX at the interface Figure 3 shows the average temperature profiles of RDX molecules at the interface during the NVE simulations. The temperature values displayed in Figure 3 are the temperature profiles from MD simulations based on the total kinetic energy of each layer, usually defined as T=2Ek/3Nk. At this point, it is relevant to mention that since a single RDX molecule contains 21 11



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 12 of 35

atoms, the total number of degrees of freedom is 63. The three translational degrees of freedom account for 4.8 % of the total degrees of freedom. This shows that the total temperature profiles presented in Figure 3 are mainly dominated by the intramolecular vibrational contributions. It can be seen that the interface temperature increases from an initial value of 550K until it reaches a plateau of the thermal equilibrium. For both the reactive and non-reactive MD simulations, the interface temperature rises faster in the single layer system than the double layer system. As a result, the temperature reaches thermal equilibrium plateau quicker in the single layer configuration than in the double layer system. Also for both types of simulations, the RDX interface gets equilibrated at lower temperatures for the double layer interface than the single layer interface. This longer time and lower temperature increase in the double layer configurations are

12


Page 13 of 35

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


(a)

(b)

Figure 3. Temperature profiles of the RDX molecules at interface. The three initial gas phase temperatures of 1000, 2000, and 3000 K are shown in the panels, (a) Non-reactive MD with SB potential, (b) Reactive MD with ReaxFF potential fully in line with the fact that there are more degrees of freedom to distribute the same amount of thermal energy from hot gas-phase molecules with respect to the single layer. Figure 3 shows that the incorporation of anharmonic potential does not have any significant effect on the equilibrium temperature.

(a)

(b) 13



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 14 of 35

Figure 4: Temperature profiles of RDX molecules from the non-reactive simulations, (a) Larger supercell with 84,672 atoms and 4x cross-sectional area, (b) Effect of number of interface layers. The entire unreacted RDX (14 layers) was allowed to interact with the hot gas-phase molecules. Figure 4a shows the size of interface cross-sectional area does not have any significant effect on the temperature profiles. The overall trend in temperature profiles for all the three temperatures is similar to the trend obtained for the basic supercell (Figures 3a). The temperature profiles for the simulations, in which the entire unreacted crystal (14 layers) was allowed to interact, are shown in Figure 4b. The observed trends are similar to single and double layer simulations. As expected, the time required for reaching the thermal equilibrium for the entire crystal is higher than the single and double layer simulations. Figure 3 also shows that, for both the single layer and double layer configurations, the interface equilibrium temperature attained in non-reactive simulations is higher than the corresponding interface equilibrium temperature of reactive simulations. This difference in equilibrium interface temperatures between the non-reactive and reactive simulations becomes more evident with the increase in temperature of the gas-phase molecules e.g. for double layer configuration, after 60 ps, the difference between the interface equilibrium temperatures for nonreactive and reactive simulations is less than 50K for 1000K case, is approximately 100K for 2000K simulations and then it increases up to 350K for 3000K simulations. We did not observe any chemical reactions in reactive MD simulations up to the gas-phase temperature of 2000K during 100 picoseconds. However, when the gas-phase temperature was set to 3000K, the uppumped energy activated some interface RDX molecules within 100 picoseconds forming intermediate radicals in the interface.

14


Page 15 of 35

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 5 shows the close-up view of the interface RDX molecules from the double layer harmonic non-reactive simulations. The harmonic SB potential used for non-reactive simulations predicts that the α-RDX undergoes solid-solid phase transition between 480 and 490K followed by melting between 500 and 510K43. As a result, the interface RDX layers melt during the NVT equilibration at 550K (Figure 5, 0ps). With the further increase in temperature due to the uppumping of energy from the hot gas-phase molecules, the interface layers continue to drift from crystalline arrangement of α-RDX forming an amorphous molten interface. Although Figure 5 shows the interface structure only for the double layer configuration of basic supercell, the similar amorphous molten structure of RDX interface is observed in the larger supercell (Figure 2) and for the basic supercell in which the entire unreacted RDX (14 layers) was allowed to interact with the gas-phase.

Figure 5: Close-up view of double layer interface during the energy up-pumping from the hot gas-phase molecules at 1000K in non-reactive MD simulations. Color scheme: red – gas phase, green - double layer interface, blue – unreacted RDX crystal In all the cases, the interface temperatures displayed in Figure 3 are high enough to initiate the chemical decomposition of RDX

44

in longer duration. The thermal excitation of

15



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 16 of 35

RDX at the interface is likely to excite key vibrational modes that could trigger chemical reactivity, as no thermochemical reaction can occur without some amount of vibrational activation. Previous experimental and computational studies45-47 have shown that excitation of NN and C-N bond stretching modes play a key role in initiating RDX by concerted ring fission and dissociation of NO2 groups. So in the next section, we focus on the excitation of various vibrational modes in RDX interface. 3.2 Excitation of Vibrational Modes In order to identify which vibrational modes were excited during the energy transfer process, we calculated the time averaged power spectrum of the interface RDX for both the nonreactive and reactive MD simulations. In addition, we also investigated the energy stored in the different bond stretching and bending modes of the interface RDX. 3.2.1 Power Spectrum of Interface RDX The power spectrum was calculated only for the double layer simulations. For calculating the power spectrum, we initially calculated velocity auto-correlation function (VACF) given as, = 〈( ) × ( + )〉. For calculating VACF, for each MD trajectory, the velocities are sampled every 2.5fs up to 50ps. For enhanced statistical averaging, the VACF was calculated using shifting register technique in which every 5th sample is assumed as new time origin. The power spectrum of the interface rdx molecules is then calculated by taking Fourier transform of

VACF to convert time domain into frequency domain, () = (−) (). Figure 6 shows the power spectrum of interface rdx molecules obtained from the non-reactive and reactive MD simulations. We did not calculate the power spectrum for the reactive

16


Page 17 of 35

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


simulations when the gas-phase temperature was set to 3000K due to the initiation of some RDX molecules in the interface. Although the reactive potential was not explicitly tuned for vibrational frequencies, the reactive spectrum shows qualitative similarities with the non-reactive spectra. For all three potentials, the intensity of the peaks increases with the increase in interface temperature while the location of the peaks remain the same. This trend shows that the similar types of vibrational modes are excited by the up pumping of thermal energy in the simulated temperature range. Figure 6 also shows that the peaks become broader with the increase in temperature. The peak around 3000 cm-1 region corresponds to the C-H bond stretching mode. The reactive spectra show a blue shift for the C-H peak. The next dominant region in power spectrum is at around 1400- 1700 cm-1. The region near 1500 cm-1 corresponds to N-N and N-O bond stretching modes whereas region near 1400 cm-1 corresponds to C-N bond stretching vibrations. The previous studies9 with

(a)

(b)

17



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 18 of 35

(d)

(c)

Figure 6: Calculated power spectrum of interface RDX molecules at different gas-phase temperatures, (a) Harmonic SB potential, basic supercell (b) Anharmonic SB potential, basic supercell (c) Harmonic SB potential, larger supercell, (d) ReaxFF potential harmonic bond stretching potential have also shown that C-H2 rocking mode contributes to the peak around 1400 cm-1. In addition to these modes, we also observe a peak around 800 cm-1 which is mostly contributed by NO2 stretching, CH2 twisting and ring breathing modes48. The power spectra of the larger supercell show excitation of similar vibrational states as the smaller system. This indicates that the cross-sectional area of the RDX interface has a nominal effect on the power spectra. Also, in all the cases, the observed intensity of phonon modes is weaker than the bond stretching and bond bending modes. This can be due to the melting of the interface layers and leading to the formation of local amorphous structure (Figure 5). Comparison between different types of potentials also shows that the anharmonic nonreactive spectra reflect the gradual transition from the non-reactive to the reactive potential. For example, the harmonic potential shows the split in the C-H peak at lower temperatures. However, such splitting of C-H mode is not observed in anharmonic non-reactive and in reactive 18


Page 19 of 35

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


ReaxFF simulations. The previous studies48 have shown that N-O bond stretching mode is responsible for the sharp peaks near 2200cm-1 region. However, similar to reactive potential, these peaks are completely absent in the power spectra obtained from anharmonic potential. Also, unlike both the non-reactive spectra which shows the single broad peak near 1400-1700 cm-1 region, the reactive spectrum shows two distinct peaks, one near 1400 cm-1 and second peak near 1750 cm-1. The power spectrum analysis clearly shows that the thermal up-pumping of energy from hot gas-phase to interface RDX excites both stretching and angular twisting vibrational modes in the interface RDX molecules. The excited bond stretching modes includes C-H bonds, N-N bonds and C-N bonds along with the ring breathing mode.

3.2.2 Nuclear Quantum Effects on Power Spectra The relative intensities of the peaks can vary if quantum effects are included in MD simulations. One of the well-known shortcoming of MD simulations is that it produces MaxwellBoltzmann statistics instead of Bose-Einstein statistics. The conventional MD simulations cannot fully capture the nuclear quantum effects. These effects can alter the power spectrum and heat capacity particularly for the materials that have large population of hydrogen atoms. So in order to investigate this quantum effect, we performed MD simulations with quantum thermal bath49 as implemented in LAMMPS. This technique accounts for quantum statistics by simply coupling the motion of an ensemble of coordinates interacting through a given force field to a fictitious

19



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 20 of 35

quantum thermal bath at a low computational cost. The MD simulations were done only for the double layer system (basic supercell) at 1000K with harmonic SB potential.

Figure 7: Effect of quantum thermal bath on the power spectrum of double layer interface RDX molecules Figure 7 shows the effect of quantum bath on power spectrum. It can be seen that incorporation of quantum effects alter the relative intensities of the different peaks and reduces the intensity C-H vibrational modes with expected broadening. However, the power spectra with quantum correction still predicts that C-N, N-N and C-H bonds along with bending modes in interface layers are excited during the energy transfer process. As expected from the Debye temperature of bonds, the C-H bond is supposed to contain negligible energy at 2000K for an isolated system. Calculation of accurate C-H bond population in a many-body system is a formidable task. In addition, the coupling to the bending mode in RDX molecule is significant. As a result, the trend indicate an up-pumping mechanism that needs further detailed investigation

20


Page 21 of 35

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


including a review of QTB method’s applicability in a complex molecular system. For proton transfer spectra studied in literature, it is seen that role of quantum correction offered by QTB is promising but the corrections are minimal at higher temperatures50. Recently it was found that QTB may be more suited for small molecular systems with low anharmonicity51 and may need further development to obtain accurate behavior. Since, only high temperature up-pumping dynamics is our main focus in the current work, we maintain that the power spectrum corrections from QTB are reasonable for the current purposed work. 3.2.3 Energy Up-pumping in Bond Stretching and Bending Modes In order to quantify the distribution of transferred energy into various vibrational modes, we further analyzed the change in each component of the potential energy in the double layer simulations. In LAMMPS simulations, one can extract the component-wise energy (bond, angle, dihedral etc) contributions for any group of atoms. Such component-based energies are used for investigating the partitioning of up-pumped energy into the various energy components. The component-based analysis of potential energy indicates that the most of the change in potential energy is due to the change in bond energy and bending energy. Figure 8 shows the change in those two energy components from both; the harmonic non-reactive and reactive MD simulations for the basic supercell at various gas-phase temperatures. The anharmonic non-reactive results are not reported here as they show similar trends as that of harmonic potential. It can be seen that the excitation of the bending vibrational modes is as dominant as bond stretching mode. The amount of energy absorbed by the bond stretching modes in both types of simulations is in good agreement with each other at all the temperatures. However, the amount of energy transferred in the bending 21



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

Non-reactive MD

Page 22 of 35

Reactive MD

Figure 8: Distribution of absorbed thermal energy into stretching and bending vibrational modes in double layer RDX interface during the non-reactive and reactive MD simulations mode is much higher than the stretching mode in non-reactive simulations at higher temperatures. For example, in non-reactive simulations, the bending mode absorbs nearly 400 kcal/mol of more energy than the bond stretching mode at 2000K. The difference in the absorbed energies between the two modes increases further to 850 kcal/mol at 3000K. For reactive simulations, the difference in the amount of energy absorbed between two modes is approximately 200 kcal/mol and it roughly remains constant at all temperatures. Although Figure 8 shows how the absorbed energy is partitioned into different vibrational modes, it does not show the relative rates at which different bond stretching modes such C-H, NN and C-N modes are excited. In order to gain further insight into the relative order at which different stretching modes are excited, we calculated the energy absorbed by C-H bonds, C-N bonds and N-N bonds as function of time for the double layer simulations at 1000K and 2000K. For evaluating the energy stored in each type of bond, we first calculated instantaneous 22


Page 23 of 35

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


interatomic distances for the atom pairs that were part of the interface RDX layers after every 2.5 fs. Then using harmonic expression that is used in SB potential, we calculated the energy stored in each bond. The energy stored in each bond type was then averaged over 250 fs to eliminate the inherent noise. This analysis was done only for 20 ps at 1000K because bond energy reaches its equilibrium in that duration. At 2000K, we extended that analysis up to 50ps. Also, this analysis was done only for non-reactive harmonic simulations because the bond energy is described by simple harmonic function. It should be noted that the incorporation of quantum nuclear effects can influence the relative ordering of stretching excitations. Figure 9 shows the average energy up-pumping rates in C-H, C-N and N-N bonds. At both the temperatures, our calculations indicate that C-N and C-H bonds store roughly twice the energy than N-N bonds during the energy transfer process. At 2000K, the C-H and C-N bonds absorb energy at the rate of 28.3 kcal/mol-ps for first 6 picoseconds and then at the rate of 2.75 kcal/mol-ps for next 40 ps. At the same temperature, the energy transfer rate in N-N bonds is 14.2 kcal/mol-ps for first 6 picoseconds and 1.3 kcal/mol-ps for next 40 ps. Figure 9 also shows that with respect to total energy transferred into the RDX interface, the thermal energy contained in N-N bonds is roughly 9% at 1000K and 5% at 2000K. This observed relative contribution of N-N bonds is lower than that estimated by Kraczek and Chung using lattice dynamics method. The authors reported that N-N bond stores at least 13% of the total thermal energy above 1000K in phonon-mediated energy transfer process. It should be noted that each RDX molecule has three N-N bonds and six C-N and C-H bonds each. As a result, on per bond basis, the amount of energy absorbed by different types of bonds is roughly the same.

23



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 24 of 35

Figure 9: The rate of up pumping of thermal energy into C-H, C-N and N-N bonds obtained from the non-reactive MD simulations 3.2.4 Effect of Temperature Gradient on Energy Up-pumping In order to investigate the effect of temperature gradient between gas phase and interface on the partitioning of adsorbed energy, we calculated the energy adsorbed by these three bonds at lower temperatures. Figure 10 shows the up pumping rate when the gas-phase temperatures were in range from 350K to 450K. It should be noted that the entire unreacted RDX crystal was allowed to interact in these calculations.

24


Page 25 of 35

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 10: Effect of temperature on the up-pumping of thermal energy into C-H, C-N and N-N bonds. The amount of energy stored in C-H and C-N bonds increases with the increase in RDX crystal temperature. At 350K, the total energy absorbed by N-N bonds is comparable to C-H and C-N bonds. Therefore, the energy absorbed per N-N bond is nearly twice that of C-H or C-N bonds on a per bond basis. However, the C-H and C-N bonds progressively store more energy with increase in temperature. At 450K, the contribution by C-H and C-N bonds is nearly twice that of N-N bonds. This energy partitioning trend at 450K is similar to the trend observed at higher temperatures for double layer systems. As discussed previously, the interface RDX layers are likely to melt around 500K. The observed trends on the energy partitioning clearly indicates that as the unreacted crystal approaches the melting temperature, the extent of excitation of C-H and C-N bonds increases until the energy stored by C-H, C-N and N-N bonds become comparable on per bond basis. Based on the power spectra and the component wise analysis of the absorbed energy at various temperature gradients between gas-phase and interface, our calculations indicate that the bending modes are as dominant as bond stretching modes and the N-N, C-N and C-H bonds show comparable levels of excitations during the energy transfer process beyond 450K. 3.3 Intermolecular Energy transfer 25



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 26 of 35

In this section, the numerical value of intermolecular energy transfer, expressed as rate (κ), from the gas-phase to the unreacted RDX molecules at the interface is computed using the following expression

κ =

∇ ∇ × ∇!

where ∇ is change in potential energy, ∇! is the corresponding change in temperature and ∇ is time required to reach the average plateau in potential energy change in the MD simulations. The potential energy profiles clearly indicate that the time required to reach the equilibrium state increases with increase in the interface layer. As a result, we have used 30ps as equilibration time for single layer configurations and 60ps for double layer configurations (Figure 11). For the systems in which the entire unreacted crystal was allowed to interact with the hot gas-phase, we have used equilibration time of 400ps for calculating energy transfer rate (Figure 4b).

26


Page 27 of 35

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


(a)

(b)

Figure 11. Potential energy change of interface RDX molecules for the basic supercell, (a) Non-

reactive SB potential, (b) Reactive ReaxFF potential The calculated values for single and double RDX layers for both the non-reactive and reactive potentials are tabulated in Table 1. For the reactive MD simulations, the potential energy did not equilibrate after 100ps when the gas-phase temperature was initialized to 3000K. This is mainly because at 3000K, approximately 28% of interface RDX molecules undergo initiation reactions and as a result, the entire interface becomes chemically active region that evolves through various endothermic and exothermic chemical reactions. Hence energy transfer rates for 3000K reactive simulations are not reported in table 1. Although the energy transfer rates for larger supercell are not included in table 1, similar energy transfer rates as basic supercell were obtained for larger supercell. 27



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 28 of 35

Table 1. Calculated energy transfer rates (κ) for basic supercell. Rate units are in J/K-mol-s. Gas Phase Initial

Non-reactive

Non-reactive

Harmonic MD

Anharmonic MD

Temperature

Single

Double

(K)

Layer

Layer

Entire Crystal, 12 Layers

Reactive MD

Single

Double

Single

Double

Layer

Layer

Layer

Layer

1000

2.85×1014 1.96×1014

2.95×1014

2.87×1014 3.32×1014

2.30×1014

2.96×1014

2000

3.13×1014 2.96×1014

3.47×1014

3.06×1014 3.45×1014

2.88×1014

3.26×1014

3000

3.18×1014 3.16×1014

3.44×1014

3.29×1014 3.42×1014

-

-

Table 1 shows that the anharmonic energy transfer rates are higher than the corresponding harmonic rates for most of the cases. This can be due to the fact that the anharmonic bond stretching function allows the atoms to move in bond dissociation region without placing large energy penalty. As a result, the equilibrium potential energy reached in anharmonic simulations is higher than the corresponding harmonic cases (Figure 11). This enables the interface atoms to rise to higher energy state by absorbing more thermal energy from the hot and confined gas-phase molecules. The table 1 also shows that the anharmonic nonreactive and the reactive MD-based energy transfer rates are in reasonably good agreement for most of the simulations. For both the anharmonic non-reactive and reactive potentials, the energy transfer rate increases from single layer to double layer. The observed energy transfer rates of the systems in which the entire crystal was allowed to interact, are higher than single and double layer interfaces at all the temperatures. The calculated energy transfer rates from the different potentials have the same order of magnitudes for all the three systems indicating that the underlying energy transfer mechanism is similar for all three interfaces.

28


Page 29 of 35

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 simulated systems in this work are different and have not reached chemical equilibrium, it is nonetheless worthwhile to compare our values to the ones reported by Dlott et al.2 In Dlott’s work, the energy transfer parameter characterizes the rate of phonon-pair conversion to vibration in naphthalene, which is a reasonable model for up-pumping for RDX and HMX.2, 4 The authors reported the energy transfer rate of 2.5×1012 J/(K-mol-s). The energy transfer rate obtained from our MD simulations is large by a factor of 100. This is mainly due to the fact that RDX molecules have different stretching modes in the form of C-N and N-N bonds and also have higher number of angular degrees of freedom in the form of bending modes that can absorb considerable amount of energy. 4. Conclusions We have performed classical MD simulations to study the process of up-pumping of thermal energy into unreacted RDX crystals at the gas-solid interface that corresponds to the conditions that may exist around a growing hot spot. The energy transfer process at such gassolid interface is not clearly understood in the condensed phase literature due to complexity of vibrational excitation dynamics and high magnitudes of temperature gradients. In order to investigate such non-equilibrium processes, both reactive and non-reactive interatomic potentials were used to study the effect of temperature gradient and interface width on the energy transfer process. The results show that the high temperature gas phase molecules in the reacted domains transfer considerable amount of thermal energy to the solid unreacted counterpart at the interface exciting key vibrational modes within the molecules. As seen in reactive MD simulations, beyond certain critical temperatures, these excited vibrational modes trigger the local chemical reactivity that can lead to the growth of hot spots. The power spectra of the interface RDX

29



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 30 of 35

molecules show that in addition to various bending modes, the N-N, C-H and C-N bond stretching modes are activated in the simulated temperature range. We investigated the nuclear quantum effects on power spectrum using QTB method. Although inclusion of the quantum effects can change the relative ordering of stretching and bending excitations, we observed excitation of similar vibrational states in the presence of quantum effects. Also, it has been demonstrated in the literature that at high temperatures, the quantum nature of the vibrational states are smeared and thus, it is less important compared to the more classical energy transfer pathways discussed in this manuscript. The observed energy transfer process is indicative of a key difference in initiation pathways compared to a doorway mode based analysis. In doorway mode based analysis, vibrational excitation of higher frequency N-N, C-H and C-N bonds happen at a later stage. In this work, a perfect RDX crystal interface exposed to hot and decomposed gas molecules resembles the energy transfer processes expected in a highly correlated but random phonon bath. This mechanism show a different order of up pumping. Moreover, the component-based analysis of potential energy indicates that the angular vibrational modes are as dominant as bond stretching modes. It was also found that the thermal excitation of different types of bonds depends upon temperature. Below 350K, the N-N bonds absorb more thermal energy than C-N and C-H bonds. However, as the interface temperature approaches the RDX melting point, then it was found that the C-N and C-H bonds of RDX rings are more likely to get excited before the N-N bonds. The energy transfer rates from both the reactive and non-reactive models were compared and we found that both the models are in reasonably good agreement with each other. Although

30


Page 31 of 35

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 calculated energy transfer rates are 100 times faster than the previous results of molecular models based on naphthalene, these results can be rationalized based on bending modes of RDX molecules which makes the system capable of absorbing thermal energy at a faster rate. The physical accuracy of simulation results is closely related to the accuracy of the potentials. The comparisons presented in this work indicate that the two force fields agree quite well with each other before the onset of reactions. Considering the differences between the force fields, the results indicate a fundamental similarity in the process of energy transfer at a molten, condensed phase RDX interface. Also, the inclusion of anharmonicity in SB potential narrows the gap in the energy transfer rates. The results also provide a strong indication that hot spot based activation models need to carefully consider the role of temperature difference between the hot spot and the unreacted energetic crystal. The sensitivity to hot spot based initiation of chemistry and the enhanced energy transfer rate can both be reconciled using the differential rates of vibrational excitation dynamics of C-N and N-N stretching and the higher number of angular degrees of freedom in such energetic crystal interfaces. Further quantum mechanical study is needed to extend these results to related temperature dependent reactions rates from a variational transition state theory or quantities such as collisional reaction RRKM rates to derive the upper limits of hot spot growth rates52-53.

31



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 32 of 35

Acknowledgments We are grateful to Professor Thomas Sewell from MU and Dana Dlott from UIUC for helpful discussions. This work is supported by the Defense Threat Reduction Agency (DTRA), Grant No. HDTRA1-15-1-0034.

References 1. Dlott, D. D., New Developments in the Physical Chemistry of Shock Compression. Annu. Rev. Phys. Chem. 2011, 62, 575. 2. Dlott, D. D.; Fayer, M. D., Shocked Molecular Solids: Vibrational up Pumping, Defect Hot Spot Formation, and the Onset of Chemistry. J. Chem. Phys. 1990, 92, 3798. 3. Tokmakoff, A.; Fayer, M. D.; Dlott, D. D., Chemical Reaction Initiation and Hot-Spot Formation in Shocked Energetic Molecular Materials. J. Phys. Chem. 1993, 97, 1901. 4. Kim, H.; Dlott, D. D., Molecular Dynamics Simulation of Nanoscale Thermal Conduction and Vibrational Cooling in a Crystalline Naphthalene Cluster. J. Chem. Phys. 1991, 94, 8203. 5. Kim, H.; Dlott, D. D., Theory of Ultrahot Molecular-Solids - Vibrational Cooling and ShockInduced Multiphonon up Pumping in Crystalline Naphthalene. J. Chem. Phys. 1990, 93, 1696-1709. 6. Kim, H.; Dlott, D. D., Molecular-Dynamics Simulation of Nanoscale Thermal Conduction and Vibrational Cooling in a Crystalline Naphthalene Cluster. J. Chem. Phys. 1991, 94, 8203-8209. 7. Fried, L. E.; Ruggiero, A. J., Energy-Transfer Rates in Primary, Secondary, and Insensitive Explosives. J. Phys. Chem. 1994, 98, 9786-9791. 8. Dawes, R.; Siavosh-Haghighi, A.; Sewell, T. D.; Thompson, D. L., Shock-Induced Melting of (100)-Oriented Nitromethane: Energy Partitioning and Vibrational Mode Heating. J. Chem. Phys. 2009, 131, 224513. 9. Boyd, S. G.; Boyd, K. J., A Computational Analysis of the Interaction of Lattice and Intramolecular Vibrational Modes in Crystalline Alpha-Rdx. J. Chem. Phys. 2008, 129, 134502. 10. Kraczek, B.; Chung, P. W., Investigation of Direct and Indirect Phonon-Mediated Bond Excitation in Alpha-Rdx. J. Chem. Phys. 2013, 138. 11. Pereverzev, A.; Sewell, T. D.; Thompson, D. L., Molecular Dynamics Study of the PressureDependent Terahertz Infrared Absorption Spectrum of Alpha- and Gamma-Rdx. J. Chem. Phys. 2013, 139. 12. Joshi, K.; Chaudhuri, S., Reactive Simulation of the Chemistry Behind the Condensed-Phase Ignition of Rdx from Hot Spots. Phys. Chem. Chem. Phys. 2015, 17, 18790-18801. 13. McGrane, S. D.; Barber, J.; Quenneville, J., Anharmonic Vibrational Properties of Explosives from Temperature-Dependent Raman. J. Phys. Chem. A 2005, 109, 9919-9927. 14. Nitzan, A.; Jortner, J., Electronic Relaxation of Small Molecules in a Dense Medium. Theor. Chim. Acta 1973, 29, 97-116. 15. Nitzan, A.; Jortner, J., Vibrational Relaxation of a Molecule in a Dense Medium. Mol. Phys. 1973, 25, 713-734.

32


Page 33 of 35

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


16. Nitzan, A.; Mukamel, S.; Jortner, J., Some Features of Vibrational-Relaxation of a Diatomic Molecule in a Dense Medium. J. Chem. Phys. 1974, 60, 3929-3934. 17. Nitzan, A.; Mukamel, S.; Jortner, J., Energy-Gap Law for Vibrational-Relaxation of a Molecule in a Dense Medium. J. Chem. Phys. 1975, 63, 200-207. 18. Mader, C., L., Numerical Modeling of Explosives and Propellants. Third ed.; CRC Press: Boca Raton, FL, 2007. 19. Smith, G. D.; Bharadwaj, R. K., Quantum Chemistry Based Force Field for Simulations of Hmx. J. Phys. Chem. B 1999, 103, 3570. 20. Wallis, E. P.; Thompson, D. L., Molecular Dynamics Simulations of Ring Inversion in Rdx. J. Chem. Phys. 1993, 99, 2661. 21. Sewell, T. D.; Thompson, D. L., Classical Dynamics Study of Unimolecular Dissociation of Hexahydro-1,3,5-Trinitro-1,3,5-Triazine (Rdx). J. Phys. Chem. 1991, 95, 6228. 22. Bedrov, D.; Chakravarthy, A.; Smith, G. D.; Sewell, T. D.; Menikoff, R.; Zaug, J., Molecular Dynamics Simulations of Hmx Crystal Polymorphs Using a Flexible Molecule Force Field. J. ComputerAided Mat. Design 2001, 8, 77. 23. Plimpton, S.; Pollock, R.; Stevens, M. In Particle-Mesh Ewald and Rrespa for Parallel Molecular Dynamics Simulations, Eighth SIAM Conference on Parallel Processing for Scientific Computing, Minneapolis, MN, Minneapolis, MN, 1997. 24. Wallis, E. P.; Thompson, D. L., Molecular-Dynamics Simulations of Ring Inversion in Rdx. J. Chem. Phys. 1993, 99, 2661-2673. 25. Plimpton, S., Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Compt. Phys. 1995, 117, 1. 26. Daw, M. S.; Baskes, M. I., Embedded-Atom Method - Derivation and Application to Impurities, Surfaces, and Other Defects in Metals. Phys. Rev. B: Condens. Matter 1984, 29, 6443-6453. 27. Martinez, L.; Andrade, R.; Birgin, E. G.; Martinez, J. M., Packmol: A Package for Building Initial Configurations for Molecular Dynamics Simulations. J. Comput. Chem. 2009, 30, 2157-2164. 28. 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. 29. Larentzos, J. P.; Rice, B. M.; Byrd, E. F. C.; Weingarten, N. S.; Lill, J. V., Parameterizing Complex Reactive Force Fields Using Multiple Objective Evolutionary Strategies (Moes). Part 1: Reaxff Models for Cyclotrimethylene Trinitramine (Rdx) and 1,1-Diamino-2,2-Dinitroethene (Fox-7). J. Chem. Theory Comput. 2015, 11, 381-391. 30. Rice, B. M.; Larentzos, J. P.; Byrd, E. F. C.; Weingarten, N. S., Parameterizing Complex Reactive Force Fields Using Multiple Objective Evolutionary Strategies (Moes): Part 2: Transferability of Reaxff Models to C-H-N-O Energetic Materials. J. Chem. Theory Comput. 2015, 11, 392-405. 31. An, Q.; Goddard, W. A.; Zybin, S. V.; Luo, S. N., Inhibition of Hotspot Formation in Polymer Bonded Explosives Using an Interface Matching Low Density Polymer Coating at the Polymer-Explosive Interface. J. Phys. Chem. C 2014, 118, 19918-19928. 32. Sergeev, O. V.; Yanilkin, A. V., Molecular Dynamics Simulation of Combustion Front Propagation in a Petn Single Crystal. Combust. Explos. Shock Waves 2014, 50, 323-332. 33. Guo, F.; Zhang, H.; Hu, H. Q.; Cheng, X. L., Effects of Hydrogen Bonds on Solid State Tatb, Rdx, and Datb under High Pressures. Chinese Physics B 2014, 23, 046501. 34. Katz, G.; Zybin, S.; Goddard, W. A.; Zeiri, Y.; Kosloff, R., Direct Md Simulations of Terahertz Absorption and 2d Spectroscopy Applied to Explosive Crystals. J. Phys. Chem. Lett. 2014, 5, 772-776. 35. Zhang, L.; Chen, L., The Effect of Pressure on Thermal Decomposition of Solid Nitromethane Via Md Simulation. Acta Phys. Sin-Ch. Ed 2013, 62, 138201.

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

Page 34 of 35

36. Liu, L. C.; Liu, Y.; Zybin, S. V.; Sun, H.; Goddard, W. A., Reaxff-/G: Correction of the Reaxff Reactive Force Field for London Dispersion, with Applications to the Equations of State for Energetic Materials. J. Phys. Chem. A 2011, 115, 11016-11022. 37. Zhang, L. Z.; Zybin, S. V.; Van Duin, A. C. T.; Goddard, W. A., Modeling High Rate Impact Sensitivity of Perfect Rdx and Hmx Crystals by Reaxff Reactive Dynamics. J. Energ. Mater. 2010, 28, 92-127. 38. Nomura, K. I.; Kalia, R. K.; Nakano, A.; Vashishta, P.; van Duin, A. C. T.; Goddard, W. A., Dynamic Transition in the Structure of an Energetic Crystal During Chemical Reactions at Shock Front Prior to Detonation. Phys. Rev. Lett. 2007, 99, 148303. 39. Strachan, A.; Kober, E. M.; van Duin, A. C. T.; Oxgaard, J.; Goddard, W. A., Thermal Decomposition of Rdx from Reactive Molecular Dynamics. J. Chem. Phys. 2005, 122, 054502. 40. Strachan, A.; van Duin, A. C. T.; Chakraborty, D.; Dasgupta, S.; Goddard, W. A., Shock Waves in High-Energy Materials: The Initial Chemical Events in Nitramine Rdx. Phys. Rev. Lett. 2003, 91, 098301. 41. Wood, M. A.; van Duin, A. C. T.; Strachan, A., Coupled Thermal and Electromagnetic Induced Decomposition in the Molecular Explosive Alpha Hmx; a Reactive Molecular Dynamics Study. J. Phys. Chem. A 2014, 118, 885-895. 42. Kamat, A. M.; van Duin, A. C. T.; Yakovlev, A., Molecular Dynamics Simulations of LaserInduced Incandescence of Soot Using an Extended Reaxff Reactive Force Field. J. Phys. Chem. A 2010, 114, 12561-12572. 43. Zheng, L.; Thompson, D. L., Molecular Dynamics Simulations of Melting of Perfect Crystalline Hexahydro-1,3,5-Trinitro-1,3,5-S-Triazine. J. Chem. Phys. 2006, 125, 084505. 44. Dreger, Z. A.; Gupta, Y. M., Phase Diagram of Hexahydro-1,3,5-Trinitro-1,3,5-Triazine Crystals at High Pressures and Temperatures. J. Phys. Chem. A 2010, 114, 8099. 45. Zhao, X.; Hintsa, E. J.; Lee, Y. T., Infrared Multiphoton Dissociation of Rdx in a Molecular Beam. J. Chem. Phys. 1988, 88, 801. 46. Chakraborty, D.; Muller, R. P.; Dasgupta, S.; W. A. Goddard, I., The Mechanism for Unimolecular Decomposition of Rdx (1,3,5-Trinitro-1,3,5-Triazine), an Ab Initio Study. J. Phys. Chem. A 2000, 104, 2261. 47. Patterson, J. E.; Dreger, Z. A.; Miao, M.; Gupta, Y. M., Shock Wave Induced Decomposition of Rdx: Time-Resolved Spectroscopy. J. Phys. Chem. A 2008, 112, 7374. 48. Izvekov, S.; Chung, P. W.; Rice, B. M., Non-Equilibrium Molecular Dynamics Simulation Study of Heat Transport in Hexahydro-1,3,5-Trinitro-S-Triazine (Rdx). Int. J. Heat Mass Transfer 2011, 54, 5623-5632. 49. Dammak, H.; Chalopin, Y.; Laroche, M.; Hayoun, M.; Greffet, J.-J., Quantum Thermal Bath for Molecular Dynamics Simulation. Phys. Rev. Lett. 2009, 103. 50. Basire, M.; Borgis, D.; Vuilleumier, R., Computing Wigner Distributions and Time Correlation Functions Using the Quantum Thermal Bath Method: Application to Proton Transfer Spectroscopy. Phys. Chem. Chem. Phys. 2013, 15, 12591-601. 51. Hernández-Rojas, J.; Calvo, F.; Noya, E. G., Applicability of Quantum Thermal Baths to Complex Many-Body Systems with Various Degrees of Anharmonicity. J. Chem. Theory Comput. 2015, 11, 861-870. 52. Losada, M.; Chaudhuri, S., Theoretical Study of Elementary Steps in the Reactions between Aluminum and Teflon Fragments under Combustive Environments. J. Phys. Chem. A 2009, 113, 593341. 53. Losada, M.; Chaudhuri, S., Finite Size Effects on Aluminum/Teflon Reaction Channels under Combustive Environment: A Rice–Ramsperger–Kassel–Marcus and Transition State Theory Study of Fluorination. J. Chem. Phys. 2010, 133, 134305.

34


Page 35 of 35

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


35


Intermolecular Energy Transfer Dynamics at a Hot-Spot Interface in

Recommend Documents