Explicitly Correlated Orbital Optimized Contracted Pair Correlation

Article pubs.acs.org/JPCA

Explicitly Correlated Orbital Optimized Contracted Pair Correlation Methods: A Short Overview Christian Lasar and Thorsten Klüner* Department of Chemistry, Carl von Ossietzky University Oldenburg, 26111 Oldenburg, Germany ABSTRACT: Designing efficient algorithms to apply highly accurate pair correlation methods to large molecules of scientific interest has remained an important field of research for a long period of time. We present a new approach toward fast algorithms, which represents an interesting alternative and extension to currently existing methods. The presented new contraction scheme saves a significant amount of memory and can be easily combined with efficient linear scaling algorithms. Additionally, the extension to orbital optimization and explicitly correlated f12-theory is demonstrated to further improve accuracy and applicablility.

■

set of variables, namely, the coefficients Aab ij , which need to be optimized and stored. The idea of local correlation methods is to introduce a basis, where many of the coefficients are negligibly small. Unfortunately, this condition is strongly system-dependent. In our new ansatz, we introduce a contraction scheme that drastically reduces the number of coefficients to be stored. The conventional coefficients Aab ij are approximated by a linear combination of integrals Aab ij(B).

INTRODUCTION Pair correlation methods are able to achieve highly accurate energies for chemical reactions.1−3 Unfortunately, their applicability is usually restricted to medium-sized molecules due to storage requirements and computational costs. These restrictions can be overcome by local correlation methods,4 which use physical and mathematical criteria to identitfy longrange interactions, which are negligible and do not need to be computed and saved. The extent to which the interactions and therefore the wave function and the electronic gradient can be compressed is strongly system-dependent. For linear molecules the compression rates will be much higher than for bulky systems such as clusters, since the number of long-range interactions differs significantly.4−9

The large number of slater determinants Φab ij is contracted by a weighted sum over the virtual space a,b, where the contraction coefficients are given by the integral expressions Rabg kl . The remaining coefficients Bklg need to be optimized to minimize ij the energy. They only depend on the occupied space i,j,k,l, which is in general much smaller than the virtual space.

CONTRACTED PAIR CORRELATION In our ansatz, we define a new way toward local correlation. The range of interactions is rigorously defined by the decay of integrals over Gaussian-type geminals and therefore can be accurately predicted by the Schwarz inequality. The compression rates are given by a formal expression that is independent of the basis set size and the spatial extension of the molecule. The basic ansatz for many pair correlation methods, except coupled cluster, can be reduced to configuration interaction singles doubles.1 We will leave out the single excitations at this point to reduce the complexity. They will be approximately included by the orbital optimization described later. With this wave function and a modified energy functional a whole class of pair correlation methods arises.1

∑

|Ψ(B)⟩ = |Φ⟩ + =

∑ R klabg |Φijab⟩ (3)

ab

The tensor elements Rabg are given by integrals over the kl g ̂ operators f . These operators are Gaussian-type two-electron 12

functions with different exponents and angular momenta. (1)

Received: April 27, 2017 Revised: June 1, 2017 Published: June 1, 2017

Here, |Φ⟩ denotes the Hartree−Fock determinant, and |Φab ij ⟩ are doubly excited determinants. This ansatz introduces a large © 2017 American Chemical Society

∑ Bijklg |Ξijklg ⟩ ijklg

|Ξijklg ⟩

Aijab|Φijab⟩

ijab

(2)

klg

■

|Ψ(A)⟩ = |Φ⟩ +

∑ Bijklg R klabg

Aijab(B) ≈

4707

DOI: 10.1021/acs.jpca.7b03960 J. Phys. Chem. A 2017, 121, 4707−4711

Article

The Journal of Physical Chemistry A g

even point of the method. However, the two main requirements leading to linear scaling algorithms are valid in our ansatz. The local molecular virtual orbitals (pair natural orbitals) are replaced by the local atomic orbitals.12 The number of atomic orbitals needed for each pair contribution is automatically restricted due to the schwarz inequality.13 Using these facts, a linear scaling algorithm will be accessible in similar ways as proposed in the literature.14 Details on this rather specific subject are beyond the scope of this paper and will be presented elsewhere.15 The most important question of course remains how accurate the contraction of the wave function is. In Figure 2

