PHYSICAL REVIEW A 72, 062502 共2005兲

Free iterative-complement-interaction calculations of the hydrogen molecule1

Yusaku Kurokawa,1 Hiroyuki Nakashima,1 and Hiroshi Nakatsuji1,2

Department of Synthetic Chemistry and Biological Chemistry, Graduate School of Engineering, Kyoto University, Katsura,
Nishikyo-ku, Kyoto 615-8510, Japan
2
Fukui Institute for Fundamental Chemistry, Kyoto University, 34-4 Takano-Nishihiraki-cho, Sakyo-ku, Kyoto 606-8103, Japan
共Received 26 June 2005; published 13 December 2005兲
The free iterative-complement-interaction 共ICI兲 method based on the scaled Schrödinger equation proposed
previously has been applied to the calculations of very accurate wave functions of the hydrogen molecule in ananalytical expansion form. All the variables were determined with the variational principle by calculating thenecessary integrals analytically. The initial wave function and the scaling function were changes to see theeffects on the convergence speed of the ICI calculations. The free ICI wave functions that were generated
automatically were different from the existing wave functions, and this difference was shown to be physicallyimportant. The best wave function reported in this paper seems to be the best worldwide in the literature fromthe variational point of view. The quality of the wave function was examined by calculating the nuclear andelectron cusps.
DOI: 10.1103/PhysRevA.72.062502

PACS number共s兲: 31.10.⫹z, 03.65.Ca, 03.65.Ge, 31.15.⫺p

I. INTRODUCTION

When the quantum chemical methodology came out, thefirst molecular application was done to hydrogen molecule.
In 1927, Heitler and London 关1兴 applied newly born quantum
mechanics to this molecule, which was a birth of quantum
chemistry. After a pioneering study of helium atom by Hylleraas 关2兴, James and Coolidge 关3兴 applied their explicitly
correlated functions to hydrogen molecule and obtained very
accurate results. Since then, Kolos, Roothaan, Wolniewitz,
Komasa, Cencek, Rychlewski, Kutzelnigg, and others 关4–12兴
have reported very accurate wave functions of this molecule.
Recently, Rychlewski and Komasa gave an extensive review
关13兴 including the explicitly correlated functions of hydrogenmolecule. We have recently developed a methodology forexactly solving the Schrödinger equation in an analytical expansion form 关14–16兴. So, in this paper, we apply our methodology to the hydrogen molecule in order to examine theefficiency of the proposed method and the quality of thecalculated energy and wave function.
The methodology is based on the study on the structure of
the exact wave function 关14–16兴. When the structure of theexact wave function is clarified, we construct such functions
and make the variables included to be optimal by using thevariational principle. However, when we formulate such atheory based on the regular Schrödinger equation for atoms
and molecules, we encounter a difficulty called the singularity problem 关15,16兴. The formulation includes the integralsof higher powers of the Hamiltonian, but such integrals diverge owing to the singularity of the Coulomb potential included in the Hamiltonian. We showed that this difficulty
could be avoided by introducing the scaled Schrödinger
equation 关15,16兴: when we reformulate our theory based on
the scaled Schrödinger equation, our theory for constructing
the exact wave function becomes free from the singularity
problem. Further, the free iterative-complement-interaction
共ICI兲 method provides a very flexible way of constructing theexact wave function in both of the choices of the initial function ␺0 and of the scaling g function.
1050-2947/2005/72共6兲/062502共11兲/$23.00

We explain in the next section our methodology applied tothe calculations of the exact wave function of the hydrogenmolecule. Detailed computational aspects are then described
and the results are summarized. We compare the present free
ICI wave functions with the existing accurate wave functionsin the field of explicitly correlated wave functions and examine the importance of the terms newly generated by thepresent free ICI method. The quality of the present wave
function is examined by calculating the nuclear and electron
cusp values 关17兴. The conclusion of this study is given in thelast section. The Appendix summarizes briefly the mathematics necessary for the present study.
II. FREE ICI METHOD APPLIED TO THE HYDROGEN
MOLECULE

We want to solve the Schrödinger equation of the hydrogen molecule,
共H − E兲␺ = 0,

共1兲with the Hamiltonian given by
H = − 21 ⵜ21 − 21 ⵜ22 − 1/r1a − 1/r1b − 1/r2a − 1/r2b + 1/r12 + 1/R,
共2兲
where 1 and 2 denote electrons, a and b two protons, and R
the internuclear distance. When we introduce the elliptic coordinates

␮i = 共ria − rib兲/R,␭i = 共ria + rib兲/R,␳ = 2r12/R,

共3兲

with i being 1 or 2, the kinetic operator and the potential
operator are written as
ⵜ2i =

冋

4
⳵
⳵
⳵
⳵
共␭2i − 1兲+
共1 − ␮2i 兲⳵ ␭i ⳵ ␮i⳵ ␮i
R2共␭2i − ␮2i 兲 ⳵ ␭i
+

再

1

1

+
␭2i − 1 1 − ␮2i

冎 册
⳵2
⳵ ␾2i

共4兲

and

062502-1©2005 The American Physical Society

PHYSICAL REVIEW A 72, 062502 共2005兲

KUROKAWA, NAKASHIMA, AND NAKATSUJI

Vne + Vee =

冉

冊

2
1
␭1
␭2
− 2
,
2 − 2
2 +
R
␭1 − ␮1 ␭2 − ␮2 ␳

共5兲

respectively, where Vne and Vee represent the nuclear attraction and electron repulsion operators, respectively.
One can see that the Hamiltonian includes the Coulombic
potential that becomes plus or minus infinity when two particles meet together. Higher powers of such a potential become strongly diverging like a well potential so that the integrals of ␺2 over such higher powers of the potential
diverge. Therefore, the ICI calculations based on the ordinary Hamiltonian become impossible since it involves suchdiverging integrals. However, this singularity problem can be
avoided by introducing the scaled Schrödinger equationg共H − E兲␺ = 0,

共6兲

where g is a scaling function that is positive and nonzero
everywhere except at the singular points r0 and even there it
satisfies
lim gV = a,
r→r0

共7兲

with a being a nonzero constant. The condition given by Eq.
共7兲 is necessary for not to eliminate the information of the
Hamiltonian at the singular points: the singularity is also an
important physics of the system.
Based on the scaled Schrödinger equation, we can formulate the ICI method that is free from the singularity problem.
The simplest ICI 共SICI兲 wave function is defined by the recursion formula as␺n+1 = 关1 + Cng共H − En兲兴␺n ,共8兲

