Bond breaking with auxiliary-field quantum Monte Carlo
Bond stretching mimics different levels of electron correlation and provides a challenging testbed for approximate many-body computational methods. Using the recently developed phaseless auxiliary-field quantum Monte Carlo (AF QMC) method, we examine bond stretching in the well-studied molecules BH and N, and in the H chain. To control the sign/phase problem, the phaseless AF QMC method constrains the paths in the auxiliary-field path integrals with an approximate phase condition that depends on a trial wave function. With single Slater determinants from unrestricted Hartree-Fock (UHF) as trial wave function, the phaseless AF QMC method generally gives better overall accuracy and a more uniform behavior than the coupled cluster CCSD(T) method in mapping the potential-energy curve. In both BH and N, we also study the use of multiple-determinant trial wave functions from multi-configuration self-consistent-field (MCSCF) calculations. The increase in computational cost versus the gain in statistical and systematic accuracy are examined. With such trial wave functions, excellent results are obtained across the entire region between equilibrium and the dissociation limit.
Quantum Monte Carlo (QMC) methods are an attractive means to treat explicitly the interacting many fermion system. Their computational cost scales favorably with system size, as a low power. The ground state wave function is obtained stochastically by Monte Carlo (MC) sampling, either in particle coordinate space  or in Slater determinant space . Except for a few special cases, however, these methods suffer from the fermion sign problem , which, if uncontrolled, causes an exponential loss of the MC signal and negates the favorable computational scaling. No formal solution has been found for this problem, but approximate methods have been developed that control it. These include the fixed-node method  in real coordinate space and constrained path methods  in Slater determinant space. The real-space fixed-node diffusion Monte Carlo (DMC) method has long been applied to a variety of solids and molecules . Recently, the phaseless auxiliary-field (AF) method was introduced which provides a framework for ab initio electronic structure calculations by QMC in Slater determinant space, within a Hilbert space defined by any single-particle basis .
The phaseless AF QMC method controls the phase problem in an approximate way by using an input trial wave function (WF) . This is a generalization of the constrained path approach  which has been applied to lattice models with Hubbard-like interactions. Compared with previous efforts  on realistic electronic systems using the standard auxiliary-field formalism , the phaseless AF QMC method overcomes the poor (exponential) scaling with system size and projection time and has statistical errors that are well-behaved.
The systematic error from the phaseless approximation has been found to be small near equilibrium geometries in a variety of systems. The method was applied using a planewave basis with pseudopotentials to several -bonded atoms, molecules, and solids  and to the transition metal molecules TiO and MnO . It has also been applied, with Gaussian basis sets, to first- and second-row atoms and molecules , to post-d elements (Ga-Br and In-I) , and to hydrogen bonded systems . The calculated all-electron total energies of first-row atoms and molecules at equilibrium geometries show typical systematic errors of no more than a few milli-hartree (m) compared to exact results. This accuracy is roughly comparable to that of CCSD(T), coupled cluster with single- and double-excited clusters plus a non-iterative correction to the energy due to triple excited clusters. In post-d systems, our results with several basis sets are in good agreement with CCSD(T) results and, for large basis sets, in excellent agreement with experiment . In almost all of these calculations, we have used as trial WF mean-field solutions from independent-electron calculations.
Bond stretching provides a difficult test for approximate correlated methods. In the dissociation limit, the unrestricted Hartree-Fock (UHF) solution gives a qualitatively correct description of the system. The intermediate region between the equilibrium and dissociated geometries represents a situation analogous to a metal-insulator transition. Due to quasi-degeneracies, there can be more than one important electronic configuration, and a single determinant often cannot adequately describe the system. Multi-configurational approaches can describe to a large degree the static correlations in the system, but often miss a large proportion of the dynamic correlations.
No general method has demonstrated the ability to consistently maintain uniformly high accuracy away from equilibrium. Coupled cluster (CC) methods , such as CCSD(T), are remarkably good in describing the equilibrium properties, but are less successful in describing systems with quasi-degeneracies such as the case in the breaking of chemical bonds . Higher order clusters have to be fully included in the iterative approach, because the perturbative corrections are based on non-degenerate perturbation theories, and usually lead to divergences for stretched nuclear geometries. Since CCSD already scales as with basis size, going to triple and higher order clusters is computationally expensive. Multi-reference CC methods could potentially solve some of these problems, but unlike the single-reference CC method, these are still not widely established . Other coupled-cluster-based approaches have been introduced recently to handle bond stretching, and this remains an active field of research; see for example Refs. .
In this paper, we test the phaseless AF QMC method away from Born-Oppenheimer equilibrium configurations. We investigate bond stretching in two well-studied molecules, BH and N, and in a hydrogen chain, H, where exact or very accurate results from full-configuration interaction (FCI) or density-matrix renormalization group (DMRG)  are available. We first use single Slater determinant trial WF’s, obtained by the unrestricted Hartree-Fock (UHF) method. It is shown that AF QMC with UHF as trial WF generally gives better overall accuracy and a more uniform behavior than CCSD(T). The use of multiple determinant trial WFs from multi-configuration self-consistent-field (MCSCF) calculations is then examined in the diatomic molecules. With these trial WFs, excellent AF QMC results are achieved across the entire potential energy surface.
The rest of the paper is organized as follows. The phaseless AF QMC method is first briefly reviewed in the next section. In Section 3, we present and discuss the potential-energy curves of the various systems. Finally, in Section 4, we conclude with a brief summary.
2The phaseless AF QMC Method
The many-body Born-Oppenheimer Hamiltonian in electronic systems can be written in second quantization, in any single-particle basis, as
where is the size of the chosen one-particle basis, and and are the corresponding creation and annihilation operators. The one-electron and two-electron matrix elements depend on the chosen basis.
The phaseless AF QMC obtains the ground state of the system by projecting from a trial WF which has a non zero overlap with the exact ground state of the system:
where is a small time-step, and is assumed to be in the form of a Slater determinant or a linear combination of Slater determinants. Using a second order Trotter decomposition, we can write . The resulting Trotter time-step error decreases with , and can be eliminated by an extrapolation to with multiple calculations.
The central idea in the AF QMC method is the use of the Hubbard-Stratonovich (HS) transformation :
Equation (Equation 3) introduces one-body operators which can be defined generally for any two-body operator by writing the latter in a quadratic form, such as , with a real number. The many-body problem as defined by is now mapped onto a linear combination of non-interacting problems defined by , interacting with external auxiliary fields. Averaging over different auxiliary-field configurations is then performed by MC techniques. Formally, this leads to a representation of as a linear combination of an ensemble of Slater determinants, . The orbitals of each are written in terms of the chosen one-particle basis and stochastically evolve in imaginary time.
Generally, the AF QMC method suffers from the sign or phase problem . The phaseless AF QMC method  used in this paper controls the phase/sign problem in an approximate manner using a trial WF, . The method recasts the imaginary-time path integral as a branching random walk in Slater-determinant space . It uses the overlap , to construct phaseless random walkers, , which are invariant under a phase gauge transformation. The resulting two-dimensional diffusion process in the complex plane of the overlap is then approximated as a diffusion process in one dimension. The ground-state energy computed with the so-called mixed estimate is approximate and not variational in the phaseless method. The error depends on the quality of , and the method becomes exact as the trial WF approaches the exact ground state of the system. This is the only error in the method that cannot be eliminated systematically.
In most applications to date , the trial WF has been a single Slater determinant taken directly from mean-field calculations. We have found  that using the UHF solution leads to better QMC energies than using the restricted Hartree-Fock (RHF) Slater determinant. This was the case even with singlets.
In this study, we will present, in addition to the single-determinant trial WF, results based on multi-determinant trial WFs obtained from MCSCF calculations. In some cases, such as bond stretching, a multi-determinant trial WF can capture some of the static correlation in the system, and thus improve the quality of the constraint in the phaseless approximation. A better trial WF will generally reduce the systematic errors of the phaseless AF QMC method.
In addition, a better trial WF will typically also lead to better statistics in the AF QMC method, for a fixed number of independent MC samples. A simple measure of the efficiency of the multi-determinant MCSCF trial WF relative to a single-determinant UHF trial WF is the following ratio, :
where is the final statistical error, and is the total number of MC samples used in the calculation. (A more precise but closely related measure is the ratio of the variances of the local energy. For its purpose here as a rough indicator, however, the difference between them is not significant.) We expect for a reasonable number of determinants in the MCSCF; in general, the better , the smaller . Since the computational cost of the phaseless AF QMC method increases linearly with the number of determinants in , the overall computational cost of the QMC calculation with respect to the single-determinant trial WF is times the number of determinants in .
3Results and discussion
To examine the performance of the phaseless AF QMC method in bond stretching, the potential-energy curves of the diatomic molecules BH and N are first studied and compared to exact FCI, near-exact DMRG, and several levels of CC methods. In addition, the symmetric and asymmetric bond stretching of an H linear chain is examined and compared to DMRG results. For the diatomic molecules, both single-determinant UHF (QMC/UHF) and multi-determinant MCSCF (QMC/MCSCF) trial WFs were used. AF QMC for the H chain used single-determinant Hartree-Fock trial wave functions.
In our calculations, all the electrons are correlated and the spherical harmonic (as opposed to Cartesian) form of the Gaussian basis functions was used. For the molecules, cc-pVDZ basis sets were used, except in the challenging case of the (triple) bond stretching of N, where calculations were also performed with the cc-pVTZ basis set . For the H chain, the minimal STO-6G basis set was used.
All of the Hartree-Fock, MCSCF, and CC calculations were carried out using NWCHEM  within C symmetry. Some of these calculations were also verified using Gaussian 98  and MOLPRO . The MCSCF energies were obtained from a complete-active-space SCF (CASSCF)  calculations. In most of the molecules, we used the RHF and UHF reference states for the CC calculations [e.g., labeled RCCSD(T) and UCCSD(T), respectively]. FCI calculations were performed using MOLPRO .
Table ? summarizes the cc-pVDZ basis set energies [in hartrees ()] obtained with a variety of methods for seven BH geometries over a range , where . The MCSCF energy was obtained by a CASSCF calculation, performed with 4 active electrons and 8 active orbitals. Figure 1 shows the potential-energy curves from selected methods.
Near the equilibrium geometry, the RCCSD(T) energies are in good agreement with the FCI energy. However, this agreement deteriorates for larger nuclear separation, and RCCSD(T) shows an unphysical dip for which increases for larger bondlength . The failure of RCCSD(T) to describe the molecule for larger bondlengths is attributed to the poor quality of the RHF WF in describing bond breaking. In the large bondlength limit, the UHF solution is better than the RHF solution. This is reflected also in CC results based on the UHF solution; the UCCSD(T) energies are in very good agreement with the FCI energy for large . The UCCSD(T) energies are in less good agreement with FCI in the intermediate region. Overall, UCCSD(T) does quite well in BH, which has a relatively small number of excitations.
As shown in Table ?, QMC/UHF energies are comparable to RCCSD(T) and in good agreement with FCI near the equilibrium geometry. As the bond is stretched, QMC/UHF energies become less accurate. The discrepancy with FCI energies is m for . In the QMC/MCSCF calculations, the multi-determinant trial WF included determinants from MCSCF with coefficient cut-offs . Thus the variational energy of our MCSCF WF is higher than the corresponding MCSCF result listed in the tables. The average value of , as defined in Eq. (Equation 4), is 0.04, and the largest value is 0.08, at the largest bondlength. The QMC/MCSCF energies are in excellent agreement with the FCI energies, to within m for all studied bondlengths.
The optimum cutoff value of determinant coefficient cut-off in the MCSCF trial WF is, of course, system-dependent. The accuracy of the QMC calculation generally improves as the cutoff is lowered, while the computational cost increases. For a small system like BH, a relatively low cutoff leads to excellent trial WFs with large efficiency gain, as the values show. Figure 2 shows the QMC errors as a function of the number of determinants included in the trial WF for three geometries of BH. For , the QMC results with MCSCF trial WFs containing determinants with coefficient cut-offs less than 0.02 (29 determinants) and less than 0.01 (52 determinants) are equivalent within statistical errors. Similarly, for , the QMC results obtained with trial WFs of 24 and 44 determinants are indistinguishable within the statistical errors. Indeed, in both of these cases, 8 determinants in the trial WF give systematic errors less than m. By contrast, for , considerably more determinants are required to achieve converged QMC systematic errors. Note that, because the MCSCF WF is in a spin restricted form, more than one determinant (many more in the case of ) is required to surpass the accuracy of QMC/UHF.
Bond stretching in N is particularly challenging, because it involves the breaking of a triple bond. As a result, N has been extensively studied . Table ? summarizes the calculated total energies, using cc-pVDZ and cc-pVTZ basis sets. Figure 3 plots a selected subset of these potential-energy curves. With the cc-pVDZ basis set, CC results based on the UHF reference state, and the near-exact DMRG energies are from Ref. . We have also verified the UCCSD and UCCSD(T) energies. For both basis sets, the CASSCF calculations are performed with 6 active electrons and 12 active orbitals.
The main features of the CC potential-energy curves of N are similar to those of BH. In contrast with the BH molecule, however, the effects beyond double excitations are substantial in N, even at the equilibrium geometry. CC results based on a RHF reference show an unphysical dip for bohr ( in Figure 3). For the cc-pVDZ basis at the larger and 4.2 bohr bond lengths, UCCSD(T) based on the UHF reference is in better agreement with DMRG than RCCSD(T). Fully including triple excitations with UCCSDT leads to a significant improvement over UCCSD(T) for all geometries except , while RCCSDT seems to be slightly worse than RCCSD(T), except at the last geometry.
QMC with an UHF trial wave function gives a better overall accuracy and a more uniform behavior than CCSD(T) in mapping the potential-energy curve in the cc-pVDZ basis. The largest difference of the QMC/UHF energies compared to DMRG is at the second to last nuclear separation, and is approximately 9 m. With QMC/MCSCF, we included in the multi-determinant trial WF all determinants with a weight larger than 0.01. This gives 65, 66, 76, 97, 82, and 58 determinants for the six bondlengths (in ascending order), respectively. As can be seen from Table ? and the inset of Figure 3, the agreement between the QMC/MCSCF and DMRG values is more uniform and the discrepancy is less than 2-3 m for all geometries.
In the QMC/MCSCF calculations for the cc-pVDZ basis set, the average value of of Eq. (Equation 4) is 0.42, and the largest value is 0.80 at the equilibrium geometry. The weight cutoff choice of 0.01 in selecting the determinants to include from the MCSCF WF was the same as in the BH calculations. This was likely too conservative as in BH. For example, with bohr, the QMC results were within statistical errors for a trial WF that included determinants with coefficient cut-offs .
The cc-pVTZ results from the various methods parallel very well the cc-pVDZ results, as can be seen from Figure 3 and Table ?. Both the QMC/UHF and QMC/MCSCF, for example, mirror each other in the two basis sets. We thus expect the accuracy of the different QMC and CC methods using the cc-pVTZ basis to be comparable to that using the cc-pVDZ basis, where DMRG results are available. For the cc-pVTZ basis set, the QMC/MCSCF calculations included determinants with coefficient cut-offs . The average value of is 0.16 and the largest value is 0.41 for bohr. Additional QMC/MCSCF calculations were performed for , , and bohr, including determinants with coefficient cut-offs , and the same energies were obtained as those in Table ? within statistical errors.
3.3Hydrogen chain: H
The hydrogen linear chain exhibits characteristic signatures of a metal-insulator transition as the interatomic distances are varied. It also provides a simple but challenging model for extended systems, where the favorable scaling of QMC will be especially valuable. Bond stretching in a linear chain of hydrogen atoms, H, was recently benchmarked with DMRG . This 50-electron system was treated using a minimal STO-6G basis set of 50 orbitals. Both symmetric and asymmetric bond stretching were considered. In the case of symmetric bond stretching, the bond between consecutive hydrogen atoms is stretched over the range bohr, and the final structure consists of 50 equidistant, nearly-independent H-atoms. In the case of asymmetric bond stretching, 25 equivalent H molecules are considered, each with a fixed bondlength of 1.4 bohr, where two consecutive hydrogen atoms belonging to two different H molecules are separated over a range of bohr, and the final structure consists of 25 equidistant, nearly-independent H molecules, each at its equilibrium bondlength. Table ? shows the results for both symmetric and asymmetric bond stretching. The RHF and UHF energies, as well as our QMC results obtained using the UHF trial wave function (or RHF when there is no UHF solution) are shown. The RCCSD(T) and DMRG energies as reported in Ref.  are also shown for comparison.
For symmetric stretching, the QMC/UHF total energies are in good agreement with the DMRG results, with the largest discrepancy being about 5 m for and bohr. As the bondlength is stretched, the correlation energy of the system increases. In view of the above results for bond stretching in diatomic molecules, it is not surprising that the discrepancy between RCCSD(T) and DMRG increases as is increased, and for , RCCSD(T) fails to converge as reported in Ref. .
For asymmetric bond stretching, the QMC energies are again in good agreement with the DMRG values. The difference between the QMC and DMRG total energies is less than m for all bond lengths. Here no distinct UHF solution was found, so the RHF Slater determinant was used as the trial WF. The RHF trial WF dissociates properly as is increased in this case, so, not surprisingly, RCCSD(T) is in good agreement with DMRG.
Bond stretching in chemistry is a non-trivial challenge for all approximate correlated methods. In this paper, we applied the recently introduced phaseless auxiliary-field QMC method to study bond-stretching in BH and N and in the H chain. The quality of the phaseless AF QMC method depends on the trial wave function that is used to control the sign/phase problem. With a single UHF Slater determinant as trial wave function, AF QMC has performed very well for molecular geometries near the equilibrium configuration, as shown by comparisons with exact values, CCSD(T) calculations, and experimental results. The results in this paper are consistent with this and extend the reach of phaseless AF QMC method beyond Born Oppenheimer equilibrium structures to bond stretching and bond breaking. For larger nuclear separations, we find that AF QMC with a single-determinant UHF solution in general gives better overall accuracy and a more uniform behavior than coupled cluster CCSD(T). In some stretched bond situations, however, QMC/UHF errors are seen to be significant. In these cases, we find that a trial wave function with a modest number of determinants usually reduces the QMC error to a few m. The QMC computational cost of the multi-determinant trial wave function scales linearly with the number of determinants, but a better trial wave function can reduce both systematic and statistical errors. Using multi-determinant trial wave functions taken directly from MCSCF calculations, the AF QMC results are in very good agreement with exact energies, and uniform behavior is seen across the entire potential-energy curve.
We thank Garnet Chan and Wirawan Purwanto for helpful discussions. This work is supported by ONR (N000140110365 and N000140510055), NSF (DMR-0535529), and ARO (48752PH) grants, and by the DOE computational materials science network (CMSN). Computations were carried out at the Center for Piezoelectrics by Design, the SciClone Cluster at the College of William and Mary, and NCSA at UIUC.
- W. M. C. Foulkes, L. Mira’s, R. J. Needs, and G. Rajagopal, Rev. Mod. Phys. 71, 33 (2001).
- J. W. Moskowitz, K. E. Schmidt, M. A. Lee, and Malvin H. Kalos, J. Chem. Phys. 77, 349 (1982); Peter J. Reynolds, David M. Ceperley, B. J. Alder, and W. A. Lester, J. Chem. Phys. 77, 5593 (1982).
- S. Zhang and H. Krakauer, Phys. Rev. Lett. 90, 136401 (2003).
- S. Zhang, J. Carlson, and J. E. Gubernatis, Phys. Rev. B 55, 7464 (1997).
- D. M. Ceperley and B. J. Alder, J. Chem. Phys. 81, 5833 (1984); J. B. Anderson in Quantum Monte Carlo: Atoms, Molecules, Clusters, Liquids and Solids, Reviews in Computational Chemistry, Vol. 13, ed. by Kenny B. Lipkowitz and Donald B. Boyd (1999); Shiwei Zhang and M. H. Kalos, Phys. Rev. Lett. 71 2159 (1993).
- Shiwei Zhang, in Theoretical Methods for Strongly Correlated Electrons, edited by D. Senechal, A.-M. Tremblay, and C. Bourbonnais (Springer, New York 2003).
- J. B. Anderson, J. Chem. Phys. 63, 1499 (1975); David M. Ceperley and B. J. Alder, Phys. Rev. Lett. 45, 566 (1980).
- N. Rom, D. M. Charutz, and D. Neuhauser, Chem. Phys. Lett. 270, 382 (1997); R. Baer, M. Head-Gordon, and D. Neuhauser, J. Chem. Phys. 109, 6219 (1998).
- P. L. Silvestrelli, S. Baroni, and R. Car, Phys. Rev. Lett. 71, 1148 (1993).
- R. Blankenbecler, D. J. Scalapino, and R. L. Sugar, Phys. Rev. D 24, 2278 (1981).
- G. Sugiyama and S. E. Koonin, Ann. Phys. (NY) 168, 1 (1986).
- S. Zhang, H. Krakauer, W. Al-Saidi, and M. Suewattana, Comp. Phys. Comm. 169, 394 (2005).
- M. Suewattana, W. Purwanto, S. Zhang, H. Krakauer, and E. J. Walter (unpublished).
- W. A. Al-Saidi, H. Krakauer, and S. Zhang, Phys. Rev. B 73, 075103 (2006).
- W. A. Al-Saidi, S. Zhang, and H. Krakauer, J. Chem. Phys. 124, 224101(2006).
- W. A. Al-Saidi, H. Krakauer, and S. Zhang, J. Chem. Phys. 125, 154110 (2006).
- W. A. Al-Saidi, H. Krakauer, and S. Zhang, J. Chem. Phys., in press (2007). Preprint available at http://arxiv.org/abs/physics/0702184.
- J. Cizek, J. Chem. Phys. 45, 4256 (1966).
- T. Crawford and H. Schaefer, Rev. Comp. Chem. 14, 33 (2000).
- P. Piecuch, V. Spirko, A. E. Kondo, and J. Paldus, J. Chem. Phys. 104, 4699 (1996).
- A. Dutta and C. D. Sherrill, J. Chem. Phys. 118, 1610 (2003).
- M. Musial and R. J. Bartlett, J. Chem. Phys. 122, 224 102 (2005).
- P. Piecuch and K. Kowalski, J. Chem. Phys. 113, 5644 (2000).
- R. J. Bartlett, Int. J. Mol. Sci. 3, 579 (2002).
- T. V. Voorhis and M. Head-Gordon, Chem. Phys. Lett. 317, 575 (2000).
- P. Piecuch, K. Kowalski, I. S. O. Pimienta, and M. J. McGuire, Int. Rev. Phys. Chem. 21, 527 (2002).
- A. D. Bochevarov, B. Temelso, and C. David Sherrill, J. Chem. Phys. 125, 4699 (2006).
- S. R. White and R. L. Martin, J. Chem. Phys. 110, 4127 (1999)
- G. K-L. Chan and M. Head-Gordon, J. Chem. Phys. 118, 8551 (2003).
- U. Schollwöck, Rev. Mod. Phys. 77, 259 (2005)
- R. L. Stratonovich, Sov. Phys. Dokl. 2, 416 (1958); J. Hubbard, Phys. Rev. Lett. 3, 77 (1959).
- T. H. Dunning, Jr., J. Chem. Phys. 90, 1007 (1989).
- T. P. Straatsma, E. Aprá, T. L. Windus, E. J. Bylaska, W. de Jong, S. Hirata, M. Valiev, M. T. Hackler, L. Pollack, R. J. Harrison, M. Dupuis, D. M. A. Smith, J. Nieplocha, V. Tipparaju, M. Krishnan, A. A. Auer, E. Brown, G. Cisneros, G. I. Fann, H. Fruchtl, J. Garza, K. Hirao, R. Kendall, J. A. Nichols, K. Tsemekhman, K. Wolinski, J. Anchell, D. Bernholdt, P. Borowski, T. Clark, D. Clerc, H. Dachsel, M. Deegan, K. Dyall, D. Elwood, E. Glendening, M. Gutowski, A. Hess, J. Jaffe, B. Johnson, J. Ju, R. Kobayashi, R. Kutteh, Z. Lin, R. Littlefield, X. Long, B. Meng, T. Nakajima, S. Niu, M. Rosing, G. Sandrone, M. Stave, H. Taylor, G. Thomas, J. van Lenthe, A. Wong, and Z. Zhang, “NWChem, A Computational Chemistry Package for Parallel Computers, Version 4.6” (2004), Pacific Northwest National Laboratory, Richland, Washington 99352-0999, USA.
- Gaussian 98, Revision A.11.4, M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, V. G. Zakrzewski, J. A. Montgomery, Jr., R. E. Stratmann, J. C. Burant, S. Dapprich, J. M. Millam, A. D. Daniels, K. N. Kudin, M. C. Strain, O. Farkas, J. Tomasi, V. Barone, M. Cossi, R. Cammi, B. Mennucci, C. Pomelli, C. Adamo, S. Clifford, J. Ochterski, G. A. Petersson, P. Y. Ayala, Q. Cui, K. Morokuma, N. Rega, P. Salvador, J. J. Dannenberg, D. K. Malick, A. D. Rabuck, K. Raghavachari, J. B. Foresman, J. Cioslowski, J. V. Ortiz, A. G. Baboul, B. B. Stefanov, G. Liu, A. Liashenko, P. Piskorz, I. Komaromi, R. Gomperts, R. L. Martin, D. J. Fox, T. Keith, M. A. Al-Laham, C. Y. Peng, A. Nanayakkara, M. Challacombe, P. M. W. Gill, B. Johnson, W. Chen, M. W. Wong, J. L. Andres, C. Gonzalez, M. Head-Gordon, E. S. Replogle, and J. A. Pople, Gaussian, Inc., Pittsburgh PA, 2002.
- MOLPRO is a package of ab initio programs written by H.-J. Werner and P. J. Knowles, with contributions from J. Almlof, R. D. Amos, M. J. O. Deegan, S. T. Elbert, C. Hampel, W. Meyer, K. A. Peterson, R. M. Pitzer, ¨ A. J. Stone, P. R. Taylor, and R. Lindh, Universitat Bielefeld, Bielefeld, Germany, University of Sussex, Falmer, Brighton, England, 1996.
- B. O. Roos, P. R. Taylor, and P. E. M. Siegbahn, Chem. Phys. 48, 157 (1980).
- P. J. Knowles and N. C. Handy, Chem. Phys. Letters 111, 315 (1984); P. J. Knowles and N. C. Handy, Comp. Phys. Commun. 54, 75 (1989).
- Garnet Kin-Lic Chan, Nihaly Kallay, and Jürgen Gauss, J. Chem. Phys. 121, 6110 (2004).
- W. D. Laidig, P. Saxe, and R. J. Bartlett, J. Chem. Phys. 86, 887 (1987).
- T. V. Voorhis and M. Head-Gordon, J. Chem. Phys. 112, 5633 (2000).
- J. Hachmann, W. Gardoen, and G. K-L. Chan, J. Chem. Phys. 125, 144101 (2006).