R klabg = ⟨ab|f12̂ |kl⟩ g

i

j

k

2 f12̂ = x12g y12g z12g exp( −αg r12 )

(4)

This ansatz represents an extension to the ideas of Höfener et al.10 who applied a similar ansatz to the coupled cluster approximation CC2. In this study, the linear combination of the Gaussian geminals was fixed, and only s-type Gaussian geminals were used. We may see later that the inclusion of p-type geminals is rather essential to reproduce accurate correlation energies in the general case. Additionally, all results were calculated using explicitly correlated f12 theory, which will also be investigated. Note that the idea of approximating expansion coefficients by a more compact expression was also investigated by Kinoshita et al. using singular value decomposition for a coupled cluster wave function.11 As already mentioned, the formal scaling of memory requirements can be drastically reduced using our novel contraction scheme (see Figure 1). While uncontracted pair

Figure 2. Percentage of the correlation energy recovered using two sand p-type Gaussian geminals with respect to the variation of exponents a and b. Both angular momenta share the same exponent in each calculation. The correlation energy is computed with the ACPF2 functional using methane and the cc-pVTZ basis set.16

we varied the exponents of two s- and p-type Gaussian geminals. All these calculations were performed using the averaged coupled pair functional (ACPF) correlation energy2 of the methane molecule in the cc-pVTZ basis set.16 The geometry was optimized with the MP2 method17 and the 631+G** basis set18 using the Gaussian 09 program package.19 Applying the exponents 0.75 and 0.1 we recover 95.6% of the conventional correlation energy of −0.219 hartree (see Figure 2). This percentage can of course be increased by including more geminals in the expansion. Note that using only s-type Gaussian geminals resulted in a maximum percentage of 94.5% recovered, even when applying 24 geminals in the expansion. This indicates that higher angular momenta are necessary to recover a high percentage of the correlation energy. Note that the deviations of the correlation energy range between 91 and 97% using the cc-pVTZ basis and a set of small molecules constructed from the atoms of the first up to the third period. Different bonding situations are contained, and two s- and ptype Gaussian geminals as defined above were applied for this test. Detailed results on this topic will be given in future publications.15 In Table 1 we present the electronic reaction energies for the model reaction C2H2 + H2→C2H4. The error

Figure 1. Formal number N of coefficients that need to be stored for the CnH2n+2 molecule using different basis sets16 compared with the contracted ansatz with two s- and p-type Gaussian geminals. 1

correlation methods formally require 2 no(no + 1)n v2 variables 1

to be stored, we only need to store 2 no(no + 1)no2ng variables (no number of occupied orbitals, nv number of virtual orbitals, ng number of geminals). As soon as the condition n v > no ng is fullfilled, a substantial amount of memory is saved independent from the specific topology of the molecule. Note that for uncontracted local correlation methods applied to linear molecules this formal scaling will be drastically improved, since a lot of pair interactions do not need to be computed and stored.4−9 Nevertheless, the same improvements also apply to our contraction scheme. Especially, for bulky systems we can improve on storage requirements compared to uncontracted local correlation methods. In this case, a lot of pairs interact strongly, and the formal scaling represents the real memory requirements more realistically. Furthermore, the memory requirement of our contraction scheme is independent from the basis set size when integral direct algorithms are used. Additionally, the number of geminals used in the expansion effects the memory consumption only as a linear, systemindependent prefactor. Considering the efficient construction of the electronic gradient, which is the time-determining step, our ansatz can be completely expressed in the very local atomic orbital basis. Currently, there is no linear scaling integral program available to compute p-type Gaussian geminals. This fact prevents us from studying the actual scaling and the break

Table 1. Electronic Reaction Energies for the Model Reaction C2H2 + H2→C2H4 Calculated in the cc-pVTZ Basis with and without Contraction