which is guaranteed to become exact at convergence and we
do not have the singularity difficulty in the course of theiteration process starting from the initial function ␺0 becauseof the existence of the scaling function g. In the present
calculation of the hydrogen molecule, we use actually thefree ICI method that is formulated from the SICI wave function as follows. We examine the right-hand side of Eq. 共8兲,collect all the independent functions, and arrange them as
兵␾k其共n兲 , k = 1 , … , M n. Here M n is the number of independentfunctions included in 兵␾k其共n兲.With this set of functions
兵␾k其共n兲, we expand our ␺n+1 asMn␺n+1 = 兺 ck,n␾k .

共9兲

k=1

The coefficients 兵ck,n其 , k = 1 , … , M n, are calculated with thevariational principle 共the ordinary Ritz variational principle
is easier to use than the variational principle for the scaled
Schrödinger equation兲.
In the free ICI method summarized above, we have two
freedoms: one is the choice of the g function and the other isthe choice of the initial function ␺0. First, we explain the
choice of the g function. Referring to the potential of thehydrogen molecule given by Eq. 共5兲, we examined two different g functions

␭21 − ␮21 ␭22 − ␮22
␳
␭1
␭2

共10兲

␭21 − ␮21 ␭22 − ␮22
−
+ ␳.
␭1
␭2

共11兲

g1 =
and
g2 = −

The function g1 is a product of the inverse of the Coulomb
potentials, g1 = 1 / Vne共1兲Vne共2兲Vee, and g2 is a sum of the inverse of each Coulomb potential, g2 = 1 / Vne共1兲 + 1 / Vne共2兲
+ 1 / Vee. In the free ICI calculation, the g2 function is more
preferable than the g1 function for the following two reasons.
First, when one of the three terms in the g1 function becomes
zero, the other terms can be arbitrary, so that the restriction
may not be complete. Second, the function g1 is more complex than the function g2, so that the choice of g1 makes the
ICI wave function more complex than the choice of g2. The
g2 function produces more plain and flexible complementsfunctions than the function g1. However, the function g1 also
produces functions that are singular, so that we have to eliminate such functions from our expansion bases of the exact
wave function: the wave function must be integrably finite
by its definition.
The initial functions ␺0 that were adopted in the present
calculations are of two kinds. One is the simplest possiblefunction given by␺共1兲
0 = exp关− ␣共r1a + r2a + r1b + r2b兲兴 = exp关− ␣共␭1 + ␭2兲兴.
共12兲
The spin part is singlet and antisymmetric, so that the spatial
part is symmetric. This is a product of the 1s orbitals centered on the two protons. The other initial function we usedin this study is given by
L

␺共L兲
0 =

exp关− ␣共␭1 + ␭2兲兴␳l−1 .
兺
l=1

共13兲

When L = 1, the initial function given by Eq. 共13兲 is identicalwith that given by Eq. 共12兲. This choice was based on thesuggestion due to Kolos and Wolniewicz 关5兴 that an inclusionof the functions of higher power in ␳ accelerates the convergence of the expansion of the wave function. Kolos and
Wolniewicz 关5兴 also noted that the inclusion of the functions
higher powers in ␭ and ␮ did not accelerate the convergence.
The orbital exponent ␣ in Eqs. 共12兲 and 共13兲 is a nonlinear
parameter and may be optimized variationally at each iteration n, but the ICI theory claims that we can get the exact
wave function with only linear expansions and so we fixed ␣
to 1.1. The internuclear distance R was also fixed to theexperimental value R = 1.4011 a . u., which is also used inother calculations to be compared with the present one. The
calculations at R = 1.4 a . u. were also carried out for comparison with the previous results in the literature.
The initial functions ␺0 given above are analytical functions, and by inserting these initial functions into Eq. 共8兲, onegets analytical functions that are necessary for constructing
the first iteration function ␺1: the kinetic operator in the
Hamiltonian is essentially a differentiation operator, and the062502-2

FREE ITERATIVE-COMPLEMENT-INTERACTION …

PHYSICAL REVIEW A 72, 062502 共2005兲

TABLE I. Ground-state energy of the hydrogen molecule calculated with the g function given by Eq.
共10兲.
Initial function 共L兲a

Number of
iterations 共n兲

Number of basis
functions 共M n兲

0
1
2
3

1
11
186
1156

1

Best value

Total energy 共a.u.兲
−0.999 780 120 198 08−1.169 313 411 900 71−1.174 470 133 104 04−1.174 475 901 628 54−1.174 475 931 397 74a

L is defined in Eq. 共13兲.

potential function and the g function are multiplicative functions. These operations are automatically done by using symbolic operation programs like MAPLE 关18兴. Then, one gets ␺1of the SICI in an analytical form. In the free ICI, one extracts
from this ␺1 all independent analytical functions, selects only
such functions that do not give divergence in the calculationsof Hamiltonian and overlap integrals, and collects them as
兵␾k其共1兲 , k = 1 , … , M 1. Then, one gives an independent variable ck,1 to each ␾k and expands one’s ␺1 of free ICI as ␺1
M1= 兺k=1ck,1␾k. The variables 兵ck,1其 are determined by the variational principle by solving the secular equation. The Hamiltonian and overlap matrices are calculated without the singularity problem for the existence of the g function. After the
diagonalization one gets the free ICI functions of the first
iteration, ␺1. The lowest solution is an approximation to theground state and the second lowest solution is an approximation to the first excited state, etc. Putting this ␺1 into Eq. 共8兲
and doing the same procedure as above, one obtains ␺2. One
repeats this iteration cycle until one gets the convergence in
one’s desired accuracy. Since the secular equation at each
iteration cycle is due to the variational principle, the energy
approaches from above the true energies, for both groundand excited states. Note that in the free ICI method, the next
iteration cycle does not require the variables 兵ck,n其 of theformer cycles, so that one can get the nth iteration free ICI
functions 兵␾k其共n兲 directly from ␺0 by applying n times theoperator part of Eq. 共8兲 to ␺0. This means that the accumulation of errors during the iteration process does not occur inthe free ICI calculations.
The Hamiltonian and overlap integrals over the basis
functions 兵␾k其共n兲 were calculated analytically by applying
and extending the method reported by James and Coolidge
关3兴. We reduced the basic integrals into the forms that can behandled with MAPLE. The details were summarized briefly in
the Appendix for convenience. As we proceed with iterations, the number of the independent functions 兵␾k其共n兲 increases and they may involve the functions whose overlap
integrals are rather close to unity. This means that the calculations must be performed in high accuracy, and we kept
60-decimal-figure accuracy throughout the calculations. The
secular equation was also solved in high precision using the
GMP 共GNU multiple precision arithmetic兲 library 关19兴.
III. RESULTS
A. Energy and wave function

