Abstract
An essentially exact groundstate calculation algorithm for fewelectron systems based on superposition of nonorthogonal Slater determinants (SDs) is described, and its convergence properties to ground states are examined. A linear combination of SDs is adopted as manyelectron wave functions, and all oneelectron wave functions are updated by employing linearly independent multiple correction vectors on the basis of the variational principle. The improvement of the convergence performance to the ground state given by the multidirection search is shown through comparisons with the conventional steepest descent method. The accuracy and applicability of the proposed scheme are also demonstrated by calculations of the potential energy curves of fewelectron molecular systems, compared with the conventional quantum chemistry calculation techniques.
Keywords:
Groundstate calculation; Nonorthogonal Slater determinants; Superposition; Fewelectron system; Multiple correction vectorBackground
In recent years, there have been many significant achievements regarding electronic structure calculations in the fields of computational physics and chemistry. However, theoretical and methodological approaches for providing practical descriptions and tractable calculation schemes of the electron–electron correlation energy with continuously controllable accuracy still remain as significant issues [115]. Although density functional theory (DFT) supplies a computationally economical and practical method, there are many unexplored problems raised by unreliable results obtained for some systems in which highly accurate electron–electron correlation energy calculations are required, since results by DFT depend significantly on the exchangecorrelation energy functional used to perform the calculation [1618].
The available quantitatively reliable methods require higher computational costs than the DFT method [18]. Although quantum Monte Carlo methods [1923] can be applied to molecular and crystal systems and show good quantitative reliability where extremely highaccuracy calculations are required, difficulties in calculating forces for optimizing atomic configurations are a considerable disadvantage and inhibit this method from becoming a standard molecular dynamics calculation technique. Configuration interaction (CI), coupled cluster, and MøllerPlesset secondorder perturbation methods, each of which use a linear combination of orthogonalized Slater determinants (SDs) as manyelectron wave functions, are standard computational techniques in quantum chemistry by which highly accurate results are obtained [24], despite suffering from basis set superposition and basis set incompleteness errors. The full CI calculation can perform an exact electron–electron correlation energy calculation in a space given by an arbitrary basis set. However, it is only applicable for small molecules with modest basis sets since the required number of SDs grows explosively on the order of the factorial of the number of basis.
The required number of SDs in order to determine groundstate energies can be drastically decreased by employing nonorthogonal SDs as a basis set. The resonating HartreeFock method proposed by Fukutome utilizes nonorthogonal SDs, and many noteworthy results have been reported [2530]. Also, Imada and coworkers [3133] and Kojo and Hirose [34,35] employed nonorthogonal SDs in path integral renormalization group calculations. Goto and coworkers developed the direct energy minimization method using nonorthogonal SDs [3639] based on the realspace finitedifference formalism [40,41]. In these previous studies, steepest descent directions and acceleration parameters are calculated to update oneelectron wave functions on the basis of a variational principle [2530,3639]. Although the steepest descent direction guarantees a secure approach to the ground state, a more effective updating process might be performed in a multidirection search.
In the present study, a calculation algorithm showing an arbitrary set of linearly independent correction vectors is employed to optimize oneelectron wave functions with Gaussian basis sets. Since the dimension of the search space depends on the number of linearly independent correction vectors, a sufficient number of correction vectors ensure effective optimization, and the iterative updating of all the oneelectron wave functions leads to smooth convergence to the ground states. The primary purpose of this article is to demonstrate the advantage of using multiple correction vectors in searching for the ground state over the conventional steepest descent search in which only one correction vector is used. As a demonstration of the accuracy and applicability of the proposed calculation algorithm, essentially exact potential energy curves of fewelectron molecular systems with long interatomic distances are described for cases where the conventional calculation methods of quantum chemistry fail.
The organization of the article is as follows. In the ‘Optimization algorithm’ section, the proposed calculation algorithm for constructing a basis set of nonorthogonal SDs by updating oneelectron wave functions with multiple correction vectors is described. The expression of the conventional steepest descent direction with a Gaussian basis set is also given for comparison. The convergence characteristics to the ground states of fewelectron systems for calculations using single and multiple correction vectors are illustrated in the ‘Applications to fewelectron molecular systems’ section. As demonstrations of the proposed calculation procedure, the convergence properties to the ground states of fewelectron atomic and molecular systems are also shown. Finally, a summary of the present study is given in the ‘Conclusions’ section.
Optimization algorithm
The calculation procedures for constructing a basis set consisting of nonorthogonal SDs for Nelectron systems using single and multiple correction vectors are described here. An Nelectron wave function ψ(r_{1}, σ_{1}, r_{2}, σ_{2},…, r_{N}, σ_{N}) is expressed by a linear combination of nonorthogonal SDs as follows:
Here, r_{i }and σ_{t }denote the position and spin index of the ith electron, respectively. L is the number of SDs, and Φ^{A}(r_{1}, σ_{1}, r_{2}, σ_{2,}…, r_{N}, σ_{N}) is the Ath SD, given by
with ϕ_{i}^{A}(r) and γ_{i}(σ_{i}) being nonorthogonal and unnormalized oneelectron basis functions and spin orbital functions, respectively. The oneelectron wave function ϕ_{i}^{A}(r) is constructed as a linear combination of Gaussian basis functions x_{s}(r) [24] as
Here, M and D_{i,s}^{A} are the number of basis functions and the sth expansion coefficient for the ith oneelectron wave function ϕ_{i}^{A}(r), respectively.
The steepest direction is implemented in the expression of the total energy functional E of the target system on the basis of the variational principle, without the constraints of orthogonality and normalization on the oneelectron wave functions. The updating procedure of the pth oneelectron wave function belongs to the Ath SD which is represented as
where a_{p}^{A} is the acceleration parameter, which is determined by the variational principle with respect to the total energy E, i.e., [28]
The component of the steepest descent vector K_{p,m}^{A} is given by
where
and
Here, denotes the element of the jth row and ith column of the matrix . The overlap integral matrix S^{AB }is defined by
where the elements S_{ij}^{AB }are the overlap integrals between the oneelectron basis functions, i.e.,
and are the matrix elements of the Hamiltonians,
and
respectively. Here, V(r) stands for an external potential.
The proposed calculation procedure employs linearly independent multiple correction vectors for updating the oneelectron wave function. The pth oneelectron wave function in the Ath SD is updated by
where C_{j}(j = 1, 2,…, L + N_{c}) and N_{c} are the expansion coefficient and the number of correction vectors, respectively. The components of the correction vectors G_{μ,m}^{A} determine N_{c} linearly independent correction functions ξ_{μ}(r) which are defined as linear combinations of Gaussian basis functions as
Since the linearly independent correction vectors can be given arbitrarily, randomly chosen values are employed in the present study. A larger number of correction vectors N_{c} realize a larger volume search space; however, the number of the linearly independent vectors N_{c} is restricted to the dimension of the space defined by the basis set used.
Thus, we have a linear combination of L + N_{c} SDs as the new Nelectron wave function
where
Figure 1 illustrates the flow of the present calculation procedure. Unrestricted HartreeFock (UHF) solutions for a target system are used for initial oneelectron wave functions. The coefficients of Equation 17 are given by solving the generalized eigenvalue equations obtained by employing the variational principle applied to the total energy, and we can have a new Nelectron wave function as a linear combination of L SDs as shown in Equation 17. Iteration of the above updating process for all the oneelectron wave functions of all SDs increasing the number of the SDs’ L leads to an essentially exact numerical solution of the ground state.
Figure 1. Flow of the present algorithm.
Applications to fewelectron molecular systems
Convergence performances for searching for the ground state of a C atom with the 631G** basis set are shown in Figure 2. The UHF solutions are adopted as initial states, and the number of employed SDs is 30. The steepest descent direction and acceleration parameter are adopted for the calculation using one correction vector (N_{c} =1), and seven randomly chosen linearly independent correction vectors are added to the steepest descent correction to create a calculation with eight correction vectors (N_{c} =8). An indispensable advantage of the multidirection search over the single steepest descent direction search is clearly demonstrated. Although the steepest descent vector gives the direction with the largest gradient, it does not necessarily point toward the global energy minimum state. On the contrary, a linear combination of multiple correction vectors can be used to obtain the minimum energy state within the given space by adopting the variation principle.
Figure 2. Effectiveness of multidirection search on total energy convergence. Effectiveness of multidirection search on total energy convergence as a function of the number of iterations for a C atom with the 631G** basis set is shown.
Figure 3 illustrates the convergence performance of the proposed method for the electron–electron correlation energy of a HF molecule with the 631G** basis set as a function of the number of employed SDs. Calculated correlation energies are shown by ratios to exact ones obtained by full CI. The convergence performance to the exact ground state is improved by increasing the number of correction vectors, since the volume of the search space for a oneelectron wave function increase with increasing N_{c}. The essentially exact groundstate energy is obtained using less than 100 nonorthogonal SDs with an error of 0.001%, compared with the exact value in which 99.5% of the electron–electron correlation energy is counted. The obtained convergence is so smooth that the accuracy of the total energy is controllable by adjusting the number of employed SDs. On the other hand, the full CI method requires over 10^{8} orthogonal SDs, and thus the reduction in the numbers of SDs is a significant advantage of adopting nonorthogonal SDs. The groundstate energy obtained by the proposed method does not depend on the components of the correction vectors; however, the rate of convergence does depend on the number of employed correction vectors N_{c}.
Figure 3. Convergence performance of the proposed method for the correlation energy. Convergence performance of the proposed method for the correlation energy of a HF molecule with the 631G** basis set as a function of the number of employed SDs is shown.
The potential energy curve calculated when a single H atom is extracted from a CH_{4} molecule as shown in Figure 4. Calculations are performed using the 631G* basis set. Although the bond lengths are close to the equilibrium one, the errors in the energies obtained by coupledcluster theory with singles and doubles (CCSD) plus perturbative triples (CCSD(T)) are a few milliHartree; at longer bond lengths, the accuracy of the results appears to deteriorate [42]. In contrast, the proposed calculation procedure ensures essentially exact ground states at all bond lengths, since no approximations are employed.
Figure 4. Potential energy curve of a CH_{4 }molecule obtained using the proposed algorithm with 631G* basis set.
Figure 5 illustrates the potential energy curve along the symmetric stretching coordinate of a H_{2}O molecule in the 321G basis set. The angle between the OH bonds is fixed at 107.6°. These results shown for the proposed calculation method, CCSD and CCSD(T) exhibit the same trends as for a CH_{4} molecule. The results for near the equilibrium bond length demonstrate comparable performance between the four methods, whereas results for long bond lengths indicate only that the proposed method has comparable performance with full CI not producing the same unphysical energy curves as CCSD and CCSD(T) around 2.3 Å [42].
Figure 5. Potential energy curve of a H_{2}O molecule obtained using the proposed algorithm with 321G basis set.
Conclusions
A reliable and tractable technique for constructing the groundstate wave function by the superposition of nonorthogonal SDs is described. Linear independent multiple correction vectors are employed in order to update oneelectron wave functions, and a conventional steepest descent method is also performed as a comparison. The dependence of convergence performance on the number of adopted correction vectors is also illustrated. The electron–electron correlation energy converges rapidly and smoothly to the ground state through the multidirection search, and an essentially exact groundstate energy is obtained with drastically fewer SDs (less than 100 SDs in the present study) compared with the number required in the full CI method. For the fewelectron molecular systems considered in the present study, essentially exact electron–electron correlation energies can be calculated even at long bond lengths for which the standard singlereference CCSD and CCSD(T) show poor results, and the practicality and applicability of the proposed calculation procedure have been clearly demonstrated. In future studies, calculations employing periodic boundary conditions and effective core potentials (ECPs) [43] will be performed. A new procedure to reduce the iteration cost should be found in order to increase the applicability of the proposed algorithm for the calculation of essentially exact groundstate energies of manyelectron systems.
Competing interests
The authors declare that they have no competing interests.
Authors’ contributions
HG conceived, planned this study, carried out the coding of the computation program, and drafted the manuscript. MK and KH participated in the discussions on the basic theory of the present method. AS performed tunings of the code and made all of calculations. All authors read and approved the final manuscript.
Acknowledgments
The present study was partially supported by a GrantinAid for the Global COE Program ‘Center of Excellence for Atomically Controlled Fabrication Technology’ (grant no. H08), a GrantinAid for Scientific Research on Innovative Areas ‘Materials Design through Computics: Complex Correlation and NonEquilibrium Dynamics’ (grant no. 22104008), a GrantinAid for Scientific Research in Priority Areas ‘Carbon Nanotube NanoElectronics’ (grant no. 19054009) and a GrantinAid for Scientific Research (B) ‘Design of Nanostructure Electrode by Electron Transport Simulation for Electrochemical Processing’ (grant no. 21360063) from the Ministry of Education, Culture, Sports, Science, and Technology (MEXT) of Japan.
References

