Effect of the Prolactin Inhibitor Cabergoline and the Gonadotropin Releasing-Hormone Agonist Deslorelin in the Suppression of Plasma Prolactin Concentrations and Egg Laying in Quail (Coturnix japonica).

Egg binding and excessive laying frequently affect avian patients, and in many cases the treatment includes suppression of egg production. Currently, for the suppression of egg production in avian patients, a gonadotropin-releasing hormone agonist, in the form of a deslorelin implant, is often used. However, the commercially available deslorelin implants have an undesired delayed onset, as well as a potential brief increase in gonadotropin secretion after administration ("flare-up" effect) that can lead to oviposition before the actual suppression of gonadotropins. The objective of this study was to investigate whether the prolactin inhibitor cabergoline suppresses ovulation and whether it could be used to bridge the time until the onset of effect by the deslorelin implant. We measured the effect of cabergoline (30 µg/kg PO q24h × 14 days), deslorelin implants (4.7 mg SC), and a combination of both on egg laying and plasma prolactin concentrations in 37 quail (Coturnix japonica) over 6 weeks. Quail were divided into 4 groups: group DesCab (deslorelin implant and cabergoline oral; n = 9); group DesPlac (deslorelin implant and placebo oral; n = 9); group PlacCab (placebo implant and cabergoline oral; n = 9); and group PlacPlac (placebo implant and placebo oral; n = 10). Regular egg laying stopped in 100% (9/9) of birds in group DesCab and 78% (7/ 9) of birds in group DesPlac within 5 days of placing the deslorelin implant. No bird ceased egg production in group PlacCab (0/9), and 10% of birds ceased egg production intermittently in group PlacPlac (1/10). Treatment with the deslorelin implant (P < .001) and with cabergoline (P = .04) had a significant (negative) influence on plasma prolactin concentrations compared with the baseline. The interaction of deslorelin and cabergoline treatment, as well as time after initiation of treatment, did not have a significant effect on plasma prolactin concentrations. These results show that daily oral cabergoline has no significant influence on egg laying and only a minor biologically nonsignificant effect on lowering the relative plasma prolactin concentrations in quail.

Identifying Collisions of Various Molecularities in Molecular Dynamics Simulations.

We present a method based on kinetic molecular theory that identifies reactions of various molecularities in molecular dynamics (MD) simulations of bulk gases. The method allows characterization of the thermodynamic conditions at which higher than bimolecular reactions are a factor in the mechanisms of complex gas-phase chemistry. Starting with Bodenstein's definition of termolecular collisions we derive analytical expressions for the frequency of higher molecularity collisions. We have developed a relationship for the ratio of the frequencies of termolecular to bimolecular collisions in terms of the temperature, density, and collision times. To demonstrate the method, we used ReaxFF in LAMMPS to carry out MD simulations for NVT ensembles of mixtures of H2-O2 over the density range 120.2-332.7 kg m-3 and temperature range 3000-5000 K. The simulations yield ReaxFF-based predictions of the relative importance of termolecular collisions O2···H2···O2 and bimolecular collisions O2···H2 in the early chemistry of hydrogen combustion.

General application of Tolman's concept of activation energy.

We present a generalization of Tolman's concept of activation energy applicable to thermal and non-thermal reactions in molecular dynamics simulations of reactions in bulk gases. To illustrate the applicability of the method, molecular dynamics calculations were carried out for the NVT ensemble to determine the activation energies of O2 + H2 → H + HO2 and 2O2 + H2 → 2HO2 from MD simulation results for [H2]/[O2] = 1 at 3000 K using the reactive force field, ReaxFF. Assuming local thermodynamic equilibrium, we define the reaction cluster local energy, the energy of the atoms participating in an individual reaction, which is conserved. The generalized Tolman activation energy (GTEa) approach is applicable to reactions of any molecularity. Although we have applied GTEa for thermal conditions, it is applicable to chemistry occurring under non-thermal conditions because it rests upon local rather than global equilibrium. We have defined the transition configurations, unique points that define a seam separating reactants and products at which the local energies of the reactants and products become equal.