We first performed the free ICI calculations of the hydrogen molecule using the scaling function g1 given by Eq. 共10兲and the simple initial function ␺共1兲0 given by Eq. 共12兲. We
summarize in Table I the calculated energies at different iteration cycles. At the bottom of the table, we gave the bestvalue obtained in the present free ICI calculations. One cansee that as the iteration proceeds, the energy approaches thebest value from above. Already at the second iteration withdimension 186, the energy is correct to five decimal figures,
which is by far beyond the chemical accuracy. At the third
iteration, we obtain the energy correct to seven decimal figures. We will show later that our wave function is different
from those existing in the literature.
We next performed the free ICI calculations using thescaling function g2 given by Eq. 共11兲 and a set of initial
functions ␺共L兲
0 共L = 1–6兲 given by Eq. 共13兲. The results are
given in Table II. When L is unity, the initial function is thesame as that used in Table I, so that the difference is due only
to the difference in the g function. In Table II the number ofthe complement functions at the fifth iteration is 832 and theenergy is −1.174 475 917 a . u., which is lower than the energy of Table I at the third iteration, −1.174 475 901 a . u.
with 1156 independent functions. This means that the g2
function given by Eq. 共11兲 is more efficient than the g1 function given by Eq. 共10兲. Actually, a multiplication of the g1
function increases the orders of all variables ␭1, ␭2, and ␳ byunity or minus unity, but the g2 function increases the orderof only one of the variables ␭1, ␭2, or ␳ by unity, so that theg2 function can produce more flexible basis functions than
the g1 function, leading to more efficient basis functions.
Table II gives a comparison of the usage of different initial functions, all with the same g2 function. As L increases,
the initial function already includes explicitly the interelectron distance ␳. If ␳ is essential, it is better to include it fromthe beginning of the calculations and this is the case as seenfrom the table: the energy with L = 1 and n = 6 is
−1.174 475 930 732 a . u. with 1788 functions, which ishigher in energy than the case of L = 4 and n = 4 with 1667
functions, −1.174 475 931 331 a . u. By using better-quality
initial function, we can get better energy with smaller number of basis functions. However, probably more important
implication of Table II is that we can always get very accurate results as we perform iterations, independent of the quality of the initial functions. When the number of the basis
functions is similar, the calculated energy is more-or-lesssimilar.
Between the differences in the g function and in the initial
function, the difference in the g function causes a larger dif-

062502-3

PHYSICAL REVIEW A 72, 062502 共2005兲

KUROKAWA, NAKASHIMA, AND NAKATSUJI

TABLE II. Ground-state energy of the hydrogen molecule calculated with the g function given by Eq.
共11兲.
La

nb

M nc

Energy 共a.u.兲

La

nb

M nc

Energy 共a.u.兲

1
1
1
1
1
1
1

0
1
2
3
4
5
6

1
5
30
114
343
832
1788

−0.999 780 120 198 080−1.164 409 776 802 471−1.172 712 604 472 602−1.174 434 056 534 598−1.174 475 391 331 891−1.174 475 917 716 333−1.174 475 930 732 940

2
2
2
2
2
2

0
1
2
3
4
5

2
9
57
215
624
1459

−1.138 078 114 832 792−1.169 819 224 009 141−1.174 397 294 989 745−1.174 474 888 511 114−1.174 475 918 013 359−1.174 475 930 608 043

3
3
3
3
3

0
1
2
3
4

3
17
114
414
1119

−1.142 973 092 050 475−1.173 036 862 957 403−1.174 471 341 851 610−1.174 475 901 329 452−1.174 475 930 609 382

4
4
4
4
4

0
1
2
3
4

4
25
174
630
1667

−1.143 082 324 090 173−1.173 448 296 529 433−1.174 475 069 845 826−1.174 475 928 234 319−1.174 475 931 331 490

5
5
5
5
5

0
1
2
3
4

5
33
260
951
2441

−1.143 084 264 123 815−1.173 489 787 488 802−1.174 475 797 933 379−1.174 475 931 085 682−1.174 475 931 391 155

6
6
6
6
6

0
1
2
3
4

6
41
346
1276
3246

−1.143 084 849 985 530−1.173 497 602 544 993−1.174 475 872 960 003−1.174 475 931 318 436−1.174 475 931 397 736

a

L of Eq. 共13兲.
Iteration number.
c
M n of Eq. 共9兲: number of basis functions at the nth iteration.
bference of the calculated results. We see from the comparisonof Tables I and II that the energy calculated with 1156 functions produced with the g1 function, −1.174 475 901 a . u. ishigher than the one calculated with only 624 functions 共abouthalf兲 produced with the g2 function, −1.174 475 918 a . u.
This example emphasizes again the importance of the choiceof the g function.
In the literature, many earlier calculations were done forthe bond length of 1.4 a.u. So we also performed the free ICI
calculations with R = 1.4 a . u. and the results are given in
Table III. The g function was due to Eq. 共11兲 and the initial
function was due to Eq. 共13兲 with L = 6. Comparing the energies for R = 1.4 and 1.4011, the latter is lower, showing thatthe minimum-energy bond length should be closer to thelatter. The experimental bond length is 1.401 12 a.u. 关20兴.
Many important studies have been reported using the explicitly correlated wave functions of the hydrogen molecule关13兴. We summarize in Table IV some of the representative

studies. As the present result, we gave our best results shownin Tables II and III for R = 1.4011 and 1.4 a.u., respectively.
First, we notice that our best results are certainly the best
worldwide in the literature for both cases of R = 1.4 and1.4011 a.u. The second best energies were reported by Cencek and Kutzelnigg 关9兴 for R = 1.4 a . u. and by Cencek andRychlewski 关10兴 for R = 1.4011 a . u. They used Gaussianfunctions.
Our free ICI wave functions are composed of the Slatertype functions and are written as␺ = 兺 ci共1 + p12兲exp关− ␣共␭1 + ␭2兲兴␭m1 i␭n2i␮1ji␮k2i␳li ,i

共14兲
where p12 is an electron exchange operator. This wave function is very simple and similar to the original wave function
due to James and Coolidge 关3兴. Our wave function and

TABLE III. Ground-state energy of hydrogen molecule at R = 1.4 a . u. calculated with the g function
given by Eq. 共10兲.
Initial function 共L兲a

Number of
iterations 共n兲

Number of basis
functions 共M n兲

Energy 共a.u.兲

6
6
6
6

0
1
2
3

6
41
346
1276

−1.143 006 074 717−1.173 494 068 035−1.174 475 655 534−1.174 475 714 138a

L of Eq. 共13兲.
062502-4

FREE ITERATIVE-COMPLEMENT-INTERACTION …

PHYSICAL REVIEW A 72, 062502 共2005兲

