Effects of Defects on Thermal Decomposition of HMX via ReaxFF

Dec 13, 2010 - Phone: +86-10-68914518. .... The first principles-based ReaxFF describes the system energy in various partial ...... octane [101347-48-...
0 downloads 0 Views 5MB Size
Download PDF
278

J. Phys. Chem. B 2011, 115, 278–287

Effects of Defects on Thermal Decomposition of HMX via ReaxFF Molecular Dynamics Simulations Ting-Ting Zhou and Feng-Lei Huang* State Key Laboratory of Explosion Science and Technology, Beijing Institute of Technology, Beijing, 100081, People’s Republic of China ReceiVed: June 23, 2010; ReVised Manuscript ReceiVed: NoVember 11, 2010

Effects of molecular vacancies on the decomposition mechanisms and reaction dynamics of condensed-phase β-HMX at various temperatures were studied using ReaxFF molecular dynamics simulations. Results show that three primary initial decomposition mechanisms, namely, N-NO2 bond dissociation, HONO elimination, and concerted ring fission, exist at both high and lower temperatures. The contribution of the three mechanisms to the initial decomposition of HMX is influenced by molecular vacancies, and the effects vary with temperature. At high temperature (2500 K), molecular vacancies remarkably promote N-N bond cleavage and concerted ring breaking but hinder HONO formation. N-N bond dissociation and HONO elimination are two primary competing reaction mechanisms, and the former is dominant in the initial decomposition. Concerted ring breaking of condensed-phase HMX is not favored at high temperature. At lower temperature (1500 K), the most preferential initial decomposition pathway is N-N bond dissociation followed by the formation of NO3 (O migration), although all three mechanisms are promoted by molecular vacancies. The promotion effect on concerted ring breaking is considerable at lower temperature. Products resulting from concerted ring breaking appear in the defective system but not in the perfect crystal. The mechanism of HONO elimination is less important at lower temperature. We also estimated the reaction rate constant and activation barriers of initial decomposition with different vacancy concentrations. Molecular vacancies accelerate the decomposition of condensed-phase HMX by increasing the reaction rate constant and reducing activation barriers. 1. Introduction As one of the most important energetic materials, HMX is widely used in many fields, such as explosives, rocket propellants, and airbag inflators, to name a few. It releases large amounts of energy through a complicated bulk decomposition process that involves unimolecular and bimolecular reactions. Understanding the underlying complex chemical processes is essential to obtain an improved model for combustion and detonation. Applying density functional theory, Lewis et al. presented four possible decomposition reaction pathways of R-HMX in the gas phase:1 N-NO2 bond dissociation, HONO elimination, C-N bond scission of the ring, and concerted ring fission. The dissociation of the N-NO2 bond is putatively the initial mechanism of nitramine decomposition in the gas phase. Lewis et al. also proposed that the HONO elimination and C-N bond scission reaction of the ring would be favorable in the condensed phase. Chakraborty suggested three distinct mechanisms for unimolecular decomposition of HMX:2 hemolytic cleavage of the N-N bond and subsequent decomposition, successive HONO elimination (energetically, the most favorable pathway), and formation of CH2O and N2O from the secondary decomposition of CH2N2O2 (methylenentriamine, or MN). Several plausible mechanisms have been suggested for the condensed-phase decomposition of HMX, and most experimental studies support the formation of CH2O, N2O, HONO, and CHN as intermediate products. These are further decomposed to gaseous products, such as CO, N2, H2O, and CO2. Brill3 proposed two competing global reactions: HMX f 4(HONO * Corresponding author. Phone: +86-10-68914518. Fax: +86-1068461702. E-mail: [email protected].

+ HCN) and HMX f 4(CH2O + N2O). Tang4 concluded that a multiple-step reaction at the condensed phase might explain experimental data more realistically. McGuire and Tarver developed a three-step decomposition model of HMX,5,6 corresponding to the unimolecular endothermic decomposition of the HMX molecule by breaking N-N and C-N bonds, weakly exothermic unimolecular decomposition of primary fragments into intermediates (such as H2CO, N2O, HONO, and HCN), and highly exothermic gas-phase reactions producing final and stable products (such as H2O, N2, CO2, and CO). Zhang7 suggested that the rupture of the N-NO2 bond accompanied by HONO elimination dominates the initial steps of HMX dissociation. Although a number of experiments and numerical simulations8,9 have been conducted to investigate the chemical and physical properties of HMX decomposition, little progress has been made in clarifying atomistic details about the complex phenomena occurring in condensed-phase chemistry under extreme conditions. To the best of our knowledge, minimal work has been carried out to investigate the influence of defects on the decomposition mechanisms of energetic materials (EM). Micromolecular structures of organic explosives are closely related to their macroscopic properties, such as sensitivities and chemical reaction. The initiation of EM is associated with the “hot spot”,10,11 presumed to be local regions of a crystal (crystal defects or deformations) in which the energy of shock compression is trapped in such regions, leading to chemical chain reactions of decomposition. Molecular vacancies, impurities, voids, pores, dislocations, and other types of defects, which play a crucial role in initial reactions, are always present in real crystals.12 However, types of defects that are responsible for the sensitivity of an explosive to the initiation of combustion

10.1021/jp105805w  2011 American Chemical Society Published on Web 12/13/2010

Effects of Defects on Thermal Decomposition of HMX and detonation have not been established. Thus, studying the physical and chemical properties of EM containing defects is important in many practical applications, representing a challenge to experimental and theoretical advancement.13–15 First-principles calculations are particularly useful in elucidating the structural and electronic properties of these materials. Using quantum chemical methods, Kuklja and Kunz conducted investigations on the effects of defects, such as vacancy complexes,11 a single vacancy,16,17 edge dislocation,18–20 and others,21–23 on the optical band gap and electronic properties of EM. Their calculations show that the compression of RDX crystal in the presence of defects appreciably reduces the optical band gap of this material. Manaa studied the effects of high pressure and molecular vacancies on the electronic structure of solid nitromethane.24 Despite these efforts, however, only a few first-principles studies are available. One of the reasons for this paucity is the complicated structures of molecular crystals and the presence of defects, as well as related deformations, removing any symmetry of the crystal structure and resulting in a considerable increase in computer memory resources required to run periodic ab initio calculations.25 Molecular dynamics (MD) simulations can reveal changes in the atomic and electronic structure that may not be captured by static calculations.26 Several theoretical studies with respect to vacancies or voids have been performed using this method.27–29 However, the complexity of the reactions in EM with numerous intermediates interacting simultaneously under extreme conditions makes it difficult to obtain a detailed atomistic understanding of the chemical processes. Classic force fields can describe the energy, structures, and vibrations of molecular systems but are unable to describe chemical reactions, limiting their range of applicability. Recently, Manaa et al. used MD with atomic interactions that were obtained using quantum mechanics (QM) to study the decomposition of HMX at a temperature and pressure close to those in CJ points.30 However, QM methods are computationally intensive, making them impractical for the study of large systems. Using the concept of partial bond order as the basis, this problem was solved with the development of the reactive force field (ReaxFF).7,31,32 ReaxFF enables MD with computational costs that are nearly as low as those incurred in simple force fields. This facilitates the computational study of large systems for many nanoseconds or longer while maintaining near QM accuracy. ReaxFF reaction dynamics is practical for studying high-temperature and highpressure MD of realistic, chemically reacting EM systems and can provide detailed information on the atomistic mechanism of the chemical reactions taking place during decomposition and subsequent reactions under extreme conditions. It has been applied successfully to the thermal decomposition of RDX and HMX.7,33,34 van Duin et al.35 studied the role of voids in the initial chemical events in RDX under shock loading using ReaxFF MD. They concluded that even small defects can play a key role in reducing the detonation threshold. We performed ReaxFF MD simulations to study the effects of molecular vacancies on the initial decomposition mechanisms, together with reaction kinetics during the decomposition of condensed-phase HMX at various temperatures. We anticipate that our simulations will assist in understanding the relationship between microstructure and macrochemical properties of solid explosives, as well as the evolution of a “hot spot”. The remainder of this paper is organized as follows: in Section 2, ReaxFF and modeling methods are briefly described. Results and analyses are presented in Section 3. Conclusions are drawn and presented in Section 4.

J. Phys. Chem. B, Vol. 115, No. 2, 2011 279 2. Modeling Method The first principles-based ReaxFF describes the system energy in various partial energy contributions,7 expressed as

Esystem ) Ebond + Elp + Eover + Eunder + Eval + Epen + Etors + Econj + EH-bond + EvdWaals + ECoulomb It can be divided into three parts: covalent interactions (bonds, angle, and torsion) based on the concept of bond order; coulomb forces that allow charge transfer between atoms; and van der Waals forces calculated between all pairs of atoms. As the most stable phase at room temperature, β-HMX was studied to investigate the role of molecular vacancies in the decomposition of HMX at various temperatures (1500, 2000, 2500, and 3000 K) using ReaxFF (version 2.0) molecular dynamics. The ReaxFF parameters have been used by Zhang et al.7 The initial structure of β-HMX was taken from the Neutron Diffraction experiment.36 Using periodic defects in the form of large supercells, vacancy concentrations of 0%, 4.27%, 8.33%, and 12.5% were modeled by the removal of 0, 2, 4, and 6 molecules from four 4 × 2 × 3 supercells containing 1344, 1288, 1232, and 1176 atoms, respectively. The removed molecules are situated far from each other to diminish the interaction among molecular vacancies, as shown in Figure 1. After optimizing the atomic positions and cell parameters to obtain the minimum energy at 0 K, these systems were heated gradually from 0 K to the target temperatures. Thereafter, isothermal-isochoric MD (NVT-MD) simulations were carried out for 25 ps at these temperatures. Temperature was controlled using a Berendsen thermostat with a damping constant of 80 fs on each system with different vacancy concentrations. The time step was set to 0.1 fs to correctly describe the decomposition process. Analysis of products formed during the simulation was performed with a 0.3 bond-order cutoff for all atom pairs to recognize molecular species. Molecular species and their compositions were recorded every 50 fs, and these data were then employed to plot time evolution figures for the products, providing detailed information on initial chemical reactions during the decomposition process. 3. Results and Analysis 3.1. Effects of Vacancy Concentrations on the Decomposition of Condensed-Phase HMX at High Temperature (2500 K). The decomposition process of perfect crystals together with the system with a vacancy concentration of 4.17% begins as an endothermic reaction of molecular breakdown. This lasts 0.2 ps to reach the maximum of the potential energy (Figure 2) until the secondary reactions are initiated in the partially decomposed solid HMX, which releases energy because of the exothermic formation of small molecules. However, the potential energy decreases from the beginning to the time scale of our simulation for vacancy concentrations of 8.33% and 12.5%, suggesting that the decomposition process begins as an exothermic reaction. The rate of potential energy decrement depends on vacancy concentration; that is, the rate reduces more quickly with elevated vacancy concentration, indicating that molecular vacancies accelerate HMX decomposition. The reason for this phenomenon can be explained by the formation of a “hot spot”, which is capable of triggering a chemical reaction. The decomposition of molecules is more likely to occur around defectssa phenomenon attributed to a shorter endothermic time and a higher HMX dissociation rate.