4708

ΔE, kJ mol−1

uncontracted

contracted

ACPF CEPA0 CCSD(T) exp.20

−215 −213 −208 −203

−220 −218

DOI: 10.1021/acs.jpca.7b03960 J. Phys. Chem. A 2017, 121, 4707−4711

Article

The Journal of Physical Chemistry A

Gaussian geminals with exponents as defined above. All the following calculations were performed using the 6-31+G** basis.18 For the hydrogen fluoride molecule, we expect a significant difference between a CAS(2,2) calculation,27 which only describes the static correlation in the dissociation limit, and the orbital optimized methods, which also include dynamic correlation. We can see in Figure 3 that the uncontracted and

introduced by the contraction scheme is smaller than 3% of the total energy difference in both functionals and therefore in the same range as the energy variation introduced by the choice of the correlation method itself. The deviation of the energy difference calculated with the respective pair correlation kJ method from the experimental reference energy of −203 mol 20 is much larger than the error introduced by contraction and mainly introduced by the insufficient basis set size. A detailed study on a larger testset of molecules reproduces similar results and will be presented elsewhere.15 To further reduce this deviation, we can include orbital optimization and f12 theory. The extension to coupled cluster methods might also be necessary in future work to reproduce experimental results as illustrated in the literature.20 Another important test for any new local correlation method is the decription of weak interactions. We examine this for the dispersion interaction of two neon atoms (see Table 2). Note that using only two s- and Table 2. ACPF Interaction Energies of the Neon Dimer at a Distance of 3.11 Å Using the 6-31+G** Basis ΔE [au]

contracted ACPF

ACPF

2s2p + diff

−0.232

−0.241

Figure 3. Bond dissociation of the hydrogen fluoride molecule.

the contracted orbital optimized configuration interaction doubles reproduce the topology of the CAS calculation while recovering a large fraction of the correlation energy compared to the coupled cluster results close to the minimum of the potential energy curve. The energies of the coupled cluster calculations in the dissociation limit diverge, since they are not variationally bound. In Table 3 we demonstrate that the orbital optimization slightly improves the electronic reaction energies for the model

p-type Gaussian geminals as above leads to convergence problems. Adding a diffuse set of one s- and p-type geminal with the exponent 0.005 solves this problem and recovers the uncontracted interaction energy quite well. Note that a comparison with experimental interaction energies will only be meaningful using large one- and two-particle basis sets and taking care of the basis set superposition error.21 Here, we only would like to show that dispersion interactions can be recovered in general.

Table 3. Orbital Optimized Electronic Reaction Energies for the Model Reaction C2H2 + H2→C2H4 Calculated in the ccpVTZ Basis with and without Contraction

■

ORBITAL OPTIMIZED CONTRACTED PAIR CORRELATION In the following, we present the extension of our contraction scheme to orbital optimized pair correlation methods.22−25 These methods can significantly improve the accuracy of molecular properties in general as shown by Bozkaya et al. for the orbital optimized coupled electron pair approximation OCEPA(0).22 The use of locality in the resulting working equations is again accessible by the atomic orbital basis. Apart from the efficient treatment of large molecules, these methods facilitate a faithful description of the dissociation of single bonds within a molecule. Orbital optimized pair correlation methods extend this ansatz by modifing the energy expression to be size-extensive. Additionally, we can apply our novel contraction scheme. Details of the implementation of the contracted orbital optimized pair correlation methods will be given elsewhere.26 To investigate the ability of orbital optimized pair correlation methods to describe static correlation, we investigate the homolytic bond breaking of the hydrogen fluoride molecule. This process needs two reference determinants for the single bond cleavage, namely, the optimized ground-state determinant and a doubly excited determinant, where the orbital describing the bond σb is replaced by the antibonding one σab * . These determinants obviously occur in our methods, since we include all doubly excited determinants in our expansion. The question is if the contraction scheme is also able to reproduce the correct dissociation behavior. Again, we use two s- and two p-type

ΔE, kJ mol−1

uncontracted

contracted