TABLE IV. History of accurate calculations of the hydrogen molecule.
Type of wave function

Reference

H-H distance 共au兲

Total energy 共au兲

Hartree-Fock
Full CI 共33␴27␲18␦16␾兲Hylleraas type
Gaussian functions
Free ICI 共extended Hylleraas type兲
Hylleraas type
Gaussian functions
Free ICI 共extended Hylleraas type兲

Sundholma
Liu-Hagstromb
James-Coolidgec
Cencek-Kutzelnigge
Present
Wolniewiczf
Cencek-Rychlewskig
Present

1.4
1.4
1.4
1.4
1.4
1.40111.40111.4011−1.133 629 573−1.1743 304 3−1.173 539d
−1.174 475 714 037−1.174 475 714 138−1.174 475 930 742−1.174 475 931 39−1.174 475 931 397 74

aReference 关21兴.
Reference 关11兴.
c
Reference 关3兴.
d
The energy and the wave function were recalculated as shown in Table VIII.
eReference 关9兴.
f
Reference 关7兴.
g
Reference 关10兴.
bJames-Coolidge wave function differ only in the powers mand n of the variables ␭1 and ␭2: m and n are always nonnegative in the James-Coolidge wave function, but they canbe even negative in the present free ICI wave function. In
other words, our ICI theory starting from the initial functions
given by Eqs. 共12兲 and 共13兲 generates not only non-negative
power terms of ␭1 and ␭2, but also negative power terms ofthem. Such functions are theoretically more relevant for describing the exact wave function of the hydrogen moleculethan the original James-Coolidge functions. We show later inthis paper some of the important roles of these negativepower functions in our wave function. The high quality ofthe present free ICI wave function would be attributed to theexistence of these negative m and n terms.
Kolos and Wolniewicz 关5–7兴 extended the JamesCoolidge wave function to a more general form like

␺ = 兺 C共1 + p12兲exp关− ␣ ␭1 − ¯␣␭2兴␭m1 ␭n2␮1j ␮k2␳l⫻ 关exp共␤ ␮1 + ¯␤␮2兲 + 共− 1兲 j+kexp共− ␤ ␮1 − ¯␤␮2兲兴,共15兲
but they always used non-negative m and n in their wave
functions. A reason for extending the James-Coolidge wave
function was that the original James-Coolidge wave function

did not have a proper asymmetric form as the interproton
distance increases up to infinity. On the other hand, Cencek
and co-workers 关8–10兴 obtained very accurate wave functions using the Gaussian-type functions22⌿ = 兺 C共1 + p12兲exp关− ␣ r1C− ¯␣r2C⬘ − ␤ r12
兴.
2

共16兲

They performed optimizations of a large number of nonlinear
parameters. Primitive questions about this type of wave function are that how well the cusp conditions 关17兴 are describedand how important these properties are in the actual wave
functions. Nevertheless, the easiness in the integral evaluations with the Gaussian bases is very important in actual
calculations of atoms and molecules.
Table IV shows also the full CI energy 关11兴 and theHartree-Fock energy 关21兴. These calculations were performed at R = 1.40 a . u., like the original James-Coolidge calculations, while the other calculations were done at R
= 1.4011 a . u. The full CI energy is the “exact” energy withina given space of the basis set. For the 共33␴27␲18␦16␾兲Gaussian-type basis, the full CI energy is −1.174 304 3 a . u.
Referring to Table II, this energy is worse than that of thecalculation, L = 2 and n = 2 , − 1.174 397 a . u., so that we estimate that this energy would be obtained with the functions
less than 57, if we use the present free ICI method.

TABLE V. First iteration nine-term wave function of the free ICI calculation shown in Table II 共L= 2 , n = 1兲a. The energy is −1.169 819 224 a . u.
No关m , n , j , k , p兴

CoefficientNo关m , n , j , k , p兴

Coefficient

1
2
3
4
5

关0, 0, 0, 0, 0兴
关0, 0, 0, 0, 1兴
关0, 0, 0, 0, 2兴
关0, 1, 0, 0, 0兴
关0, 1, 0, 0, 1兴

1.000 000 000−0.478 943 233−0.205 371 707
0.558 786 404
1.537 900 860

6
7
8
9

关−1 , 0 , 0 , 0 , 0兴关−1 , 0 , 0 , 0 , 1兴关−1 , 0 , 2 , 0 , 0兴关−1 , 0 , 2 , 0 , 1兴

−0.169 894 927
0.310 672 513
0.064 715 165
0.234 383 618

aEach basis function is normalized to unity.
062502-5

PHYSICAL REVIEW A 72, 062502 共2005兲

KUROKAWA, NAKASHIMA, AND NAKATSUJI

TABLE VI. Best five-terma extended James-Coolidge-type wave function of the hydrogen molecule
within −1 艋 mi, ni 艋 1, 兩mi − ni兩 艋 1, 0 艋 ji, ki 艋 2, and 0 艋 li 艋 1. R = 1.4 a . u. and ␣ = 0.95. The calculated
energy was −1.172 276 534 a . u.
No关m , n , j , k , l兴

CoefficientNo关m , n , j , k , l兴

Coefficient12
3

关0, 0, 0, 0, 0兴
关0, 0, 0, 0, 1兴
关0, 0, 0, 2, 0兴

1.000 000 000
0.718 124 997
0.198 062 046

45

关0, 0, 1, 1, 0兴关0 , −1 , 0 , 2 , 1兴

−0.055 207 501
0.081 476 902

aEach basis function is normalized to unity.

B. Importance of the terms with negative powers of ␭i

The free ICI method produces automatically the analyticalbasis functions that are necessary to expand the exact wave
function under a choice of initial function ␺0 and the g function. In the present choice of the initial function, the ICI
method always produced terms of ␭i that have positive, zero,and negative powers. Since these terms were produced theoretically by the ICI method, they are considered to be necessary for constructing the exact wave function. However, suchterms were not included in the previous calculations byJames and Coolidge and Kolos and Wolniewitz. So we examine here the importance of such terms in constructing theaccurate wave functions of the hydrogen molecule. Note thatfor the helium atom the importance of similar terms withnegative powers has been pointed out by Kinoshita 关22兴.
First, the importance of the terms with negative powers of
␭i is understood from the high performance of the present
free ICI calculations as presented in the preceding section. In
Table V, we show the expansion coefficients of the ninebasis functions obtained at the first iteration of the free ICI
calculation shown in Table II 共L = 2 and n = 1兲. This calculation is based on the g2 function given by Eq. 共11兲 and theinitial function ␺共2兲0 given by Eq. 共13兲. The coefficients given
in Table V are those for the terms of Eq. 共14兲 designated by关m , n , j , k , l兴, which are normalized to unity, and the firstiteration free ICI wave function ␺共2兲
1 was normalized suchthat the coefficient of the initial function, 关0,0,0,0,0兴, is unity.
One can see from Table V that the first five terms, Nos. 1–5,
are conventional ones, but the other four terms are terms
whose m are negative. You see that terms 6, 7, and 9 havelarge coefficients, showing the importance of the negativepower terms of ␭i.
Kolos and Roothaan 关4兴 calculated five-term JamesCoolidge wave functions using only positive powers and obtained an energy of −1.171 619 a . u. We have examined thebest possible five-term extended James-Coolidge-type wave

