Improving the Performance of Hybrid Functional-Based Molecular Dynamics Simulation through Screening of Hartree-Fock Exchange Forces

Speeding-up Hybrid Functional-Based Ab Initio Molecular Dynamics Using Multiple Time-stepping and Resonance-Free Thermostat

Ab initio molecular dynamics (AIMD) based on density functional theory (DFT) has become a workhorse for studying the structure, dynamics, and reactions in condensed matter systems. Currently, AIMD simulations are primarily carried out at the level of generalized gradient approximation (GGA), which is at the second rung of DFT functionals in terms of accuracy. Hybrid DFT functionals, which form the fourth rung in the accuracy ladder, are not commonly used in AIMD simulations as the computational cost involved is 100 times or higher. To facilitate AIMD simulations with hybrid functionals, we propose here an approach using multiple time stepping with adaptively compressed exchange operator and resonance-free thermostat, that could speed up the calculations by ∼30 times or more for systems with a few hundred of atoms. We demonstrate that by achieving this significant speed up and making the compute time of hybrid functional-based AIMD simulations at par with that of GGA functionals, we are able to study several complex condensed matter systems and model chemical reactions in solution with hybrid functionals that were earlier unthinkable to be performed.

Improving the scaling and performance of multiple time stepping-based molecular dynamics with hybrid density functionals

Density functionals at the level of the generalized gradient approximation (GGA) and a plane-wave basis set are widely used today to perform ab initio molecular dynamics (AIMD) simulations. Going up in the ladder of accuracy of density functionals from GGA (second rung) to hybrid density functionals (fourth rung) is much desired pertaining to the accuracy of the latter in describing structure, dynamics, and energetics of molecular and condensed matter systems. On the other hand, hybrid density functional based AIMD simulations are about two orders of magnitude slower than GGA based AIMD for systems containing ~100 atoms using ~100 compute cores. Two methods, namely MTACE and s-MTACE, based on a multiple time step integrator and adaptively compressed exchange operator formalism are able to provide a speed-up of about 7-9 in performing hybrid density functional based AIMD. In this work, we report an implementation of these methods using a task-group based parallelization within the CPMD program package, with the intention to take advantage of the large number of compute cores available on modern high-performance computing platforms. We present here the boost in performance achieved through this algorithm. This work also identifies the computational bottleneck in the s-MTACE method and proposes a way to overcome it.

Arbitrarily accurate quantum alchemy

Doping compounds can be considered a perturbation to the nuclear charges in a molecular Hamiltonian. Expansions of this perturbation in a Taylor series, i.e., quantum alchemy, have been used in the literature to assess millions of derivative compounds at once rather than enumerating them in costly quantum chemistry calculations. So far, it was unclear whether this series even converges for small molecules, whether it can be used for geometry relaxation, and how strong this perturbation may be to still obtain convergent numbers. This work provides numerical evidence that this expansion converges and recovers the self-consistent energy of Hartree-Fock calculations. The convergence radius of this expansion is quantified for dimer examples and systematically evaluated for different basis sets, allowing for estimates of the chemical space that can be covered by perturbing one reference calculation alone. Besides electronic energy, convergence is shown for density matrix elements, molecular orbital energies, and density profiles, even for large changes in electronic structure, e.g., transforming He₃ into H₆. Subsequently, mixed alchemical and spatial derivatives are used to relax H₂ from the electronic structure of He alone, highlighting a path to spatially relaxed quantum alchemy. Finally, the underlying code that allows for arbitrarily accurate evaluation of restricted Hartree-Fock energies and arbitrary order derivatives is made available to support future method development.

Hidden Hemibonding in the Aqueous Hydroxyl Radical

The existence of a two-center, three-electron hemibond in the first solvation shell of ^•OH(aq) has long been a matter of debate. The hemibond manifests in ab initio molecular dynamics simulations as a small-r feature in the oxygen radial distribution function (RDF) for H₂O···^•OH, but that feature disappears when semilocal density functionals are replaced with hybrids, suggesting a self-interaction artifact. Using periodic simulations at the PBE0+D3 level, we demonstrate that the hemibond is actually still present (as evidenced by delocalization of the spin density) but is obscured by the hydrogen-bonded feature in the RDF due to a slight elongation of the hemibond. Computed electronic spectra for ^•OH(aq) are in excellent agreement with experiment and confirm that hemibond-like configurations play an outsized role in the spectroscopy due to an intense charge-transfer transition that is strongly attenuated in hydrogen-bonded configurations. Apparently, 25% exact exchange (as in PBE0) is insufficient to eliminate delocalization of unpaired spins.

Achieving an Order of Magnitude Speedup in Hybrid-Functional- and Plane-Wave-Based Ab Initio Molecular Dynamics: Applications to Proton-Transfer Reactions in Enzymes and in Solution

Ab initio molecular dynamics (MD) with hybrid density functionals and a plane wave basis is computationally expensive due to the high computational cost of exact exchange energy evaluation. Recently, we proposed a strategy to combine adaptively compressed exchange (ACE) operator formulation and a multiple time step integration scheme to reduce the computational cost significantly [J. Chem. Phys. 2019, 151, 151102 ]. However, it was found that the construction of the ACE operator, which has to be done at least once in every MD time step, is computationally expensive. In the present work, systematic improvements are introduced to further speed up by employing localized orbitals for the construction of the ACE operator. By this, we could achieve a computational speedup of an order of magnitude for a periodic system containing 32 water molecules. Benchmark calculations were carried out to show the accuracy and efficiency of the method in predicting the structural and dynamical properties of bulk water. To demonstrate the applicability, computationally intensive free-energy computations at the level of hybrid density functional theory were performed to investigate (a) methyl formate hydrolysis reaction in neutral aqueous media and (b) proton-transfer reaction within the active-site residues of the class C β-lactamase enzyme.

Role of hemibonding in the structure and ultraviolet spectroscopy of the aqueous hydroxyl radical

The presence of a two-center, three-electron hemibond in the solvation structure of the aqueous hydroxl radical has long been debated, as its appearance can be sensitive to self-interaction error in density functional theory.