280

J. Phys. Chem. B, Vol. 115, No. 2, 2011

Zhou and Huang

Figure 1. 4 × 2 × 3 supercells of β-HMX with 0, 2, 4, and 6 molecular vacancies, identified by open circles, corresponding to vacancy concentrations of 0% (a); 4.17% (b); 8.33% (c); and 12.5% (d).

Figure 2. Evolution of potential energy per molecule in NVT-MD simulations for vacancy concentrations of 0%, 4.17%, 8.33%, and 12.5% at 2500 K; for clarity, the inset shows data at the early stage of decomposition.

To evaluate the degree of HMX decomposition for vacancy concentrations of 0%, 4.17%, 8.33%, and 12.5% at 2500 K, the remaining HMX and main products are plotted in Figure 3. Only species with a relatively high abundance are shown in this figure. These simulations reveal that the rupture of the N-NO2 bond accompanied by HONO elimination dominates the initial stage of HMX dissociation. The initial decomposition products NO2 and HONO are formed faster and in larger quantities than other products at the early stage. Previous studies on HMX decomposition pathways2,37 have elucidated that the first dissociated fragment from solid HMX is NO2 for its lower reaction barrier. Oyumi and Brill38–40 pointed out that the rupture

of the N-NO2 bond is the key step in the decomposition process of HMX. Yuji Kohno et al.41 found that only the N-NO2 bond is lengthened considerably by the electronic correlation effect, ascribed to the antibonding character of the N-NO2 bond in the lowest unoccupied molecular orbital. After reaching the maximum, NO2 begins to reduce dramatically because of secondary reactions that consume NO2. As the dominant final products of HMX decomposition, H2O and N2 appear shortly after the initial dissociation of NO2. In the time scale of the simulation, the amount of H2O and N2 rises monotonically with increasing time, and the rate of H2O formation is remarkably faster than that of N2. In addition, NO, HNO3, and CO2, for example, are observed in the decomposition process. We compared the evolution of HMX, dissociated NO2, HONO, CH2O, N2O, CH2N, and primary final products H2O and N2 for different vacancy concentrations (Figure 4(a)-(h)). The rate of HMX decomposition increases rapidly with the rise in vacancy concentration. Complete HMX decomposition takes about 1, 0.85, 0.80, and 0.75 ps for vacancy concentrations of 0%, 4.33%, 8.17%, and 12.5%, respectively. Figure 4(b) reveals the evolution of NO2, which can be divided into two stages: the dramatic growth from the beginning of the simulation and the reduction until an equilibrated value is reached. Although change tendency is not affected by vacancy concentration, the entire evolution is accelerated with the increase in vacancy concentration. It takes almost the same time (about 1 ps) for NO2 to reach the maximum; however, the maximum values differ by about 1.75, 1.75, 1.9, and 2.1 for vacancy concentrations of 0%, 4.33%, 8.17%, and 12.5%, respectively, implying that molecular vacancies accelerate the rupture of the N-NO2 bond. After reaching the maximum, NO2

