Influence of Silicon on the Detonation Performance of Energetic

Oct 2, 2018 - ... in civilian and military fields, but there is no clear understanding about how the specific atom type and molecular structure contro...
0 downloads 0 Views 934KB Size
Subscriber access provided by Kaohsiung Medical University

C: Physical Processes in Nanomaterials and Nanostructures

The Influence of Silicon to the Detonation Performance of Energetic Materials from First Principles Molecular Dynamics Simulations Dezhou Guo, Dezhao Guo, Fenglei Huang, and Qi An J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.8b08305 • Publication Date (Web): 02 Oct 2018 Downloaded from http://pubs.acs.org on October 4, 2018

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

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

Page 1 of 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

The Journal of Physical Chemistry

The Influence of Silicon to the Detonation Performance of Energetic Materials from First Principles Molecular Dynamics Simulations Dezhou Guo1,2, Dezhao Guo3, Fenglei Huang1*, and Qi An2,* 1

State Key Laboratory of Explosion Science and Technology, Beijing Institute of Technology,

Beijing 100081, People’s Republic of China. 2

Department of Chemical and Materials Engineering, University of Nevada, Reno, Nevada

89557, United States 3

School of Mechatronical Engineering, Beijing Institute of Technology, Beijing 100081,

People’s Republic of China.

1 Environment ACS Paragon Plus

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

ABSTRACT: The development of new generation energetic materials (EMs) with improved detonation performance is of considerable importance for applications in civilian and military fields, but there is no clear understanding about how the specific atom type and molecular structure control detonating properties. To gain an atomistic-level understanding of the influence of Si element on the detonation properties of EMs, we carried out RxMD(cQM) procedure, combining ReaxFF reactive molecular dynamics with quantum mechanics molecular dynamics, to predict the thermodynamics parameters of Chapman Jouguet (CJ) state as measurements of detonation performance. We find that the detonation temperature for tetrakis(nitratomethyl)silane (Si-PETN) is higher than that for pentaerythritol tetranitrate (PETN) because of high energy release while forming Si products. However, lower detonation pressure and detonation velocity for Si-PETN system than those for PETN system were found because Si atoms attract nearby oxygen atoms from other molecules or fragments resulting in cluster products and leading to less gas products formation. Our results indicate that silicon based energetic materials can keep active in oxygen deficient conditions. This study uncovers how the specific atoms determine the detonation properties of EMs from atomic perspective, providing useful information for designing EMs with improved properties.

2 Environment ACS Paragon Plus

Page 2 of 23

Page 3 of 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

The Journal of Physical Chemistry

INTRODUCTION Energetic materials (EMs) are specific materials which release large amount of energy by fast exothermic decomposing chemical reactions1. Mixing EMs with metal or metalloid elements is a practical method to increase the performance of an explosive formulation2. Usually, silicon is considered as an important element for applications at decoy flairs and irradiation agents, in propulsion engines, and in hydrocarbon fuels.3 In addition, silicon-based polymers were investigated for polymer binders because of such excellent properties as high flexibility at low temperature, low rubber-glass transition, and low glass transition temperature. This character is essential for preventing temperature induced fractures in explosive charges, which can lead to incomplete detonation4. Moreover, although silicon-based materials have such advantages as excellent chemical stability, thermal stability, and biocompatibility, pure Si based compounds are rarely applied in the EMs. The first silicon containing EM, Si-PETN (tetrakis(nitratomethyl)silane), was synthesized recently by Klapotke’s group. Si-PETN is a silicon analogue of the commonly used EM PETN (pentaerythretyl tetranitrate). The Si-PETN can be synthesized by reacting the nitration of hydroxymethylsilane, Si(CH2OH)4, with nitric acid5. The molecular structure of Si-PETN is almost identical to PETN in which the central C atom is replaced by Si atom. Surprisingly, SiPETN exhibits an easy decomposition character and an unexpected extreme sensitivity than PETN. Because of the chemical affinity of silicon and carbon, it is interesting to compare properties of Si-PETN with PETN. Quantum mechanics (QM)5-8 and molecular dynamics (MD) simulations9 were carried out to reveal the molecular and structural factors for dramatically

3 Environment ACS Paragon Plus

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

