Accurate Methods for Large Molecular Systems

Reformulation of All ONIOM-Type Molecular Fragmentation Approaches and Many-Body Theories Using Graph-Theory-Based Projection Operators: Applications to Dynamics, Molecular Potential Surfaces, Machine Learning, and Quantum Computing

We present a graph-theory-based reformulation of all ONIOM-based molecular fragmentation methods. We discuss applications to (a) accurate post-Hartree-Fock AIMD that can be conducted at DFT cost for medium-sized systems, (b) hybrid DFT condensed-phase studies at the cost of pure density functionals, (c) reduced cost on-the-fly large basis gas-phase AIMD and condensed-phase studies, (d) post-Hartree-Fock-level potential surfaces at DFT cost to obtain quantum nuclear effects, and (e) novel transfer machine learning protocols derived from these measures. Additionally, in previous work, the unifying strategy discussed here has been used to construct new quantum computing algorithms. Thus, we conclude that this reformulation is robust and accurate.

Intermolecular interactions in clusters of ethylammonium nitrate and 1-amino-1,2,3-triazole

The intermolecular interaction energies, including hydrogen bonds (H-bonds), of clusters of the ionic liquid ethylammonium nitrate (EAN) and 1-amino-1,2,3-triazole (1-AT) based deep eutectic propellants (DeEP) are examined. 1-AT is introduced as a neutral hydrogen bond donor (HBD) to EAN in order to form a eutectic mixture. The effective fragment potential (EFP) is used to examine the bonding interactions in the DeEP clusters. The resolution of the Identity (RI) approximated second order Møller-Plesset perturbation theory (RI-MP2) and coupled cluster theory (RI-CCSD(T)) are used to validate the EFP results. The EFP method predicts that there are significant polarization and charge transfer effects in the EAN:1-AT complexes, along with Coulombic, dispersion and exchange repulsion interactions. The EFP interaction energies are in good agreement with the RI-MP2 and RI-CCSD(T) results. The quasi-atomic orbital (QUAO) bonding and kinetic bond order (KBO) analyses are additionally used to develop a conceptual and semi-quantitative understanding of the H-bonding interactions as a function of the size of the system. The QUAO and KBO analyses suggest that the H-bonds in the examined clusters follow the characteristic hydrogen bonding three-center four electron interactions. The strongest H-bonding interactions between the (EAN)1:(1-AT)n and (EAN)2:(1-AT)n (n = 1-5) complexes are observed internally within EAN; that is, between the ethylammonium cation [EA]+ and the nitrate anion ([NO3]-). The weakest H-bonding interactions occur between [NO3]- and 1-AT. Consequently, the average strengths of the H-bonds within a given (EAN)x:(1-AT)n complex decrease as more 1-AT molecules are introduced into the EAN monomer and EAN dimer. The QUAO bonding analysis suggests that 1-AT in (EAN)x:(1-AT)n can act as both a HBD and a hydrogen bond acceptor simultaneously. It is observed that two 1-AT molecules can form H-bonds to each other. Although the KBOs that correspond to H-bonding interactions in [EA]+:1-AT, [NO3]-:1-AT and between two 1-AT molecules are weaker than the H-bonds in EAN, those weak H-bond networks with 1-AT could be important to form a stable DeEP.

Graph-|Q⟩⟨C|: A Quantum Algorithm with Reduced Quantum Circuit Depth for Electronic Structure

The accurate determination of chemical properties is known to have a critical impact on multiple fundamental chemical problems but is deeply hindered by the steep algebraic scaling of electron correlation calculations and the exponential scaling of quantum nuclear dynamics. With the advent of new quantum computing hardware and associated developments in creating new paradigms for quantum software, this avenue has been recognized as perhaps one way to address exponentially complex challenges in quantum chemistry and molecular dynamics. In this paper, we discuss a new approach to drastically reduce the quantum circuit depth (by several orders of magnitude) and help improve the accuracy in the quantum computation of electron correlation energies for large molecular systems. The method is derived from a graph-theoretic approach to molecular fragmentation and enables us to create a family of projection operators that decompose quantum circuits into separate unitary processes. Some of these processes can be treated on quantum hardware and others on classical hardware in a completely asynchronous and parallel fashion. Numerical benchmarks are provided through the computation of unitary coupled-cluster singles and doubles (UCCSD) energies for medium-sized protonated and neutral water clusters using the new quantum algorithms presented here.

Accurate Calculation of Isomerization and Conformational Energies of Larger Molecules Using Explicitly Correlated Local Coupled Cluster Methods in Molpro and ORCA