Anisotropy in surface-initiated melting of the triclinic molecular crystal 1,3,5-triamino-2,4,6-trinitrobenzene: A molecular dynamics study.

Surface-initiated melting of 1,3,5-triamino-2,4,6-trinitrobenzene (TATB), a triclinic molecular crystal, was investigated using molecular dynamics simulations. Simulations were performed for the three principal crystallographic planes exposed to vacuum, with the normal vectors to the planes given by b × c, c × a, and a × b (where a, b, and c define the edge vectors of the unit cell), denoted as (100), (010), and (001), respectively. The best estimate of the normal melting temperature for TATB is 851 ± 5 K. The nature and extent of disordering of the crystal-vacuum interface depend on the exposed crystallographic face, with the (001) face exhibiting incomplete melting and superheating. This is attributed to the anisotropy of the inter-molecular hydrogen bonding and the propensity of the crystal to form stacking faults in directions approximately perpendicular to the (100) and (010) faces. For all three crystal orientations, formation of molecular vacancies in the lattice at the crystal-vacuum (or crystal-quasi-liquid layer) interface precedes the complete loss of order at the interface.

Pressure effects on the relaxation of an excited nitromethane molecule in an argon bath.

Classical molecular dynamics simulations were performed to study the relaxation of nitromethane in an Ar bath (of 1000 atoms) at 300 K and pressures 10, 50, 75, 100, 125, 150, 300, and 400 atm. The molecule was instantaneously excited by statistically distributing 50 kcal/mol among the internal degrees of freedom. At each pressure, 1000 trajectories were integrated for 1000 ps, except for 10 atm, for which the integration time was 5000 ps. The computed ensemble-averaged rotational energy decay is ∼100 times faster than the vibrational energy decay. Both rotational and vibrational decay curves can be satisfactorily fit with the Lendvay-Schatz function, which involves two parameters: one for the initial rate and one for the curvature of the decay curve. The decay curves for all pressures exhibit positive curvature implying the rate slows as the molecule loses energy. The initial rotational relaxation rate is directly proportional to density over the interval of simulated densities, but the initial vibrational relaxation rate decreases with increasing density relative to the extrapolation of the limiting low-pressure proportionality to density. The initial vibrational relaxation rate and curvature are fit as functions of density. For the initial vibrational relaxation rate, the functional form of the fit arises from a combinatorial model for the frequency of nitromethane "simultaneously" colliding with multiple Ar atoms. Roll-off of the initial rate from its low-density extrapolation occurs because the cross section for collision events with L Ar atoms increases with L more slowly than L times the cross section for collision events with one Ar atom. The resulting density-dependent functions of the initial rate and curvature represent, reasonably well, all the vibrational decay curves except at the lowest density for which the functions overestimate the rate of decay. The decay over all gas phase densities is predicted by extrapolating the fits to condensed-phase densities.

Calculation of anharmonic couplings and THz linewidths in crystalline PETN.

We have developed a method for calculating the cubic anharmonic couplings in molecular crystals for normal modes with the zero wave vector in the framework of classical mechanics, and have applied it, combined with perturbation theory, to obtain the linewidths of all infrared absorption lines of crystalline pentaerythritol tetranitrate in the terahertz region (<100 cm(-1)). Contributions of the up- and down-conversion processes to the total linewidth were calculated. The computed linewidths are in qualitative agreement with experimental data and the results of molecular dynamics simulations. Quantum corrections to the linewidths in the terahertz region are shown to be negligible.

Molecular dynamics simulations of shock waves in hydroxyl-terminated polybutadiene melts: mechanical and structural responses.