function allowing the powers in Eq. 共14兲 to change within−1 艋 mi, ni 艋 1, 兩mi − ni兩 艋 1, 0 艋 ji, ki 艋 2, and 0 艋 li 艋 1. The
bond length R and the orbital exponent ␣ are the same as
those used by Kolos and Roothaan 关4兴, R = 1.4 a . u. and ␣= 0.95. Table VI shows the calculated best wave function. It
includes one term that has negative n and its coefficient is thesecondsmallest.
The
calculated
energywas−1.172 276 534 a . u. We also performed a similar calculation
within the positive range of 0 艋 mi, ni 艋 1, 0 艋 ji, ki 艋 2, and0 艋 li 艋 1 and obtained the best wave function as shown in
Table VII. The calculated energy was −1.171 998 568 a . u.
The last term 关1,1,0,0,0兴 was not included in the original
5-term Kolos-Roothaan wave function, though it was included in their 12-term wave function. They included instead
the term 关0,1,0,0,0兴. The other four terms were common in
all of these three wave functions. It is concluded that within
5-term James-Coolidge wave functions of hydrogen molecule, the best result is obtained by allowing negative-powerterms of ␭i.
Another piece of evidence was obtained by performing
calculations similar to the original James-Coolidge one. We
recalculated the James-Coolidge 13-term wave function andshow the result in Table VIII. The calculated energy was−1.173 539 685 a . u., which was slightly different from theoriginal one, E = 1.173 501 a . u., reported by James andCoolidge later 关3兴. 共For the 5- and 11-term wave functions,the present recalculations gave the same energies as thosegiven in the original paper 关3兴.兲 Then, we replaced the smallest two terms, term Nos. 5 and 7, of the original 13 termswith the two terms from −1 艋 m, n 艋 1, 兩m − n兩 艋 1, 0 艋 j, k艋 2, and 0 艋 l 艋 1. The lowest-energy result is given in Table
IX that has the energy −1.173 962 233 a . u., which is substantially lower than the original James-Coolidge energy
given in Table VIII. In Table IX, term Nos. 12 and 13 are thereplaced ones and term No. 12 has positive m and n, but term
No. 13 has negative m and n. Again, this shows the importance of the negative-power terms in the calculations of thehydrogen molecule.

TABLE VII. Best five-terma James-Coolidge-type wave function of the hydrogen molecule within thepositive range of 0 艋 mi, ni 艋 1, 0 艋 ji, ki 艋 2, and 0 艋 li 艋 1. R = 1.4 a . u. and ␣ = 0.95. The calculated energy
was −1.171 998 568 a . u.
No关m , n , j , k , l兴

CoefficientNo关m , n , j , k , l兴

Coefficient12
3

关0, 0, 0, 0, 0兴
关0, 0, 0, 0, 1兴
关0, 0, 0, 2, 0兴

1.000 000 000
0.802 238 155
0.282 362 668

4
5

关0, 0, 1, 1, 0兴
关1, 1, 0, 0, 0兴−0.067 134 733−0.087 098 459

aEach basis function is normalized to unity.
062502-6

FREE ITERATIVE-COMPLEMENT-INTERACTION …

PHYSICAL REVIEW A 72, 062502 共2005兲

TABLE VIII. Recalculated 13-terma wave function of James and Coolidge with R = 1.4 a . u. and ␣= 3 / 4. The calculated energy was −1.173 539 685 a . u.
No.

关m , n , j , k , l兴

CoefficientNo.

关m , n , j , k , l兴

Coefficient

1
2
3
4
5
6
7

关0, 0, 0, 0, 0兴
关0, 0, 0, 0, 1兴
关0, 0, 0, 0, 2兴
关0, 0, 0, 2, 0兴
关0, 0, 0, 2, 1兴
关0, 0, 1, 1, 0兴
关0, 0, 1, 1, 1兴

1.000 000 000
0.553 586 789−0.099 140 955
0.195 335 087
0.044 417 387−0.062 129 657−0.006 994 087

8
9
10
11
12
13

关1, 0, 0, 0, 0兴
关1, 0, 0, 2, 0兴
关1, 0, 1, 1, 0兴
关1, 0, 2, 0, 0兴
关2, 0, 0, 0, 0兴
关1, 0, 0, 0, 1兴

−0.933 408 702−0.052 106 263
0.044 627 723−0.074 756 346
0.315 198 529−0.188 359 478

aBasis functions are normalized to unity.

We also have some intuitive explanations about the importance of the negative-power terms of ␭i. Figure 1 shows
the illustrations of the functions ␾ = exp共−␭兲␭m␮2 with positive and negative m. Figure 1共a兲 is for m = 1,2,3 and Fig. 1共b兲for m = 0 , −1 , −2, and the two protons were assumed to be atz = 0 and 1. Note that the variables ␭i and ␮i describe thewave function outside and inside, respectively, of the twonuclei. We see that the functions with positive m do not show
a nice cusplike behavior near the positions of the nuclei, but
those with negative m show reasonable cusplike behaviors at
the nuclear positions. This also explains the importance ofthe functions with negative powers of ␭i.

and electron cusp condition. The nuclear cusp condition isexpressed for the hydrogen molecule as 关4兴
FN共␭2, ␮2兲 =

−1
⳵⌿
⳵⌿
␭1
− ␮1
R␭1⌿
⳵ ␭1
⳵ ␮1

冊册

␭1=1,␮1=−1= ZA = 1,
共17兲

which is the condition when electron 1 is at proton a. The
left-hand side 共LHS兲 of Eq. 共17兲 is a function of the electron
coordinate 2, so that Eq. 共17兲 must be satisfied at any coordinate of electron 2, except at the additional singular point␭2 = 1 and ␮2 = −1, where electron 2 also collides with protona. The electron cusp condition is given by

C. Cusp properties

