Marc Cawkwell

Homogeneous dislocation nucleation in cyclotrimethylene trinitramine under shock loading

The propagation of shock waves normal to (111) in the energetic molecular crystal cyclotrimethylene trinitramine (RDX) has been studied using large-scale molecular dynamics simulations. Partial dislocation loops with Burgers vector 0.16[010] are nucleated homogeneously on (001) at Rankine–Hugoniot shock pressures greater than 1.3 GPa. Calculations of the [010] cross-section of the (001) generalized stacking fault energy surface as a function of applied pressure along [001] reveals that the stacking fault enclosed by the partial dislocation loops is rendered metastable by a stress-induced change in molecular conformation. Furthermore, large-scale molecular dynamics simulations performed on quasi-two-dimensional (111)-oriented single crystals show a two-wave elastic-plastic response with a “galloping” plastic wave. We propose that the onset of homogeneous dislocation nucleation accounts for the abrupt change in the elastic-plastic response of macroscopic (111)-oriented RDX single crystals observed in recent...

Energy conserving, linear scaling Born-Oppenheimer molecular dynamics

Born-Oppenheimer molecular dynamics simulations with long-term conservation of the total energy and a computational cost that scales linearly with system size have been obtained simultaneously. Linear scaling with a low pre-factor is achieved using density matrix purification with sparse matrix algebra and a numerical threshold on matrix elements. The extended Lagrangian Born-Oppenheimer molecular dynamics formalism [A. M. N. Niklasson, Phys. Rev. Lett. 100, 123004 (2008)] yields microcanonical trajectories with the approximate forces obtained from the linear scaling method that exhibit no systematic drift over hundreds of picoseconds and which are indistinguishable from trajectories computed using exact forces.

Anomalous hardening under shock compression in (021)-oriented cyclotrimethylene trinitramine single crystals

We recently proposed that the change observed in the elastic-plastic response of (111)-oriented cyclotrimethylene trinitramine (RDX) crystals under shock compression is caused by an anomalous hardening that is mediated by the homogeneous nucleation of partial dislocation loops with Burgers vector 0.16[010] on (001) {Cawkwell et al., [J. Appl. Phys. 107, 063512 (2010)]}. The orientation dependencies of the (001)[010] slip system suggested that (021)-oriented RDX crystals should also display an anomalous hardening. Molecular dynamics simulations of (021)-oriented RDX crystals confirm that this slip system is activated at a shock pressure 1.34<P≤1.54 GPa. Plate impact experiments on (021)-oriented RDX single crystals show a two-wave elastic-plastic response at 1.0 GPa and an almost overdriven response at 2.25 GPa that is entirely consistent with the theoretical prediction.

Generalized extended Lagrangian Born-Oppenheimer molecular dynamics

Extended Lagrangian Born-Oppenheimer molecular dynamics based on Kohn-Sham density functional theory is generalized in the limit of vanishing self-consistent field optimization prior to the force evaluations. The equations of motion are derived directly from the extended Lagrangian under the condition of an adiabatic separation between the nuclear and the electronic degrees of freedom. We show how this separation is automatically fulfilled and system independent. The generalized equations of motion require only one diagonalization per time step and are applicable to a broader range of materials with improved accuracy and stability compared to previous formulations.

Computing the Density Matrix in Electronic Structure Theory on Graphics Processing Units.

The self-consistent solution of a Schrödinger-like equation for the density matrix is a critical and computationally demanding step in quantum-based models of interatomic bonding. This step was tackled historically via the diagonalization of the Hamiltonian. We have investigated the performance and accuracy of the second-order spectral projection (SP2) algorithm for the computation of the density matrix via a recursive expansion of the Fermi operator in a series of generalized matrix-matrix multiplications. We demonstrate that owing to its simplicity, the SP2 algorithm [Niklasson, A. M. N. Phys. Rev. B2002, 66, 155115] is exceptionally well suited to implementation on graphics processing units (GPUs). The performance in double and single precision arithmetic of a hybrid GPU/central processing unit (CPU) and full GPU implementation of the SP2 algorithm exceed those of a CPU-only implementation of the SP2 algorithm and traditional matrix diagonalization when the dimensions of the matrices exceed about 2000 × 2000. Padding schemes for arrays allocated in the GPU memory that optimize the performance of the CUBLAS implementations of the level 3 BLAS DGEMM and SGEMM subroutines for generalized matrix-matrix multiplications are described in detail. The analysis of the relative performance of the hybrid CPU/GPU and full GPU implementations indicate that the transfer of arrays between the GPU and CPU constitutes only a small fraction of the total computation time. The errors measured in the self-consistent density matrices computed using the SP2 algorithm are generally smaller than those measured in matrices computed via diagonalization. Furthermore, the errors in the density matrices computed using the SP2 algorithm do not exhibit any dependence of system size, whereas the errors increase linearly with the number of orbitals when diagonalization is employed.