Effects of Defects on Thermal Decomposition of HMX

J. Phys. Chem. B, Vol. 115, No. 2, 2011 281

Figure 3. Evolution of products per molecule in NVT-MD simulations for vacancy concentrations of 0% (a), 4.17% (b), 8.33% (c), and 12.5% (d) at 2500 K.

reduction begins. Reduction takes nearly the same time to decrease from the maximum to 0.50, indicating that the quantity of NO2 molecules reduces much more rapidly under higher concentration. This phenomenon can be explained by the increase in secondary reactions at higher vacancy concentration, resulting in a rapid consumption of NO2 molecules and a faster formation of intermediates and stable products. Notably, the evolution of HONO is delayed by molecular vacancies, which is contrary to that of NO2. Figure 4(c) shows that HONO begins to appear at about 0.1 ps in all the systems with different vacancy concentrations. The amount of HONO reaches almost the same maximum (about 0.55) at about 1.8, 2.5, 3.3, and 4 ps for vacancy concentrations of 0%, 4.17%, 8.33%, and 12.5%, respectively, suggesting that molecular vacancies effectively slow down the evolution of HONO. This process involves the breakage of both N-NO2 and C-H bonds caused by the short nonbonded N-O · · · H-C distance.2 We propose that intermolecular reaction can partially contribute to the formation of HONO. If O and H come from the same molecule, the likelihood of O in the nitro group depriving H in the methyl group cannot be influenced by molecular vacancies because these do not affect the intramolecular nonbonded N-O · · · H-C distance. To elucidate this speculation, we performed 50 ps NVT-MD simulation on a unit cell at 2000 K, which shows intermolecular transfer of H from -CH2 to -NO2. At 11.5 ps, the population of HONO in the four systems is reduced to an equilibrated number (about 0.25), maintained for 8.5 ps, and then decreased at a very slow speed. This indicates that the dissociation of HONO is a reversible process. HONO

