^{1}, Thomas D. Sewell

^{1,a)}and Donald L. Thompson

^{1,b)}

### Abstract

The structural relaxation of crystalline nitromethane initially at *T* = 200 K subjected to moderate (∼15 GPa) supported shocks on the (100), (010), and (001) crystal planes has been studied using microcanonical molecular dynamics with the nonreactive Sorescu–Rice–Thompson force field [D. C. Sorescu, B. M. Rice, and D. L. Thompson, J. Phys. Chem. B **104**, 8406 (2000)]. The responses to the shocks were determined by monitoring the mass density, the intermolecular, intramolecular, and total temperatures (average kinetic energies), the partitioning of total kinetic energy among Cartesian directions, the radial distribution functions for directions perpendicular to those of shock propagation, the mean-square displacements in directions perpendicular to those of shock propagation, and the time dependence of molecular rotational relaxation as a function of time. The results show that the mechanical response of crystalline nitromethane strongly depends on the orientation of the shock wave. Shocks propagating along [100] and [001] result in translational disordering in some crystal planes but not in others, a phenomenon that we refer to as *plane-specific disordering*; whereas for [010] the shock-induced stresses are relieved by a complicated structural rearrangement that leads to a paracrystalline structure. The plane-specific translational disordering is more complete by the end of the simulations (∼6 ps) for shock propagation along [001] than along [100]. Transient excitation of the intermolecular degrees of freedom occurs in the immediate vicinity of the shock front for all three orientations; the effect is most pronounced for the [010] shock. In all three cases excitation of molecular vibrations occurs more slowly than the intermolecular excitation. The intermolecular and intramolecular temperatures are nearly equal by the end of the simulations, with 400–500 K of net shock heating. Results for two-dimensional mean-square molecular center-of-mass displacements, calculated as a function of time since shock wave passage in planes perpendicular to the direction of shock propagation, show that the molecular translational mobility in the picoseconds following shock wave passage is greatest for [001] and least for the [010] case. In all cases the root-mean-square center-of-mass displacement is small compared to the molecular diameter of nitromethane on the time scale of the simulations. The calculated time scales for the approach to thermal equilibrium are generally consistent with the predictions of a recent theoretical analysis due to Hooper [J. Chem. Phys. **132**, 014507 (2010)].

We are grateful to Dr. Ali Siavosh-Haghighi for many fruitful discussions and assistance with some of the technical aspects of the calculations. This work was supported by a DOD MURI grant managed by the Army Research Office under Grant Number W911NF-05–1-0265 and by the Los Alamos National Security, LLC Laboratory Directed Research and Development program.

I. INTRODUCTION

II. METHODS

A. Force field

B. Computational details

C. Analysis of simulation results

III. RESULTS

A. Shock velocity and shock pressure

B. Qualitative features of the shocked states

C. Density

D. Temperature

E. 2D molecular translational disordering

F. 2D molecular mass transport

G. Molecular rotational relaxation

IV. SUMMARY AND CONCLUSIONS

### Key Topics

- Shock waves
- 65.0
- Shock wave effects
- 8.0
- Crystal structure
- 7.0
- Single crystals
- 7.0
- Lattice constants
- 6.0

## Figures

Schematic diagram of the approach used to generate shock waves. The leftmost (tan) region is the stationary piston comprised of rigid molecules. A slab containing flexible molecules (blue and green regions) equilibrated to 200 K impacts with speed *u* _{ p } from the right onto the stationary piston. The interface between the postshock material (blue region) and the preshock material (green region) is the shock front position (solid red line) that advances at shock velocity *u* _{ s } from left to right.

Schematic diagram of the approach used to generate shock waves. The leftmost (tan) region is the stationary piston comprised of rigid molecules. A slab containing flexible molecules (blue and green regions) equilibrated to 200 K impacts with speed *u* _{ p } from the right onto the stationary piston. The interface between the postshock material (blue region) and the preshock material (green region) is the shock front position (solid red line) that advances at shock velocity *u* _{ s } from left to right.

Nitromethane crystal structure and the definition of molecular layers used in the analysis of the MD results. (a) Initial configuration of the piston for the [100] shock simulation; the 96 molecules between the two dashed lines comprise one molecular layer. (b) Initial configuration of the piston for the [010] shock simulation; only one-half of the vertical size of the 20*a* × 65*b* × 12*c* system studied is shown; the 240 molecules between the two dashed lines comprise one molecular layer (480 molecules in the full-size sample). (c) Initial configuration of the piston for the [001] shock simulation; the 80 molecules between the two dashed lines comprise one molecular layer.

