The isothermal second-order elastic stiffness tensor and isotropic moduli of β-1,3,5,7- tetranitro-1,3,5,7-tetrazoctane (β-HMX) were calculated, using the P21/n space group convention, from molecular ...dynamics for hydrostatic pressures ranging from 10−4 to 30 GPa and temperatures ranging from 300 to 1100 K using a validated all-atom flexible-molecule force field. The elastic stiffness tensor components were calculated as derivatives of the Cauchy stress tensor components with respect to linear strain components. These derivatives were evaluated numerically by imposing small, prescribed finite strains on the equilibrated β-HMX crystal at a given pressure and temperature and using the equilibrium stress tensors of the strained cells to obtain the derivatives of stress with respect to strain. For a fixed temperature, the elastic coefficients increase substantially with increasing pressure, whereas, for a fixed pressure, the elastic coefficients decrease as temperature increases, in accordance with physical expectations. Comparisons to previous experimental and computational results are provided where possible.
•Heat current in solids can be filtered without affecting thermal conductivity.•Filtering reduces oscillations in the current correlation function substantially.•Filtering reduces oscillations in the ...time-dependent thermal conductivity.•New fitting procedure for the time-dependent thermal conductivity is developed.
The closely related Green–Kubo and Helfand moment approaches are applied to obtain the thermal conductivity tensor of β-1,3,5,7-tetranitro-1,3,5,7-tetrazoctane (β-HMX) at T=300 K and P=1 atm from equilibrium molecular dynamics (MD) simulations. Direct application of the Green–Kubo formula exhibits slow convergence of the integrated thermal conductivity values even for long (120 ns) simulation times. To partially mitigate this slow convergence we developed a numerical procedure that involves filtering of the MD-calculated heat current. The filtering is accomplished by physically justified removal of the heat-current component which is given by a linear function of atomic velocities. A double-exponential function is fitted to the integrated time-dependent thermal conductivity, calculated using the filtered current, to obtain the asymptotic values for the thermal conductivity. In the Helfand moment approach the thermal conductivity is obtained from the rates of change of the averaged squared Helfand moments. Both methods are applied to periodic β-HMX supercells of four different sizes and estimates for the thermal conductivity of the infinitely large crystal are obtained using Matthiessen’s rule. Both approaches yield similar although not identical thermal conductivity values. These predictions are compared to experimental and other theoretically determined values for the thermal conductivity of β-HMX.
We apply the Green–Kubo (G–K) approach to obtain the thermal conductivity tensor of β‐1,3,5,7‐tetranitro‐1,3,5,7‐tetrazoctane (β‐HMX) as a function of pressure and temperature from equilibrium ...molecular dynamics (MD) simulations. Direct application of the G–K formula exhibits slow convergence of the integrated thermal conductivity values even for long (120 ns) accumulated simulation times. To partially mitigate this slow convergence, we apply a recently implemented numerical procedure that involves physically justified filtering of the MD‐calculated G–K heat current and fitting the integrated time‐dependent thermal conductivity to a physically motivated double‐exponential function. The thermal conductivity tensor was determined for pressures 0.1 MPa
≤P≤
${\le P\le }$
30 GPa and temperatures 300 K
≤T≤
${\le T\le }$
900 K. The thermal conductivity
καβT,P
${{\kappa }^{\alpha \beta }\left(T,P\right)}$
increases with increasing pressure, by approximately an order of magnitude over the interval considered, and decreases with increasing temperature. The MD predictions are compared to experimental and other theoretically determined values for the thermal conductivity of β‐HMX. A simple, semi‐empirical fitting form is proposed that captures the behavior of
κααT,P
${{\kappa }^{\alpha \alpha }\left(T,P\right)}$
over the pressure and temperature intervals studied.
Using shock‐induced quasi‐2D pore collapse in β‐HMX as a specific case study, we address three practical questions that arise when designing and interpreting molecular dynamics (MD) simulations of ...explicit shock wave propagation in crystals. How sensitive are the overall results to sample thickness perpendicular to the (quasi) plane of the problem? For impacts on a given crystal surface, how much does the transverse orientation of the sample matter? And, for a given sample size and orientation, how much run‐to‐run variability exists among results for independent but statistically equivalent realizations of the shock? The first and second questions are interrelated but pertain individually to assessing the roles of finite‐size and crystal anisotropy effects, respectively, on the simulated collapse mechanisms and associated local thermo‐mechanical states in the sample during and after collapse. The third addresses the confidence with which the results of individual simulations can be regarded as representative. All three questions become particularly important if the MD predictions are intended to serve as “ground truth” for validation of continuum mesoscale models. Here, quasi‐2D samples of (010)‐oriented single‐crystal β‐HMX (β‐1,3,5,7‐tetranitro‐1,3,5,7‐tetrazocane) containing a cylindrical pore were subjected to reverse‐ballistic impacts, resulting in explicit, supported shock propagation initially along 010, for impact speeds up
=1.0 and 0.5 km s−1. The individual samples differed in either the thickness perpendicular to the sample plane or the transverse crystal orientation normal to 010 in the quasi‐2D cell. Three independent realizations of the shock were performed for a selected case to assess run‐to‐run variability. Comparisons of qualitative features of the collapse, temperature and pressure distributions in the samples, and time scales for pore collapse suggest that there is little sensitivity to sample thickness for the same crystal orientation and moderate sensitivity to transverse crystal orientation for samples of the same thickness. Run‐to‐run variability is evident to the eye in side‐by‐side system observations. However, overall mechanisms of collapse, distributions for temperature and pressure in the samples, and time scales for collapse are in near‐quantitative agreement among the realizations.
All‐atom molecular dynamics (MD) and Eulerian continuum simulations, performed using the LAMMPS and SCIMITAR3D codes, respectively, were used to study thermo‐mechanical aspects of the shock‐induced ...collapse of an initially cylindrical 50 nm diameter pore in single crystals of 1,3,5‐triamino‐2,4,6‐trinitrobenzene (TATB). Three impact speeds, 0.5 km s−1, 1.0 km s−1 and 2.0 km s−1, were used to generate the shocks. These impact conditions are expected to yield collapse mechanisms ranging from predominantly visco‐plastic to hydrodynamic. For the MD studies, three crystal orientations (relative to shock‐propagation direction) were examined that span the limiting cases with respect to the crystalline anisotropy in TATB. An isotropic constitutive model was used for the continuum simulations, thus crystal anisotropy is absent. The evolution of spatiotemporally resolved quantities during collapse is reported including local pressure, temperature, pore size and shape, and material flow. Multiple models for the melting curve and specific heat were studied. Within the isotropic elastic/perfectly plastic continuum framework and for the range of impact conditions studied, the specific heat and melting curve sub‐models are shown to have modest effects on the continuum hotspot predictions in the present inert calculation. Treating the MD predictions as ‘ground truth’, albeit with a classical rather than quantum‐like heat capacity, it is clear that – at a minimum – an extension of the constitutive model to account for crystal plasticity and anisotropic strength will be required for a high‐fidelity continuum description.
All-atom molecular dynamics (MD) simulations were used to study shock wave loading in oriented single crystals of the highly anisotropic triclinic molecular crystal ...1,3,5-triamino-2,4,6-trinitrobenzene (TATB). The crystal structure consists of planar hydrogen-bonded sheets of individually planar TATB molecules that stack into graphitic-like layers. Shocks were studied for seven systematically prepared crystal orientations with limiting cases that correspond to shock propagation exactly perpendicular and exactly parallel to the graphitic-like layers. The simulations were performed for initially defect-free crystals using a reverse-ballistic configuration that generates explicit, supported shocks. Final longitudinal stress components are between ≈8.5 and ≈10.5 GPa for the 1.0 km s–1 impact speed studied. Orientation-dependent properties are reported including shock speeds, stresses, temperatures, compression ratios, and local material strain rates. Spatiotemporal maps of the temperature, stress tensor, material flow, and molecular orientations reveal complicated processes that arise for specific shock directions. The results indicate that TATB shock response is highly sensitive to crystal orientation, with significant qualitative differences for the time evolution of the stress tensor and temperature, elastic/inelastic compression response, defect formation and growth, critical von Mises stress, and strain rates during shock rise that span nearly an order of magnitude. A variety of inelastic deformation mechanisms are identified, ranging from crumpling of graphitic-like layers to dislocation-mediated plasticity to intense shear strain localization. To our knowledge, these are the first systematic MD simulations and analysis of explicit shock wave propagation along nontrivial crystal directions in a triclinic molecular crystal.
Heterogeneous energetic materials (EMs) subjected to mechanical shock loading exhibit complex thermo-mechanical processes which are driven by the high temperature, pressure, and strain rate behind ...the shock. These lead to spatial energy localization in the microstructure, colloquially known as “hotspots” where chemistry may commence and possibly culminating in detonation. Shock-induced pore collapse is one of the dominant mechanisms by which localization occurs. In order to physically predict the shock sensitivity of energetic materials under these extreme conditions, we formulate a multiplicative crystal plasticity model with key features inferred from molecular dynamics (MD) simulations. Within the framework of thermodynamics, we incorporate the pressure dependence of both monoclinic elasticity and critical resolved shear stress into the crystal plasticity formulation. Other fundamental mechanisms, such as strain hardening and pressure-dependent melting curves, are all inferred from atomic-scale computations performed across relevant intervals of pressure and temperature. To handle the extremely large deformation and the evolving geometry of the self-contact due to pore collapse, we leverage the capabilities of the Material Point Method (MPM) to track the interface via the Lagrangian motion of material points and the Eulerian residual update to avoid the mesh distortion issue. This combination of features enables us to simulate the shock-induced pore collapse and associated hotspot evolution with a more comprehensive physical underpinning, which we apply to the monoclinic crystal β-HMX. Treating MD predictions of the pore collapse as ground truth,head-to-head validation comparisons between MD and MPM predictions are made for samples with identical sample geometry and similar boundary conditions for reverse-ballistic impact speeds ranging from 0.5kms−1 to 2.0kms−1. Comparative studies are performed to reveal the importance of incorporating a frictional contact algorithm, pressure-dependent elastic stiffness, and non-Schmid critical resolved shear stress in the mesoscale model.
The second‐order elastic stiffness tensor and isotropic moduli of β‐octahydro‐1,3,5,7‐tetranitro‐1,3,5,7‐tetrazocine (β‐HMX, P21/n space group setting) on the 0 K isotherm are presented for ...hydrostatic pressures between 10−4 GPa and 30 GPa. The results were obtained from molecular statics simulations using a validated all‐atom flexible‐molecule force field. Comparisons to previous experimental and computational determinations are provided.
We present a machine learning framework to train and validate neural networks to predict the anisotropic elastic response of a monoclinic organic molecular crystal known as β$$ \beta $$‐HMX in the ...geometrical nonlinear regime. A filtered molecular dynamic (MD) simulations database is used to train neural networks with a Sobolev norm that uses the stress measure and a reference configuration to deduce the elastic stored free energy functional. To improve the accuracy of the elasticity tangent predictions originating from the learned stored free energy, a transfer learning technique is used to introduce additional tangential constraints from the data while necessary conditions (e.g., strong ellipticity, crystallographic symmetry) for the correctness of the model are either introduced as additional physical constraints or incorporated in the validation tests. Assessment of the neural networks is based on (1) the accuracy with which they reproduce the bottom‐line constitutive responses predicted by MD, (2) the robustness of the models measured by detailed examination of their stability and uniqueness, and (3) the admissibility of the predicted responses with respect to mechanics principles in the finite‐deformation regime. We compare the training efficiency of the neural networks under different Sobolev constraints and assess the accuracy and robustness of the models against MD benchmarks for β$$ \beta $$‐HMX.
PDF