An overview of the approximations in the explicitly correlated local coupled cluster methods PNO-LCCSD(T)-F12 in Molpro and DLPNO-CCSD(T)F12 in ORCA is given. Options to select the domains of projected atomic orbitals (PAOs), pair natural orbitals (PNOs), and triples natural orbitals (TNOs) in both programs are described and compared in detail. The two programs are applied to compute isomerization and conformational energies of the ISOL24 and ACONFL test sets, where the former is part of the GMTKN55 benchmark suite. Thorough studies of basis set effects are presented for selected systems. These revealed large intramolecular basis set superposition effects that make it practically impossible to reliably determine the complete basis set (CBS) limits without including explicitly correlated terms. The latter strongly reduce the basis set dependence and at the same time also errors caused by the local domain approximations. On the basis of these studies, the PNO-LCCSD(T)-F12 method is applied to determine new reference energies for the above-mentioned benchmark sets. We are confident that our results should agree within a few tenths of a kcal mol-1 with the (unknown) CCSD(T)/CBS values, which therefore allowed us to define computational settings for accurate explicitly correlated local coupled cluster methods with moderate computational effort. With these protocols, especially PNO-LCCSD(T)-F12b/AVTZ', reliable reference values for comprehensive benchmark sets can be generated efficiently. This can significantly advance the development and evaluation of the performance of approximate electronic structure methods, especially improved density functional approximations or machine learning approaches.

Treating Semiempirical Hamiltonians as Flexible Machine Learning Models Yields Accurate and Interpretable Results

Quantum chemistry provides chemists with invaluable information, but the high computational cost limits the size and type of systems that can be studied. Machine learning (ML) has emerged as a means to dramatically lower the cost while maintaining high accuracy. However, ML models often sacrifice interpretability by using components such as the artificial neural networks of deep learning that function as black boxes. These components impart the flexibility needed to learn from large volumes of data but make it difficult to gain insight into the physical or chemical basis for the predictions. Here, we demonstrate that semiempirical quantum chemical (SEQC) models can learn from large volumes of data without sacrificing interpretability. The SEQC model is that of density-functional-based tight binding (DFTB) with fixed atomic orbital energies and interactions that are one-dimensional functions of the interatomic distance. This model is trained to ab initio data in a manner that is analogous to that used to train deep learning models. Using benchmarks that reflect the accuracy of the training data, we show that the resulting model maintains a physically reasonable functional form while achieving an accuracy, relative to coupled cluster energies with a complete basis set extrapolation (CCSD(T)*/CBS), that is comparable to that of density functional theory (DFT). This suggests that trained SEQC models can achieve a low computational cost and high accuracy without sacrificing interpretability. Use of a physically motivated model form also substantially reduces the amount of ab initio data needed to train the model compared to that required for deep learning models.

Low-Data Deep Quantum Chemical Learning for Accurate MP2 and Coupled-Cluster Correlations

Accurate ab initio prediction of electronic energies is very expensive for macromolecules by explicitly solving post-Hartree-Fock equations. We here exploit the physically justified local correlation feature in a compact basis of small molecules and construct an expressive low-data deep neural network (dNN) model to obtain machine-learned electron correlation energies on par with MP2 and CCSD levels of theory for more complex molecules and different datasets that are not represented in the training set. We show that our dNN-powered model is data efficient and makes highly transferable predictions across alkanes of various lengths, organic molecules with non-covalent and biomolecular interactions, as well as water clusters of different sizes and morphologies. In particular, by training 800 (H2O)8 clusters with the local correlation descriptors, accurate MP2/cc-pVTZ correlation energies up to (H2O)128 can be predicted with a small random error within chemical accuracy from exact values, while a majority of prediction deviations are attributed to an intrinsically systematic error. Our results reveal that an extremely compact local correlation feature set, which is poor for any direct post-Hartree-Fock calculations, has however a prominent advantage in reserving important electron correlation patterns for making accurate transferable predictions across distinct molecular compositions, bond types, and geometries.

ONIOM meets : efficient, accurate, and robust multi-layer simulations across the periodic table

The computational treatment of large molecular structures is of increasing interest in fields of modern chemistry. Accordingly, efficient quantum chemical approaches are needed to perform sophisticated investigations on such systems. This engaged the development of the well-established "Our own N-layered integrated molecular orbital and molecular mechanics" (ONIOM) multi-layer scheme [L. W. Chung et al., Chem. Rev., 2015, 115, 5678-5796]. In this work, we present the specific implementation of the ONIOM scheme into the xtb semi-empirical extended tight-binding program package and its application to challenging transition-metal complexes. The efficient and broadly applicable GFNn-xTB and -FF methods are applied in the ONIOM framework to elucidate reaction energies, geometry optimizations, and explicit solvation effects for metal-organic systems with up to several hundreds of atoms. It is shown that an ONIOM-based combination of density functional theory, semi-empirical, and force-field methods can be used to drastically reduce the computational costs and thus enable the investigation of huge systems at almost no significant loss in accuracy.

Toward Accurate Prediction of Ion Mobility in Organic Semiconductors by Atomistic Simulation

A multiscale scheme (MLMS: Multi-Level Multi-Scale) to predict the ion mobility (μ) of amorphous organic semiconductors is proposed, which was successfully applied to the hole mobility predictions of 14 organic systems. An inverse relationship between μ and reorganization energy is observed due to local polaronic distortions. Another moderate inverse correlation between μ and distribution of site energy change exists, representing the effects of geometric flexibility. The former and the latter represent the intramolecular and intermolecular geometric effects, respectively. In addition, a linear correlation between transfer coupling and μ is observed, showing the importance of orbital overlaps between monomers. Especially, the highest hole mobility of C6-2TTN is due to its large transfer coupling. On the other hand, another high hole mobility of CBP turned out to come from the high first neighbor density (ρ_FND) of its first self-solvation, emphasizing the proper description of amorphous structural configurations with a sufficiently large number of monomers. In general, systems with either unusually high transfer coupling or high first neighbor density can potentially have high μ regardless of geometric effects. Especially, the newly suggested design parameter, ρ_FND, is pointing to a new direction as opposed to the traditional π-conjugation strategy.