The exact wave function must satisfy some necessary
conditions like the variational condition, virial theorem, cusp
condition, etc. In the free ICI wave function given by Eq.
共14兲, the coefficients 兵ci其 were calculated with the variational
principle. We did not optimize the exponent ␣, since thetheorem of ICI guarantees that we can approach the exact
result with only linear variables. The cusp property is known
to be a rather sensitive property, reflecting the quality of thecalculated wave functions. Below, we examine how well thepresent ICI wave functions satisfy the cusp conditions.
Kato 关17兴 presented two cusp conditions for the exact
many-electron wave functions: the nuclear cusp condition冋 冉

Fe共␭i, ␮i兲 =

冋 册

1
1 ⳵⌿
= ,
⌿ ⳵ r12 r12=0 2

共18兲

which is the condition when two electrons meet together. The
LHS is also the function of the variables other than ␳ = r12
and so it must be 1 / 2 everywhere except for the additional
singular points.
Figure 2 shows the plots of the function FN along theiteration process of the free ICI calculations using the g2
function and the initial function ␺共6兲
0 : at the fourth iteration ofthis calculation we have obtained the worldwide best valueof the energy of the hydrogen molecule shown in Tables II
and IV. Figure 2共a兲 is along the z coordinate that connects

TABLE IX. The wave function constructed from the 11 terms of James and Coolidge 13 terms plus 2
terms from −1 艋 m, n 艋 1, 兩m − n兩 艋 1, 0 艋 j, and k 艋 2, and 0 艋 l 艋 1 for the H2 molecule with R = 1.4 a . u. and␣ = 3 / 4.a,b The calculated energy was −1.173 962 233 a . u.
No.

关m , n , j , k , l兴

CoefficientNo.

关m , n , j , k , l兴

Coefficient

1
2
3
4
5
6
7

关0, 0, 0, 0, 0兴
关0, 0, 0, 0, 1兴
关0, 0, 0, 0, 2兴
关0, 0, 0, 2, 0兴
关0, 0, 1, 1, 0兴
关1, 0, 0, 0, 0兴
关1, 0, 0, 2, 0兴

1.000 000 000
0.650 858 318−0.059 439 543
0.138 956 703−0.041 545 078−0.933 425 960−0.018 305 773

8
9
10
11
12
13

关1, 0, 1, 1, 0兴
关1, 0, 2, 0, 0兴
关1, 0, 0, 0, 1兴
关2, 0, 0, 0, 0兴
关1, 1, 0, 0, 1兴关−1 , −1 , 0 , 2 , 1兴

0.015 613 455−0.033 975 753−0.404 443 002
0.337 386 255
0.071 193 197
0.031 686 115aBasis functions are normalized to unity.
The last two terms, Nos. 12 and 13, are the selected terms.

b062502-7

PHYSICAL REVIEW A 72, 062502 共2005兲

KUROKAWA, NAKASHIMA, AND NAKATSUJI

FIG. 1. 共Color兲 共a兲 Plot of the function ␾ = exp共−␭兲␭m␮2 with
m = 1 共green兲, m = 2 共red兲, and m = 3 共blue兲. The two protons are atz = 0 and z = 1 a . u. 共b兲 Plot of the function ␾ = exp共−␭兲␭m␮2 withm = −2 共green兲, m = −1 共red兲, and m = 0 共blue兲. The two protons areat z = 0 and z = 1 a . u.

two protons lying at z = 0 and z = 1.4011 a . u., and Fig. 2共b兲 isalong an axis lying parallel to the z coordinate and distantfrom it by 1 a.u. The figures show the local values of thecusp function FN when electron 2 moves along these axes.
Though the cusp value of the free ICI wave function is very
different from unity at the first iteration 共n = 1兲, it becomescloser and closer to unity as the iteration proceeds. At n = 4,
the cusp values are unity everywhere except at the positionof nucleus a 关the origin of Fig. 2共a兲兴, where an additional
singularity occurs. There, not only does nucleus a collidewith the two electrons, but also the two electrons collide witheach other: this is a new singularity point not expressed by
Eq. 共15兲. It is remarkable that at z = 1.4011 a . u. of Fig. 2共a兲,the value of FN becomes unity at n = 4: there, the two cusp
conditions for the two pairs—electron 1, nucleus a, and electron 2, nucleus b—are satisfied independently.
Figure 3 shows plots of the function Fe along the iteration
process of the same free ICI wave functions given in Fig. 2
above. Figure 3共a兲 shows the value of the cusp function Fe
when the colliding two electrons move along the z axis andFig. 3共b兲 along the axis separated by 1 a.u. from the z axis.
Again, the cusp value is improved as the iteration proceeds:

FIG. 2. 共Color兲 共a兲 Plots of the nuclear-electron cusp function
FN along the z coordinate that connects two protons lying at z = 0
and z = 1.4011 a . u. at different iteration processes of the free ICI
calculations using the g2 function and the initial function ␺共6兲
0 : iteration number n = 1 共blue兲, n = 2 共green兲, n = 3 共violet兲, and n = 4
共red兲. The wave function at n = 4 is our best wave function reported
in this paper. 共b兲 Plots of the nuclear-electron cusp function FN
along the axis lying parallel to the z coordinate that connects two
protons lying at z = 0 and z = 1.4011 a . u. at different iteration processes of the free ICI calculations using the g2 function and theinitial function ␺共6兲
0 : iteration number n = 1 共blue兲, n = 2 共green兲, n= 3 共violet兲, and n = 4 共red兲. The wave function at n = 4 is our best
wave function reported in this paper.

the cusp value of the wave function at n = 1 is far from onehalf, but it comes closer and closer to one-half as n increases.
At the nuclear positions a and b along the z axis, the cusp
value is different from one half even for the best wave function since there an additional singularity occurs: there, twocolliding electrons collide farther with the nucleus. Except atthese two points, the cusp function Fe of our best free ICI
wave function 共n = 4兲 is flat and takes one-half everywhere
along the z axis and along the axis apart 1 a.u. from the z
axis.
IV. CONCLUSION

The free ICI method developed in earlier papers has been
applied to the hydrogen molecule to calculate its very accu-

062502-8

FREE ITERATIVE-COMPLEMENT-INTERACTION …

PHYSICAL REVIEW A 72, 062502 共2005兲

However, it was also true that any initial functions havegiven the converging series of energies to the same bestvalue. This was also true for the two g functions examined
here. We have also shown that the cusp values, rather sensitive properties of the wave function, were substantially improved as the iteration proceeds. Recently, we have extended
the present methodology to the relativistic case 关23兴. Accurate and systematic descriptions of the relativistic electronic
structures of many-electron systems are a very important
subject that must be explored in the near future.
ACKNOWLEDGMENT

This research was supported by a grant for Creative Scientific Research from the Ministry of Education, Science,
Culture, and Sports of Japan.
APPENDIX