increased sensitivity of Si-PETN. The initial DFT calculations5 indicated that the tendency to form Si-O bonds at initial decomposition step increases the sensitivity of Si-PETN. More extensive studies by Liu6 concluded that the most likely first step for decomposing Si-PETN is that O atom attacks the central Si atom leading to the break of Si-C bond and to the formation of Si-O bond. This rearrangement reaction is quite high exothermic with low energy barrier (32 kcal/mol), leading to the colossal sensitivity of Si-PETN. Furthermore, Murray7 pointed out that the lower energy barrier of Si-PETN (32 kcal/mol) than PETN (80 kcal/mol) arises from the initial stages of the decomposition process, in which the positive σ-hole on the silicon interacts with the negative linking oxygen, and then the central silicon is easier than carbon to temporarily expand its coordination sphere. In addition, Zhou et al.8,9 suggested that the thermal decomposition mechanisms of Si-PETN arise from the Si–O bond formation including the interaction between Si and O atoms from the decomposed fragments, the Si–CH2–O realignment in the sectionally dissociated fragments, along with the molecular interactions between Si and NO2, and that the intramolecular rearrangement of Si-O bonds which releases large amount of energy to promote further decomposition. As a result, fast temperature increase and reaction zone expansion might happen in Si-PETN. However, there is no clear understanding about the influence of Si element on the detonation properties of EMs. Based on the similar molecular structures with parallel contacts in two crystals but forming dramatic different sensitivities, we expect that there may be big differences of detonation properties between Si-PETN and PETN. The RxMD(cQM)10, combined RxMD with QMD procedure, is a first-principles-based approach to in silico predict the thermochemical parameters at the CJ state as a measure of performance for new materials. To provide valuable information for understanding the nature of detonation properties of Si-PETN and PETN, we first used ReaxFF9 to predict the structure and

4 Environment ACS Paragon Plus

Page 4 of 23

Page 5 of 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

The Journal of Physical Chemistry

thermodynamical properties at high T and P. ReaxFF was optimized by comparison to QM results to depict the decomposition reactions in a fast way. ReaxFF based reactive molecular dynamics (RMD) is nearly as accurate as density functional theory (DFT), but with much lower cost which is like non-reactive force fields11-17. Therefore, it makes the RMD a practical way to investigate molecular systems for as long as nanoseconds, which enables to simulate complicated processes at high temperature and high pressure conditions. Then we follow the RMD trajectories with quantum mechanics dynamics (QMD) simulations for accurately describing the equilibrated fully reacted states to predict the CJ state and properties. Since the simulations can directly output the detailed detonation products at the CJ state without any pre-assumptions of chemical distributions, we believe that this method provides critical information to help understand the difference of detonation mechanisms between condensed Si-PETN and PETN. 2 METHODS AND PROCEDURES 2.1 Simulation Models. The systems examined herein are the molecular crystals of Si-PETN and PETN. The initial unit cell structure of PETN was taken from the Cambridge Structural Database available at the Cambridge; the structure of Si-PETN was made according to the experiment and DFT calculations by substituting the central carbon atom by silicon6. The unit cell was optimized using conjugate gradients method in VASP package. The dispersion effect was considered using the D3 approach18. The optimized cell parameters for PETN are a = b = 9.33 Å, c = 6.63 Å (at 0 K), which leads to a density of 1.81 g/cm3, agreeing very well with the X-ray experimental values of a = b = 9.28 Å, c= 6.61 Å, leading to a density of 1.84 g/cm3 at 100 K19. Thus, the DFT-D3 methods can accurately describe the molecular crystal structures. We predict that the

5 Environment ACS Paragon Plus

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

Page 6 of 23

cell parameters for Si-PETN at 0 K are a = b = 9.57 Å, c = 6.86 Å, leading to a density of 1.76 g/cm3. Then we built a 1×1×2 supercell of each crystal containing 4 molecules (116 atoms) for our further calculations, as shown in Figure 1.

Figure 1. Supercell structure replicated twice alone c direction of (a) PETN containing 4 molecules. (b) Si-PETN containing 4 molecules. The C, N, H, O, and Si atoms are represented by brown, light

blue, white, red, and dark blue balls, respectively. 2.2 Method of RxMD(cQM) to Predict Chapman-Jouguet State The thermodynamic quantities of initial inert state and final reacted states across the detonation wave are connected by the ZND detonation model through the below equations of mass, momentum, and energy conservation:

(1), (2), (3),

6 Environment ACS Paragon Plus

Page 7 of 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

The Journal of Physical Chemistry

where ρ is the density, D is the detonation wave velocity, u is the products velocity behind the detonation wave, e is the specific internal energy, p is the pressure, and v is the specific volume. Here, the “specific” means the quantity per unit mass and the subscript “0” refers to the unshocked state. Based on eqn (1) and (2), the Hugoniot function of eqn (3) can be written as

