J. Chem. Theory Comput. 2009, 5, 1573–1584
1573
Electrostatically Embedded Many-Body Approximation for Systems of Water, Ammonia, and Sulfuric Acid and the Dependence of Its Performance on Embedding Charges Hannah R. Leverentz and Donald G. Truhlar* Department of Chemistry and Supercomputing Institute, UniVersity of Minnesota, Minneapolis, Minnesota 55455-0431 Received February 24, 2009
Abstract: This work tests the capability of the electrostatically embedded many-body (EE-MB) method to calculate accurate (relative to conventional calculations carried out at the same level of electronic structure theory and with the same basis set) binding energies of mixed clusters (as large as 9-mers) consisting of water, ammonia, sulfuric acid, and ammonium and bisulfate ions. This work also investigates the dependence of the accuracy of the EE-MB approximation on the type and origin of the charges used for electrostatically embedding these clusters. The conclusions reached are that for all of the clusters and sets of embedding charges studied in this work, the electrostatically embedded three-body (EE-3B) approximation is capable of consistently yielding relative errors of less than 1% and an average relative absolute error of only 0.3%, and that the performance of the EE-MB approximation does not depend strongly on the specific set of embedding charges used. The electrostatically embedded pairwise approximation has errors about an order of magnitude larger than EE-3B. This study also explores the question of why the accuracy of the EE-MB approximation shows such little dependence on the types of embedding charges employed.
1. Introduction To compute properties of a chemical system often requires one to find a balance between computational cost and accuracy. A variety of relatively low-cost classical mechanical and semiempirical quantum mechanical methods allow one to calculate the properties of large (hundreds to thousands of atoms) systems quickly (sometimes within a fraction of a second), but, without problematic parametrization against experimental data, these methods are often incapable of providing more than qualitative accuracy for properties derived from a potential energy surface (PES). At the other extreme, calculations based on the first principles of quantum mechanics [such as coupled cluster1 (CC) or configuration interaction2 (CI) theory] have been developed that in principle could be carried to nearly arbitrary levels of quantitative accuracy3 but that in practice may be used to calculate the energies only of systems containing a few atoms because of the methods’ high computational cost. Thus, much
effort has been expended in order to find a broadly applicable method that can accurately calculate the energy of a large system at a cost that would be reasonable for use in either molecular dynamics (MD) or Monte Carlo (MC) simulations. Fragment-based approaches4-13 are one class of methods that attempt to accomplish this goal. These methods involve breaking the large system into subsystems (which will be called fragments) that are small enough to be treated at some desired level of electronic structure theory. Often, an attempt is made to polarize each fragment by representing the “missing” fragments as point charges or continuous charge density distributions, and the large system’s total energy is then calculated as some linear combination of the fragments’ energies and sometimes of the energies of pairs and trimers of the fragments as well. The electrostatically embedded many-body (EE-MB) method,13-17 which will be described in greater detail in Section 2, is a relatively simple fragment-based method that is computationally inexpensive because it does not involve
10.1021/ct900095d CCC: $40.75 2009 American Chemical Society Published on Web 05/01/2009
1574
J. Chem. Theory Comput., Vol. 5, No. 6, 2009
the self-consistent determination of embedding point charges or charge distributions. In the formal EE-MB approximation, each fragment (or monomer), pair of fragments (dimer), and sometimes group of three or more fragments (trimer or higher oligomer), is embedded in a predetermined set of point charges (called embedding charges or background charges) that represents the fragments that are not explicitly included in the electronic structure calculation of a given monomer, dimer, or trimer. When tested on water clusters and on mixed clusters of water and ammonia, the EE-MB approximation showed itself to be a cost-effective way to accurately calculate the total energy of a system of noncovalently interacting molecules at virtually any desired level of electronic structure theory.13-17 The present work continues to explore the EE-MB approximation by looking at two additional aspects of the EE-MB calculations, as described in the next two paragraphs. First, the present study applies the EE-MB approximation to more complicated mixed systems than any on which it has yet been tested; the largest clusters considered in this article are formed from six water molecules, one ammonia molecule, and two sulfuric acid molecules. Clusters of this type were selected because these molecules are thought to be the fundamental components of clusters formed during the early stages of atmospheric nucleation processes.18 In addition, these clusters test the EE-MB approximation’s ability to predict accurate energies (compared to the “full” quantum mechanical calculation by the same electronic structure method) for systems involving both large and small fragments (the large fragment being sulfuric acid with five heavy atoms and the small fragments being water and ammonia with only one heavy atom each) as well as ions or charge transfer complexes because several of the configurations considered in this article correspond to clusters of ammonia, sulfuric acid, bisulfate ion, ammonium ion, and water rather than clusters of only ammonia, sulfuric acid, and water. Second, the present study compares various ways to obtain the embedding charges and tests how sensitively the accuracy of the EE-MB approximation depends on the resulting sets of embedding charges. Typically the sets of background charges that represent the “missing” monomers are determined by performing some kind of population analysis or charge analysis on the electron density matrices of the isolated and optimized gas-phase monomers. Using these predetermined sets of background charges has several advantages relative to using charges that depend on the configuration under consideration: (1) it lowers the cost of the EE-MB calculation by precluding the need to perform additional self-consistent field calculations to determine the “best” background charges for each configuration, and (2) it maintains the straightforward availability of analytic gradients and Hessians (if they are already available for a given method of electronic structure theory) by removing the embedding charges’ dependence on the specific geometry of the system. However, one might argue that using such an inflexible set of embedding charges may not adequately polarize each fragment and could potentially compromise the accuracy of the EE-MB approximation. Therefore, in the
Leverentz and Truhlar
present study we also test some inexpensive ways to obtain embedding charges that do depend on the specific geometry of each system being studied, and we compare the EE-MB results from those geometry-dependent (GD) charges with those from the geometry-independent (GI) charges that would be used in the formal EE-MB approximation. One should note that the formal EE-MB approximation would be more easily applied to dynamical simulations17 that require fast calculations of PES gradients, but that either the formal EEMB approximation or one that uses geometry-dependent background charges would be convenient for Monte Carlo simulations, where the calculation of PES gradients is not required. The outline of the rest of this paper is as follows: Section 2 briefly reviews the theoretical underpinnings of the EEMB approximation, Section 3 describes the computational methods used to perform the tests in this study and also gives the details of how the various sets of background charges were obtained, Section 4 presents the results and discusses their significance, and Section 5 summarizes our conclusions.
2. Theory The EE-MB approximation, like several other fragment-based methods, is based on the many-body expansion of a system’s total energy. Once a system has been fragmented into N monomers, the many-body expansion expresses the system’s total energy as a sum of the energetic contributions of the one-body (i.e., individual monomer) interactions (V1), the two-body interactions (V2), the three-body interactions (V3), and so on up to the N-body term, as shown in eq 1. E ) V1 + V2 + V3 + ... + VN
(1)
If one denotes the energy of one of the monomers as though it had the geometry it has in the cluster but were alone in a vacuum as Ei (where i runs over the arbitrary labels given to the monomers), the energy of dimer as Eij, and the energy of a trimer as Eijk, then the first three terms on the right-hand side of eq 1 are defined in eqs 2 through 4; the definitions of the remaining terms can be inferred from these equations. N
∑ Ei
(2)
∑ (Eij - Ei - Ej)
(3)
V1 )
i)1
N
V2 )
i