We briefly summarize here the method of analytical calculations of the integrals that appear in the variational calculations of the free ICI wave function given in the form of Eq.
共14兲. The formulation is a small generalization of the one
presented by James and Cooligde in their original paper 关3兴.
Ruedenberg 关24兴 also gave formulas related to the present
case.
The integrals we have to calculate are the Hamiltonianand overlap integrals
Hij =FIG. 3. 共Color兲 共a兲 Plots of the electron-electron cusp function
Fe along the z coordinate that connects two protons lying at z = 0
and z = 1.4011 a . u. at different iteration processes of the free ICI
calculations using the g2 function and the initial function ␺共6兲
0 : iteration number n = 1 共blue兲, n = 2 共green兲, n = 3 共violet兲, and n = 4
共red兲. The wave function at n = 4 is our best wave function reported
in this paper. 共b兲 Plots of the electron-electron cusp function Fe
along the axis lying parallel to the z coordinate that connects two
protons lying at z = 0 and z = 1.4011 a . u. at different iteration processes of the free ICI calculations using the g2 function and theinitial function ␺共6兲
0 : iteration number n = 1 共blue兲, n = 2 共green兲, n= 3 共violet兲, and n = 4 共red兲. The wave function at n = 4 is our best
wave function reported in this paper.

rate wave functions. We could obtain the best variational
wave functions in the literature for both R = 1.4011 and 1.4
a.u. with a rather simple Hylleraas-James-Coolidge form.
The difference of the present free ICI wave function from theprevious James-Coolidge ones was the existence of thenegative-power terms of ␭i. These terms were generated automatically by the free ICI formalism, so that these functions
are believed to be essential for effective descriptions of theexact wave function of the hydrogen molecule. We haveshown some results showing the importance of thesenegative-power terms.
We have examined two different g functions and several
initial functions. A clever choice of the g function was more
important than the different choices of the initial function.

冕冕
冕冕

Sij =

␾ iH ␾ j d ␶ 1d ␶ 2 ,共A1兲␾ i␾ j d ␶ 1d ␶ 2 ,where the functions ␾i have the form

␾i = 共1 + p12兲exp关− ␣共␭1 + ␭2兲兴␭m1 i␭n2i␮1ji␮k2i␳li ,共A2兲

where m and n are positive and negative integers and zeroand j, k, and l are non-negative integers. Applying the
Hamiltonian on ␾i, we obtainj nj jj kj ljH␾i = 兺 d j共1 + p12兲exp关− ␣共␭1 + ␭2兲兴␭m
1 ␭2 ␮1 ␮2 ␳ ,j

共A3兲which are similar in form to Eq. 共14兲, and this process can beeasily done with the use of MAPLE or other similar program.
The Jacobian for the elliptic coordinate is given byd ␶ 1d ␶ 2 =R6 2
共␭ − ␮21兲共␭22 − ␮22兲d␭1d␭2d␮1d␮2d␸1d␸2 .
64 1
共A4兲

Then, all the integrals reduce to the sums of the integrals
I=062502-9

冕

n j k lexp关− 2␣共␭1 + ␭2兲兴␭m
1 ␭2␮1␮2␳ d␭1d␭2d␮1d␮2d␸1d␸2 .

共A5兲PHYSICAL REVIEW A 72, 062502 共2005兲

KUROKAWA, NAKASHIMA, AND NAKATSUJI

Now, we summarize the necessary formulas to calculatethese integrals. First, the variables of ␳ are transformed byusing the relation

␳2= ␭21 + ␭22 + ␮21 + ␮22 − 2 − 2 ␭1␭2␮1␮2− 2关共␭21 − 1兲共␭22 − 1兲共1 − ␮21兲共1 − ␮22兲兴1/2cos共␸1 − ␸2兲,共A6兲

D␶N = 2共2␶ + 1兲⬁⬁

␳ =兺 兺

␶=0 N=0D␶N P␶N冉冊 冉冊

Z共␯,m,n, j,k,l兲 =2

共N ⬎ 0兲.

冕

exp关− 2␣共␭1 + ␭2兲兴n j k l ␯⫻␭m
1 ␭2␮1␮2␳ M d␭1d␭2d␮1d␮2d␸1d␸2 ,where
M = 关共␭21 − 1兲共␭22 − 1兲共1 − ␮21兲共1 − ␮22兲兴1/2cos共␸1 − ␸2兲.

共A7兲

共A10兲
Then, the integral I defined by Eq. 共A5兲 is written with this Z
function as

which is Von Neuman’s expansion, whereD␶0 = 2␶ + 1,

册

共A9兲

␭1 N ␭2 N
Q
P 共␮1兲P␶N共␮2兲␭2 ␶ ␭1 ␶⫻cos N共␸1 − ␸2兲,共␶ − N兲!
共␶ + N兲!

P and Q are the unnormalized complex associated Legendre
functions of the first and second kinds, and we take uppervariable when ␭2 艌 ␭1 and lower variable otherwise.
We introduce the function Z defined bywhich is obtained from the second cosine formula and the
expansion−1冋

I = Z共0,m,n, j,k,l兲.
共A8兲

共A11兲

From Eq. 共A6兲, we obtain the recursion formula

Z共␯,m,n, j,k,l兲 = Z共␯,m + 2,n, j,k,l − 2兲 + Z共␯,m,n + 2, j,k,l − 2兲 + Z共␯,m,n, j + 2,k,l − 2兲 + Z共␯,m,n, j,k + 2,l − 2兲
− 2Z共␯,m,n, j,k,l − 2兲 − 2Z共␯,m + 1,n + 1, j + 1,k + 1,l − 2兲 − 2Z共␯ + 1,m,n, j,k,l − 2兲,case 共a兲:

Z共␯,m,n, j,k,0兲

共l,even兲,case 共b兲:

Z共␯,m,n, j,k,− 1兲

共l,odd兲.

⬁⬁which implies that the value of the function Z reduces finallyto the sum of the following terms, depending on the parity of
the integer l:

Z共␯,m,n, j,k,− 1兲 = 兺 兺 D共␶,N兲R共␶, ␯,N, j兲R共␶, ␯,N,k兲
␶=0 N=0⫻H共␶, ␯,N,m,n, ␣兲f共␯,N兲,
共A13兲

D共␶,N兲 = D␶N ,Z共␯,m,n, j,k,0兲 = L共␯,m, ␣兲L共␯,n, ␣兲U共␯, j兲U共␯,k兲f 0共␯兲,
共A14兲

R共␶, ␯,N, j兲 =

冕

1

−1

L共␯,m, ␣兲 =

冕⬁exp共− 2␣␭兲共␭2 − 1兲␯/2␭md␭,