Palmer IJ, Brown WB, Hillier IH: Simulation of the charge transfer absorption of the H_{2}O/O_{2 }van der Waals complex using high level ab initio calculations.
J Chem Phys 1996, 104:3198. Publisher Full Text

Kowalski K, Piecuch P: The method of moments of coupledcluster equations and the renormalized CCSD[T], CCSD(T), CCSD(TQ), and CCSDT(Q) approaches.
J Chem Phys 2000, 113:18. Publisher Full Text

Gwaltney SR, Sherrill CD, HeadGordon M: Secondorder perturbation corrections to singles and doubles coupledcluster methods: General theory and application to the valence optimized doubles model.
J Chem Phys 2000, 113:3548. Publisher Full Text

Tsuzuki S, Honda K, Uchimaru T, Mikami M, Tanabe K: Origin of attraction and directionality of the π/π interaction: model chemistry calculations of benzene dimer interaction.
J Am Chem Soc 2002, 124:104. PubMed Abstract  Publisher Full Text

Dutta A, Sherrill CD: Full configuration interaction potential energy curves for breaking bonds to hydrogen: an assessment of singlereference correlation methods.
J Chem Phys 2003, 118:1610. Publisher Full Text

Abrams ML, Sherrill CD: Full configuration interaction potential energy curves for the X 1Σg+, B 1Δg, and B’ 1Σg+ states of C2: a challenge for approximate methods.
J Chem Phys 2004, 121:9211. PubMed Abstract  Publisher Full Text