1 H = e − e0 − ( p + p0 )(v0 − v ) (4). 2

Here, if H = 0, the energy before and after the shock wave are conservative. Therefore, the Hugoniot states of both inert material and fully reacted material are determined by a set of (p,v,T) parameters at H = 0. In our simulations, for the initial un-shocked state, the equilibrium properties are derived by averaging the 10 ps constant volume and constant temperature (NVT) QMD simulation at T= 300 K. For the final reacted state, reaching equilibrium needs hundreds of picoseconds NVT simulations at specific temperatures and densities, which require several orders of magnitude beyond the capability of current QMD simulations. The combination of initial RxMD simulation until the final equilibrated state with an additional 30 ps QMD simulation is a practical way to accurately achieve detonated state10. Thus, we first carried out the NVT-RMD simulation with ReaxFF for 270 ps for PETN system and 150 ps for Si-PENT system until reaching fully decomposition states. Then we extracted the final atoms positions and velocities from ReaxFF as the input to the QMD simulation for another 30 ps to examine the final equilibrium state and get an accurate first-principle description. The internal energy e and the pressure P of the final reacted states are determined by averaging the last 10 ps QMD results. So far, the Hugoniot value is determined using eqn (4).

7 Environment ACS Paragon Plus

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

In order to accurately determine the Hugoniot curve for the detonation products, five sets of simulations of different temperatures are required. At each specific temperature, five Hugoniot values are obtained through changing material’s density, bracketing both H > 0 and H < 0 cases. Then the relationship between H and V/V0 is fitted using spline function to find the volume compression ratio V/V0 satisfying H = 0. The locations of the intersections between Hugoniot value at certain temperatures and H = 0 axis give the exact compression volume ratios for the fully reacted Hugoniot curve. We determine five H = 0 points to describe the curve accurately. Then the CJ point is located as the point of tangency of the Rayleigh line and the completely decomposed Hugoniot curve. Rayleigh line is a straight line connecting points corresponding to the initial and final states on a equation of state (P-V) plot for a substance subjected to a shock wave. The Hugoniot curve can be described by a quadratic polynomial and the Rayleigh line can be expressed as a line:

P = a0 + a1n + a2 n 2 (5)

P = an + b

where

n = V / V0

(6)

.

The Rayleigh line begins from point (1,0) because the initial pressure of the unshocked material can be neglected compared with the huge pressure after the shock front. Thus, we can obtain the CJ point volume-compression-ratio as:

nCJ = 1 − 1 +

a1 + a0 a2

(7).

8 Environment ACS Paragon Plus

Page 8 of 23

Page 9 of 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

The Journal of Physical Chemistry

The detonation velocity DCJ can be calculated by putting the PCJ and nCJ into the Hugoniot relations:

DCJ =

PCJ − P0 ρ0 (1 − nCJ )

(10).

The temperature at CJ state (TCJ) is then obtained using VCJ/V0 in the functional form of the T~V/V0 curve. We performed the RMD simulations using the Large-Scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) code20. We set up a time step as 0.1 fs and periodic boundary condition along three dimensions. For QMD simulations, we applied the Born-Oppenheimer molecular dynamics (BO-MD) approach implemented in VASP package with plane wave basis set21,22. The Perdew-Burke-Ernzerh (PBE) functional and the projector augmented-wave method (PAW) are applied accounting for the exchange-correlation interaction and the core-valence interaction, respectively23,24. The PAW potentials for C, H, N, O and Si considered the 2s22p2 electrons, 1s1 electron, 2s22p3 electrons, 2s22p4 electrons and 3s23p2 electrons as valence states, respectively. We used the tetrahedron method with Blöchl corrections to determine the electron partial occupancies25. The kinetic energy cut-off of 500 eV was used for the plane wave expansions. We set the convergence criteria to a 1 × 10-5 eV and 1 × 10-2 eV Å-1 for the electronic wave function and geometry optimization, respectively. To perform the Brillouin zone integration, we used the Γ-centered symmetry-reduced Monkhorst−Pack meshes with 1×1×1 k-point grid mesh in the Brillouin zone. We considered van der Waals interactions using DFTD3 method with Becke-Jonson damping approach26.

9 Environment ACS Paragon Plus

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

Page 10 of 23

The detailed products information is obtained by our molecular analytical procedure which uses the connectivity matrix based on bond orders calculated from ReaxFF with 0.1 ps intervals10. Our algorithm recognizes the molecule by comparing their bond orders with cutoff values and then assigns them with unique ID numbers to track their reaction routes. The fake bonds arising from the momentary fluctuations of bond orders are avoided by using a time window of 1.0 ps.