Graph-Theoretic Molecular Fragmentation for Potential Surfaces Leads Naturally to a Tensor Network Form and Allows Accurate and Efficient Quantum Nuclear Dynamics

Molecular fragmentation methods have revolutionized quantum chemistry. Here, we use a graph-theoretically generated molecular fragmentation method, to obtain accurate and efficient representations for multidimensional potential energy surfaces and the quantum time-evolution operator, which plays a critical role in quantum chemical dynamics. In doing so, we find that the graph-theoretic fragmentation approach naturally reduces the potential portion of the time-evolution operator into a tensor network that contains a stream of coupled lower-dimensional propagation steps to potentially achieve quantum dynamics with reduced complexity. Furthermore, the fragmentation approach used here has previously been shown to allow accurate and efficient computation of post-Hartree-Fock electronic potential energy surfaces, which in many cases has been shown to be at density functional theory cost. Thus, by combining the advantages of molecular fragmentation with the tensor network formalism, the approach yields an on-the-fly quantum dynamics scheme where both the electronic potential calculation and nuclear propagation portion are enormously simplified through a single stroke. The method is demonstrated by computing approximations to the propagator and to potential surfaces for a set of coupled nuclear dimensions within a protonated water wire problem exhibiting the Grotthuss mechanism of proton transport. In all cases, our approach has been shown to reduce the complexity of representing the quantum propagator, and by extension action of the propagator on an initial wavepacket, by several orders, with minimal loss in accuracy.

Graph Theoretic Molecular Fragmentation for Multidimensional Potential Energy Surfaces Yield an Adaptive and General Transfer Machine Learning Protocol

Over a series of publications we have introduced a graph-theoretic description for molecular fragmentation. Here, a system is divided into a set of nodes, or vertices, that are then connected through edges, faces, and higher-order simplexes to represent a collection of spatially overlapping and locally interacting subsystems. Each such subsystem is treated at two levels of electronic structure theory, and the result is used to construct many-body expansions that are then embedded within an ONIOM-scheme. These expansions converge rapidly with many-body order (or graphical rank) of subsystems and have been previously used for ab initio molecular dynamics (AIMD) calculations and for computing multidimensional potential energy surfaces. Specifically, in all these cases we have shown that CCSD and MP2 level AIMD trajectories and potential surfaces may be obtained at density functional theory cost. The approach has been demonstrated for gas-phase studies, for condensed phase electronic structure, and also for basis set extrapolation-based AIMD. Recently, this approach has also been used to derive new quantum-computing algorithms that enormously reduce the quantum circuit depth in a circuit-based computation of correlated electronic structure. In this publication, we introduce (a) a family of neural networks that act in parallel to represent, efficiently, the post-Hartree-Fock electronic structure energy contributions for all simplexes (fragments), and (b) a new k-means-based tessellation strategy to glean training data for high-dimensional molecular spaces and minimize the extent of training needed to construct this family of neural networks. The approach is particularly useful when coupled cluster accuracy is desired and when fragment sizes grow in order to capture nonlocal interactions accurately. The unique multidimensional k-means tessellation/clustering algorithm used to determine our training data for all fragments is shown to be extremely efficient and reduces the needed training to only 10% of data for all fragments to obtain accurate neural networks for each fragment. These fully connected dense neural networks are then used to extrapolate the potential energy surface for all molecular fragments, and these are then combined as per our graph-theoretic procedure to transfer the learning process to a full system energy for the entire AIMD trajectory at less than one-tenth the cost as compared to a regular fragmentation-based AIMD calculation.

A double exponential coupled cluster theory in the fragment molecular orbital framework

Fragmentation-based methods enable electronic structure calculations for large chemical systems through partitioning them into smaller fragments. Here, we have developed and benchmarked a dual exponential operator-based coupled cluster theory to account for high-rank electronic correlation of large chemical systems within the fragment molecular orbital (FMO) framework. Upon partitioning the molecular system into several fragments, the zeroth order reference determinants for each fragment and fragment pair are constructed in a self-consistent manner with two-body FMO expansion. The dynamical correlation is induced through a dual exponential ansatz with a set of fragment-specific rank-one and rank-two operators that act on the individual reference determinants. While the single and double excitations for each fragment are included through the conventional rank-one and rank-two cluster operators, the triple excitation space is spanned via the contraction between the cluster operators and a set of rank-two scattering operators over a few optimized fragment-specific occupied and virtual orbitals. Thus, the high-rank dynamical correlation effects within the FMO framework are computed with rank-one and rank-two parametrization of the wave operator, leading to significant reduction in the number of variables and associated computational scaling over the conventional methods. Through a series of pilot numerical applications on various covalent and non-covalently bonded systems, we have shown the quantitative accuracy of the proposed methodology compared to canonical, as well as FMO-based coupled-cluster single double triple. The accuracy of the proposed method is shown to be systematically improvable upon increasing the number of contractible occupied and virtual molecular orbitals employed to simulate triple excitations.