Juhasz T, Mazziotti DA: Perturbation theory corrections to the twoparticle reduced density matrix variational method.
J Chem Phys 2004, 121:1201. PubMed Abstract  Publisher Full Text

RochaRinza T, Vico LD, Veryazov V, Roos BO: A theoretical study of singlet lowenergy excited states of the benzene dimer.
Chem Phys Lett 2006, 426:268. Publisher Full Text

Du S, Francisco JS: The OH radicalH_{2}O molecular interaction potential.
J Chem Phys 2006, 124:224318. PubMed Abstract  Publisher Full Text

Benedek NA, Snook IK: Quantum Monte Carlo calculations of the dissociation energy of the water dimer.
J Chem Phys 2006, 125:104302. PubMed Abstract  Publisher Full Text

Bonfanti M, Martinazzo R, Tantardini GF, Ponti A: Physisorption and diffusion of hydrogen atoms on graphite from correlated calculations on the Hcoronene model system.
J Phys Chem C 2007, 111:5825. Publisher Full Text

Beaudet TD, Casula M, Kim J, Sorella S, Martin RM: Molecular hydrogen adsorbed on benzene: insights from a quantum Monte Carlo study.
J Chem Phys 2008, 129:164711. PubMed Abstract  Publisher Full Text

Ma J, Michaelides A, Alfe D: Binding of hydrogen on benzene, coronene, and graphene from quantum Monte Carlo calculations.
J Chem Phys 2011, 134:134701. PubMed Abstract  Publisher Full Text