Table 1. Bond order cut-off values for different atom pairs. C

H

O

N

Si

C

0.55 0.40 0.80 0.30 0.30

H

0.55 0.40 0.55 0.30

O

0.65 0.55 0.30

N

0.45 0.30

Si

0.30

3 RESUTLS AND DISCUSSION In order to predict the CJ state of PETN and Si-PETN, we firstly determine the Hugoniot states of completely reacted states by computing from a long cook-off simulation holding the system at prescribed temperature and volume compression. Detonated state should be in the appropriate temperatures and volume compressions range where the final Hugoniot values are close to zero. In these ranges, we carried out NVT-MD simulation (constant volume, constant temperature and constant number of atoms) using ReaxFF for 25 cases (namely, five sets of temperatures are set, and for each temperature five different densities are considered). Then, in each case, we performed 270 ps - RMD for PETN system and 150 ps - RMD for Si-PETN system to reach the final equilibrium steady state, based on the energy release tendency. The

10 Environment ACS Paragon Plus

Page 11 of 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

The Journal of Physical Chemistry

total energy reveals the dependence of degree of thermal decomposition: for PETN, it increases initially because of the endothermic O–NO2 bond breaking; but for Si-PETN, it decreases directly for Si–O bond formation, as shown in Figure 2a and 2c. Then the secondary reactions are initiated, in which molecules break into small pieces along with releasing energy dramatically. For example, Figure 2a and 2c shows the time evolutions of total energy at V/V0 = 0.65 with T= 3000 K for PETN, and at V/V0 = 0.65 with T= 3600 K for Si-PETN, respectively. For PETN system, the total energy increases at the first 1 ps and decreases exponentially during the following 25 picoseconds, indicating the fast energy release from violent exothermic chemical reactions. Comparing with PETN system, since the Si-PETN molecule is easier to decompose and its system was heated to a higher temperature, Si-PETN system has a faster total energy decreasing process, which is about 15 ps. During the next 100 ps for PETN system and 50 ps for Si-PETN system, the energy decreases slowly indicating the system gradually reaches an equilibrium steady state. Then we performed an additional 100 ps to verify the unchanging of the detonated properties. Finally, we extracted the atoms positions and velocities at the end of RMD simulations (270 ps for PETN and 150 ps for Si-PETN) and use these data as input for QMD simulations for another 30 ps to obtain an accurate first-principles description. The system has reached equilibrium since the total energy remains constant for the last 20 ps and the stable products (H2O, CO2, N2, CO etc.) keep unchanged for the last 10 ps in the QMD simulations (Fig. 2b and Fig. 2d). The properties of fully decomposed states, such as pressure and the total energy are determined by averaging the last 10 ps QMD simulations.

11 Environment ACS Paragon Plus

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

Figure 2. Time evolution of the total energy per unit mass of NVT simulation. The initial energy is set to zero as a reference. (a) ReaxFF-MD for the first 270 ps (inserted for the initial 4 ps total energy evolution) and (b) QMD for the last 30 ps (inserted for the last 10 ps products distribution) at T = 3100 K and V/V0 = 0.65 for PETN system, (c) ReaxFF-MD for the first 150 ps (inserted for the initial 4 ps total energy evolution) and (d) QMD for the last 30 ps (inserted for the last 10 ps products distribution) at T = 3600 K and V/V0 = 0.65 for Si-PETN system.

For a specific temperature, a family of equilibrated reacted state (Hugoniot values) is obtained from eqn (4) by changing the volume-compression-ratio and using spline fitting. Five isotherms curves at various temperatures are plotted in order to describe the Hugoniot curve accurately. The intersections between Hugoniot function and the Hg = 0 axis (the dashed line) provides the locations of the final Hugoniot curve. Therefore, five different compression ratios are required to bracket the Hg = 0 case, as shown in Figure 3.

12 Environment ACS Paragon Plus

Page 12 of 23

Page 13 of 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

The Journal of Physical Chemistry

Figure 3. The spline fitted relationship of the Hugoniot value (Hg) and the volume-compression-ratio (V/V0) at various temperatures for PETN and Si-PETN. The fully reacted Hugoniot states are determined by the intersections of Hg = 0 line (the dashed line) and the spline fitted curves.