OACPF OCEPA0 CCSD(T) exp20

−214 −212 −208 −203

−218 −217

reaction C2H2 + H2→C2H4. Note that the frozen core approximation is used during our calculations, which may introduce small errors. Using orbital-optimized coupled cluster may further improve on these energy differences compared to experimental data.23−25 The error introduced by the contraction scheme is again always smaller than 2% of the total energy difference. Further improvements on energy differences can be achieved using explicitly correlated f12 theory to modify the ansatz for the wave function.

■

EXPLICITLY CORRELATED CONTRACTED PAIR CORRELATION Explicitly correlated wave functions are part of extensive research, since they offer a way to drastically improve the otherwise poor basis set convergence of correlation methods.28−31 The contraction approach shown above is the basis for explicitly correlated wave functions. Since we are also using p-type Gaussian geminals, the working equations become slightly more complicated, and details of the implementation will be shown elsewhere.32 While in conventional explicitly correlated theory the linear combination of Gaussian geminals 4709

DOI: 10.1021/acs.jpca.7b03960 J. Phys. Chem. A 2017, 121, 4707−4711

Article

The Journal of Physical Chemistry A is kept fixed, we allow its optimization by the dependence of Bijklg on the geminal index g and therefore increase the variational freedom of the wave function. The ansatz for the wave function is modified in the following way: |Ξijklg ⟩ =

∑ Rklabg|Φijab⟩ + ab

1 (∑ R klPQg |ΦijPQ ⟩ − 2 PQ

properties,22 the application of these methods to rather large molecules is highly interesting. Additionally, we show first results including explicitly correlated f12 theory. Improved convergence toward the experimental data is demonstrated. The extension of our ansatz to coupled cluster methods will be an important future development. A combination with the very accurate pair natural orbital methods may also be a promising future research field.14

∑ Rklpqg|Φijpq⟩) pq

(5)

■

The indices p,q correspond to the current basis set, while the capital P,Q assume a complete one-particle basis and will introduce the explicit correlation. Certain explicitly correlated integral types cannot be solved analytically in an efficient way. In these cases, we assume that the computational basis is complete {p} → {P} and therefore return to the conventional working equations. We currently do not employ an auxilary basis set33 for these integrals, since the virtual space is contracted anyway. This extension might be given in future work. If the used basis set is very small the assumption that it is complete may introduce large errors. Thus, it is recommended to use sufficiently large basis sets of at least tiple zeta quality. We have shown before that it is necessary to include p-type Gaussian geminals into the wave function. These functions do not improve on the description of the Coulomb cusp, and the corresponding explicitly correlated integrals would be expensive and complicated to calculate. We therefore do not include them into the explicitly correlated part of our calculations. In Table 4 we show reaction energies of explicitly correlated calculations. Note that without an auxiliary basis a rather large

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We like to thank Prof. Dr. Klopper, Prof. Dr. Staemmler, and R. Röhse for their time and very helpful discussions. Additionally, we like to thank the DFG Major Research Instrumentation Programme (INST 184/108-1 FUGG) for funding the HPC Cluster HERO located at the University of Oldenburg.

■

uncontracted

contracted

expl. corr. contracted

ACPF CEPA0 exp20

−215 −213 −203

−220 −218

−217 −215

REFERENCES