Breaking covalent bonds in the context of the many-body expansion (MBE). I. The purported "first row anomaly" in XHn (X = C, Si, Ge, Sn; n = 1-4)

We present a new, novel implementation of the Many-Body Expansion (MBE) to account for the breaking of covalent bonds, thus extending the range of applications from its previous popular usage in the breaking of hydrogen bonds in clusters to molecules. A central concept of the new implementation is the in situ atomic electronic state of an atom in a molecule that casts the one-body term as the energy required to promote it to that state from its ground state. The rest of the terms correspond to the individual diatomic, triatomic, etc., fragments. Its application to the atomization energies of the XHn series, X = C, Si, Ge, Sn and n = 1-4, suggests that the (negative, stabilizing) 2-B is by far the largest term in the MBE with the higher order terms oscillating between positive and negative values and decreasing dramatically in size with increasing rank of the expansion. The analysis offers an alternative explanation for the purported "first row anomaly" in the incremental Hn-1X-H bond energies seen when these energies are evaluated with respect to the lowest energy among the states of the XHn molecules. Due to the "flipping" of the ground/first excited state between CH2 (3B1 ground state, 1A1 first excited state) and XH2, X = Si, Ge, Sn (1A1 ground state, 3B1 first excited state), the overall picture does not exhibit a "first row anomaly" when the incremental bond energies are evaluated with respect to the molecular states having the same in situ atomic states.

Adsorption Kinetics of Imidacloprid, Acetamiprid and Methomyl Pesticides in Aqueous Solution onto Eucalyptus Woodchip Derived Biochar

This work reports the application of a biochar (BC) derived from eucalyptus wood chips to remove pesticides (imidacloprid, acetamiprid and methomyl) from water. The pseudo-second order kinetic adsorption model is the best fit describing the adsorption of pesticides on BC. Furthermore, the Langmuir model correlated well with the adsorption isotherm data for acetamiprid and methomyl, while the Freundlich model was selected to explain the adsorption of imidacloprid on BC. The maximum adsorption capacities for methomyl, imidacloprid and acetamiprid on the BC material are 32.42, 14.75 and 4.87 mg g−1, respectively. The highest adsorption capacity of methomyl on the BC surface could be the result of multilayer adsorption suggested by the adsorption isotherm studies, with imidacloprid (or acetamiprid) monolayer being adsorbed on the BC surface. The structure, functional groups of pesticides, including their polarity, all played an important role contributing to the performance of biochar sorbent. Preferable interactions between the studied pesticides and the BC surface may include π-π interactions and hydrogen bonding. The steric aromatic entity in adsorbed imidacloprid and acetamiprid on the BC surface may hinder the possibility of other pesticide molecules approaching the available sorption sites on the surface.

Graph-Theory-Based Molecular Fragmentation for Efficient and Accurate Potential Surface Calculations in Multiple Dimensions

We present a multitopology molecular fragmentation approach, based on graph theory, to calculate multidimensional potential energy surfaces in agreement with post-Hartree-Fock levels of theory but at the density functional theory cost. A molecular assembly is coarse-grained into a set of graph-theoretic nodes that are then connected with edges to represent a collection of locally interacting subsystems up to an arbitrary order. Each of the subsystems is treated at two levels of electronic structure theory, the result being used to construct many-body expansions that are embedded within an ONIOM scheme. These expansions converge rapidly with the many-body order (or graphical rank) of subsystems and capture many-body interactions accurately and efficiently. However, multiple graphs, and hence multiple fragmentation topologies, may be defined in molecular configuration space that may arise during conformational sampling or from reactive, bond breaking and bond formation, events. Obtaining the resultant potential surfaces is an exponential scaling proposition, given the number of electronic structure computations needed. We utilize a family of graph-theoretic representations within a variational scheme to obtain multidimensional potential surfaces at a reduced cost. The fast convergence of the graph-theoretic expansion with increasing order of many-body interactions alleviates the exponential scaling cost for computing potential surfaces, with the need to only use molecular fragments that contain a fewer number of quantum nuclear degrees of freedom compared to the full system. This is because the dimensionality of the conformational space sampled by the fragment subsystems is much smaller than the full molecular configurational space. Additionally, we also introduce a multidimensional clustering algorithm, based on physically defined criteria, to reduce the number of energy calculations by orders of magnitude. The molecular systems benchmarked include coupled proton motion in protonated water wires. The potential energy surfaces and multidimensional nuclear eigenstates obtained are shown to be in very good agreement with those from explicit post-Hartree-Fock calculations that become prohibitive as the number of quantum nuclear dimensions grows. The developments here provide a rigorous and efficient alternative to this important chemical physics problem.

On the Accuracy of QM/MM Models: A Systematic Study of Intramolecular Proton Transfer Reactions of Amino Acids in Water