The mechanical and structural responses of hydroxyl-terminated cis-1,4-polybutadiene melts to shock waves were investigated by means of all-atom non-reactive molecular dynamics simulations. The simulations were performed using the OPLS-AA force field but with the standard 12-6 Lennard-Jones potential replaced by the Buckingham exponential-6 potential to better represent the interactions at high compression. Monodisperse systems containing 64, 128, and 256 backbone carbon atoms were studied. Supported shock waves were generated by impacting the samples onto stationary pistons at impact velocities of 1.0, 1.5, 2.0, and 2.5 km s(-1), yielding shock pressures between approximately 2.8 GPa and 12.5 GPa. Single-molecule structural properties (squared radii of gyration, asphericity parameters, and orientational order parameters) and mechanical properties (density, shock pressure, shock temperature, and shear stress) were analyzed using a geometric binning scheme to obtain spatio-temporal resolution in the reference frame centered on the shock front. Our results indicate that while shear stress behind the shock front is relieved on a ∼0.5 ps time scale, a shock-induced transition to a glass-like state occurs with a concomitant increase of structural relaxation times by several orders of magnitude.

A classical trajectory study of the dissociation and isomerization of C2H5.

Motivated by photodissociation experiments in which non-RRKM nanosecond lifetimes of the ethyl radical were reported, we have performed a classical trajectory study of the dissociation and isomerization of C2H5 over the energy range 100-150 kcal/mol. We used a customized version of the AIREBO semiempirical potential (Stuart, S. J.; et al. J. Chem. Phys. 2000, 112, 6472-6486) to more accurately describe the gas-phase decomposition of C2H5. This study constitutes one of the first gas-phase applications of this potential form. At each energy, 10,000 trajectories were run and all underwent dissociation in less than 100 ps. The calculated dissociation rate constants are consistent with RRKM models; no evidence was found for nanosecond lifetimes. An analytic kinetics model of isomerization/dissociation competition was developed that incorporated incomplete mode mixing through a postulated divided phase space. The fits of the model to the trajectory data are good and represent the trajectory results in detail through repeated isomerizations at all energies. The model correctly displays single exponential decay at lower energies, but at higher energies, multiexponential decay due to incomplete mode mixing becomes more apparent. At both ends of the energy range, we carried out similar trajectory studies on CD2CH3 to examine isotopic scrambling. The results largely support the assumption that a H or a D atom is equally likely to dissociate from the mixed-isotope methyl end of the molecule. The calculated fraction of products that have the D atom dissociation is ∼20%, twice the experimental value available at one energy within our range. The calculated degree of isotopic scrambling is non-monotonic with respect to energy due to a non-monotonic ratio of the isomerization to dissociation rate constants.

Molecular dynamics study of the pressure-dependent terahertz infrared absorption spectrum of α- and γ-RDX.

Terahertz infrared absorption spectra of the α and γ polymorphs of 1,3,5-trinitro-1,3,5-triazacyclohexane (RDX) were predicted using two different theoretical approaches based on molecular dynamics simulations. The thermodynamic conditions studied were T = 298 K and hydrostatic pressures P = 0.0, 1.0, and 2.0 GPa for α-RDX and P = 3.0, 5.2, and 7.0 GPa for γ-RDX. The spectra obtained using the two methods are similar but not identical. In the case of α-RDX for pressure P = 0.0 GPa both spectra agree reasonably well with experimental data. The predicted spectra for α-RDX exhibit red-shifting (mode softening) of the main absorption peak with increasing pressure while for γ-RDX the spectra exhibit overall blue-shifting with increasing pressure.

Post-shock relaxation in crystalline nitromethane.

