Article pubs.acs.org/JPCA
Cooperative Roles of Charge Transfer and Dispersion Terms in Hydrogen-Bonded Networks of (H2O)n, n = 6, 11, and 16 Suehiro Iwata,*,† Pradipta Bandyopadhyay,‡ and Sotiris S. Xantheas§ †
Department of Chemistry, Faculty of Science and Technology, Keio University, Kohoku, Yokohama 223-8522, Japan, and Toyota Physical and Chemical Research Institute, Nagakute, Aichi 480-1192, Japan ‡ School of Computational and Integrative Sciences, Jawaharlal Nehru University, New Delhi, India § Physical Sciences Division, Pacific Northwestern National Laboratory, 902 Battelle Boulevard, P.O. Box 999, MS K1-83, Richland, Washington 99352, United States S Supporting Information *
ABSTRACT: The perturbation expansion based on the locally-projected molecular orbital (LPMO PT) was applied to the study of the hydrogenbonded networks of water clusters with up to 16 molecules. Utilizing the local nature of the occupied and excited MOs on each monomer, the chargetransfer and dispersion terms are evaluated for every pair of molecules. The two terms are strongly correlated with each other for the hydrogen-bonded pairs. The strength of the hydrogen bonds in the clusters is further classified by the types of the hydrogen donor and acceptor water molecules. The relative energies evaluated with the LPMO PT among the isomers of (H2O)6, (H2O)11, and (H2O)16 agree very well with those obtained from CCSD(T) calculations with large basis sets. The binding energy of the LPMO PT is approximately free of the basis set superposition errors caused both by the orbital basis inconsistency and by the configuration basis inconsistency.
■
INTRODUCTION The properties of hydrogen bonds in water cluster networks have extensively been studied using both wave function based theories and density functional theories (DFT).1−13 In the series of previous studies, one of the authors (S.S.X.) reported the most stable isomers of water clusters (H2O)n (n = 2−21), using high levels of correlated wave function theories with large basis sets.14−18 These earlier results emphasized the necessity of using correlated levels of theory with large basis sets to correctly evaluate both the magnitude of the absolute as well as the relative binding energies among isomers of the same cluster size to account properly for the dispersion energy. One of the reasons for requiring extensive basis sets in the study of these hydrogen-bonded molecular clusters is to minimize the basis set superposition error (BSSE), which is inherent to the expansion methods in approximating the electronic wave function and density of the monomers and of the assembled clusters. The approximations have to be well-balanced both for the monomers and the interacting systems. In the wave function based theories the expansion methods are utilized in two steps; the first step is to expand the one-electron functions, called the molecular orbitals (MO), and the next step is to expand the many-electron function, usually as a sum of the Slater determinants of the electron configurations. The BSSE is caused by the unbalanced approximations in those two expansions and results from the orbital basis inconsistency (OBI) and from the configuration basis inconsistency (CBI) in © 2013 American Chemical Society
the wave functions of the monomers and of the interacting systems.19,20 The best way to avoid the BSSE caused by OBI is to use extensive basis sets (such as the aug-cc-pVQZ) or to rely on an extrapolation scheme leading to the complete basis set (CBS) limit. Because the Hartree−Fock (HF) wave function for a closed shell system can be described by a single determinant (single electron configuration), there is no error due to CBI at this level. As a result, the binding energy evaluated by a HF wave function, in particular after the counterpoise correction (CP) procedure,21,22 smoothly and quickly converges to a limit with the cardinal number x of aug-cc-pVxZ.23 The binding energies with x = T and Q are already close to the CBS limit determined by the extrapolation.24,25 The BSSE for the correlation energy, which is due to CBI, on the other hand, behaves differently in the basis set dependence. The counterpoise corrected correlation energy evaluated by MP2 and CCSD(T) converges slowly to the CBS limit, as has been shown in several previous studies.16,24,25 Halkier et al. recommended the separate extrapolation with different formulas for the HF energy and for the correlation energy.25 The results of these earlier studies suggested that the uncorrected correlation energies plus the CP corrected HF Received: April 17, 2013 Revised: June 25, 2013 Published: June 27, 2013 6641
dx.doi.org/10.1021/jp403837z | J. Phys. Chem. A 2013, 117, 6641−6651
The Journal of Physical Chemistry A
Article
added after the supermolecule calculation in DFT.33,34 One of the examples is the “dispersionless DFT + Dispersion” by Chalasinski and co-workers,35,36 which uses the Pauli blockade (PB) procedure to determine the Kohn−Sham orbitals. As was clearly stated by Gutowski and Pile in their original paper of the PB procedure,37 BSSE has to be removed by the counterpoise procedure. This is because the converged energy should recover the Hartree−Fock energy of the supermolecule wave function. The PB procedure was developed only to interpret the Hartree−Fock interaction energy between closed-shell systems. Another interesting development on this subject is the threeparameter doubly hybrid density functionals (DHDF, named XYG3) by Xin and co-workers,38−40 which has been successfully applied for dispersion-dominant dimers. Because they are based on the sum formula of the MP2 at the second hybrid step, they are not immune of the CBI as well as of the OBI. Recent progress in this field is reviewed by Sherrill and coworkers33 and by Klimes and Michaelides.34 In this respect, our approach can be considered as the wave function version of the “+D” method without any parameters. The purpose of the present study is as follows: 1. To examine the applicability and the restriction of the LPMO 3SPT + Dispersion for the description of water clusters. 2. To demonstrate the correlation of the charge-transfer and dispersion terms for the hydrogen-bonded pairs of water dimers. 3. To categorize the strength of the hydrogen bonds by the participating water molecules that are classified either as donors or acceptors of hydrogen bonds. 4. To determine the factors that are responsible for the energetic stabilization of the most stable isomers of water clusters. The present analysis thus aims at the understanding of the stability of water clusters and their hydrogen-bonded networks.
binding energies are closer to the CBS limit binding energy than the fully CP corrected ones. The convergence behavior of the correlation energy with the basis set may be systemdependent, and so care should be exercised in the application of the CP procedure to the correlation energy. One of the present authors (S.I.) demonstrated that, after the removal of the BSSE in the HF binding energy by the CP procedure or by the method described below, the binding energy evaluated by the CP uncorrected correlation methods converges quickly and is close to the CBS limit.26 Because the dispersion interaction plays a key role in noncovalent interactions, such as the hydrogen bond, electron correlation is required to properly describe those systems. However, in the supermolecule calculations the correlated wave functions are expressed as linear combinations of electronic configurations, and the separation of the intra- and intermolecular contributions to the correlation energy is difficult. Therefore, special care has to be taken to avoid the CBI. One approach to differentiate between the intra- and intermolecular correlation is the use of localized orbitals. The local MP2 formulation of Pulay is a viable approach,27 but the occupied MOs are canonical MOs though they are localized around atoms in molecules, and the OBI at the Hartree−Fock level cannot be removed. Werner and co-workers have developed efficient implementations of that method and applied it for various systems.28 An alternative approach is to variationally determine the local MOs for every fragment of the cluster. The local MOs are expanded in terms of the basis sets centered on the atomic centers in each molecule. It is called the locally projected MO (LPMO), because of the form of the effective Hartree−Fock equations.29 The energy calculated with the Slater determinant for the whole cluster, constructed from the nonorthogonal LPMOs, is free of OBI, but it extremely underestimates the binding energy because of the strict restriction of electron delocalization over the other molecules.30 After defining the excited MOs that are local to each molecule, the Slater determinants of the excited electronic configurations have to be added to the original Slater determinant, and the MP2 type perturbation theory (LPMO PT) has been developed. Several test calculations demonstrated that the third-order correction within the single excitations (LPMO 3SPT) approximates the CP corrected HF binding energy when the augmented basis functions are used.20,31 If small basis sets such as the cc-pvdz are used, the LPMO 3SPT approach partly introduces BSSE, and a simple, approximate procedure for removing the error was previously introduced.31 To avoid the CBI for the electron correlation case, only the pair type excitations for the intramolecular correlation and the dispersion type excitations for the intermolecular correlation are chosen.26 The further approximation to only include the dispersion types (LPMO 3SPT + Dispersion) was previously successfully tested for the weak electron−donor−acceptor complexes and small water clusters.32 The present work represents an extension of that previous work. The binding energy evaluated by DFT also contains the BSSE due to OBI, unless the complete basis set is used to determine the electron density. As for the configuration basis inconsistency (CBI), because the electron correlation is evaluated through the functional, the parameters in the used functionals could correct for it, if the supermolecule approach is adopted. But the analysis may not be simple. There are numerous recent papers, in which the dispersion correction is
■
THEORETICAL AND COMPUTATIONAL PROCEDURE Energy Analysis. Over the last 10 years, one of the authors (S.I.) has been developing the perturbation expansion theory based on the locally projected molecular orbital (LPMO) approach. The methods and procedures behind the main ideas of the present study are already documented in previous papers.20,26,31,32,41,42 Here we describe only a few salient features that pertain to the current discussion. The zero-order wave function ΨLPMO for a molecular cluster is a single Slater determinant constructed from sets of local MOs variationally determined by the coupled Hartree−Fock equations with the strong constraint that the occupied MOs for each molecule are expanded only in terms of the basis sets on that molecule. Because the MOs of the different molecules are not orthogonal, the Slater determinant has to be properly normalized. The energy evaluated with this Slater determinant includes the exchange−repulsion, electrostatic interaction as well as the induction (polarization) terms caused by the other molecules, but it can be proved that the charge-transfer terms between the molecules are excluded. One of the important aspects is that the energy is free of the basis set superposition error (BSSE) caused by the orbital basis inconsistency (OBI). LPMOs are canonical for each molecule but not for the whole cluster, and therefore, the first-order wave function starts with the single excitations as 6642
dx.doi.org/10.1021/jp403837z | J. Phys. Chem. A 2013, 117, 6641−6651
The Journal of Physical Chemistry A
Article
Table 1. Relative Energies and Energy Components (kJ mol−1) of Various Water Cluster Isomers isomers
geom.
MP2
CCSD(T)
EDisp
ECT+LE
ELPMO BindE
0.0 0.9 1.9 3.5
−73.4 −73.7 −69.2 −64.0
−56.7 −59.0 −62.7 −63.5
−59.1 −55.6 −55.4 −58.1
0.0 −2.1 −1.4 5.9 −2.9 6.7 10.6 −0.7 9.2 10.7
−166.7 −167.9 −164.0 −165.6 −164.1 −164.1 −163.1 −164.7 −163.8 −168.0
−147.6 −140.7 −149.4 −139.5 −147.9 −146.2 −142.5 −142.8 −137.2 −135.3
−103.8 −111.6 −106.0 −107.1 −109.0 −101.1 −101.9 −111.3 −106.4 −105.6
0.0 2.3 2.6 1.5 2.8
−264.6 −258.4 −257.9 −256.8 −268.7
−207.1 −223.1 −222.5 −219.5 −211.6
−192.4 −180.3 −181.0 −186.4 −181.0
3SPT+Disp (H2O)6
prism cage book cyclic
a a a a
0.0 1.1 3.1 8.6
(0.0a) (0.3) (1.0) (4.2)
0.0 1.9 4.8 11.3
X434 X443 X515 P3189 X551 P1226 P2476 X4412 P4205 P2877
b b b c b c c b c c
0.0 1.6 1.3 2.6 2.6 3.6 6.0 6.1 7.5 7.2
[0.0d] [1.8] [−0.5] [2.9] [0.2] [2.6] [5.5] [4.3] [7.6] [8.9]
0.0 1.3 2.6 2.7 3.9 4.8 6.7 6.9 7.1 7.3
4444-a boat-a boat-b antiboat 4444-b
f f f f f
0.0f −1.6 −0.9 −0.2 1.9
(H2O)11 [0.0e] [1.5] [0.8] [3.1] [1.5] [3.8] [6.2] [5.1] [7.6] [8.5] (H2O)16
a
0.0f 1.1 1.8 2.1 2.3
Reference 15. The MP2/CBS limit is given in parentheses.15 bReference 17. MP2: aug-cc-pVDZ, CCS(T): aug-cc-pVDZ. cPresent work. See text. LPMO 3SPT/aug-cc-pVTZ + MP2/aug-cc-pVDZ. eLPMO 3SPT/aug-cc-pVTZ + CCSD(T)/aug-cc-pVDZ. fRefererence 18. aug-cc-pVTZ.
d
X≠Y
Ψ1ST =
∑
|LE X ⟩ +
Mol = X
∑
the change of the intramolecular electron correlation is canceled out when evaluating the binding energy. The calculated binding energy in this approximation can be written as
|CTX → Y ⟩
Mol = X,Y
X