This work presents a systematic assessment of QM/QM' and QM/MM models with respect to direct QM calculations for the tautomerization (neutral to zwitterion) reactions of amino acids (glycine, alanine, valine, aspartate, and neutral and protonated histidine) solvated in a 160 water cluster. The effect of varying QM region size and choice of embedding potentials, including fixed-charge and polarizable molecular mechanics force fields (TIP3P and EFP) and various semiempirical QM methods (PM7, GFN2-xTB, DFTBA, DFTB3, HF-3c, and PBEh-3c), on the accuracy of the models was examined. A surprising finding was that molecular mechanics force fields outperformed many of the semiempirical methods. Generally, the errors in the QM/QM' and QM/MM models converge slowly with respect to the QM region size, requiring 50 or more waters to be included in the QM region before the error in the model falls below 1 kcal mol^-1 of its pure QM result. Different QM region selection schemes were also compared, and it was found that selection based on Natural Population Analysis (NPA) atomic charges significantly reduced the error in the QM/QM' and QM/MM models particularly if a low-quality embedding potential was used. It is envisaged that these results will be useful for the development of future hybrid QM models.

Effective Fragment Potential-Based Molecular Dynamics Studies of Diffusion in Acetone and Hexane

To facilitate more reliable descriptions of transport properties in liquids, molecular dynamics (MD) simulations are performed based on the effective fragment potential (EFP) method derived from first-principles quantum mechanics (in contrast to MD based upon empirically fitted potentials). The EFP method describes molecular interactions in terms of Coulomb, polarization/induction, dispersion, exchange-repulsion, and charge-transfer interactions. The EFP MD simulations described in this paper, performed on hexane and acetone, are able to track the mean-square displacement of molecules for sufficient time to reliably extract translational diffusion coefficients. The results reported here are in reasonable agreement with experiment.

Weighted-Graph-Theoretic Methods for Many-Body Corrections within ONIOM: Smooth AIMD and the Role of High-Order Many-Body Terms

We present a weighted-graph-theoretic approach to adaptively compute contributions from many-body approximations for smooth and accurate post-Hartree-Fock (pHF) ab initio molecular dynamics (AIMD) of highly fluxional chemical systems. This approach is ONIOM-like, where the full system is treated at a computationally feasible quality of treatment (density functional theory (DFT) for the size of systems considered in this publication), which is then improved through a perturbative correction that captures local many-body interactions up to a certain order within a higher level of theory (post-Hartree-Fock in this publication) described through graph-theoretic techniques. Due to the fluxional and dynamical nature of the systems studied here, these graphical representations evolve during dynamics. As a result, energetic "hops" appear as the graphical representation deforms with the evolution of the chemical and physical properties of the system. In this paper, we introduce dynamically weighted, linear combinations of graphs, where the transition between graphical representations is smoothly achieved by considering a range of neighboring graphical representations at a given instant during dynamics. We compare these trajectories with those obtained from a set of trajectories where the range of local many-body interactions considered is increased, sometimes to the maximum available limit, which yields conservative trajectories as the order of interactions is increased. The weighted-graph approach presents improved dynamics trajectories while only using lower-order many-body interaction terms. The methods are compared by computing dynamical properties through time-correlation functions and structural distribution functions. In all cases, the weighted-graph approach provides accurate results at a lower cost.

Scalable Electron Correlation Methods. 8. Explicitly Correlated Open-Shell Coupled-Cluster with Pair Natural Orbitals PNO-RCCSD(T)-F12 and PNO-UCCSD(T)-F12