Molecular dynamics simulations of shocked (100)-oriented crystalline nitromethane were carried out to determine the rates of relaxation behind the shock wave. The forces were described by the fully flexible non-reactive Sorescu-Rice-Thompson force field [D. C. Sorescu, B. M. Rice, and D. L. Thompson, J. Phys. Chem. B 104, 8406 (2000)]. The time scales for local and overall thermal equilibration in the shocked crystal were determined. The molecular center-of-mass and atomic kinetic energy distributions rapidly reach substantially different local temperatures. Several picoseconds are required for the two distributions to converge, corresponding to establishment of thermal equilibrium in the shocked crystal. The decrease of the molecular center-of-mass temperature and the increase of the atomic temperature behind the shock front exhibit essentially exponential dependence on time. Analysis of covalent bond distance distributions ahead of, immediately behind, and well behind the shock front showed that the effective bond stretching potentials are essentially harmonic. Effective force constants for the C-N, C-H, and N-O bonds immediately behind the shock front are larger by factors of 1.6, 2.5, and 2.0, respectively, than in the unshocked crystal; and by factors of 1.2, 2.2, and 1.7, respectively, compared to material sufficiently far behind the shock front to be essentially at thermal equilibrium.

A classical trajectory study of the intramolecular dynamics, isomerization, and unimolecular dissociation of HO2.

The classical dynamics and rates of isomerization and dissociation of HO2 have been studied using two potential energy surfaces (PESs) based on interpolative fittings of ab initio data: An interpolative moving least-squares (IMLS) surface [A. Li, D. Xie, R. Dawes, A. W. Jasper, J. Ma, and H. Guo, J. Chem. Phys. 133, 144306 (2010)] and the cubic-spline-fitted PES reported by Xu, Xie, Zhang, Lin, and Guo (XXZLG) [J. Chem. Phys. 127, 024304 (2007)]. Both PESs are based on similar, though not identical, internally contracted multi-reference configuration interaction with Davidson correction (icMRCI+Q) electronic structure calculations; the IMLS PES includes complete basis set (CBS) extrapolation. The coordinate range of the IMLS PES is limited to non-reactive processes. Surfaces-of-section show similar generally regular phase space structures for the IMLS and XXZLG PESs with increasing energy. The intramolecular vibrational energy redistribution (IVR) at energies above and below the threshold of isomerization is slow, especially for O-O stretch excitations, consistent with the regularity in the surfaces-of-section. The slow IVR rates lead to mode-specific effects that are prominent for isomerization (on both the IMLS and XXZLG) and modest for unimolecular dissociation to H + O2 (accessible only on the XXZLG PES). Even with statistical distributions of initial energy, slow IVR rates result in double exponential decay for isomerization, with the slower rate correlated with slow IVR rates for O-O vibrational excitation. The IVR and isomerization rates computed for the IMLS and XXZLG PESs are quantitatively, but not qualitatively, different from one another with the largest differences ascribed to the ~2 kcal/mol difference in the isomerization barrier heights. The IMLS and XXZLG results are compared with those obtained using the global, semi-empirical double-many-body expansion DMBE-IV PES [M. R. Pastrana, L. A. M. Quintales, J. Brandão, and A. J. C. Varandas, J. Chem. Phys. 94, 8073 (1990)], for which the surfaces-of-section display more irregular phase space structure, much faster IVR rates, and significantly less mode-specific effects in isomerization and unimolecular dissociation. The calculated IVR results for all three PESs are reasonably well represented by an analytic, coupled three-mode energy transfer model.

Molecular dynamics simulations of shock waves in oriented nitromethane single crystals: plane-specific effects.

