Jeremy Schofield

Using a classical potential as an efficient importance function for sampling from an ab initio potential

In this paper the ab initio potential of mean force for the formic acid–water system is calculated in a Monte Carlo simulation using a classical fluctuating charge molecular mechanics potential to guide Monte Carlo updates. The ab initio energies in the simulation are calculated using density-functional theory (DFT) methods recently developed by Salahub et al. [J. Chem. Phys. 107, 6770 (1997)] to describe hydrogen-bonded systems. Importance sampling methods are used to investigate structural changes and it is demonstrated that using a molecular mechanics importance function can improve the efficiency of a DFT simulation by several orders of magnitude. Monte Carlo simulation of the system in a canonical ensemble at T=300 K reveals two chemical processes at intermediate time scales: The rotation of the H2O bonded to HCOOH, which takes place on a time scale of 3 ps, and the dissociation of the complex which occurs in 24 ps. It is shown that these are the only important structural “reactions” in the formic ac...

Mixed quantum-classical molecular dynamics: Aspects of the multithreads algorithm

The mixed quantum-classical Liouville equation is derived from a semiclassical perspective starting from the full quantum Schrodinger equation. An asymptotic numerical scheme for solving the equation is discussed which relies on propagating swarms of interacting “threads” which represent the density matrix or other observable. It is demonstrated that this “multithreads” method performs extremely well on simple one-dimensional model systems designed to test nonadiabatic molecular dynamic methods, yielding essentially exact results for a variety of models.

Solutions of mixed quantum-classical dynamics in multiple dimensions using classical trajectories

The multithreads algorithm for solving the mixed quantum-classical Liouville equation is extended to systems in which multiple classical degrees of freedom couple explicitly to a quantum subsystem. The method involves evolving a discrete set of matrices representing operators positioned at classical phase space coordinates according to precise dynamical rules dictated by evolution equations. The propagation scheme is based on the Trotter expansion of the time evolution operator and involves trajectory (thread) branching and pruning operations at each time step. The method is tested against exact numerical solution of the quantum dynamics for two models in which the nonadiabatic evolution of two heavy coordinates (nuclei) induces changes in population in two electronic states. It is demonstrated that the multithreads algorithm provides a good quantitative as well as qualitative description of the dynamics for branching ratios and populations as a function of time. Critical performance issues such as the computational demand of the method, energy conservation, and how the scheme scales with the number of classical degrees of freedom coupled to the quantum subsystem are discussed.

Mapping quantum-classical Liouville equation: Projectors and trajectories

The evolution of a mixed quantum-classical system is expressed in the mapping formalism where discrete quantum states are mapped onto oscillator states, resulting in a phase space description of the quantum degrees of freedom. By defining projection operators onto the mapping states corresponding to the physical quantum states, it is shown that the mapping quantum-classical Liouville operator commutes with the projection operator so that the dynamics is confined to the physical space. It is also shown that a trajectory-based solution of this equation can be constructed that requires the simulation of an ensemble of entangled trajectories. An approximation to this evolution equation which retains only the Poisson bracket contribution to the evolution operator does admit a solution in an ensemble of independent trajectories but it is shown that this operator does not commute with the projection operators and the dynamics may take the system outside the physical space. The dynamical instabilities, utility, and domain of validity of this approximate dynamics are discussed. The effects are illustrated by simulations on several quantum systems.

Exact and asymptotic solutions of the mixed quantum-classical Liouville equation

In this article, an exact surface-hopping procedure and an approximate asymptotic method for performing molecular dynamics based on a mixed quantum-classical Liouville equation [J. Chem. Phys. 110, 8919 (1999)] for partially Wigner transformed dynamical variables of a coupled quantum subsystem and classical bath are elaborated. The methods are based upon writing the equations of motion in a basis set in which quantum transitions do not alter the classical trajectory, and therefore avoid ad-hoc momentum jump approximations and are free of singular kernels associated with sampling momenta. Results obtained utilizing the new trajectory methods are presented for a model two-level system bilinearly coupled to a classical harmonic oscillator. These results are compared to results obtained from standard methods of performing mixed quantum-classical dynamics. The new methods perform well for the model system over a wide range of initial kinetic energies.