We present explicitly correlated open-shell pair natural orbital local coupled-cluster methods, PNO-RCCSD(T)-F12 and PNO-UCCSD(T)-F12. The methods are extensions of our previously reported PNO-R/UCCSD methods (J. Chem. Theory Comput., 2020, 16, 3135-3151, https://pubs.acs.org/doi/10.1021/acs.jctc.0c00192) with additions of explicit correlation and perturbative triples corrections. The explicit correlation treatment follows the spin-orbital CCSD-F12b theory using Ansatz 3*A, which is found to yield comparable or better basis set convergence than the more rigorous Ansatz 3C in computed ionization potentials and reaction energies using double- to quaduple-ζ basis sets. The perturbative triples correction is adapted from the spin-orbital (T) theory to use triples natural orbitals (TNOs). To address the coupling due to off-diagonal Fock matrix elements, the local triples amplitudes are iteratively solved using small domains of TNOs, and a semicanonical (T0) domain correction with larger domains is applied to reduce the domain errors. The performance of the methods is demonstrated through benchmark calculations on ionization potentials, radical stabilization energies, reaction energies of fragmentations and rearrangements in radical cations, and spin-state energy differences of iron complexes. For a few test sets where canonical calculations are feasible, PNO-RCCSD(T)-F12 results agree with the canonical ones to within 0.4 kcal mol^-1, and this maximum error is reduced to below 0.2 kcal mol^-1 when large local domains are used. For larger systems, results using different thresholds for the local approximations are compared to demonstrate that 1 kcal mol^-1 level of accuracy can be achieved using our default settings. For a couple of difficult cases, it is demonstrated that the errors from individual approximations are only a fraction of 1 kcal mol^-1, and the overall accuracy of the method does not rely on error compensations. In contrast to canonical calculations, the use of spin-orbitals does not lead to a significant increase of computational time and memory usage in the most expensive steps of PNO-R/UCCSD(T)-F12 calculations. The only exception is the iterative solution of the (T) amplitudes, which can be avoided without significant errors by using a perturbative treatment of the off-diagonal coupling, known as (T1) approximation. For most systems, even the semicanonical approximation (T0) leads only to small errors in relative energies. Our program is well parallelized and capable of computing accurate correlation energies for molecules with 100-200 atoms using augmented triple-ζ basis sets in less than a day of elapsed time on a small computer cluster.

Dividing the Periodic Box into Subdivisions with Their Surroundings to Accelerate Molecular Dynamics Simulation with High Accuracy

A major field of current research in chemistry and biology is the development of the tools that enable in situ analysis of complex systems. However, the long-time dynamics simulation for an extremely large system in solution is almost impossible by an all-atom force field combined with an explicit solvent model. The results show that the larger the periodic box is, the closer the properties of the system are to the experimental values. Therefore, how can we carry out simulations for systems that are fast, accurate, and large enough? A method of dividing the periodic box into subdivisions with their surroundings (DBSS) is presented here, and it clearly increases the computation speed without losing accuracy and enables the simulation of extremely large systems by strongly decreasing the dimension of the charge matrix. The DBSS method divides a single periodic box or unit in an extremely large system into several subdivisions with a suitable choice according to atomic coordinates. This method ensures that these subdivisions are always changing and allows the atoms to communicate with each other. Intermolecular communication is important for molecular properties and functions but is not possible with other fragment methods. The partial charges are calculated in each subdivision with an overlapping surrounding used to take hydrogen bond interaction between the subdivisions into account. This is an iterative process because the charge population will be recalculated at intervals during the dynamics simulations. Taking a water system as an example, each subdivision is extended by 4 Å as the surrounding. The computation time scales almost linearly with the size of the system, and the slope is small. MD simulations for several properties have been performed by the ABEEM-DBSS method. The results indicate that the ABEEM-DBSS method can accurately simulate the properties of water system, and the accuracy can reach or approach that of the experimental data or of other water potentials. Interestingly, the properties become closer to the experimental data as the sizes of the periodic box increase, further validating the need for the simulation of a large system and demonstrating the value of the DBSS method.

Recent developments in the general atomic and molecular electronic structure system

A discussion of many of the recently implemented features of GAMESS (General Atomic and Molecular Electronic Structure System) and LibCChem (the C++ CPU/GPU library associated with GAMESS) is presented. These features include fragmentation methods such as the fragment molecular orbital, effective fragment potential and effective fragment molecular orbital methods, hybrid MPI/OpenMP approaches to Hartree-Fock, and resolution of the identity second order perturbation theory. Many new coupled cluster theory methods have been implemented in GAMESS, as have multiple levels of density functional/tight binding theory. The role of accelerators, especially graphical processing units, is discussed in the context of the new features of LibCChem, as it is the associated problem of power consumption as the power of computers increases dramatically. The process by which a complex program suite such as GAMESS is maintained and developed is considered. Future developments are briefly summarized.

How accurate are approximate quantum chemical methods at modelling solute-solvent interactions in solvated clusters?

In this paper, the performance of a wide range of DFT methods is assessed for the calculation of interaction energies of thermal clusters of a solute in water. Three different charge states (neutral, proton transfer transition state and zwitterion) of glycine were solvated by 1 to 40 water molecules as sampled from molecular dynamics simulations. While some ab initio composite methods that employ insufficiently large basis sets incurred significant errors even for a cluster containing only 5 water molecules relative to the W1X-2 benchmark, the DLPNO-CCSD(T)/CBS and DSD-PBEP86 (triple zeta basis set) levels of theory predicted very accurate interaction energies. These levels of theory were used to benchmark the performance of 16 density functionals from different rungs of Jacob's Ladder. Of the Rung 4 functionals examined, the ωB97M-V and ωB97X-V functionals stood out for predicting absolute interaction energies in 40-water clusters with mean absolute deviations (MAD) ∼4 kJ mol-1. The B3LYP-D3(BJ) functional performed exceptionally well with a MAD ∼1.7 kJ mol-1 and is the overall best performing method. Calculations of relative interaction energies allow for cancellation of systematic errors, including basis set truncation and superposition errors, and the ωB97M-V and B3LYP-D3(BJ) double zeta basis set calculations yielded relative interaction energies that are within ∼3 kJ mol-1 of the benchmark. The ONIOM approximation provides another strategy for accelerating the calculation of accurate absolute interaction energies provided that the calculations have converged with respect to the size of the "high-level-layer".

Scalable Electron Correlation Methods. 6. Local Spin-Restricted Open-Shell Second-Order Møller-Plesset Perturbation Theory Using Pair Natural Orbitals: PNO-RMP2

We present a (near) linear scaling implementation of high-spin open-shell Møller-Plesset perturbation theory using pair natural orbitals (PNO-RMP2). The theory is based on a new variant of open-shell MP2 which is fully spin-adapted and uses a single set of spin-free amplitudes, as in closed-shell MP2. This method, denoted SROMP2, is invariant to unitary orbital transformations within the closed, open, and virtual orbital subspaces. Accordingly, only a single set of PNOs per spatial orbital pair is needed, and the efficiency is similar to closed-shell calculations. The PNOs are obtained using a semicanonical approximation with large domains of projected atomic orbitals (PAOs). Linear scaling is achieved provided that the open-shell orbitals are local, and distant pairs are treated by multipole approximations. The method is efficiently parallelized. The convergence of ionization and reaction energies as a function of the PAO and PNO domain sizes is demonstrated and found to be very similar as for closed-shell calculations. The suitability of the PNOs for explicitly correlated PNO-RCCSD-F12 calculations is also tested. So far, this method is only simulated using a conventional program with appropriate projections to the PAO and PNO subspaces. It is demonstrated for radical stabilization energies as well as ionization potentials that the errors caused by the local domain approximations with our default thresholds are negligible.

Do Better Quality Embedding Potentials Accelerate the Convergence of QM/MM Models? The Case of Solvated Acid Clusters

This study examines whether the use of more accurate embedding potentials improves the convergence of quantum mechanics/molecular mechanics (QM/MM) models with respect to the size of the QM region. In conjunction with density functional theory calculations using the ωB97X-D functional, various embedding potentials including the TIP3P water model, the effective fragment potential (EFP), and semi-empirical methods (PM6, PM7, and DFTB) were used to simulate the deprotonation energies of solvated acid clusters. The calculations were performed on solvated neutral (HA) and cationic (HB⁺) acids clusters containing 160 and 480 water molecules using configurations sampled from molecular dynamics simulations. Consistently, the ωB97X-D/EFP model performed the best when using a minimal QM region size. The performance for the other potentials appears to be highly sensitive to the charge character of the acid/base pair. Neutral acids display the expected trend that semi-empirical methods generally perform better than TIP3P; however, an opposite trend was observed for the cationic acids. Additionally, electronic embedding provided an improvement over mechanical embedding for the cationic systems, but not the neutral acids. For the best performing ωB97X-D/EFP model, a QM region containing about 6% of the total number of solvent molecules is needed to approach within 10 kJ mol^-1 of the pure QM result if the QM region was chosen based on the distance from the reaction centre.

Benchmarking the Effective Fragment Potential Dispersion Correction on the S22 Test Set

The usual modeling of dispersion interactions in density functional theory (DFT) is often limited by the use of empirically fitted parameters. In this study, the accuracies of the popular empirical dispersion corrections and the first-principles derived effective fragment potential (EFP) dispersion correction are compared by computing the DFT-D and HF-D equilibria interaction energies and intermolecular distances of the S22 test set dimers. Functionals based on the local density approximation (LDA) and generalized gradient approximation (GGA), as well as hybrid functionals, are compared for the DFT-D calculations using coupled cluster CCSD(T) at the complete basis set (CBS) limit as the reference method. In general, the HF-D(EFP) method provides accurate equilibrium dimerization energies and intermolecular distances for hydrogen-bonded systems compared to the CCSD(T)/CBS reference data without using any empirical parameters. For dispersion-dominant and mixed systems, the structures and interaction energies obtained with the B3LYP-D(EFP) method are similar to or better than those obtained with the other DFT-D and HF-D methods. Overall, the first-principles derived -D(EFP) correction presents a robust alternative to the empirical -D corrections when used with the B3LYP functional for dispersion-dominant and mixed systems or with Hartree-Fock for hydrogen-bonded systems.

Accurate prediction of the structure and vibrational spectra of ionic liquid clusters with the generalized energy-based fragmentation approach: critical role of ion-pair-based fragmentation

The GEBF method with the ion-pair-based fragmentation has been developed to facilitate ab initio calculations of general ionic liquid clusters.

Quantum mechanical force fields for condensed phase molecular simulations

Molecular simulations are powerful tools for providing atomic-level details into complex chemical and physical processes that occur in the condensed phase. For strongly interacting systems where quantum many-body effects are known to play an important role, density-functional methods are often used to provide the model with the potential energy used to drive dynamics. These methods, however, suffer from two major drawbacks. First, they are often too computationally intensive to practically apply to large systems over long time scales, limiting their scope of application. Second, there remain challenges for these models to obtain the necessary level of accuracy for weak non-bonded interactions to obtain quantitative accuracy for a wide range of condensed phase properties. Quantum mechanical force fields (QMFFs) provide a potential solution to both of these limitations. In this review, we address recent advances in the development of QMFFs for condensed phase simulations. In particular, we examine the development of QMFF models using both approximate and ab initio density-functional models, the treatment of short-ranged non-bonded and long-ranged electrostatic interactions, and stability issues in molecular dynamics calculations. Example calculations are provided for crystalline systems, liquid water, and ionic liquids. We conclude with a perspective for emerging challenges and future research directions.

Quantum mechanical calculation of nanomaterial-ligand interaction energies by molecular fractionation with conjugated caps method

Molecular fractionation with conjugate caps (MFCC) method is introduced for the efficient estimation of quantum mechanical (QM) interaction energies between nanomaterial (carbon nanotube, fullerene, and graphene surface) and ligand (charged and neutral). In the calculations, nanomaterials are partitioned into small fragments and conjugated caps that are properly capped, and the interaction energies can be obtained through the summation of QM calculations of the fragments from which the contribution of the conjugated caps is removed. All the calculations were performed by density functional theory (DFT) and dispersion contributions for the attractive interactions were investigated by dispersion corrected DFT method. The predicted interaction energies by MFCC at each computational level are found to give excellent agreement with full system (FS) calculations with the mean energy deviation just a fractional kcal/mol. The accurate determination of nanomaterial-ligand interaction energies by MFCC suggests that it is an effective method for performing QM calculations on nanomaterial-ligand systems.

Quantum Chemical Methods for the Prediction of Energetic, Physical, and Spectroscopic Properties of Ionic Liquids

The accurate prediction of physicochemical properties of condensed systems is a longstanding goal of theoretical (quantum) chemistry. Ionic liquids comprising entirely of ions provide a unique challenge in this respect due to the diverse chemical nature of available ions and the complex interplay of intermolecular interactions among them, thus resulting in the wide variability of physicochemical properties, such as thermodynamic, transport, and spectroscopic properties. It is well understood that intermolecular forces are directly linked to physicochemical properties of condensed systems, and therefore, an understanding of this relationship would greatly aid in the design and synthesis of functionalized materials with tailored properties for an application at hand. This review aims to give an overview of how electronic structure properties obtained from quantum chemical methods such as interaction/binding energy and its fundamental components, dipole moment, polarizability, and orbital energies, can help shed light on the energetic, physical, and spectroscopic properties of semi-Coulomb systems such as ionic liquids. Particular emphasis is given to the prediction of their thermodynamic, transport, spectroscopic, and solubilizing properties.

Trends in Two- and Three-Body Effects in Multiscale Clusters of Ionic Liquids

Applications of higher correlated levels of ab initio theory to condensed systems require a significant amount of computational resources. The recent development of the fragment molecular orbital (FMO) approach alleviates this issue by splitting the system into individual fragments and achieves the accuracy of the method by accounting for all possible two-body and three-body interactions. In this work a comprehensive application of the FMO approach in combination with a second order of Møller-Plesset perturbation theory method, MP2, is presented for multiscale clusters of ionic liquids such as [C₁mim]X, [C₁mpyr]X, [C₂py]X, and [NMe₄]X, where X = chloride and tetrafluoroborates, BF₄^-, with the clusters varying in size from 4, 8, 16, to 32 ion pairs. Reliable cutoff criteria for the inclusion of two-body and three-body interactions are identified for both HF energy and MP2 correlation energy to achieve the desired accuracy of 1 kJ mol^-1. The importance of two-body and three-body interactions in ionic liquids is also discussed.

Atomic Spectral Methods for Ab Initio Molecular Electronic Energy Surfaces: Transitioning From Small-Molecule to Biomolecular-Suitable Approaches

Continuing attention has addressed incorportation of the electronically dynamical attributes of biomolecules in the largely static first-generation molecular-mechanical force fields commonly employed in molecular-dynamics simulations. We describe here a universal quantum-mechanical approach to calculations of the electronic energy surfaces of both small molecules and large aggregates on a common basis which can include such electronic attributes, and which also seems well-suited to adaptation in ab initio molecular-dynamics applications. In contrast to the more familiar orbital-product-based methodologies employed in traditional small-molecule computational quantum chemistry, the present approach is based on an "ex-post-facto" method in which Hamiltonian matrices are evaluated prior to wave function antisymmetrization, implemented here in the support of a Hilbert space of orthonormal products of many-electron atomic spectral eigenstates familiar from the van der Waals theory of long-range interactions. The general theory in its various forms incorporates the early semiempirical atoms- and diatomics-in-molecules approaches of Moffitt, Ellison, Tully, Kuntz, and others in a comprehensive mathematical setting, and generalizes the developments of Eisenschitz, London, Claverie, and others addressing electron permutation symmetry adaptation issues, completing these early attempts to treat van der Waals and chemical forces on a common basis. Exact expressions are obtained for molecular Hamiltonian matrices and for associated energy eigenvalues as sums of separate atomic and interaction-energy terms, similar in this respect to the forms of classical force fields. The latter representation is seen to also provide a long-missing general definition of the energies of individual atoms and of their interactions within molecules and matter free from subjective additional constraints. A computer code suite is described for calculations of the many-electron atomic eigenspectra and the pairwise-atomic Hamiltonian matrices required for practical applications. These matrices can be retained as functions of scalar atomic-pair separations and employed in assembling aggregate Hamiltonian matrices, with Wigner rotation matrices providing analytical representations of their angular degrees of freedom. In this way, ab initio potential energy surfaces are obtained in the complete absence of repeated evaluations and transformations of the one- and two-electron integrals at different molecular geometries required in most ab inito molecular calculations, with large Hamiltonian matrix assembly simplified and explicit diagonalizations avoided employing partitioning and Brillouin-Wigner or Rayleigh-Schrödinger perturbation theory. Illustrative applications of the important components of the formalism, selected aspects of the scaling of the approach, and aspects of "on-the-fly" interfaces with Monte Carlo and molecular-dynamics methods are described in anticipation of subsequent applications to biomolecules and other large aggregates.