can easily decompose into HO, which contributes to the formation of H2O. As a result of the delayed evolution of HONO, the formation of H2O (Figure 4(d)) slows down, and the population at higher vacancy concentrations is smaller. Concerted ring breaking of condensed-phase HMX is promoted by molecular vacancies. A minimal amount of CH2N2O2 (MN) is observed in the decomposition process of perfect HMX crystals. However, it occurs at about 0.7, 0.35, and 0.15 ps for vacancy concentrations of 4.17%, 8.33%, and 12.5%, respectively. The population of MN increases with the rise in concentration, showing that molecular vacancies advance concerted ring breaking of HMX. Suryanarayana42 proposed a concerted decomposition mechanism of HMX to four MN molecules, which can further decompose into CH2O and N2O. Morgan and Bayer,43 in their pyrolysis study using electron spin resonance, detected CH2N and NO2 radicals as decomposition products from MN intermediates. In our simulations, we also observed CH2O, N2O, and CH2N, which are more likely to appear at higher vacancy concentrations (Figure 4(e)-(g)). This agreement reveals that defects exist in real crystals. CH2N2O also appears at the early stage of the HMX decomposition process and increases with the rise in vacancy concentration. The formation of N2 is acceleratedsa result ascribed to the promoted concerted ring breaking of HMX by molecular vacancies (Figure 4(h)). In conclusion, we propose that at high temperature, the effects of molecular vacancies on the three mechanisms differ despite the existence of the three primary initial decomposition mechanisms of the condensed-phase HMX, namely, N-NO2 bond

282

J. Phys. Chem. B, Vol. 115, No. 2, 2011

Zhou and Huang

Figure 4. Evolution of HMX (a), dissociated NO2 (b), HONO (c), H2O (d), CH2O (e), N2O (f), CH2N (g), and N2 (h) per molecule for vacancy concentrations of 0%, 4.17%, 8.33%, and 12.5% at 2500 K.

dissociation, HONO elimination, and concerted ring fission followed by C-N bond scission. Molecular vacancies notably promote N-N bond cleavage and concerted ring breaking, whereas they hinder the formation of HONO. The population

of MN is rather small compared with that of NO2 and HONO, indicating that the ring breaking of condensed-phase HMX is not preferential in spite of the promotion by molecular vacancies. Thus, we propose that N-N bond dissociation and HONO

Effects of Defects on Thermal Decomposition of HMX

J. Phys. Chem. B, Vol. 115, No. 2, 2011 283

Figure 5. Evolution of potential energy (a) and the number of remaining HMX (b) per molecule in NVT-MD simulations for vacancy concentrations of 0% and 12.5% at 1500 K.

elimination are two primary competing reaction mechanisms, with the former playing a leading role. 3.2. Effects of Vacancy Concentrations on the Decomposition of Condensed-Phase HMX at Lower Temperature (1500 K). To elucidate the influence of molecular vacancy on the decomposition of HMX at lower temperature, we performed two simulations at 1500 K. The potential energy (Figure 5(a)) initially increases in the endothermic process of molecule breakdown until the secondary reactions are initiated in partially decomposed solid HMX. At this point, HMX begins to release energy because of the exothermic formation of small fragments. The endothermic stage lasts about 16 and 7.5 ps for vacancy concentrations of 0% and 12.5%, respectively, demonstrating that molecular vacancies dramatically accelerate the decomposition of HMX. The evolution of the remaining HMX molecules (Figure 5(b)) is strong evidence. It takes 18.5 ps for HMX at a vacancy concentration of 12.5% to completely decompose; however, some undissociated HMX molecules are still observed in the perfect crystal at the end of the simulation. Figure 6 displays the time evolution of the total population of fragments per single HMX molecule produced in the two systems at 1500 K. The populations of H2O, N2, OH, and NO are extremely small; H, CH2O, as well as CO2 are not found in any of the systems, and the amounts of NO2, HONO, and HNO3 are less than those at 2500 K, demonstrating the dramatic promoting effect of temperature on the decomposition of HMX. CH2O2N2, CH2N, or CHN was not observed in the decomposition process of perfect HMX in the time scale of our simulations, but these species appear in the system with the vacancy concentration of 12.5%, elucidating that molecular vacancies strongly promote the ring breaking of HMX. Under molecular vacancies, the population of CH2O2N2 together with CH2N is comparable to that of HONO, indicating that concerted ring breaking competes with HONO elimination at 1500 K. This is quite different from that at 2500 K. Figure 6 shows that NO2 dominates the initial decomposition process and highly surpasses other products, suggesting that the rupture of N-NO2 is the most preferential pathway. The evolution of the dissociated NO2 molecules for vacancy concentrations of 0% and 12.5% (Figure 7(a)) shows that the rate of NO2 formation in the defective crystal rises more quickly than that in the perfect crystal, revealing that molecular vacancies accelerate N-N bond dissociation. After reaching the maximum (about 1 at 19 ps), the population of NO2 begins to