Numerical implementation of the exact dynamics of free rigid bodies

In this paper the exact analytical solution of the motion of a rigid body with arbitrary mass distribution is derived in the absence of forces or torques. The resulting expressions are cast into a form where the dependence of the motion on initial conditions is explicit and the equations governing the orientation of the body involve only real numbers. Based on these results, an efficient method to calculate the location and orientation of the rigid body at arbitrary times is presented. This implementation can be used to verify the accuracy of numerical integration schemes for rigid bodies, to serve as a building block for event-driven discontinuous molecular dynamics simulations of general rigid bodies, and for constructing symplectic integrators for rigid body dynamics.

Molecular dynamics simulation of a graphite-supported copper nanocluster: thermodynamic properties and gas adsorption

Molecular dynamics simulations were conducted for a cubic Cu cluster supported on a graphite bilayer. The Sutten–Chen and Lennard–Jones potentials were used for metal–metal and metal–graphite interactions, respectively. Heating and cooling processes were performed by NVT simulations at different temperatures in the range 200 to 1800 K. The melting point was identified on the basis of caloric and heat capacity curves. The calculated melting point was 770 K, far below the bulk melting point of crystalline copper. Several phenomena such as the appearance of a hysteresis (irreversibility) in caloric curves, surface melting, and cluster-induced surface wetting were justified from the results. The simulation of cluster in the presence of gas atmosphere showed that the CO gas is adsorbed more than H2 and it has a greater impact on the clusters structure.

Discontinuous molecular dynamics for semiflexible and rigid bodies

A general framework for performing event-driven simulations of systems with semiflexible or rigid bodies interacting under impulsive forces is outlined. The method consists of specifying a means of computing the free evolution of constrained motion, evaluating the times at which interactions occur, and determining the consequences of interactions on subsequent motion. Algorithms for computing the times of interaction events and carrying out efficient event-driven simulations are discussed. The semiflexible case and the rigid case differ qualitatively in that the free motion of a rigid body can be computed analytically and need not be integrated numerically.

Mode coupling and generalized hydrodynamics

In this paper a general hydrodynamic mode coupling theory of equilibrium fluctuations in simple liquids is developed from molecular considerations. The approach developed here avoids the shortcomings of previous mode coupling theories by adopting a complete hierarchy of equations for the slow modes of a hydrodynamic system and solving it formally through the N ordering approximation scheme developed previously by Machta and Oppenheim. A series is obtained from the hierarchy of equations which allows the generalized transport coefficients to be obtained exactly in the thermodynamic limit up to arbitrary order in the wavevector, frequency and mode coupling parameters. A self-consistent equation for the transport coefficients is formulated up to Burnett order in the wavevector and all orders in the mode coupling parameter. The results are similar to phenomenological models used by previous researchers.

Efficient ab initio sampling methods in rate constant calculations for proton-transfer reactions

In this article, the classical potential based importance Monte Carlo sampling method of Iftimie et al. [J. Chem. Phys. 113, 4852, (2000)] is applied to an ab initio simulation of the proton transfer tautomerization reaction of malonaldehyde in an aprotic, nonpolar solvent. It is demonstrated that ad hoc bond-energy bond-order relations derived from bond evolution theory combined with Pauling’s valence bond ideas can be used to construct a molecular mechanics guidance potential for the simulation of the proton transfer reaction which improves the statistics of the calculation by three orders of magnitude. The sampling method is extended to simulations in which quantum effects are treated using the imaginary time path-integral representation. A new algorithm based on multiple Markov chain theory is introduced by which it is possible to obtain very short integrated correlation lengths in calculations of quantum static correlation functions.