Molecular dynamics simulations of supported shock waves (shock pressure P(s) ∼ 15 GPa) propagating along the [110], [011], [101], and [111] directions in crystalline nitromethane initially at T = 200 K were performed using the nonreactive Sorescu-Rice-Thompson force field [D. C. Sorescu, B. M. Rice, and D. L. Thompson, J. Phys. Chem. B 104, 8406 (2000)]. These simulations, combined with those from a preceding study of shocks propagating along [100], [010], and [001] directions in nitromethane for similar conditions of temperature and shock pressure [L. He, T. D. Sewell, and D. L. Thompson, J. Chem. Phys. 134, 124506 (2011)], have been used to study the post-shock relaxation phenomena. Shocks along [010] and [101] lead to a crystal-crystal structure transformation. Shocks propagating along [011], [110], [111], [100], and [001] exhibit plane-specific disordering, which was characterized by calculating as functions of time the 1D mean square displacement (MSD), 2D radial distribution function (RDF), and 2D orientation order parameter P(2)(θ) in orthogonal planes mutually perpendicular to the shock plane; and by calculating as functions of distance behind the shock front the Cartesian components of intermolecular, intramolecular, and total kinetic energies. The 2D RDF results show that the structural disordering for shocks along [100], [110], and [111] is strongly plane-specific; whereas for shocks along [001] and [011], the loss of crystal structural order is almost equivalent in the orthogonal planes perpendicular to the shock plane. Based on the entire set of simulations, there is a trend for the most extensive disordering to occur in the (010) and (110) planes, less extensive disordering to occur in the (100) plane, and essentially no disordering to occur in the (001) plane. The 2D P(2)(θ) and 1D MSD profiles show, respectively, that the orientational and translational disordering is plane-specific, which results in the plane-specific structural disordering observed in the 2D RDF. By contrast, the kinetic energy partitioning and redistribution do not exhibit plane specificity, as shown by the similarity of spatial profiles of the Cartesian components of the intermolecular, intramolecular, and total kinetic energies in orthogonal planes perpendicular to the shock plane.

Molecular dynamics simulations of shock waves in oriented nitromethane single crystals.

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)].

Molecular dynamics study of the crystallization of nitromethane from the melt.

The crystallization of nitromethane, CH(3)NO(2), from the melt on the (100), (010), (001), and (110) crystal surfaces at 170, 180, 190, 200, 210, and 220 K has been investigated using constant-volume and -temperature (NVT) molecular dynamics simulations with a realistic, fully flexible force field [D. C. Sorescu, B. M. Rice, and D. L. Thompson, J. Phys. Chem. B 104, 8406 (2000)]. The crystallization process and the nature of the solid-liquid interface have been investigated by computing the molecular orientations, density, and radial distribution functions as functions of time and location in the simulation cell. During crystallization the translational motion of the molecules ceases first, after which molecular rotation ceases as the molecules assume proper orientations in the crystal lattice. The methyl groups are hindered rotors in the liquid; hindrance to rotation is reduced upon crystallization. The width of the solid-liquid interface varies between 6 and 13 Å (about two to five molecular layers) depending on which crystal surface is exposed to the melt and which order parameter is used to define the interface. The maximum rate of crystallization varies from 0.08 molecules ns(-1) Å(-2) for the (010) surface at 190 K to 0.41 molecules ns(-1) Å(-2) for the (001) surface at 220 K.

Effects of Epinephrine, Detomidine, and Butorphanol on Assessments of Insulin Sensitivity in Mares.

Sympathoadrenal stimulation may perturb results of endocrine tests performed on fractious horses. Sedation may be beneficial; however, perturbation of results may preclude useful information. Four experiments were designed to 1) determine the effects of epinephrine on insulin response to glucose (IR2G), 2) assess the effects of detomidine (DET), alone or combined with butorphanol (DET/BUT), on IR2G and glucose response to insulin (GR2I), and 3) assess the effects of BUT alone on IR2G. In Experiment 1, mares were administered saline or epinephrine (5 µg/kg BW) immediately before infusion of glucose (100 mg/kg BW). Glucose stimulated (P < .05) insulin release in controls at 5 minutes that persisted through 30 minutes; insulin was suppressed (P < .05) by epinephrine from 5 to 15 minutes, rising gradually through 30 minutes. Experiments 2 (IR2G) and 3 (GR2I) were conducted as triplicated 3 × 3 Latin squares with the following treatments: saline (SAL), DET, and DET/BUT (all administered at .01 mg/kg BW). Glucose stimulated (P < .05) insulin release that persisted through 30 minutes in SAL mares; DET and DET/BUT severely suppressed (P < .0001) the IR2G. Sedation did not affect resting glucose and had inconsistent effects on the GR2I when mares were treated with 50 mIU/kg BW recombinant human insulin. Butorphanol had no effect on IR2G. In conclusion, adrenergic agonists severely suppress the IR2G and cannot be used for sedation for this test. The use of DET did not alter the GR2I, and therefore may be useful for conducting this test in fractious horses.