Figure 6. Evolution of products per molecule in NVT-MD simulations for vacancy concentrations of 0% (a) and 12.5% (b) at 1500 K.

decrease in the defective crystal, whereas the number of NO2 in the perfect crystal continues to increase at the end of the simulation. Figure 7(b) shows that HONO begins to appear at about 11 and 0.5 ps for vacancy concentrations of 0% and 12.5%, respectively. HONO increases more quickly, and the amount is greater under molecular vacancies than that in the perfect crystal; thus, molecular vacancies promote the formation of HONO at 1500 K, which disagrees with the results at 2500 K. The population of HONO remains less than that of NO3 at 1500 K, revealing that at lower temperatures, it is difficult for H in the methyl group to migrate to O in the nitro group, resulting in less HONO formation.

284

J. Phys. Chem. B, Vol. 115, No. 2, 2011

Zhou and Huang

Figure 7. Evolution of dissociated NO2 (a) and HONO (b) per molecule for vacancy concentrations of 0% and 12.5% at 1500 K.

The effects of molecular vacancies on the initial decomposition of HMX are more obvious at lower temperature. Figure 8 shows that the acceleration effect of molecular vacancies on the evolutions of HMX, NO2, as well as HONO is more obvious at 1500 K than at 2500 K. Apart from molecular vacancies, temperature also plays a crucial role in the decomposition process of HMX. HMX decomposes more rapidly, and the population of NO2 together with HONO is considerably larger at 2500 K than at 1500 K (Figure 8). These phenomena show that the effects of molecular vacancies and temperature are two competing factors that strongly influence the decomposition process of HMX. The former wield a more important influence at lower temperature. At lower temperature, three initial decomposition mechanisms appear, with N-N bond dissociation as the most preferential initial decomposition pathway followed by the formation of NO3 (O migration) instead of HONO elimination (H migration). The effects of molecular vacancies at lower and high temperatures differ. Molecular vacancies not only promote N-N bond dissociation and concerted ring breaking but also accelerate HONO elimination. However, the population of HONO is still small compared with those of NO2 and NO3, indicating that the mechanism of HONO elimination is less important at lower temperature. In addition, the promotion effect of molecular vacancies on concerted ring breaking is quite evident. Fragments resulting from concerted ring breaking are not observed in the perfect crystal but appear in the defective one in a quantity comparable with that of HONO. This shows that the mechanism of concerted ring breaking becomes more important at lower temperature, in contrast to that at high temperature. 3.3. Effect of Molecular Vacancies on the Reaction Kinetics of HMX Decomposition. To look into the effect of vacancies on the kinetic parameters of the decomposition process of HMX, chemical reactive kinetics was employed. The rate of reaction is often described by the classical Arrhenius equation

dR E ) f(R)k(T) ) f(R)A exp dt RT

( )

(1)

where R is the reaction progress; t is time; f(R) is the reaction model; k(T) is the temperature-dependent rate constant; A is the pre-exponential factor; E is the activation energy; R is the universal gas constant; and T is temperature. The Arrhenius law (linearized form) relates the rate constant to the activation barrier and temperature

ln(k) ) ln(A) -

E RT

(2)

The calculated time evolution profiles for the decomposition of HMX generally follow exponential fitting first-order kinetics. Thus, the first-order expression is applied to evaluate the rate constant

R(t) ) 1 - exp(-kt)

(3)

dR ) k exp(-kt) ) k(1 - R) ) kf(R) dt

(4)

Thus