The completely detonated Hugoniot curve and the evolution of pressure as a function of volume are determined by five Hg = 0 points with their volume-compression-ratios and corresponding pressures, as shown in Figure 4a for PETN and 4c for Si-PETN. This is fitted by the quadratic polynomial. We determine the CJ point as the tangential point between the Rayleigh line and the Hugoniot curve, represented as red dots in Figure 4a and 4c. Thus, the CJ points are located at V/V0 = 0.72 for PETN and V/V0 = 0.70 for Si-PETN. The relation of temperature versus V/V0 is also fitted by a quadratic polynomial function, as shown in Figure 4b for PETN and 4d for Si-PETN, from which the temperature at the CJ state is derived.

13 Environment ACS Paragon Plus

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

Figure 4. The Hugoniot curve describing the completely detonated states, and the CJ points for (a) PETN and (c) Si-PETN, and the curve of detonated temperature and compressed volume proportion satisfying Hg = 0 for (b) PETN and (d) Si-PETN. The CJ points (red dot) were derived as the tangent point of the Rayleigh line and the fully reacted Hugoniot curve. The CJ volume-compression-ratio and the CJ pressure are V/V0 = 0.72 and PCJ = 32.23േ2.15 GPa for PETN and are V/V0 = 0.7 and PCJ = 27.78േ2.01 GPa for Si-PETN. The quadratic polynomial fitted temperature at the CJ points (the red dot) are TCJ = 2960 K for PETN and TCJ = 3638 K for Si-PETN.

The detonation properties can be well described at the Chapman-Jouguet (CJ) state because the reaction products at the end of reaction zone reach chemical equilibrium in CJ condition. As shown in Table 2, the predicted detonation velocity from our simulation is DCJ = 7.974േ0.262 km/s at ρ0 = 1.81 g/cm3 for PETN, in good agreement with experimental data of DCJ = 7.975

14 Environment ACS Paragon Plus

Page 14 of 23

Page 15 of 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

The Journal of Physical Chemistry

km/s at ρ0 = 1.67 g/cm3 and 8.26 km/s at ρ0 = 1.76 g/cm3, validating our results. No data available for Si-PETN, we predict the detonation velocity DCJ = 7.254േ0.257 km/s and the CJ temperature TCJ = 3638േ112 K at ρ0 = 1.76 g/cm3 for Si-PETN. Table 2. Detonation properties at the CJ state for PETN and Si-PETN27

PETN

Si-PETN

MD&QMD

Exp1

Exp2

MD&QMD

Density (g/cm3)

1.81

1.67

1.76

1.76

PCJ (GPa)

32.23±2.15

27.78േ2.01

TCJ (K)

2960 ± 101

3638 ± 112

VCJ(cm3/g)

0.398±0.017

0.398±0.019

DCJ(m/s)

7974±262

7975

8260

7254±257

For CJ temperature, the Si-PETN system is TCJ = 3638 ± 112 K, 22.9% higher than those of PETN system which is TCJ = 2960 ± 101 K. This is because the decomposition mechanism of SiPETN includes molecule rearrangement and the formation of Si-O bonds which release large amount of energy, promoting more complicated chemical reactions with huge energy releasing. The high energy releasing as heat or radiation makes the silicon usually used as irradiation agents and decoy flairs28. For CJ pressure, the Si-PETN system is PCJ = 27.78േ2.01 Gpa, which is 13.8% lower than PCJ = 32.23േ2.15 Gpa of PETN. As a result, lower detonation velocity is achieved as well: 7254 േ 257 m/s for Si-PETN vs 7974 േ 262 km/s for PETN. This high temperature with low pressure phenomena were also observed in our previous reactive molecular dynamics study of shock and thermal expansion simulations8,16.

15 Environment ACS Paragon Plus

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

Page 16 of 23

Si-PETN has a molecular structure nearly identical but only a central carbon atom difference to its carbon analogue PETN, but Si-PETN shows a high temperature with low pressure detonated properties. We devote this to the existence of Si element in Si-PETN leading to different types of detonation products distribution compared with PETN system. Thus, we analyze the main detonated products at the CJ states with the fragment analysis program by averaging the last 10 ps QMD simulations, as shown in Table 3. Table 3. Detonation Products predicted at the CJ state for PETN and Si-PETN, compared with experimental data29 PETN

Density (g/cm3)

Si-PETN

MD&QMD

Experiment confined

MD&QMD

1.81

1.73-1.74

1.76

Main Products

N2

2.00±0.00

2.00

1.50±0.00

(mol/mol)

H2O

1.80±0.09

3.50

0.51±0.00

HO

1.60±0.15

CO

1.00±0.12

1.69