Application of interpolating moving least squares fitting to hypervelocity collision dynamics: O(3P) + HCl.

We use an automated interpolating moving least-squares (IMLS) algorithm, which generates a fitted ab initio surface for systems of arbitrary topology, to construct a global OHCl ((3)A'') surface at the UB3LYP/aug-cc-pVTZ level of theory. This analytic PES includes all reaction channels and OHCl geometries with energies up to 144 kcal/mol (6.25 eV) above the O + HCl asymptote. The fitted surface was combined with the quasiclassical trajectory method to study the dynamics of the O((3)P) + HCl reaction at hyperthermal collision energies. The fitted PES greatly improves energy conservation during trajectory integration and eliminates problems with ab initio convergence, which are often encountered during direct dynamics studies. The more extensive trajectory calculations yield new insight into the title reaction and agree well with previous experimental studies and direct dynamics results.

Ab initio wavenumber accurate spectroscopy: 1CH2 and HCN vibrational levels on automatically generated IMLS potential energy surfaces.

We report here calculated J = 0 vibrational frequencies for (1)CH(2) and HCN with root-mean-square error relative to available measurements of 2.0 cm(-1) and 3.2 cm(-1), respectively. These results are obtained with DVR calculations with a dense grid on ab initio potential energy surfaces (PESs). The ab initio electronic structure calculations employed are Davidson-corrected MRCI calculations with double-, triple-, and quadruple-zeta basis sets extrapolated to the complete basis set (CBS) limit. In the (1)CH(2) case, Full CI tests of the Davidson correction at small basis set levels lead to a scaling of the correction with the bend angle that can be profitably applied at the CBS limit. Core-valence corrections are added derived from CCSD(T) calculations with and without frozen cores. Relativistic and non-Born-Oppenheimer corrections are available for HCN and were applied. CBS limit CCSD(T) and CASPT2 calculations with the same basis sets were also tried for HCN. The CCSD(T) results are noticeably less accurate than the MRCI results while the CASPT2 results are much poorer. The PESs were generated automatically using the local interpolative moving least-squares method (L-IMLS). A general triatomic code is described where the L-IMLS method is interfaced with several common electronic structure packages. All PESs were computed with this code running in parallel on eight processors. The L-IMLS method provides global and local fitting error measures important in automatically growing the PES from initial ab initio seed points. The reliability of this approach was tested for (1)CH(2) by comparing DVR-calculated vibrational levels on an L-IMLS ab initio surface with levels generated by an explicit ab initio calculation at each DVR grid point. For all levels ( approximately 200) below 20 000 cm(-1), the mean unsigned difference between the levels of these two calculations was 0.1 cm(-1), consistent with the L-IMLS estimated mean unsigned fitting error of 0.3 cm(-1). All L-IMLS PESs used in this work have comparable mean unsigned fitting errors, implying that fitting errors have a negligible role in the final errors of the computed vibrational levels with experiment. Less than 500 ab initio calculations of the energy and gradients are required to achieve this level of accuracy.

Shock-induced melting of (100)-oriented nitromethane: Energy partitioning and vibrational mode heating.

