Daniel S. Lambrecht

Current Status of the AMOEBA Polarizable Force Field

Molecular force fields have been approaching a generational transition over the past several years, moving away from well-established and well-tuned, but intrinsically limited, fixed point charge models toward more intricate and expensive polarizable models that should allow more accurate description of molecular properties. The recently introduced AMOEBA force field is a leading publicly available example of this next generation of theoretical model, but to date, it has only received relatively limited validation, which we address here. We show that the AMOEBA force field is in fact a significant improvement over fixed charge models for small molecule structural and thermodynamic observables in particular, although further fine-tuning is necessary to describe solvation free energies of drug-like small molecules, dynamical properties away from ambient conditions, and possible improvements in aromatic interactions. State of the art electronic structure calculations reveal generally very good agreement with AMOEBA for demanding problems such as relative conformational energies of the alanine tetrapeptide and isomers of water sulfate complexes. AMOEBA is shown to be especially successful on protein-ligand binding and computational X-ray crystallography where polarization and accurate electrostatics are critical.

Linear-scaling atomic orbital-based second-order Møller–Plesset perturbation theory by rigorous integral screening criteria

A Laplace-transformed second-order Moller-Plesset perturbation theory (MP2) method is presented, which allows to achieve linear scaling of the computational effort with molecular size for electronically local structures. Also for systems with a delocalized electronic structure, a cubic or even quadratic scaling behavior is achieved. Numerically significant contributions to the atomic orbital (AO)-MP2 energy are preselected using the so-called multipole-based integral estimates (MBIE) introduced earlier by us [J. Chem. Phys. 123, 184102 (2005)]. Since MBIE provides rigorous upper bounds, numerical accuracy is fully controlled and the exact MP2 result is attained. While the choice of thresholds for a specific accuracy is only weakly dependent upon the molecular system, our AO-MP2 scheme offers the possibility for incremental thresholding: for only little additional computational expense, the numerical accuracy can be systematically converged. We illustrate this dependence upon numerical thresholds for the calculation of intermolecular interaction energies for the S22 test set. The efficiency and accuracy of our AO-MP2 method is demonstrated for linear alkanes, stacked DNA base pairs, and carbon nanotubes: e.g., for DNA systems the crossover toward conventional MP2 schemes occurs between one and two base pairs. In this way, it is for the first time possible to compute wave function-based correlation energies for systems containing more than 1000 atoms with 10 000 basis functions as illustrated for a 16 base pair DNA system on a single-core computer, where no empirical restrictions are introduced and numerical accuracy is fully preserved.

Rigorous integral screening for electron correlation methods

We derive rigorous multipole-based integral estimates (MBIE) in order to account for the distance dependence occurring in atomic-orbital (AO) formulations of electron correlation theory, where our focus is on AO-MP2 theory within a Laplace scheme. We find for the exact transformed integral products an extremely early onset of a linear-scaling behavior and a very small number of significant products. To preselect the significant integral products we adapt our MBIE method as rigorous upper bound. In this way it is possible to exploit the favorable scaling behavior observed and to reduce the scaling of estimated products asymptotically to linear, without sacrificing accuracy or reliability. By separating Coulomb- and exchange-type contractions only half-transformed integrals need to be computed. Furthermore, our scheme of rigorously preselecting transformed integral products via MBIE seems to offer particularly interesting perspectives for a direct formation of half- or fully transformed integrals by using multipole expansions and auxiliary basis sets.

Multipole-based integral estimates for the rigorous description of distance dependence in two-electron integrals

We derive multipole-based integral estimates (MBIE) as rigorous and tight upper bounds to four-center two-electron integrals in order to account for the 1/R distance decay between the charge distributions, which is missing in the Schwarz screening commonly used in ab initio methods. Our screening criteria are valid for all angular momenta and can be formulated for any order of multipoles. We have found the expansion limited to dipoles to be sufficiently tight for estimating the integrals in Hartree-Fock and density-functional theories, while the screening effort is negligible. For, e.g., a DNA fragment with 1052 atoms and 10,674 basis functions (6-31G*) the exchange part is faster by a factor of 2.1 as compared to the Schwarz screening both within our linear exchange scheme, whereas a smaller factor of 1.3 is gained for the Coulomb part within the continuous fast multipole method. Most importantly, our new MBIE screening is perfectly suited to exploit the strong distance decay of electron-correlation effects of at least 1/R4 in atomic-orbital-based formulations of correlation methods.

