Search

For electronic structure calculations on large systems, solving the self-consistent-field (SCF) equations is one of the key bottlenecks. Density matrix search methods provide an efficient linear scaling approach for circumventing the traditional O(N3) diagonalization procedure for solving the SCF equations. The conjugate gradient density matrix search (CG-DMS) method is a successful implementation of this approach. An alternative density matrix search method, QN–DMS, employs direct inversion in the iterative subspace using a quasi-Newton (QN) step to estimate the error vector. For linear polyglycine chains of 10–100 residues, the present approach scales linearly and is 30% faster than CG-DMS. For clusters of up to 300 water molecules, this method shows smoother and efficient convergence, and displays nearly linear scaling.

In the context of the recently developed Atom-centered Density Matrix Propagation (ADMP) approach to ab initio molecular dynamics, a formal analysis of the deviations from the Born—Oppenheimer surface is conducted. These deviations depend on the fictitious mass and on the magnitude of the commutator of the Fock and density matrices. These quantities are found to be closely interrelated and the choice of the fictitious mass provides a lower bound on the deviations from the Born—Oppenheimer surface. The relations are illustrated with an example calculation for the Cl−(H2O)25 cluster. We also show that there exists a direct one-to-one correspondence between approximate Born—Oppenheimer dynamics, where SCF convergence is restricted by a chosen threshold value for the commutator of the Fock and density matrices, and extended Lagrangian dynamics performed using a finite value for the fictitious mass. The analysis is extended to the nuclear forces used in the ADMP approximation. The forces are shown to be more general than those standardly used in Born—Oppenheimer dynamics, with the addition terms in the nuclear forces depending on the commutator mentioned above.

In a recently developed approach to ab initio molecular dynamics (ADMP), we used an extended Lagrangian to propagate the density matrix in a basis of atom centered Gaussian functions. Results of trajectory calculations obtained by this method are compared with the Born–Oppenheimer approach (BO), in which the density is converged at each step rather than propagated. For NaCl, the vibrational frequency with ADMP is found to be independent of the fictitious electronic mass and to be equal to the BO trajectory result. For the photodissociation of formaldehyde, H2CO→H2+CO, and the three body dissociation of glyoxal, C2H2O2→H2+2CO, very good agreement is found between the Born–Oppenheimer trajectories and the extended Lagrangian approach in terms of the rotational and vibrational energy distributions of the products. A 1.2 ps simulation of the dynamics of chloride ion in a cluster of 25 water molecules was used as a third test case. The Fourier transform of the velocity–velocity autocorrelation function showed ...

The photodissociation of s-tetrazine via a three-body fragmentation channel (C2N2H2→2HCN+N2) has been studied by ab initio direct classical trajectory calculations using Hartree–Fock and density functional methods with split valence and polarized basis sets [HF/3-21G, HF/6-31G(d) and B3LYP/6-31G(d)]. The calculated transition state is planar. At our most reliable method (CBS-APNO), the heat of reaction and barrier height are −53.0 kcal mol−1 and 41.1 kcal mol−1, respectively. To simulate the experimental photolysis of s-tetrazine, trajectories were started from a microcanonical ensemble at the transition state with 12 kcal mol−1 excess energy distributed among the vibrational modes and the transition vector. At all levels of theory, the HCN product has a very broad rotational distribution, ranging up to J=64, and has extensive excitation of the bending vibration. By contrast, N2 is produced with low J and with only a small amount of vibrational excitation in agreement with experiments. In accord with the experiment, the relative translation motion of the products receives about 80–83% of the available energy.

Photodissociation of glyoxal via the H2CO+CO channel has been investigated by ab initio classical trajectory calculations using Becke’s three-parameter hybrid functional method with split valence and polarized basis set [B3LYP/6-311G(d,p)]. To model the experimental conditions, trajectories were started from a microcanonical ensemble at the transition state with 8.5 kcal/mol excess energy distributed among the vibrational modes and the transition vector. The CO product was produced with a broad rotational distribution but with almost no vibration excitation. When combined with the results from the H2+2CO channel, the calculated vibrational and rotational distributions of CO are in excellent agreement with the experimental observations. The rotational distribution of H2CO was very broad ranging up to J=85. The H2CO product has significant vibrational excitation in the out-of-plane bending, CH2 rocking, CH2 scissoring, and CO stretching modes. For both the H2+2CO and the CO+H2CO channels, the majority of av...

We propose and implement an alternative approach to the original Car–Parrinello method where the density matrix elements (instead of the molecular orbitals) are propagated together with the nuclear degrees of freedom. Our new approach has the advantage of leading to an O(N) computational scheme in the large system limit. Our implementation is based on atom-centered Gaussian orbitals, which are especially suited to deal effectively with general molecular systems. The methodology is illustrated by applications to the three-body dissociation of triazine and to the dynamics of a cluster of a chloride ion with 25 water molecules.