Using eq 3, we obtain the reaction rate constants at 1500, 2000, 2500, and 3000 K of the initial decomposition of HMX with different vacancy concentrations (Table 1). The reaction rate constant becomes larger with the increase in vacancy concentration, indicating that molecular vacancies can accelerate the decomposition of HMX, especially at lower temperatures. The reaction rate constant is nearly 3-4 times larger for the vacancy concentration of 12.5% than for the perfect crystal when temperature is below 2000 K. The solid line in Figure 9 represents the linear fit based on eq 2, exhibiting typical Arrhenius behavior of initial decomposition of HMX, which allows estimation of the activation barrier and logarithm of the pre-exponential factor. The comparable estimated values of E and ln(A) are listed in Table 2. In the perfect crystal, the calculated ln(A) value is between the values obtained by Gibbs and Henson, and the value of E is smaller than both of these authors’ results. The discrepancy may be explained by the following reasons. The values derived by Gibbs from calorimetry experiments and by Henson from a compilation of ignition times are global, integrating several chemical steps in the decomposition processes, whereas the data obtained from our ReaxFF simulations are from the first step of HMX decomposition (N-N bond dissociation) at 1500-3000 K. Calorimetric experiments are always performed below 600 K, requiring extrapolation of measured rates to higher temperatures, which trigger additional reaction pathways that are not taken into account in these experiments. Moreover, products obtained from calorimetric experiments undergo isentropic expansion, which is prohibited in the constant-volume simulations of the current study. Henson’s results are obtained by fitting the data gathered from confined thermal explosion, fast pyrolysis, laser

Effects of Defects on Thermal Decomposition of HMX

J. Phys. Chem. B, Vol. 115, No. 2, 2011 285

Figure 8. Comparison of evolutions of HMX, NO2, and HONO at different temperatures for vacancy concentrations of 0% and 12.5%.

TABLE 1: Reaction Rate Constant k (ps-1) of HMX Decomposition with Different Vacancy Concentrations at Different Temperatures concentration T (K)

0% (perfect)

1500 2000 2500 3000

0.1158 1.5919 13.3367 25.8911

4.17%

17.0066

8.33%

17.0999

12.5% 0.3588 6.1875 19.1477 31.1487

ignition, impact-induced shear, and frictional heating. Thus, the overall rates may be affected by hot-spot formation and shockinduced heating, which do not manifest in our simulations for the perfect crystal.

The typical Arrhenius behavior of initial decomposition of HMX indicates that the thermal decomposition of solid HMX is controlled by nucleation and growth of reaction sites. Comparing the values between perfect and defective crystals, the activation energy for the vacancy concentration of 12.5% is 5.8 kcal/mol less than that of the perfect crystal, demonstrating that molecular vacancies can remarkably reduce activation barriers and promote chemical reactions. This can be understood from the perspective of the hot-spot formation mechanism, which suggests that it is easier to form hot spots around defects, further triggering chemical reaction and decomposition of molecules. In addition, the above-mentioned analysis shows that molecular vacancies strongly promote N-N bond dissociation, the reaction barrier of which is lower than that of HONO

286

J. Phys. Chem. B, Vol. 115, No. 2, 2011

Zhou and Huang

Figure 9. Logarithm of HMX initial decomposition rate constant (ln(k)) vs inverse temperature (1/T) in the 1500-3000 K range. The square and solid line represent calculated values and linear fitting for the perfect crystal, and the triangle and dashed line represent calculated values and linear fitting for the vacancy concentration of 12.5%. The linear equation with correlation coefficient is also depicted.

TABLE 2: Chemical Reaction Kinetic Parameters for HMX Thermal Decomposition ln A (s-1)

model

perfect crystal

Gibbs44 Henson45 this study

45.36 29.35 36.49

E (kcal/mol) vacancy concentration ) 12.5%

perfect crystal

vacancy concentration ) 12.5%

35.89

52.7 35.6 32.7

26.9

elimination and concerted ring breaking. This may be another possible explanation for the decrement of the activation barrier of condensed-phase HMX decomposition. 4. Conclusions We investigated the effects of molecular vacancies on the decomposition mechanisms along with reaction dynamics of HMX at various temperatures using the ReaxFF MD method. Several conclusions are drawn. (1) Molecular vacancies accelerate the decomposition of HMX, which can be explained by the formation of a hot spot. As one type of hot spot, vacancies are capable of triggering a chemical reaction. The decomposition of molecules is more likely to occur around defects. Real crystals have high concentrations of imperfections even at low temperatures and are more sensitive to initiation than that of perfect crystals.46 Thus, we propose that at room temperature the influence of defects on the storage, safety, and aging of HMX should not be ignored. (2) Although the three initial decomposition mechanisms exist at high temperature, the effects of molecular vacancies on these mechanisms differ. They notably promote N-N bond cleavage and concerted ring breaking but hinder the formation of HONO. N-N bond dissociation and HONO elimination are the primary competing reaction mechanisms, and the former is the dominant mechanism. The population of products resulting from concerted ring breaking is rather small compared with those of NO2 and HONO, indicating that the ring breaking of condensed-phase HMX is not a preferential mechanism. (3) At lower temperature, the most preferential initial decomposition pathway is N-N bond dissociation followed by the formation of NO3 (O migration) instead of HONO elimination (H migration), although all three decomposition mechanisms are promoted by molecular vacancies. Another important