The Performance of Density Functionals for Sulfate–Water Clusters

The performance of 24 density functionals, Hartree-Fock, and MP2 is assessed with respect to the CCSD(T)/CBS* energetics of 49 sulfate-water clusters with three to six water molecules. Included among the density functionals are GGA, meta-GGA, hybrid GGA, hybrid meta-GGA, and double hybrid density functionals, as well as the LDA. Three types of dispersion corrections (VV10, XDM, and -D) are tested in conjunction with these functionals. The 26 methods are compared using the relative and binding energies of the sulfate-water clusters as the main criteria. It was discovered that a majority of the tested density functionals are unable to simultaneously capture the physics necessary to describe both the relative and binding energies of the anionic solvation clusters. The three density functionals with the best overall performance are XYG3, ωB97X-2, and XYGJ-OS. The only other density functional that performs comparably to these three double hybrids is M11. A majority of the density functionals that contain a fraction of exact exchange tend to perform well only for the relative energies, while functionals lacking exact exchange generally perform poorly with respect to both criteria. However, the meta-GGA functional, M11-L, stands out due to its superior performance for the relative energies. While dispersion correction functionals cannot replace the accuracy provided by MP2 correlation, it is shown that the proper combination of a hybrid GGA functional (LC-ωPBE) with a dispersion correction functional (VV10) can lead to drastic improvements in the binding energies of the parent functional, while preserving its performance with respect to the relative energies. Ultimately, however, MP2 has the best overall performance out of the 26 benchmarked methods.

Ab Initio Simulations Reveal that Reaction Dynamics Strongly Affect Product Selectivity for the Cracking of Alkanes over H-MFI

Product selectivity of alkane cracking catalysis in the H-MFI zeolite is investigated using both static and dynamic first-principles quantum mechanics/molecular mechanics simulations. These simulations account for the electrostatic- and shape-selective interactions in the zeolite and provide enthalpic barriers that are closely comparable to experiment. Cracking transition states for n-pentane lead to a metastable intermediate (a local minimum with relatively small barriers to escape to deeper minima) where the proton is shared between two hydrocarbon fragments. The zeolite strongly stabilizes these carbocations compared to the gas phase, and the conversion of this intermediate to more stable species determines the product selectivity. Static reaction pathways on the potential energy surface starting from the metastable intermediate include a variety of possible conversions into more stable products. One-picosecond quasiclassical trajectory simulations performed at 773 K indicate that dynamic paths are substantially more diverse than the potential energy paths. Vibrational motion that is dynamically sampled after the cracking transition state causes spilling of the metastable intermediate into a variety of different products. A nearly 10-fold change in the branching ratio between C2/C3 cracking channels is found upon inclusion of post-transition-state dynamics, relative to static electronic structure calculations. Agreement with experiment is improved by the same factor. Because dynamical effects occur soon after passing through the rate-limiting transition state, it is the dynamics, and not only the potential energy barriers, that determine the catalytic selectivity. This study suggests that selectivity in zeolite catalysis is determined by high temperature pathways that differ significantly from 0 K potential surfaces.

Distance-dependent Schwarz-based integral estimates for two-electron integrals: Reliable tightness vs. rigorous upper bounds

A new integral estimate for four-center two-electron integrals is introduced that accounts for distance information between the bra- and ket-charge distributions describing the two electrons. The screening is denoted as QQR and combines the most important features of the conventional Schwarz screening by Häser and Ahlrichs published in 1989 [J. Comput. Chem. 10, 104 (1989)] and our multipole-based integral estimates (MBIE) introduced in 2005 [D. S. Lambrecht and C. Ochsenfeld, J. Chem. Phys. 123, 184101 (2005)]. At the same time the estimates are not only tighter but also much easier to implement, so that we recommend them instead of our MBIE bounds introduced first for accounting for charge-distance information. The inclusion of distance dependence between charge distributions is not only useful at the SCF level but is particularly important for describing electron-correlation effects, e.g., within AO-MP2 theory, where the decay behavior is at least 1/R(4) or even 1/R(6). In our present work, we focus on studying the efficiency of our QQR estimates within SCF theory and demonstrate the performance for a benchmark set of 44 medium to large molecules, where savings of up to a factor of 2 for exchange integrals are observed for larger systems. Based on the results of the benchmark set we show that reliable tightness of integral estimates is more important for the screening performance than rigorous upper bound properties.