(1) Pulay, P.; Saebo, S.; Meyer, W. An Efficient Reformulation of the Closed-Shell Self-Consistent Electron Pair Theory. J. Chem. Phys. 1984, 81, 1901−1905. (2) Gdanitz, R. J.; Ahlrichs, R. The Averaged Coupled-Pair Functional (ACPF): A Size-Extensive Modification of MR CI(SD). Chem. Phys. Lett. 1988, 143, 413−420. (3) Meyer, W. PNO-CI and CEPA Studies of Electron Correlation Effects. Theoretica chimica acta 1974, 35, 277−292. (4) Saebø, S.; Pulay, P. Local Treatment of Electron Correlation. Annu. Rev. Phys. Chem. 1993, 44, 213−236. (5) Schütz, M.; Hetzer, G.; Werner, H.-J. Low-Order Scaling Local Electron Correlation Methods. I. Linear Scaling Local MP2. J. Chem. Phys. 1999, 111, 5691−5705. (6) Hetzer, G.; Schütz, M.; Stoll, H.; Werner, H.-J. Low-Order Scaling Local Correlation Methods II: Splitting the Coulomb Operator in Linear Scaling Local Second-Order Møller-Plesset Perturbation Theory. J. Chem. Phys. 2000, 113, 9443−9455. (7) Schütz, M. Low-Order Scaling Local Electron Correlation Methods. III. Linear Scaling Local Perturbative Triples Correction (T). J. Chem. Phys. 2000, 113, 9986−10001. (8) Schütz, M.; Werner, H.-J. Low-Order Scaling Local Electron Correlation Methods. IV. Linear Scaling Local Coupled-Cluster (LCCSD). J. Chem. Phys. 2001, 114, 661−681. (9) Schütz, M. Low-Order Scaling Local Electron Correlation Methods. V. Connected Triples Beyond (T): Linear Scaling Local CCSDT-1b. J. Chem. Phys. 2002, 116, 8772−8785. (10) Höfener, S.; Tew, D. P.; Klopper, W.; Helgaker, T. The Geminal Basis in Explicitly Correlated Wave Functions. Chem. Phys. 2009, 356, 25−30. Moving Frontiers in Quantum Chemistry: Electron Correlation, Molecular Properties and Relativity. (11) Kinoshita, T.; Hino, O.; Bartlett, R. J. Singular Value Decomposition Approach for the Approximate Coupled-Cluster Method. J. Chem. Phys. 2003, 119, 7756−7762. (12) Meyer, W. PNO-CI Studies of Electron Correlation Effects. I. Configuration Expansion by Means of Nonorthogonal Orbitals, and Application to the Ground State and Ionized States of Methane. J. Chem. Phys. 1973, 58, 1017−1035. (13) Saebø, S.; Pulay, P. Local Configuration Interaction: An Efficient Approach for Larger Molecules. Chem. Phys. Lett. 1985, 113, 13−18. (14) Riplinger, C.; Neese, F. An Efficient and Near Linear Scaling Pair Natural Orbital Based Local Coupled Cluster Method. J. Chem. Phys. 2013, 138, 034106.

Table 4. Electronic Reaction Energies for the Model Reaction C2H2+H2 →C2H4 Calculated in the cc-pVTZ Basis with and without Explicit Correlation ΔE, kJ mol−1

AUTHOR INFORMATION

basis set of triple-ζ quality is necessary to make the resolution of the identity approximation valid.33 The energies are slightly kJ closer to the experimental value of −203 mol compared to the contracted ansatz.20 The larger deviation compared to the uncontracted ansatz is introduced by the insufficient accuracy of the contraction scheme for the triple-ζ basis set. Recent results indicate that d- and maybe even f-type Gaussian geminals are needed for this high-accuracy energy regime. Further studies including these Gaussian geminals will have to be investigated in future publications.32 Note that the application of the contraction scheme to coupled cluster wave functions may also improve the results.20

■

SUMMARY AND OUTLOOK In conclusion, we have presented a new way for approaching and extending local correlation methods in general. The presented contraction scheme drastically reduces the storage requirements and recovers a large percentage of the correlation energy of a given pair correlation method. Additionally, we show the application of contraction to orbital optimized pair correlation methods. This includes the use of locality in an easy way for these methods and may lead to the development of efficient computer codes in future work. Since orbital optimized pair correlation methods lead to very accurate molecular 4710

DOI: 10.1021/acs.jpca.7b03960 J. Phys. Chem. A 2017, 121, 4707−4711

Article