Graph-based linear scaling electronic structure theory

We show how graph theory can be combined with quantum theory to calculate the electronic structure of large complex systems. The graph formalism is general and applicable to a broad range of electronic structure methods and materials, including challenging systems such as biomolecules. The methodology combines well-controlled accuracy, low computational cost, and natural low-communication parallelism. This combination addresses substantial shortcomings of linear scaling electronic structure theory, in particular with respect to quantum-based molecular dynamics simulations.

Extended Lagrangian Born-Oppenheimer molecular dynamics simulations of the shock-induced chemistry of phenylacetylene

The initial chemical events that occur during the shock compression of liquid phenylacetylene have been investigated using self-consistent tight binding molecular dynamics simulations. The extended Lagrangian Born-Oppenheimer molecular dynamics formalism enabled us to compute microcanonical trajectories with precise conservation of the total energy. Our simulations revealed that the first density-increasing step under shock compression arises from the polymerization of phenylacetylene molecules at the acetylene moiety. The application of electronic structure-based molecular dynamics with long-term conservation of the total energy enabled us to identify electronic signatures of reactivity via monitoring changes in the HOMO-LUMO gap, and to capture directly adiabatic shock heating, transient non-equilibrium states, and changes in temperature arising from exothermic chemistry in classical molecular dynamics trajectories.

Atomistic study of ordinary 1/2<110] screw dislocations in single-phase and lamellar γ-TiAl

Computer simulation of the core structure and glide of ordinary screw dislocations in single-phase L10 TiAl and in two lamellae forming a twin γ/γ-interface has been performed using recently constructed Bond-Order Potentials (BOPs). BOPs represent a semi-empirical, numerically efficient scheme that works within the orthogonal tight-binding approximation and is able to capture the directionality of bonding. We have studied dislocation glide in perfect L10 TiAl and along a twin interface, transmission of an ordinary screw dislocation between lamellae, and the core structure, mobility and detachment of an interfacial screw dislocation from a twin boundary under applied shear stresses in directions parallel and perpendicular to a (111) plane. Our results show that the glide of ordinary straight screw dislocations under applied stresses in L10 TiAl is characterized by zigzag movement on two conjugated {111} planes. The non-planar core of the screw dislocation is distorted asymmetrically when the elastic centre of the dislocation is close to a twin γ/γ-interface and the dislocation moves on one of the (111) planes, depending on the magnitude of the corresponding Schmid factor. Ordinary dislocations become ordinary interfacial dislocations when they reach the interface. With increasing applied stress they can glide into the adjacent lamella, leaving no remnant interfacial dislocation.

Geometric integration in Born-Oppenheimer molecular dynamics

Geometric integration schemes for extended Lagrangian self-consistent Born-Oppenheimer molecular dynamics, including a weak dissipation to remove numerical noise, are developed and analyzed. The extended Lagrangian framework enables the geometric integration of both the nuclear and electronic degrees of freedom. This provides highly efficient simulations that are stable and energy conserving even under incomplete and approximate self-consistent field (SCF) convergence. We investigate three different geometric integration schemes: (1) regular time reversible Verlet, (2) second order optimal symplectic, and (3) third order optimal symplectic. We look at energy conservation, accuracy, and stability as a function of dissipation, integration time step, and SCF convergence. We find that the inclusion of dissipation in the symplectic integration methods gives an efficient damping of numerical noise or perturbations that otherwise may accumulate from finite arithmetics in a perfect reversible dynamics.

Efficient parallel linear scaling construction of the density matrix for Born-Oppenheimer molecular dynamics.

We present an algorithm for the calculation of the density matrix that for insulators scales linearly with system size and parallelizes efficiently on multicore, shared memory platforms with small and controllable numerical errors. The algorithm is based on an implementation of the second-order spectral projection (SP2) algorithm [ Niklasson, A. M. N. Phys. Rev. B 2002 , 66 , 155115 ] in sparse matrix algebra with the ELLPACK-R data format. We illustrate the performance of the algorithm within self-consistent tight binding theory by total energy calculations of gas phase poly(ethylene) molecules and periodic liquid water systems containing up to 15,000 atoms on up to 16 CPU cores. We consider algorithm-specific performance aspects, such as local vs nonlocal memory access and the degree of matrix sparsity. Comparisons to sparse matrix algebra implementations using off-the-shelf libraries on multicore CPUs, graphics processing units (GPUs), and the Intel many integrated core (MIC) architecture are also presented. The accuracy and stability of the algorithm are illustrated with long duration Born-Oppenheimer molecular dynamics simulations of 1000 water molecules and a 303 atom Trp cage protein solvated by 2682 water molecules.