mechanism is concerted ring breaking, which is considerably advanced by molecular vacancies. We did not find concerted ring breaking in the perfect crystal but observed it in the defective crystal in abundance. The mechanism of HONO elimination becomes less crucial at lower temperature. (4) Molecular vacancies increase the reaction rate constant and reduce the activation barrier, accelerating the decomposition process of HMX, especially at lower temperatures. Reaction rate constants for the vacancy concentration of 12.5% are nearly 3-4 times larger than that of the perfect crystal at temperatures below 2000 K. The activation energy for the vacancy concentration of 12.5% is 5.8 kcal/mol less than that of the perfect crystal. The formation of final products is not fully completed within 25 ps of the simulation time, especially for the low-temperature cases. Our study reveals that ReaxFF MD is a prospective tool for the investigation of the effect of defects on EM decomposition. Thus, with suitable advanced computational facilities, chemical processes at the atomistic level of more complicated systems for longer simulation times can be studied. The interactions between defects and defect deformation were disregarded in this study. A previous study revealed that vacancies attract one another and high vacancy associations may also be energetically favored.11 In addition, displacements of the nearest neighbor molecules around the vacancy would lower the total energy of the defective crystal. To model the real defective solids and obtain more accurate results, constructing considerably larger supercells is needed, although this incurs exceedingly higher computational costs. Acknowledgment. This work was supported by the National Natural Science Foundation of China (Grant No. 10832003). We thank Prof. Goddard in Caltech for generously providing the ReaxFF code and helpful discussions. We also appreciate the valuable ideas and suggestions of Dr. Song Huajie and Dr. Shi Yiding. Note Added after ASAP Publication. This paper was published ASAP on December 13, 2010. The Supporting Information description and refs 12, 28, 29, and 35 were revised. The updated paper was reposted on December 20, 2010. Supporting Information Available: Initial structures of β-HMX with vacancy concentrations of 0%, 4.17%, 8.13%, and 12.5% (perfect-original.cif, 4.17-original.cif, 8.33-original.cif, and 12.5-original.cif). This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Lewis, J. P.; Glaesemann, K. R.; VanOpdorp, K.; Voth, G. A. J. Phys. Chem. A 2000, 104, 11384–11389. (2) Chakraborty, D.; Muller, R. P.; Dasgupta, S.; Goddard, W. A., III. J. Phys. Chem. A 2001, 105, 1302–1314. (3) Brill, T. B. J. Propul. Power 1995, 11, 740. (4) Tang, C.-J.; Lee, Y. J.; Litzinger, T. A. J. Propul. Power 1999, 15, 296. (5) McGuire, R. P.; Tarver, C. M. SeVenth Symposium (International) on Detonation; Naval Surface Weapons Center NSWC MP 82-334, Annapolis, MD, 1981; p 56. (6) Tarver, C. M.; Chidester, S. K.; Nichols, A. L. J. Phys. Chem. 1996, 100, 5794–5799. (7) Zhang, L. Z.; Zybin, S. V.; van Duin, A. C. T.; Dasgupta, S.; Goddard, W. A., III. J. Phys. Chem. A 2009, 113, 10619–10640. (8) Zel’dovich, Ya. B.; Kompaneets, A. S. Theory of Detonation; Academic Press: New York, 1960 and references therein. (9) (a) Fickett, W.; Davis, W. C. Detonation; University of Califonia Press: Berkeley, 1979. (b) Fickett, W. Introduction to Detonation Theory; University of Califonia Press: Berkeley, 1985 and references therein.

Effects of Defects on Thermal Decomposition of HMX (10) Bowden, F. P.; Yoffe, Y. D. Internation and Growth of Explosion in Liquids and Solids; Cambridge University Press: London, 1952; pp 6465. (11) Kuklja, M. M.; Kunz, A. B. J. Phys. Chem. B 1999, 103, 8427– 8431. (12) Sharma, J.; Coffey, C. S.; Ramaswamy, A. L.; Armstrong, R. W. In Decomposition, Combustion and Detonation Chemistry of Energetic Materials; Materials Research Society Symposium Proceedings, Vol. 418, 1996; pp 257-264. (13) Stoneham, A. M. Theory of Defects in Solids; Clarendon: Oxford, 1975. (14) Itoh, N. Defect Processes Induced by Electronic Excitation of Insulations; World Scientific Press: Singerpore, 1989; p 278. (15) Ashok, S.; Chevallier, J.; Sopori, B. L.; Tabe, M.; Kiesel, P. Symposium E: Semiconductor Defect Engineering-Materials, Synthetic Structures, and DeVices; Materials Research Society: Warrendale, PA, 2005. (16) Kuklja, M. M.; Kunz, A. B. J. Phys. Chem. Solids 2000, 61, 35– 44. (17) Kuklja, M. M.; Kunz, A. B. J. Chem. Phys. 2000, 112, 3417–3423. (18) Kuklja, M. M.; Kunz, A. B. J. Appl. Phys. 2000, 87, 2215–2218. (19) Kuklja, M. M.; Aduev, B. P.; Aluker, E. D.; Krasheninin, V. L.; Krechetov, A. G.; Mitrofanov, A. Y. J. Appl. Phys. 2001, 89, 4156–4166. (20) Kuklja, M. M.; Kunz, A. B. J. Appl. Phys. 2001, 89, 4962–4970. (21) Kuklja, M. M.; Kunz, A. B. J. Appl. Phys. 1999, 86, 4428–4434. (22) Kunz, A. B.; Kuklja, M. M.; Botcher, T. R.; Russell, T. P. Thermochim. Acta 2002, 384, 279–284. (23) Kuklja, M. M.; Rashkeev, S. N.; Zerilli, F. J. Appl. Phys. Lett 2006, 89, 071904. (24) Elstner, M.; Frauenhelm, T.; Manaa, M. R. J. Chem. Phys. 2002, 117, 788–799. (25) Kuklja, M. M.; Rashkeeev, S. N. Phys. ReV. B 2007, 75, 104111. (26) Reed, E. J.; Joannopoulos, J. D.; Fried, L. E. Phys. ReV. B 2000, 62, 16500–16509. (27) White, C. T.; Barrettm, J. J. C.; Mintmire, J. W.; Elert, M. L.; Robertson, D. H. In Decomposition, Combustion and Detonation Chemistry of Energetic Materials; MRS Symposium Proceedings, Brill, T. B., Russel, T. P., Tao, W. C., Wardle, R. B., Eds.: Materials Research Society: Pittsburgh, PA, 1996; Vol. 418, pp 277-280. (28) (a) Tsai, D. H. In Structure of crystal defects in damaged RDX as reVealed by an AFM; Shock Compression of Condensed Matter, 1997; Schmidt, S. C., Dandekar, D. P., Forbes, J. W., Eds.: American Institute of Physics: Woodbury, NY, 1998, pp 281-286. (b) Tsai, D. H. In Structure and Properties of Energetic Materials; MRS Symposium Proceedings,