0.50±0.06

CO2

1.50±0.28

3.39

1.50±0.12

0.77±0.08

Other molecules (mol/mol) carbon clusters

2.67±0.14

1.84±0.12

composition

C0.94H1.04O1.72

C0.68H1.78N0.54O1.60

silicon clusters

0.75±0.10

composition

C1.00H3.79N0.33O5.57Si1.33

For PETN system, N2, H2O, HO, CO and CO2 are the dominant detonation products. Compared with the calorimeter bomb test results29, the calculated total number of N2 molecules agrees well with experiments. However, the experimental data of 3.50 mol/mol of H2O is almost twice of the simulation results of 1.80 mol/mol. This is because our constant-volume simulations

16 Environment ACS Paragon Plus

Page 17 of 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

The Journal of Physical Chemistry

measure the product composition with temperature of ~3000 K, and many water molecules were not stable and dissociated into small fragments, such as HO, which is also found in the firstprinciples molecular dynamics study13, 30, 31. Since the number of HO fragments is 1.60 mol/mol in the system, we believe that our results are comparable to the experiment data. Our prediction for the number of CO2 and CO molecules from PETN is also lower than the experiment. The reason is that it takes seconds for materials in calorimetric experiment turning into detonation products, and the isentropic expansion along with cooling down to the ambient conditions is inevitable. However, we obtain the hot dense products at the CJ condition in the RMD and QMD simulations for hundreds of picoseconds. Thus, some carbon, hydrogen and oxygen atoms remain aggregated in clusters. We regarded that these clusters are in gas phase due to the simple overall formulation of C0.94H1.04O1.72, which will decompose and produce additional CO2 and CO products after expansion12. Thus, during the detonation process of PETN, the high transition percentage of mass from solid to gases devoted to high pressure and high detonation velocity. For Si-PETN, we found that N2, H2O, HO, CO and CO2 are the dominant detonation products as well. However, compared with PETN system, significant decreased amounts of H2O were found in Si-PETN. This can be explained from two factors: (1) the CJ temperature of SiPETN is even higher than PETN, so more water molecules are likely to exist; (2) silicon atoms attract surrounded oxygen atoms of water molecules, leading to their further dissociation. This explains why silicon based energetic materials can keep active in oxygen deficient conditions, such as in water, and can prevent incomplete detonation. In addition, although the amounts of CO2 found in Si-PETN are equal to those in PETN, but CO is fewer, indicating that some of the C atoms are trapped into aggregates. We divided the clusters into two groups based on with or without silicon atoms. For carbon clusters, 1.84 mol/mol were found with an averaged

17 Environment ACS Paragon Plus

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

composition of C0.68H1.78N0.54O1.60; for silicon clusters, 0.75 mol/mol were found with a total formula of C1.00H3.79N0.33O5.57Si1.33. These silicon aggregates containing 18.75% of the carbon atoms, 35.53% of the hydrogen atoms, 6.19% of the nitrogen atoms, and 34.81% of the oxygen atoms, arises from the fact that Si atoms attract surrounded oxygen atoms from HO, CO, and CO2, suppressing the generation of gaseous products, which makes a lower detonation pressure while detonating, leading to a lower detonation velocity. 4 CONCLUSIONS In summary, we carried out RxMD(cQM) simulations, a combination method of ReaxFFRMD and QMD simulations, to examine the influence of silicon element to detonation properties for EMs, by comparing PETN and its silicon analogue Si-PETN. We found that compared with PETN, the pressure from the same volume compression rate of completely decomposed Si-PETN system is lower, but the temperature is higher. For CJ temperature, the Si-PETN system is 22.9% higher than PETN system, but for CJ pressure, the SiPETN system is 13.8% lower than PETN. As a result, lower detonation velocity is calculated for Si-PETN than for PETN. The detonation products in PETN system are in gas phase including only few simple carbon clusters, leading to a high detonation pressure and a high detonation velocity. But with Si element involved in EMs, firstly, a large amount of heat is produced during Si products formation, resulting in a high detonating temperature; secondly, Si atoms attract surrounded oxygen atoms to form silicon aggregates, reducing the generation of gaseous products, which makes a lower detonation pressure and a lower detonation velocity. In addition, since silicon

18 Environment ACS Paragon Plus

Page 18 of 23

Page 19 of 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

The Journal of Physical Chemistry