Nitromethane crystal structure and the definition of molecular layers used in the analysis of the MD results. (a) Initial configuration of the piston for the [100] shock simulation; the 96 molecules between the two dashed lines comprise one molecular layer. (b) Initial configuration of the piston for the [010] shock simulation; only one-half of the vertical size of the 20*a* × 65*b* × 12*c* system studied is shown; the 240 molecules between the two dashed lines comprise one molecular layer (480 molecules in the full-size sample). (c) Initial configuration of the piston for the [001] shock simulation; the 80 molecules between the two dashed lines comprise one molecular layer.

(a) Center-of-mass position for each layer along the shock direction [010] (CMP_{ y }) as a function of the time during shock propagation (solid lines). The dashed line connecting shock front positions as time evolves tracks the shock trajectory in the *y–t* plane; the CMP_{ y } velocity is *u* _{ c }. The numbers 21, 31, 51, 111, and 121 denote specific layer numbers. (b) Layers with the maximum gradient of kinetic energy [(*K* ^{′} _{ y })_{max}] plotted as a function of time. The calculation of the maximum gradient of the kinetic energy for each layer is the same as described in Ref. 31.

(a) Center-of-mass position for each layer along the shock direction [010] (CMP_{ y }) as a function of the time during shock propagation (solid lines). The dashed line connecting shock front positions as time evolves tracks the shock trajectory in the *y–t* plane; the CMP_{ y } velocity is *u* _{ c }. The numbers 21, 31, 51, 111, and 121 denote specific layer numbers. (b) Layers with the maximum gradient of kinetic energy [(*K* ^{′} _{ y })_{max}] plotted as a function of time. The calculation of the maximum gradient of the kinetic energy for each layer is the same as described in Ref. 31.

Schematic diagram of the approach used to define fixed-volume bins for trajectory analysis. In the laboratory frame the spatial origin for the bins translates linearly to the right with time, so as to remain stationary in the reference frame centered on the shock front. The width of the bins Δλ is constant. Bins for postshock and preshock materials are labeled by negative and positive integers, respectively; the bin that contains the shock front at any given instant in time is denoted *bin 1*. The number of postshock [*N* _{ o }(*t*)] and preshock [*N* _{r}(*t*)] bins increases and decreases, respectively, with time *t*. The total number of bins *N* _{ b }(*t*) = *N* _{ o }(*t*) + *N* _{ r }(*t*) decreases linearly in time due to compression of the flexible slab as the shock wave advances through the material. Bins, or partial bins, at least 14 Å from either end of the flexible slab (denoted by crosses in the left-most and right-most bins in the figure) are excluded from the analysis of properties.

Schematic diagram of the approach used to define fixed-volume bins for trajectory analysis. In the laboratory frame the spatial origin for the bins translates linearly to the right with time, so as to remain stationary in the reference frame centered on the shock front. The width of the bins Δλ is constant. Bins for postshock and preshock materials are labeled by negative and positive integers, respectively; the bin that contains the shock front at any given instant in time is denoted *bin 1*. The number of postshock [*N* _{ o }(*t*)] and preshock [*N* _{r}(*t*)] bins increases and decreases, respectively, with time *t*. The total number of bins *N* _{ b }(*t*) = *N* _{ o }(*t*) + *N* _{ r }(*t*) decreases linearly in time due to compression of the flexible slab as the shock wave advances through the material. Bins, or partial bins, at least 14 Å from either end of the flexible slab (denoted by crosses in the left-most and right-most bins in the figure) are excluded from the analysis of properties.

Projections of snapshots of molecular centers of mass for the three shock simulations studied at the respective instants of maximum compression. Projections along (a) [010] for the [100] shock; (b) [001] for the [100] shock; (c) [001] for the [010] shock; (d) [100] for the [010] shock; (e) [100] for the [001] shock; and (f) [100] followed by a counterclockwise rotation of 50° around the direction [001] for the [001] shock. The bottom abscissa is the distance from the shock front and the top abscissa is the time since the passage of the shock front. The regions to the left of the yellow lines are the fixed pistons; the red solid lines denote the instantaneous shock front positions at the ends of the simulations. The black dashed lines mark the approximate postshock times at which the 2D RDFs shown in Fig. 9 converge to essentially stable values for the [100] and [001] shocks.