J. Phys. Chem. B, Vol. 115, No. 2, 2011 287 Liedenberg, D. H., Armstrong, R. W., Gilman, J. J., Eds.; Materials Research Society: Pittsburgh, PA, 1993; Vol. 296, pp 113-120. (29) Phillips, L.; Oran, E. S. In Shock Compression of Condensed Matter; Proceedings of the American Physical Society Topical Conference; Shmidt, S. C., Dick, R. D., Forbes, J. W., Tasker, D. G., Eds.; American Institute of Physics: North-Holland, Amsterdam, 1992; Vol. 1082, pp 143-146. (30) Manaa, M. R.; Fried, L. E.; Melius, C. F.; Elstner, M.; Frauenheim, T. J. Phys. Chem. A 2002, 106, 9024–9029. (31) van Duin, A. C. T.; Dasgupta, S.; Lorant, F.; Goddard, W. A., III. J. Phys. Chem. A 2001, 105, 9396–9409. (32) Strachan, A.; van Duin, A. C. T.; Chakraborty, D.; Dasgupta, S.; Goddard, W. A., III. Phys. ReV. Lett. 2003, 91, 098301. (33) Strachan, A.; Kober, E.; van Duin, A. C. T.; Oxgaard, J.; Goddard, W. A., III. J. Chem. Phys. 2005, 122, 054502. (34) Nomura, K.; Kalia, R. K.; Nakano, A.; Vashishta, P.; van Duin, A. C. T.; Goddard, W. A., III. Phys. ReV. Lett. 2007, 99, 148303. (35) Strachan, A.; van Duin, A. C. T.; Goddard, W. A., III. Initial Chemical Events in the Energetic Material Rdx under Shock Loading: Role of Defects. Shock Compression of Condensed Matter, 2003; Furnish, M. D.; Forbes, J. W.; Gupta, Y. M., Eds.; American Institute of Physics: Portland, Oregon, 2004; Vol. 706, pp 895-898. (36) Choi, C. S.; Boutin, H. P. Acta Crystallogr. 1970, B26, 1235. (37) Chakraborty, D.; Muller, R. P.; Dasgupta, S.; Goddard, W. A., III. J. Comput.-Aided Mater. Des. 2001, 8, 203. (38) Oyumi, Y.; Brill, T. B. Combust. Flame 1985, 62, 213–224. (39) Brill, T. B.; Oyumi, Y. J. Phys. Chem. 1986, 90, 2679–2682. (40) Oyumi, Y.; Brill, T. B. Propellents, Explos., Pyrotech. 1988, 13 (3), 69–73. (41) Yuji, K.; Kazuyoshi, U.; Akira, I. J. Phys. Chem. 1996, 100, 4701– 4712. (42) Suryanarayana, B.; Graybush, R. J.; Autera, J. R. Chem. Ind. London 1967, 52, 2177. (43) Morgan, C. U.; Bayer, R. A. Combust. Flame 1979, 36, 99–101. (44) Gibbs, T. R.; Popolato, A., Eds. LASA ExplosiVe Property Data; University of California Press: Berkeley, CA, 1980. (45) Henson, B. F.; Smilowitz, L.; Asay, B. W.; Dickson, P M.; Howe, P. M. 12th Symposium (International) on Detonation; ONR 333-05-2, 2002; p 987. (46) Kuklja, M. M. Appl. Phys. A: Mater. Sci. Process. 2003, 76, 359–366.

JP105805W