The Journal of Physical Chemistry A (15) Lasar, C.; Klüner, T. Contracted Pair Correlation Methods. to be published. (16) Dunning, T. H., Jr. Gaussian Basis Sets for Use in Correlated Molecular Calculations. I. The Atoms Boron through Neon and Hydrogen. J. Chem. Phys. 1989, 90, 1007−1023. (17) Møller, C.; Plesset, M. S. Note on an Approximation Treatment for Many-Electron Systems. Phys. Rev. 1934, 46, 618−622. (18) Hehre, W. J.; Ditchfield, R.; Pople, J. A. Self-Consistent Molecular Orbital Methods. XII. Further Extensions of Gaussian-Type Basis Sets for Use in Molecular Orbital Studies of Organic Molecules. J. Chem. Phys. 1972, 56, 2257−2261. (19) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E. Gaussian 09, Revision B.01, Gaussian Inc.: Wallingford, CT, 2009. (20) Klopper, W.; Bak, K. L.; Jørgensen, P.; Olsen, J.; Helgaker, T. Highly Accurate Calculations of Molecular Electronic Structure. J. Phys. B: At., Mol. Opt. Phys. 1999, 32, R103. (21) Tew, D. P.; Klopper, W. A Comparison of Linear and Nonlinear Correlation Factors for Basis Set Limit Møller-Plesset Second Order Binding Energies and Structures of He2, Be2, and Ne2. J. Chem. Phys. 2006, 125, 094302. (22) Bozkaya, U.; Sherrill, C. D. Orbital-Optimized CoupledElectron Pair Theory and its Analytic Gradients: Accurate Equilibrium Geometries, Harmonic Vibrational Frequencies, and Hydrogen Transfer Reactions. J. Chem. Phys. 2013, 139, 054104. (23) Scuseria, G. E., III; Schaefer, H. F. The Optimization of Molecular Orbitals for Coupled Cluster Wavefunctions. Chem. Phys. Lett. 1987, 142, 354−358. (24) Sherrill, C. D.; Krylov, A. I.; Byrd, E. F. C.; Head-Gordon, M. Energies and Analytic Gradients for a Coupled-Cluster Doubles Model using Variational Brueckner Orbitals: Application to Symmetry Breaking in O4+. J. Chem. Phys. 1998, 109, 4171−4181. (25) Krylov, A. I.; Sherrill, C. D.; Byrd, E. F. C.; Head-Gordon, M. Size-Consistent Wave Functions for Nondynamical Correlation Energy: The Valence Active Space Optimized Orbital Coupled-Cluster Doubles Model. J. Chem. Phys. 1998, 109, 10669−10678. (26) Lasar, C.; Klüner, T. Orbital Optimized Contracted Pair Correlation Methods. to be published. (27) Roos, B. O.; Taylor, P. R.; Siegbahn, P. E. A Complete Active Space {SCF} Method (CASSCF) Using a Density Matrix Formulated Super-CI Approach. Chem. Phys. 1980, 48, 157−173. (28) Kutzelnigg, W. r12-Dependent Terms in the Wavefunction as Closed Sums of Partial Wave Amplitudes for Large l. Theoretica chimica acta 1985, 68, 445−469. (29) Hättig, C.; Klopper, W.; Köhn, A.; Tew, D. P. Explicitly Correlated Electrons in Molecules. Chem. Rev. 2012, 112, 4−74. PMID: 22206503. (30) Kong, L.; Bischoff, F. A.; Valeev, E. F. Explicitly Correlated R12/F12 Methods for Electronic Structure. Chem. Rev. 2012, 112, 75− 107. PMID: 22176553. (31) Ten-no, S.; Noga, J. Explicitly Correlated Electronic Structure Theory from R12/F12 Ansätze. Wiley Interdisciplinary Reviews: Computational Molecular Science 2012, 2, 114−125. (32) Lasar, C.; Klüner, T. Explicitly Correlated Orbital Optimized Contracted Pair Correlation Methods. to be published. (33) Klopper, W.; Samson, C. C. M. Explicitly Correlated SecondOrder Møller-Plesset Methods with Auxiliary Basis Sets. J. Chem. Phys. 2002, 116, 6397−6410.

4711

DOI: 10.1021/acs.jpca.7b03960 J. Phys. Chem. A 2017, 121, 4707−4711