Article pubs.acs.org/JPCB
Anisotropic Responses and Initial Decomposition of CondensedPhase β‑HMX under Shock Loadings via Molecular Dynamics Simulations in Conjunction with Multiscale Shock Technique Ni-Na Ge,†,‡ Yong-Kai Wei,‡ Zhen-Fei Song,† Xiang-Rong Chen,*,‡ Guang-Fu Ji,*,† Feng Zhao,† and Dong-Qing Wei*,§ †
National Key Laboratory of Shock Wave and Detonation Physics, Institute of Fluid Physics, Chinese Academy of Engineering Physics, Mianyang 621999, China ‡ Key Laboratory of High Energy Density Physics and Technology of Ministry of Education, College of Physical Science and Technology, Sichuan University, Chengdu 610064, China § State Key Laboratory of Explosion Science and Technology, Beijing Institute of Technology, Beijing 00081, China ABSTRACT: Molecular dynamics simulations in conjunction with multiscale shock technique (MSST) are performed to study the initial chemical processes and the anisotropy of shock sensitivity of the condensed-phase HMX under shock loadings applied along the a, b, and c lattice vectors. A self-consistent charge density-functional tight-binding (SCC-DFTB) method was employed. Our results show that there is a difference between lattice vector a (or c) and lattice vector b in the response to a shock wave velocity of 11 km/s, which is investigated through reaction temperature and relative sliding rate between adjacent slipping planes. The response along lattice vectors a and c are similar to each other, whose reaction temperature is up to 7000 K, but quite different along lattice vector b, whose reaction temperature is only up to 4000 K. When compared with shock wave propagation along the lattice vectors a (18 Å/ps) and c (21 Å/ps), the relative sliding rate between adjacent slipping planes along lattice vector b is only 0.2 Å/ ps. Thus, the small relative sliding rate between adjacent slipping planes results in the temperature and energy under shock loading increasing at a slower rate, which is the main reason leading to less sensitivity under shock wave compression along lattice vector b. In addition, the C−H bond dissociation is the primary pathway for HMX decomposition in early stages under high shock loading from various directions. Compared with the observation for shock velocities Vimp = 10 and 11 km/s, the homolytic cleavage of N−NO2 bond was obviously suppressed with increasing pressure.
1. INTRODUCTION Octahydro-1,3,5,7-tetranitro-1,3,5,7-tetrazocine (HMX) has been widely used in weaponry and civilian areas, and its initial chemical processes have been investigated for many years. Glascoe et al.,1 Henson et al.,2 and Brill and Karpowicz3 found that HMX experienced a polymorphic transition from the β to δ phases during thermal decomposition. Ji et al.4 and Sharia and Kuklja5−7 found that the bonds of N−NO2 were the weakest, which could trigger explosion. Prasad et al.8 and Zhao et al.9 as well as Mcguire and Tarver10 suggested that C−N bond scission played an important role in the decomposition reactions for the condensed phase. Recently, Cui et al.11 indicated that N−NO2 and C−N were easily compressed with increasing pressure and became relatively strong and stable. In contrast, C−H bond became weaker under high pressures, which is also supported by other theoretical studies.12 The theoretical studies of HMX have made significant progress; however, the chemical processes upon mechanical loading, such as shock waves, are rare to see. Over the years, except for the understanding of chemical processes of energetic materials HMX, there was another hot © 2014 American Chemical Society
research, which is the sensitivity of explosives upon impact. The anisotropic response of high-energy material HMX has widely been investigated in both theoretical and experimental studies.13−22 In theoretical studies, Conroy et al.,16 Xiao et al.,17 and Zhou et al.18 investigated the anisotropy of sensitive in condensed-phase HMX through first-principles density functional theory (DFT), classical molecular dynamics (MD), and compressive shear reactive dynamics (CS-RD) model, respectively. In experimental studies, most reseachers19−21 indicated that the anisotropy of shock sensitivity exhibited by HMX is significant. For example, Menikoff et al.20 showed that the calculated effective yield strengths were 0.18 and 0.31 GPa for the directions normal to the (0 1 1) and (0 1 0) planes. Dick et al.21 indicated that the (0 1 0) orientation has larger elastic precursors and does not have the regular plastic deformation. However, Hooks and Ramos22 found that the onset of shockinduced reactions occurred nearly at the same pressure for Received: March 11, 2014 Revised: June 25, 2014 Published: June 25, 2014 8691
dx.doi.org/10.1021/jp502432g | J. Phys. Chem. B 2014, 118, 8691−8699
The Journal of Physical Chemistry B
Article
2. COMPUTATIONAL DETAILS The MD in conjunction with multiscale shock technique (MSST) is based on the Navier−Stokes equations for compressible flow.30 Instead of simulating a shock wave within a large computational cell with many atoms (the direct approach), the computational cell of the multiscale technique follows a Lagrangian point through the shock wave, which could be simulated with fewer atoms and lower computational cost.31 The multiscale shock technique simulates steady shock waves (that is, constant shock speed) by the time-evolving equations of motion for atoms and volume of the computational cell to constrain the shock-propagation-vector stress to the Rayleigh line and energy to the Hugoniot relation. For a specified shock speed, the Hugoniot relation E − E0 = 1/2(P + P0)(V0 − V) could be obtained by the conservation of mass, momentum, and energy across the shock front.32,33 Here E is the energy, P is the negative of the diagonal component of the stress tensor in the vector of the shock, and V is the volume. A subscript 0 refers to the preshocked state, while quantities without subscripts refer to the postshocked state. The thermodynamic path connecting the initial state of the system to its final (Hugoniot) state could be described by the Rayleigh line P − P0 = U2ρ0(1 − ρ0/ρ) (where U is the shock velocity and ρ is the density).26 These two relations describe a steady planar shock wave within continuum theory, which has to be obeyed for steady shock waves. Our electronic-structure calculations are based on the selfconsistent charge density functional tight binding (SCCDFTB) scheme,34 which is based on a second-order expansion of the Kohn−Sham total energy in density functional theory (DFT) with respect to charge density fluctuations. This method allows description of total energies, atomic forces, and charge transfer in a self-consistent manner. It has been successfully tested on organic and bioorganic systems35,36 and shown to accurately predict reaction energies.37−39 It was applied to the study of nitromethane and TATB under pressure33,38 which has yielded comparable results to those obtained from the other density functional calculations. Furthermore, it is also applied to study the reactivity of shocked perfect crystal TATB.39 In this study, the simulation supercell was taken from the experimentally determined X-ray crystal structure.40 It is wellknown that β-HMX is a monoclinic molecular crystals. a, b, and c represent the three monoclinic crystallographic axes of HMX in the space group P21/n. The calculated equilibrium structure of the unit-cell HMX is showed in Table 1. The conjugate gradient was used for the geometry optimization. The convergence criterion for the maximum geometry change
shocks delivered in the (0 1 1) and (0 1 0) vectors. Consequently, even though the theoretical and experimental studies of the HMX have made significant progress, the coupling between the mechanical loading and anisotropy of shock sensitivity is not yet clear. It is generally known that the coupling between shock loading induced mechanical response and complicated chemistry makes the study of high-energy materials very challenging. Understanding the chemical processes is required to build an improved model for combustion or detonation under shock. Recently, nonequilibrium MD simulations based on the quantum mechanical (QM) methods (quantum molecular dynamics simulation) are capable of studying the chemical processes in the initial stages of shock detonation in condensed systems. Strachan et al.23 developed the ReaxFF that well describes the chemical and mechanical properties of RDX at various collision velocities. They found the primary reactions involving NO2, OH, NO, and N2, NO2 being the main products for low impact velocities (