H共␶, ␯,N,m,n, ␣兲 =1

U共␯, j兲 =

冕

1

2 ␯/2

共1 − ␮ 兲 ␮ d␮ ,
j

冕冕
␲

−␲

␲

−␲

共1 − ␮2兲␯/2 P␶N共␮兲␮ jd␮ ,

冕冕
⬁

⬁

1

1

共A15兲

⫻P␶N

冉冊 冉冊

␭1 N ␭2Qd␭1d␭2␭2 ␶ ␭1冕 冕⬁

cos␯共␸1 − ␸2兲d␸1d␸2 .

=062502-10

␭1

d␭1

1

When ␯ is odd, f 0共␯兲 is zero, so that Z共␯ , m , n , j , k , 0兲 = 0.
For case 共b兲, we substitute Eq. 共A7兲 into Eq. 共A9兲 andobtain the Z value written in the form

exp关− 2␣共␭1 + ␭2兲兴

2
2
n
␯/2
⫻␭m
1 ␭2关共␭1 − 1兲共␭2 − 1兲兴

−1

f 0共 ␯ 兲 =

共A16兲

where

For case 共a兲, the Z value is written as

where

共A12兲

d␭2 exp关− 2␣共␭1 + ␭2兲兴

1

2
2
n
n m
␯/2
⫻共␭m
1 ␭2 + ␭1␭2 兲关共␭1 − 1兲共␭2 − 1兲兴

⫻P␶N共␭2兲Q␶N共␭1兲

FREE ITERATIVE-COMPLEMENT-INTERACTION …

f共␯,N兲 =

冕冕
␲

␲

−␲

−␲

PHYSICAL REVIEW A 72, 062502 共2005兲
␯ min共j,k兲+␯

Z共␯,m,n, j,k,− 1兲 = 兺cos␯共␸1 − ␸2兲cos N共␸1 − ␸2兲d␸1d␸2 ,

N=0

共A17兲
and D␶N is given by Eq. 共A8兲. In Eq. 共A16兲, we have two兺

␶=ND共␶,N兲R共␶, ␯,N, j兲R共␶, ␯,N,k兲

⫻H共␶, ␯,N,m,n, ␣兲f共␯,N兲,

共A18兲

infinite summations, but they actually reduce to finite sums
because of the following reasons. It is nonzero only when 共1兲␶ 艌 N, from the properties of the Legendre P functions, 共2兲
N 艋 ␯ from the property of the function f共␯ , N兲, and 共3兲 ␶
艋 j + ␯ from the property of R共␶ , ␯ , N , j兲. Further since thesevariables are all non-negative integers, Eq. 共A16兲 is rewrittenaswhere min共j , k兲 means the smaller value of j and k.
The functions appearing in Eqs. 共A14兲 and 共A18兲 are all
calculated with the MAPLE program 关18兴, for example, so thatthe Hamiltonian integral Hij and the overlap integral Sij areexpressed by a linear combination of Z共0 , m , n , j , k , l兲. Thus,
after the diagonalization of the secular equation, HC= ESC,
one obtains the energy E and the coefficient vector C.

关1兴 W. Heitler and F. London, Z. Phys. 44, 455 共1927兲.
关2兴 E. A. Hylleraas, Z. Phys. 54, 347 共1929兲.
关3兴 H. M. James and A. S. Coolidge, J. Chem. Phys. 1, 825
共1933兲; 3, 129 共1935兲.
关4兴 W. Kolos and C. C. Roothaan, Rev. Mod. Phys. 32, 205
共1960兲; 32, 219 共1960兲.
关5兴 W. Kolos and L. Wolniewicz, J. Chem. Phys. 41, 3663 共1964兲;43, 2429 共1965兲; W. Kolos and L. Wolniewicz, Chem. Phys.
Lett. 24, 457 共1974兲; W. Kolos, J. Chem. Phys. 101, 1330
共1994兲.
关6兴 L. Wolniewicz, J. Chem. Phys. 99, 1851 共1993兲.
关7兴 L. Wolniewicz, J. Chem. Phys. 103, 1792 共1995兲.
关8兴 J. Komasa, W. Cencek, and J. Rychlewski, Phys. Rev. A 46,
2351 共1992兲; J. Rychlewski, W. Cencek, and J. Komasa,
Chem. Phys. Lett. 229, 657 共1994兲.
关9兴 W. Cencek and W. Kutzelnigg, J. Chem. Phys. 105, 5878
共1996兲.
关10兴 W. Cencek and J. Rychlewski, cited in Ref. 关13兴.
关11兴 J. W. Liu and S. Hagstrom, Phys. Rev. A 48, 166 共1993兲.
关12兴 S. A. Alexander and R. L. Coldwell, J. Chem. Phys. 121,
11557 共2004兲.
关13兴 J. Rychlewski, and J. Komasa, in Explicitly Correlated Wave
Functions in Chemistry and Physics—Theory and Applica-

tions, edited by J. Rychlewski 共Kluwer Academic, Dordrecht,
2003兲, pp. 91–147.
关14兴 H. Nakatsuji, J. Chem. Phys. 113, 2949 共2000兲; H. Nakatsuji
and E. R. Davidson, ibid. 115, 2000 共2001兲; H. Nakatsuji, ibid.
115, 2465 共2001兲; 116, 1811 共2002兲; H. Nakatsuji and M.
Ehara, ibid. 117, 9 共2002兲; 122, 194108 共2005兲; H. Nakatsuji,
Bull. Chem. Soc. Jpn. 78, 1705 共2005兲.
关15兴 H. Nakatsuji, Phys. Rev. Lett. 93, 030403 共2004兲.
关16兴 H. Nakatsuji, this issue, Phys. Rev. A 72, 062110 共2005兲.
关17兴 T. Kato, Commun. Pure Appl. Math. 10, 151 共1957兲.
关18兴 Computer code MAPLE, Waterloo Maple Inc., Waterloo, Ontario, Canada; see http://www.maplesoft.com
关19兴 About the GMP library, see http://swox.com/gmp/ and http://
www.cs.nyu.edu/exact/core/gmp/
关20兴 K. P. Huber and G. Herzberg, Constants of Diatomic Molecules, Molecular Spectra and Molecular Structures, Vol. 4
共Van Nostrand Reinhold, New York, 1979兲.
关21兴 D. Sundholm, Chem. Phys. Lett. 149, 251 共1988兲.
关22兴 T. Kinoshita, Phys. Rev. 105, 1490 共1957兲.
关23兴 H. Nakatsuji and H. Nakashima, Phys. Rev. Lett. 95, 050407
共2005兲.
关24兴 K. Ruedenberg, J. Chem. Phys. 19, 1459 共1951兲.

062502-11