Booth GH, Cleland D, Thom AJW, Alavi A: Breaking the carbon dimer: the challenges of multiple bond dissociation with full configuration interaction quantum Monte Carlo methods.
J Chem Phys 2011, 135:084104. PubMed Abstract  Publisher Full Text

Robinson JB, Knowles P: Approximate variational coupled cluster theory.
J Chem Phys 2011, 135:044113. PubMed Abstract  Publisher Full Text

Feibelman PJ, Hammer B, Norskov JK, Wagner F, Scheffler M, Stumpf R, Watwe R, Dumesic J: The CO/Pt(111) puzzle.
J Phys Chem B 2001, 105:4018. Publisher Full Text

Hu QM, Reuter K, Scheffler M: Towards an exact treatment of exchange and correlation in materials: application to the “CO adsorption puzzle” and other systems.
Phys Rev Lett 2007, 98:176103. PubMed Abstract  Publisher Full Text

Foulkes WMC, Mitas L, Needs RJ, Rajagopal G: Quantum Monte Carlo simulations of solids.
Rev Mod Phys 2001, 73:33. Publisher Full Text

Silverstrelli PL, Baroni S, Car R: Auxiliaryfield quantum Monte Carlo calculations for systems with longrange repulsive interactions.
Phys Rev Lett 1993, 71:1148. PubMed Abstract  Publisher Full Text