Projections of snapshots of molecular centers of mass for the three shock simulations studied at the respective instants of maximum compression. Projections along (a) [010] for the [100] shock; (b) [001] for the [100] shock; (c) [001] for the [010] shock; (d) [100] for the [010] shock; (e) [100] for the [001] shock; and (f) [100] followed by a counterclockwise rotation of 50° around the direction [001] for the [001] shock. The bottom abscissa is the distance from the shock front and the top abscissa is the time since the passage of the shock front. The regions to the left of the yellow lines are the fixed pistons; the red solid lines denote the instantaneous shock front positions at the ends of the simulations. The black dashed lines mark the approximate postshock times at which the 2D RDFs shown in Fig. 9 converge to essentially stable values for the [100] and [001] shocks.

The mass density as a function of distance from the shock front (bottom abscissa) and time since shock passage (top abscissa). (a) Results for bin widths Δλ set equal to the unit cell parameters (*a*, *b*, and *c* for the [100], [010] and [001] shock directions, respectively); and (b) as in (a) except that Δλ = 2.0 Å for all three cases. Shocked material is associated with negative values for the position and positive values for time; the time is calculated based on the shock velocity in the [100] case, which is only slightly different from those for the [010] and [001] cases since the shock velocities for all three cases are very similar (see Table II).

The mass density as a function of distance from the shock front (bottom abscissa) and time since shock passage (top abscissa). (a) Results for bin widths Δλ set equal to the unit cell parameters (*a*, *b*, and *c* for the [100], [010] and [001] shock directions, respectively); and (b) as in (a) except that Δλ = 2.0 Å for all three cases. Shocked material is associated with negative values for the position and positive values for time; the time is calculated based on the shock velocity in the [100] case, which is only slightly different from those for the [010] and [001] cases since the shock velocities for all three cases are very similar (see Table II).

The intermolecular (lattice), intramolecular (vibrational–rotational), and total temperatures as functions of distance (bottom abscissa) and time (top abscissa) from the shock front. (a) Dashed curves: [100] shock; solid curves: [001] shock. The postshock time in (a) is determined as described in the caption for Fig. 6. (b) [010] shock; the postshock time is determined using the [010] shock velocity.

The intermolecular (lattice), intramolecular (vibrational–rotational), and total temperatures as functions of distance (bottom abscissa) and time (top abscissa) from the shock front. (a) Dashed curves: [100] shock; solid curves: [001] shock. The postshock time in (a) is determined as described in the caption for Fig. 6. (b) [010] shock; the postshock time is determined using the [010] shock velocity.

As in Fig. 7 except the Cartesian components of total kinetic energy are shown.

As in Fig. 7 except the Cartesian components of total kinetic energy are shown.

2D center-of-mass RDFs for shocks along (a) [100], (b) [010], and (c) [001].

2D center-of-mass RDFs for shocks along (a) [100], (b) [010], and (c) [001].

2D MSDs of molecular centers of mass vs time.

2D MSDs of molecular centers of mass vs time.

Orientational order parameter *P* _{2}(*θ*) for C–N bond vectors as functions of time.

Orientational order parameter *P* _{2}(*θ*) for C–N bond vectors as functions of time.

Molecular rotational distribution functions. (a) Angle distribution functions for five evenly spaced times between 0.2 and 1.0 ps behind the shock front. The results for the [001] shock are shown by the dashed curves and those for [100] by the solid curves. (b) Angle distribution function for five evenly spaced times between 1.0 and 5.0 ps behind the shock front. The results for the [001] shock are shown by the dashed curves and those for [100] by the solid curves. [(c) and (d)]: Same as (a) and (b) except the results are for the [010] shock.

Molecular rotational distribution functions. (a) Angle distribution functions for five evenly spaced times between 0.2 and 1.0 ps behind the shock front. The results for the [001] shock are shown by the dashed curves and those for [100] by the solid curves. (b) Angle distribution function for five evenly spaced times between 1.0 and 5.0 ps behind the shock front. The results for the [001] shock are shown by the dashed curves and those for [100] by the solid curves. [(c) and (d)]: Same as (a) and (b) except the results are for the [010] shock.

## Tables

Lattice parameters of nitromethane at *T* = 200 K and *P* = 1 atm.

Lattice parameters of nitromethane at *T* = 200 K and *P* = 1 atm.

Shock velocity, pressure, and temperature for shocks along [100], [010], and [001].^{a}

Shock velocity, pressure, and temperature for shocks along [100], [010], and [001].^{a}

Article metrics loading...

Full text loading...

Commenting has been disabled for this content