atoms are able to attract oxygen atoms from other molecules, silicon-based EMs can keep active in oxygen deficient conditions, which is useful in preventing incomplete detonation. Our findings suggest that mixing Si-PETN with other EMs is an excellent approach to increase their adaptation to some particular conditions. For example, it can be used under low temperature conditions because of the huge energy release of Si-PETN. It could also be effectively applied in oxygen deficient conditions to prevent incomplete detonation due to the attraction of Si atoms to surrounded oxygen atoms.

AUTHOR INFORMATION

Corresponding Author *E-mail: [email protected], [email protected]

Notes The authors declare no competing financial interests. ACKNOWLEDGMENT This work is supported by the American Chemical Society Petroleum Research Fund (PRF# 58754-DNI6). REFERNECES (1) Becuwe, A.; Delclos. A. Low-sensitivity Explosive Compounds for Low Vulnerability Warheads. Prop., Explos., Pyrotech. 1993, 18, 1–10. (2) Urbański, T. Chemistry and Technology of Explosives; Pergamon Press: Oxford, 1984.

19 Environment ACS Paragon Plus

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

(3) Lewis, M. J.; Chang, J. S. Joint Jet-a/silane/hydrogen Reaction Mechanism. J. Propul. Power

2000, 16, 365–367. (4) dos Anjos, D. S. C.; Revoredo, E. C. V.; Galembeck, A. Silicone-polyacrylate Chemical Compatibilization with Organosilanes. Polym. Eng. Sci. 2010, 50, 606–612. (5) Klapotke, T. M.; Krumm, B.; Rainer llg; Troegel, D.; Tacke R. The Sila-explosives Si(CH2N3)4 and Si(CH2ONO2)4:  Silicon Analogues of the Common Explosives Pentaerythrityl Tetraazide, C(CH2N3)4, and Pentaerythritol Tetranitrate, C(CH2ONO2)4. J. Am. Chem. Soc.

2007, 129, 6908–6915. (6) Liu, W. G.; Zybin, S. V.; Dasgupta, S.; Klapotke, T. M.; Goddard, W. A. III. Explanation of the Colossal Detonation Sensitivity of Silicon Pentaerythritol Tetranitrate (Si-PETN) explosive. J. Am. Chem. Soc. 2009, 131, 7490–7491.

(7) Murray, J. S.; Lane, P.; Nieder, A.; Klapotke, T. M. and Politzer, P. Enhanced Detonation Sensitivities of Silicon Analogs of PETN: Reaction Force Analysis and The Role of σ–hole Interactions. Theor. Chem. Acc. 2010, 127, 345–354. (8) Zhou, T.; Liu, L.; Goddard, W. A. III; Zybin, S. V. and Huang, F. ReaxFF Reactive Molecular Dynamics on Silicon Pentaerythritol Tetranitrate Crystal Validates the Mechanism for the Colossal Sensitivity. Phys. Chem. Chem. Phys. 2014, 16, 23779-23791. (9) Zhou, T.; Cheng, T.; Zybin, S. V.; Goddard, W. A. III and Huang, F. Reaction Mechanisms and Sensitivity of Silicon Nitrocarbamate and Related Systems from Quantum Mechanics Reaction Dynamics. J. Mater. Chem. A 2018, 6, 5082-5097. (10) Zhou, T.; Zybin, S. V.; Goddard, W. A.; Cheng, T.; Naserifar, S.; Jaramillo-Botero, A.; Huang, F. Predicted Detonation Properties at the Chapman-Jouguet State for Proposed Energetic Materials (MTO and MTO3N) from Combined Reaxff and Quantum Mechanics Reactive Dynamics. Phys. Chem. Chem. Phys. 2018, 20, 3953-3969. (11) van Duin, A. C. T.; Dasgupta, S.; Lorant, F. and Goddard, W. A. III. ReaxFF:  A Reactive Force Field for Hydrocarbons. J. Chem. Phys. A 2001, 105, 9396–9409.

20 Environment ACS Paragon Plus

Page 20 of 23

Page 21 of 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

The Journal of Physical Chemistry

(12) Guo, D.; Zybin, S. V.; An, Q.; Goddard, W. A. III and Huang, F. Prediction of the Chapman–Jouguet Chemical Equilibrium State in a Detonation Wave from First Principles Based Reactive Molecular Dynamics. Phys. Chem. Chem. Phys. 2016, 18, 2015-2022. (13) Guo, D.; An, Q.; Zybin, S. V.; Goddard, W. A. III; Huang, F. and Tang, B. The Co-crystal of TNT/CL-20 Leads to Decreased Sensitivity toward Thermal Decomposition from First Principles based Reactive Molecular Dynamics. J. Mater. Chem. A 2015, 3, 5409–5419. (14) Guo, D.; An, Q.; Goddard, W. A. III; Zybin, S. V. and Huang, F. Compressive Shear Reactive Molecular Dynamics Studies Indicating that Cocrystals of TNT/CL-20 Decrease Sensitivity. J. Phys. Chem. C, 2014, 118, 30202–30208. (15) Liu, L.; Liu, Y.; Zybin, S.V.; Sun, H.; Goddard, W. A. III. ReaxFF-lg: 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. (16) An, Q.; Goddard, W. A. III; Zybin, S. V.; Jaramillo-Botero, A.; and Zhou, T. Highly Shocked Polymer Bonded Explosives at a Nonplanar Interface: Hot-spot Formation Leading to Detonation. J. Phys. Chem. C 2013, 117, 26551−26561. (17) An, Q.; Zybin, S. V.; Goddard, W. A. III; Jaramillo-Botero, A.; Blanco M.; Luo, S. N. Elucidation of the Dynamics for Hot-spot Initiation at Nonuniform Interfaces of Highly Shocked Materials. Phys. Rev. B 2011, 84, 220101. (18) Grimme, S.; Antony, J.; Ehrlich, S.; and Krieg, S. A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (DFT-D) for the 94 Elements H-Pu, J. Chem. Phys. 2010, 132, 154104.

(19) Zhurova, E. A.; Stash, A. I.; Tsirelson, V. G.; Zhurov, V. V.; Bartashevich, E. V.; Potemkin, V. A.; Pinkerton, A. A. Atoms-in-molecules Study of Intra-and Intermolecular Bonding in the Pentaerythritol Tetranitrate Crystal. J. Am. Chem. Soc. 2006, 128, 14728-14734. (20) Plimpton, S. Fast Parallel Algorithms for Short-range Molecular Dynamics. J. Comput. Phys. 1995, 117, 1-19.

21 Environment ACS Paragon Plus

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

(21) Kresse, G. Ab Initio Molecular Dynamics for Liquid Metals. J. Non-Cryst. Solids 1995, 192-193, 222-229.

(22) Kresse, G.; Furthmüller, J. Efficiency of Ab-Initio Total Energy Calculations for Metals and Semiconductors Using a Plane-Wave Basis Set. Comput. Mater. Sci. 1996, 6, 15-50. (23) Paier, J.; Hirschl, R.; Marsman, M.; Kresse, G. The Perdew–Burke–Ernzerhof ExchangeCorrelation Functional Applied to the G2-1 Test Set Using a Plane-Wave Basis Set. J. Chem. Phys. 2005, 122, 234102.

(24) Kresse, G.; Joubert, D. From Ultrasoft Pseudopotentials to the Projector Augmented-Wave Method. Phys. Rev. B 1999, 59, 1758-1775. (25) Blöchl, P. E.; Jepsen, O. Andersen, O. K. Improved Tetrahedron Method for Brillouin-Zone Integrations. Phys. Rev. B 1994, 49, 16223-16233. (26) Grimme, S.; Ehrlich, S.; Goerigk, L. Effect of the Damping Function in Dispersion Corrected Density Functional Theory. J. Comput. Chem. 2011, 32, 1456-1465. (27) Davis, W. C.; Craig, B. G.; Ramsay, J. B. Failure of the Chapman-Jouguet Theory for Liquid and Solid Explosives. Phys. Fluids, 1965, 8, 2169–2182. (28) Simone, D.; Bruno C.; Hidding, B. Silanes as Fuels for Scramjets and Their Applications. J. Propulsion 2006, 22, 1006–1011.

(29) Ornellas, D. L. Calorimetric Determinations of the Heat and Products of Detonation for Explosives: October 1961 to April 1982; Lawrence Livermore Laboratory, University of

California, Livermore, CA, 1982. (30) Wu, C. J.; Fried, L. E.; Yang, L. H.; Goldman, N. and Bastea, S. Catalytic behaviour of dense hot water. Nat. Chem. 2009, 1, 57–62. (31) Liu, Y.; Zybin, S. V.; Guo, J.; van Duin, A. C. T. and Goddard, W. A. III. Reactive Dynamics Study of Hypergolic Bipropellants: Monomethylhydrazine and Dinitrogen Tetroxide. J. Phys. Chem. B 2012, 116, 14136–14145.

22 Environment ACS Paragon Plus

Page 22 of 23

Page 23 of 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

The Journal of Physical Chemistry

TOC Graphic

23 Environment ACS Paragon Plus