Efficient distance-including integral screening in linear-scaling Møller-Plesset perturbation theory.

Efficient estimates for the preselection of two-electron integrals in atomic-orbital based Møller-Plesset perturbation theory (AO-MP2) theory are presented, which allow for evaluating the AO-MP2 energy with computational effort that scales linear with molecular size for systems with a significant HOMO-LUMO gap. The estimates are based on our recently introduced QQR approach [S. A. Maurer, D. S. Lambrecht, D. Flaig, and C. Ochsenfeld, J. Chem. Phys. 136, 144107 (2012)], which exploits the asympotic decay of the integral values with increasing bra-ket separation as deduced from the multipole expansion and combines this decay behavior with the common Schwarz bound to a tight and simple estimate. We demonstrate on a diverse selection of benchmark systems that our AO-MP2 method in combination with the QQR-type estimates produces reliable results for systems with both localized and delocalized electronic structure, while in the latter case the screening essentially reverts to the common Schwarz screening. For systems with localized electronic structure, our AO-MP2 method shows an early onset of linear scaling as demonstrated on DNA systems. The favorable scaling behavior allows to compute systems with more than 1000 atoms and 10,000 basis functions on a single core that are clearly not accessible with conventional MP2 methods. Furthermore, our AO-MP2 method is particularly suited for parallelization and we present benchmark calculations on a protein-DNA repair complex comprising 2025 atoms and 20,371 basis functions.

Charge-transfer and the hydrogen bond: Spectroscopic and structural implications from electronic structure calculations

The absolutely localized molecular orbital (ALMO) model is a fully variational approach which permits polarization of molecules interacting in a cluster while prohibiting charge-transfer (or dative interactions) between individual molecules. The ALMO model can be applied within any density functional theory calculation--the B3LYP functional is employed in this work. ALMO DFT calculations of observables such as optimized geometry, vibrational frequencies and their intensities, and vertical detachment energies are performed for the water dimer, the chloride-water complex and the cyanide-water complex. The vibrational spectra are obtained both within the harmonic approximation and by quasiclassical trajectory simulations. By comparing these ALMO DFT calculations with full DFT calculations using precisely the same functional and basis, the role of charge-transfer on observables in these model hydrogen bonding systems can be assessed. The results can be further interpreted using ALMO-based energy decomposition analysis, which help to reveal the origin of sensitivity or insensitivity of observables to dative interactions. Analysis of the results also suggests that the B3LYP functional, while qualitatively adequate, appears to somewhat overestimate charge-transfer effects.

A Linear-Scaling MP2 Method for Large Molecules by Rigorous Integral-Screening Criteria

Abstract A brief review of our linear-scaling method for atomic-orbital (AO) second-order Møller-Plesset perturbation theory (MP2) is given. The key feature of our method is the rigorous preselection of numerically significant four-center two-electron integrals based on multipole-based integral estimates (MBIE) that do not only account for the exponential coupling between the Gaussian-type basis functions forming charge distributions, but also for the 1/R coupling between the charge distributions. This coupling turns for the required integral products in AO-MP2 into at least a 1/R4 or even a 1/R6 decay behavior. Using MBIE we attain linear scaling, which is illustrated for DNA fragments with up to 1052 atoms and 10 674 basis functions as computed on a single processor. The largest molecule calculated in our present work at the scaled-opposite spin (SOS-) AO-MP2 level is an RNA system comprising 1664 atoms and 19 182 basis functions. Furthermore, we present results for the use of Cholesky-decomposed pseudo-density matrices in Laplace-based MP2, that offers the advantage of exploiting occupied/virtual blocking both with and without auxiliary basis sets.