A study of the structural relaxation of nitromethane subsequent to shock loading normal to the (100) crystal plane performed using molecular dynamics and a nonreactive potential was reported recently [J. Chem. Phys. 131, 064503 (2009)]. Starting from initial temperatures of T(0)=50 and 200 K, shocks were simulated using impact velocities U(p) ranging from 0.5 to 3.0 km s(-1); clear evidence of melting was obtained for shocks initiated with impacts of 2.0 km s(-1) and higher. Here, we report the results of analyses of those simulation data using a method based on the Eckart frame normal-mode analysis that allows partitioning of the kinetic energy among the molecular degrees of freedom. A description of the energy transfer is obtained in terms of average translational and rotational kinetic energies in addition to the rates of individual vibrational mode heating. Generally, at early times postshock a large superheating of the translational and rotational degrees of freedom (corresponding to phonon modes of the crystal) is observed. The lowest frequency vibrations (gateway modes) are rapidly excited and also exhibit superheating. Excitation of the remaining vibrational modes occurs more slowly. A rapid, early excitation of the symmetric C-H stretch mode was observed for the shock conditions T(0)=50 K and U(p)=2.0 km s(-1) due to a combination of favorable alignment of molecular orientation with the shock direction and frequency matching between the vibration and shock velocity.

Shock-induced melting of (100)-oriented nitromethane: structural relaxation.

Molecules subjected to shock waves will, in general, undergo significant intramolecular distortion and exhibit large amplitude orientational and translational displacements relative to the unshocked material. The analysis of molecular dynamics simulations of strongly perturbed materials is complicated, particularly when the goal is to express time-dependent molecular-scale properties in terms of structural or geometric descriptors/properties defined for molecules in the equilibrium geometry. We illustrate the use of the Eckart-Sayvetz condition in a molecular dynamics study of the response of crystalline nitromethane subjected to supported shock waves propagating normal to (100). The simulations were performed with the nonreactive but vibrationally accurate force field due to Sorescu et al. [J. Phys. Chem. B 104, 8406 (2000)]. Shocks were initiated with impact velocities of U(p)=0.5, 1.0, 2.0, and 3.0 km s(-1) in crystals at initial temperatures of T0=50 and 200 K. Statistical precision in the analysis was enhanced through the use of a spatiotemporal reference frame centered on the advancing shock front, which was located as a function of time using the gradient of the kinetic energy along the shock direction. The Eckart-Sayvetz condition provides a rigorous approach by which the alignment can be obtained between a coordinate frame for a perturbed molecule and one in a convenient reference frame (e.g., one based on the equilibrium crystal structure) for analyses of the molecules in the material as the system evolves toward equilibrium. Structural and dynamic properties of the material corresponding to orientation in the lattice, translational symmetry, and mass transport (orientational order parameters, two dimensional radial distribution functions, and self-diffusion coefficients, respectively) were computed as functions of time with 4 fs resolution. The results provide clear evidence of melting for shocks initiated by impacts of at least U(p)=2.0 km s(-1) and provide insights into the evolution of changes at the molecular-mode level associated with the onset of the melting instability in shocked crystal.

Interpolating moving least-squares methods for fitting potential energy surfaces: using classical trajectories to explore configuration space.

We develop two approaches for growing a fitted potential energy surface (PES) by the interpolating moving least-squares (IMLS) technique using classical trajectories. We illustrate both approaches by calculating nitrous acid (HONO) cis-->trans isomerization trajectories under the control of ab initio forces from low-level HF/cc-pVDZ electronic structure calculations. In this illustrative example, as few as 300 ab initio energy/gradient calculations are required to converge the isomerization rate constant at a fixed energy to approximately 10%. Neither approach requires any preliminary electronic structure calculations or initial approximate representation of the PES (beyond information required for trajectory initial conditions). Hessians are not required. Both approaches rely on the fitting error estimation properties of IMLS fits. The first approach, called IMLS-accelerated direct dynamics, propagates individual trajectories directly with no preliminary exploratory trajectories. The PES is grown "on the fly" with the computation of new ab initio data only when a fitting error estimate exceeds a prescribed tight tolerance. The second approach, called dynamics-driven IMLS fitting, uses relatively inexpensive exploratory trajectories to both determine and fit the dynamically accessible configuration space. Once exploratory trajectories no longer find configurations with fitting error estimates higher than the designated accuracy, the IMLS fit is considered to be complete and usable in classical trajectory calculations or other applications.