For the integration of the classical equations of motion in the Born–Oppenheimer approach, each time the energy and gradient of the potential energy surface are needed, a properly converged wave function is calculated. If Hessians (second derivatives) can be calculated, significantly larger steps can be taken in the numerical integration of the equations of motion without loss of accuracy. Even larger steps can be taken with a Hessian-based predictor–corrector algorithm. Since updated Hessians are used successfully in quasi-Newton methods for geometry optimization, it should be possible to improve the performance of trajectory calculations using updated Hessians. The Murtagh–Sargent (MS) update, the Powell-symmetric–Broyden (PSB) update and Bofill’s update (a weighted combination of MS and PSB) were tested, and Bofill’s update was found to be the best. Slightly smaller step sizes were needed with Hessian updating to maintain good conservation of the energy, but this was more than compensated by the reduction in total computational cost. An overall factor of 3 in speed-up was obtained for trajectories of systems containing 4 to 6 heavy atoms computed at the HF/3-21G level.

Classical trajectories can be computed directly from electronic structure calculations without constructing a global potential-energy surface. When the potential energy and its derivatives are needed during the integration of the classical equations of motion, they are calculated by electronic structure methods. In the Born–Oppenheimer approach the wave function is converged rather than propagated to generate a more accurate potential-energy surface. If analytic second derivatives (Hessians) can be computed, steps of moderate size can be taken by integrating the equations of motion on a local quadratic approximation to the surface (a second-order algorithm). A more accurate integration method is described that uses a second-order predictor step on a local quadratic surface, followed by a corrector step on a better local surface fitted to the energies, gradients, and Hessians computed at the beginning and end points of the predictor step. The electronic structure work per step is the same as the second-ord...

Advances in the computation of the Coulomb, exchange, and correlation contributions to Gaussian-based Hartree–Fock and density functional theory Hamiltonians have demonstrated near-linear scaling with molecular size for these steps. These advances leave the O(N3) diagonalization bottleneck as the rate determining step for very large systems. In this work, a conjugate gradient density matrix search (CG-DMS) method has been successfully extended and computationally implemented for use with first principles calculations. A Cholesky decomposition of the overlap matrix and its inverse is used to transform to and back from an orthonormal basis, which can be formed in near-linear time for sparse systems. Linear scaling of CPU time for the density matrix search and crossover of CPU time with diagonalization is demonstrated for polyglycine chains containing up to 493 atoms and water clusters up to 900 atoms.

Ion chromatography experiments on carbon clusters have provided a powerful tool for characterizing the products of the laser ablation of graphite. Using this technique, several families of carbon clusters have been observed, and their role in a plausible fullerene formation process has been hypothesized. In this work, we have examined the experimental mobility results from a theoretical perspective. Our most interesting finding is the existence of a family of three-dimensional 2 + 4 cycloaddition products whose members match the experimental mobilities of the so-called `ring III` family over a range of cluster sizes, whereas previous studies have asserted that the `ring III` clusters are planar. In agreement with previous research, we find that the `ring I` and `ring II` families consist of monocyclic and bicycle rings, respectively. However, these families should be broadly defined so as to include ring structures with carbon branches, because short carbon branches have only a negligible effect on cluster mobility. 28 refs., 6 figs., 6 tabs.

The metastable $^{4}\mathrm{\ensuremath{\Pi}}_{\mathit{u}}$ state of the ${\mathrm{Ne}}_{2}^{\mathrm{\ensuremath{-}}}$ anion is examined with the use of the numerical Hartree-Fock method. The results suggest that due to its similarity to the experimentally discovered $^{4}\mathrm{\ensuremath{\Pi}}_{\mathit{g}}$ state of ${\mathrm{He}}_{2}^{\mathrm{\ensuremath{-}}}$, this state should also be long lived.

The curvy-extended-Lagrangian molecular-dynamics (ELMD) approach [J. M. Herbert and M. Head-Gordon, J. Chem. Phys. 121, 11542 (2004)] is similar to atom-centered density-matrix propagation (ADMP) [H. B. Schlegel, J. M. Millam, S. S. Iyengar, G. A. Voth, A. D. Daniels, G. E. Scuseria, and M. J. Frisch, J. Chem. Phys. 114, 9758 (2001); S. S. Iyengar, H.B. Schlegel, J.M. Millam, G.A. Voth, G.E. Scuseria, and M.J. Frisch, ibid.115, 10291 (2001); H.B. Schlegel, S.S. Iyengar, X. Li, J.M. Millam, G.A. Voth, G.E. Scuseria, and M.J. Frisch, ibid. 117, 8694 (2002); S.S. Iyengar, H.B. Schlegel, G.A. Voth, J.M. Millam, G.E. Scuseria, and M.J. Frisch, Israel J. Chem. 42, 191 (2002)] and based on Car-Parrinello [Phys. Rev. Lett. 55, 2471 (1985)] extended-Lagrangian [H.C. Andersen, J. Chem. Phys. 72, 2384 (1980)] molecular dynamics. Similarities between curvy-ELMD and ADMP arise from using unconverged electronic single-particle density matrices within Gaussian basis functions as dynamical variables. Curvy-ELMD differs f...

A generalization is presented here for a newly developed approach to ab initio molecular dynamics, where the density matrix is propagated with Gaussian orbitals. Including a tensorial fictitious mass facilitates the use of larger time steps for the dynamics process. A rigorous analysis of energy conservation is presented and used to control the deviation of the fictitious dynamics trajectory from the corresponding Born–Oppenheimer dynamics trajectory. These generalizations are tested for the case of the Cl−(H2O)25 cluster. It is found that, even with hydrogen atoms present in the system, no thermostats are necessary to control the exchange of energy between the nuclear and the fictitious electronic degrees of freedom.