Zhang S, Krakauer H, Zhang S: Quantum Monte Carlo method using phasefree random walks with Slater determinants.
Phys Rev Lett 2003, 90:136401. PubMed Abstract  Publisher Full Text

AlSaidi WA, Krakauer H, Zhang S: Auxiliaryfield quantum Monte Carlo study of TiO and MnO molecules.

Suewattana M, Purwanto W, Zhang S, Krakauer H, Walter E: Phaseless auxiliaryfield quantum Monte Carlo calculations with plane waves and pseudopotentials: applications to atoms and molecules.

Purwanto W, Krakauer H, Zhang S: Pressureinduced diamond to βtin transition in bulk silicon: A quantum Monte Carlo study.

Szabo A, Ostlund NS: Modern Quantum Chemistry: Introduction to Advanced Electronic Structure Theory. London: Macmillan; 1982.

Fukutome H: Theory of resonating quantum fluctuations in a fermion system—resonating HartreeFock approximation—.
Prog Theor Phys 1988, 80:417. Publisher Full Text

Ikawa A, Yamamoto S, Fukutome H: Orbital optimization in the resonating HartreeFock approximation and its application to the one dimensional Hubbard model.
J Phys Soc Jpn 1993, 62:1653. Publisher Full Text

Igawa A: A method of calculation of the matrix elements between the spinprojected nonorthogonal Slater determinants.
Int J Quantum Chem 1995, 54:235. Publisher Full Text

Tomita N, Tenno S, Yanimura Y: Ab initio molecular orbital calculations by the resonating HartreeFock approach: superposition of nonorthogonal Slater determinants.
Chem Phys Lett 1996, 263:687. Publisher Full Text

Tenno S: Superposition of nonorthogonal Slater determinants towards electron correlation problems.

Okunishi T, Negishi Y, Muraguchi M, Takeda K: Resonating Hartree–Fock approach for electrons confined in two dimensional square quantum dots.
Jpn J Appl Phys 2009, 48:125002. Publisher Full Text

Imada M, Kashima T: Pathintegral renormalization group method for numerical study of strongly correlated electron systems.
J Phys Soc Jpn 2000, 69:2723. Publisher Full Text

Kashima T, Imada M: Pathintegral renormalization group method for numerical study on ground states of strongly correlated electronic systems.
J Phys Soc Jpn 2001, 70:2287. Publisher Full Text

Noda Y, Imada M: Quantum phase transitions to chargeordered and Wignercrystal states under the interplay of lattice commensurability and longrange Coulomb interactions.
Phys Rev Lett 2002, 89:176803. PubMed Abstract  Publisher Full Text

Kojo M, Hirose K: Pathintegral renormalization group treatments for manyelectron systems with longrange repulsive interactions.
Surf Interface Anal 2008, 40:1071. Publisher Full Text

Kojo M, Hirose K: Firstprinciples pathintegral renormalizationgroup method for Coulombic manybody systems.

Goto H, Hirose K: Totalenergy minimization of fewbody electron systems in the realspace finitedifference scheme.
J Phys: Condens Matter 2009, 21:064231. Publisher Full Text

Goto H, Yamashiki T, Saito S, Hirose K: Direct minimization of energy functional for fewbody electron systems.
J Comput Theor Nanosci 2009, 6:2576. Publisher Full Text

Goto H, Hirose K: Electron–electron correlations in squarewell quantum dots: direct energy minimization approach.
J Nanosci Nanotechnol 2011, 11:2997. PubMed Abstract  Publisher Full Text

Sasaki A, Kojo M, Hirose K, Goto H: Realspace finitedifference approach for multibody systems: pathintegral renormalization group method and direct energy minimization method.
J Phys: Condens Matter 2011, 23:434001. Publisher Full Text

Chelikowsky JR, Troullier N, Saad Y: Finitedifferencepseudopotential method: electronic structure calculations without a basis.
Phys Rev Lett 1994, 72:1240. PubMed Abstract  Publisher Full Text

Hirose K, Ono T, Fujimoto Y, Tsukamoto S: FirstPrinciples Calculations in RealSpace Formalism. London: Imperial College Press; 2005.

Knowles PJ, Cooper B: A linked electron pair functional.
J Chem Phys 2010, 133:224106. PubMed Abstract  Publisher Full Text

Trail JR, Needs RJ: Smooth relativistic Hartree–Fock pseudopotentials for H to Ba and Lu to Hg.
J Chem Phys 2005, 122:174109. PubMed Abstract  Publisher Full Text