Explicitly correlated Helium wave function and hyperspherical coordinates

TL;DR

This paper introduces a new class of explicitly correlated wave functions for Helium-like atoms using hyperspherical coordinates, achieving highly accurate ground state energies and satisfying cusp conditions.

Contribution

It proposes a novel functional form of wave functions explicitly depending on interelectronic and hyperspherical coordinates, improving accuracy for Helium atom calculations.

Findings

Achieved the best ground state energy for Helium with hyperspherical coordinates.

Wave functions satisfy Kato cusp conditions within 6.10^{-4} accuracy.

Included negative powers of R in wave functions as a necessary condition for Schrödinger equation solutions.

Abstract

Wave functions of a new functional kind have been proposed for Helium-like atoms in this work . These functions explicitly depend on interelectronic and hyperspherical coordinates. The best ground state energy for the Helium atom $ -2.903724376677 a.u.$ has been calculated with variational method with basis set of simple functions with a single exponential parameter. To the author's knowledge, this is the best result with use of hyperspherical coordinates so far. Comparable result has been obtained for the hydrogen anion. For Helium atom, our best wave functions matched the Kato cusp conditions within the accuracy below $6.10^{-4} $. An important feature of proposed wave functions is the inclusion of negative powers of $R=\sqrt(r^{2}_{1}+r^{2}_{2})$ in combination with positive powers of $r_{12}$ into the wave function. We showed that this is necessary condition for proposed wave function to be a formal solution of Schr\"odinger equation.

Tables

Table 1. Table 1: Ground state energies of He and H - (at.units.) using wave functions of the type ( 10 ). Number of basis functions after elimination of linear dependencies N 𝑁 N , infinite summation truncations I m , J m , K m subscript 𝐼 𝑚 subscript 𝐽 𝑚 subscript 𝐾 𝑚 I_{m},J_{m},K_{m} and optimised exponents ζ 𝜁 \zeta for each basis set are also shown.

He
Basis set	$N$	$I_m$	$J_m$	$K_m$	$\zeta$	Energy
A160	100	3	3	5	2.429563	$-2.90372198$
A438	231	5	5	5	2.460426	$-2.90372427$
A438a	236	6	5	5	2.680000	$-2.90372427$
A860	456	7	7	5	2.249919	$-2.90372433$
Exacta						$-2.90372438$
H-
A860	530	7	7	5	0.851281	$-0.52775094$
Exacta						$-0.52775102$
aRef.[11]

Table 2. Table 2: Ground state energy (at.units) of He atom and H - ion using ( 12 ).

He
Basis set	$N$	$\zeta$	Energy
B288	137	2.300232	$-2.9037240795$
C288	129	2.063932	$-2.9037242527$
B595	222	2.197048	$-2.9037243113$
C595	193	2.182277	$-2.9037243717$
B1068	404	2.750647	$-2.9037243755$
C1068	390	2.342957	$-2.9037243767$
C944a	351	2.484067	$-2.9037243750$
C950b	366	2.339823	$-2.9037243766$
Exacta			$-2.9037243770$
H-
C1068	358	0.738130	$-0.5277509907$
Exacta			$-0.5277510165$
aRef.[11]

Table 3. Table 3: Cusp parameters [see Eqs.( 18 - 19 )]

Basis set	${\nu }_1$	${\nu }_{12}$	$⟨\delta \left({𝐫}_1\right)⟩$	$⟨\delta \left({𝐫}_{12}\right)⟩$
A438a	-1.974804	0.491942	1.809316	0.106425
A438	-1.976622	0.492143	1.809494	0.106435
A860	-1.985193	0.494084	1.809791	0.106400
B1068	-1.999224	0.498531	1.810405	0.106349
C1086	-2.000097	0.499408	1.810432	0.106348
C944a	-2.000034	0.496380	1.810429	0.106372
C950b	-2.000014	0.498969	1.810429	0.106351
Ref.values	$-2{.}^a$	${0.5}^a$	${1.810429}^b$	${0.106345}^b$
aExact, see [24]
bRef.[25]

Equations

$$H=$$

$$R=\sqrt{r_1^2+r_2^2}$$

$$t=\frac{r_2^2-r_1^2}{r_2^2+r_1^2}$$

$$r_1=\frac{R\sqrt{1-t}}{\sqrt{2}}$$

$$r_2=\frac{R\sqrt{1+t}}{\sqrt{2}}$$

$$H=$$

$$V_{nuc}=\frac{Z\sqrt{2}}{R}(\frac{1}{\sqrt{1-t}}+\frac{1}{\sqrt{1+t}})$$

$$V_{nuc}=-\frac{2\sqrt{2}Z}{R\sqrt{1-t^2}}(1-\sum _{k=1}^{\infty }\frac{(4k-3)!!}{4k!!}{t}^{2k})$$

$$V_{nuc}=-\frac{2Z}{R\sqrt{1-t^{2}}}\Big{(}1+\frac{1}{2}\sqrt{1-t^{2}}+\sum_{k=2}^{\infty}\frac{(-1)^{k+1}(2k-3)!!}{2k!!}\sqrt{1-t^{2}}^{k}\Big{)}$$

$${\Psi }_1=\sum _{i,j}^{\infty }\sum _{k=0}^{\infty }\sum _{l=0}^1{C}_{ijkl}{R}^i{r}_{12}^j{t}^{2k}{\sqrt{1-t^2}}^l{e}^{-\zeta R}\frac{\alpha (1)\beta (2)-\alpha (2)\beta (1)}{\sqrt{2}}.$$

$${\Psi }_2=\sum _{i,j}^{\infty }\sum _{k=0}^{\infty }{C}_{ijk}{R}^i{r}_{12}^j{\sqrt{1-t^2}}^k{e}^{-\zeta R}\frac{\alpha (1)\beta (2)-\alpha (2)\beta (1)}{\sqrt{2}}.$$

$$\displaystyle\Psi_{3}=\Big{(}\sum_{i}^{\infty}\sum_{j=0}^{\infty}R^{i}r^{j}_{12}e^{-\zeta R}\big{(}\sum_{k=0}^{\infty}A_{i,j,k}t^{2k}(\sqrt{1-t}+\sqrt{1+t})+$$

$$\displaystyle\sum_{k=1}^{\infty}B_{i,j,k}t^{2k-1}(\sqrt{1-t}-\sqrt{1+t})\big{)}+$$

$$\displaystyle\sum_{i}^{\infty}\sum_{j=0}^{\infty}\sum_{k=0}^{\infty}R^{i}r^{j}_{12}e^{-\zeta R}\sum_{l=0}^{1}C_{ijkl}t^{2k}\sqrt{1-t^{2}}^{l}\Big{)}\frac{\alpha(1)\beta(2)-\alpha(2)\beta(1)}{\sqrt{2}}.$$

$$\mathbf{Hc}=ϵ\mathbf{Sc}$$

$$8{\pi }^2\int {\Phi }_{i,j,k,l}{r}_1{r}_2{r}_{12}d{r}_{1}d{r}_{2}d{r}_{12}=$$

$$=8{\pi }^2\int {\Phi }_{i,j,k,l}(R,t,r_{12})J_{\varphi (R,t)}{R}^2\frac{\sqrt{1-t^2}}{2}dRdtd{r}_{12}$$

$$\mathcal{I}_{1}=2\pi^{2}\int\limits_{0}^{\infty}R^{3}dR\int\limits_{0}^{1}dt\int\limits_{\big{|}R\sqrt{\frac{1-t}{2}}-R\sqrt{\frac{1+t}{2}}\big{|}}^{R\sqrt{\frac{1-t}{2}}+R\sqrt{\frac{1+t}{2}}}r_{12}f(R,t,r_{12})dr_{12}\qquad r_{2}>r_{1}$$

$$\mathcal{I}_{2}=2\pi^{2}\int\limits_{0}^{\infty}R^{3}dR\int\limits_{-1}^{0}dt\int\limits_{\big{|}R\sqrt{\frac{1-t}{2}}-R\sqrt{\frac{1+t}{2}}\big{|}}^{R\sqrt{\frac{1-t}{2}}+R\sqrt{\frac{1+t}{2}}}r_{12}f(R,t,r_{12})dr_{12}\qquad r_{1}>r_{2}$$

$$f(R,t,r_{12})=R^I{r}_{12}^J{t}^K{\sqrt{1-t^2}}^L{\sqrt{1-t}}^M{\sqrt{1+t}}^N\exp 2(-\zeta R)$$

$$\displaystyle\int\limits_{\big{|}R\sqrt{\frac{1-t}{2}}-R\sqrt{\frac{1+t}{2}}\big{|}}^{R\sqrt{\frac{1-t}{2}}+R\sqrt{\frac{1+t}{2}}}r_{12}^{J+1}dr_{12}=\left\{\begin{array}[]{lll}\frac{2}{J+2}(\frac{R}{\sqrt{2}})^{J+2}\sum_{i=1}^{J^{\prime}}{J+2\choose i}(1-t)^{\frac{J+2-i}{2}}(1+t)^{\frac{i}{2}}&\textrm{if $r_{2}>r_{1}$}\\ \\ \frac{2}{J+2}(\frac{R}{\sqrt{2}})^{J+2}\sum_{i=1}^{J^{\prime}}{J+2\choose i}(1+t)^{\frac{J+2-i}{2}}(1-t)^{\frac{i}{2}}&\textrm{if $r_{1}>r_{2}$}\\ \end{array}\right.$$

$$\mathcal{I}=\underset{0}{\overset{1}{\int }}{t}^K{\sqrt{1-t^2}}^L{\sqrt{1-t}}^M{\sqrt{1+t}}^{N}dt,$$

$$\mathcal{I}=\underset{0}{\overset{1}{\int }}\frac{(t^K{\sqrt{1-t^2}}^L{\sqrt{1-t}}^M{\sqrt{1+t}}^N}{(\sqrt{1-t}+\sqrt{1+t}{)}^O}dt.$$

$${\nu }_1$$

$${\nu }_{12}$$

$$I_p(i,\zeta )={\int }_0^{\infty }dR{R}^i(\log R{)}^p\exp -\zeta .R$$

$$\Psi =\sum _{j=0}^{\infty }{\Phi }_j(R,t)r_{12}^j$$

$$V_{nuc}=-\frac{2Z\sqrt{2}}{R}\big{(}1+\sum_{k=1}^{\infty}\frac{(4k-1)!!}{(4k)!!}t^{2k}\big{)}.$$

$$\displaystyle\Big{(}(j+2)(j+3)+2(j+2)t\frac{\partial}{\partial t}\Big{)}\Phi_{j+2}=\Big{(}H_{j}-E\Big{)}\Phi_{j}+\Phi_{j+1}$$

$$H_j=-\frac{1}{2}\frac{{\partial }^2}{\partial R^2}-\frac{5}{2R}\frac{\partial }{\partial R}-\frac{2(1-t^2)}{R^2}\frac{{\partial }^2}{\partial t^2}+\frac{6t}{R^2}\frac{\partial }{\partial t}$$

$$-\frac{j}{R}\frac{\partial }{\partial R}+\frac{2jt}{R^2}\frac{\partial }{\partial t}-\frac{Z\sqrt{2}}{R\sqrt{1-t}}-\frac{Z\sqrt{2}}{R\sqrt{1+t}}.$$

$$(2+2t\frac{\partial }{\partial t}){\Phi }_1$$

$$(6+4t\frac{\partial }{\partial t}){\Phi }_2$$

$$(1+\frac{2}{3}t\frac{\partial }{\partial t}){\Phi }_2=\frac{1}{3}(H_0-E)(1+t\frac{\partial }{\partial t}){\Phi }_1+\frac{1}{6}{\Phi }_1$$

Full text

1 Introduction

It is well known from the time of Hylleraas [1, 2], that one of the necessary conditions of a relevant convergence towards the exact nonrelativistic ground state Helium energy is the inclusion of $r_{12}$ terms into the wave function. Until now many different methods have been suggested where the $r_{12}$ function is used to construct the wave function. These methods could be divided to variational [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12], correlated-function hyperspherical-harmonic (HH) methods [13, 14, 15] (nice overview of HH methods is in the paper of Krivec [15]) and ICI method [16, 17, 18]. Generally we can say that variational methods converge to globally optimized solutions while hyperspherical-harmonic methods converge pointwise. However it is known (Bartlett et. all [19]), that wave function created from only the Slater functions and powers of $r_{1}$,$r_{2}$ and $r_{12}$ is not the exact one. It was Bartlett [20] and Fock [21] who first pointed out that exact wave function must depend on coordinate $R=\sqrt{r_{1}^{2}+r_{2}^{2}}$. Fock’s wave function depends on the powers of $R^{2}$ multiplied by powers of the logarithmic term $\ln{R^{2}}$. Due to the mathematical difficulty of direct implementation of Fock’s approach it was probably first implemented (in some approximate way) in work of Frankowski and Pekeris [5] in 1966. Their results were improved by Freund, Huxtable and Morgan [6], each group used wave functions that were some mixture of original Hylleraas coordinates and powers of the logarithmic term of the coordinate $s=r_{1}+r_{2}$.

Our effort is concentrated on finding a competitive approach with use of simple functions as coordinate system that could be generalized to more than two electron atoms. We proposed the approach where the wave function depends besides the powers of $r_{12}$ on the hyperradial coordinate $R=\sqrt{r_{1}^{2}+r_{2}^{2}}$ and on the coordinate $t=(r_{2}^{2}-r_{1}^{2})/(r_{2}^{2}+r_{1}^{2})$. The function $t$ has only one shortage, it is not well defined at the point where $r_{1}=0$ and $r_{2}=0$, so it has not derivative in this point. But in the rest area it is well defined, so it is not a significant defect. This function is more simple and easier to handle than that Fock’s one. Probably the first attempt to use these coordinates directly in variational approach was done in [22] , but authors did not use explicitly correlated functions.

As a first test of quality of the basis set all calculations were done with the same exponential scale factor $\zeta$ in $\exp{(-\zeta R)}$. An important feature of this approach is inclusion of negative powers of $R$ with combination of positive powers of $r_{12}$. Incorporation of these terms is – similar to the approach of Kinoshita $[3]$ – necessary for the wave function to be a formal solution of Schrödinger equation. As a first step the variational method has been used. We want to show that this (or similar) proposed basis sets open new possibilities in finding of the proper wave function that satisfies all cusp conditions. We believe that this knowledge will increase a chance to construct general few electron (with number of electrons more than two) atomic (or molecular) wave functions, because of a pairwise character of electron-electron and nuclear-electron interactions.

2 Hamiltonian transformation and basis set construction

The Hamiltonian for Helium atom in S basic state in coordinates $r_{1}$, $r_{2}$ and $r_{12}$ reads

[TABLE]

As it was already mentioned in the introduction, a function created from only the powers of $r_{1}$, $r_{2}$ and $r_{12}$ is not an exact eigenfunction of the Hamiltonian (1) due to the cross terms $\frac{r^{2}_{2}}{r_{1}r_{12}}\frac{\partial^{2}}{\partial r_{1}\partial r_{12}}$ and $\frac{r^{2}_{1}}{r_{2}r_{12}}\frac{\partial^{2}}{\partial r_{2}\partial r_{12}}$. Let us make the following transformation of the coordinate system:

[TABLE]

so the $r_{1}$ and $r_{2}$ coordinates can be expressed as

[TABLE]

After somehow lengthy manipulations we obtain the following transformed Hamiltonian

[TABLE]

We see that the cross terms that caused problems in (1) are now transformed to somewhat more convenient expressions, namely, the differential terms now always contain at most a single singularity factor [compared with double singularity factors $(r_{1}r_{12})^{-1}$, $(r_{2}r_{12})^{-1}$ in (1)] with all the unconvenient double singularity factors occuring now only in potential terms. It seems rather natural to include some combinations of the powers of $R$, $r_{12}$ and $t$ into the wave function. Moreover, we must carefully treat the singularities in potentials and the antisymmetry of the wave function. If we could expand the potential term

[TABLE]

in a Taylor series around $t=0$ in the whole range of $t$, the solution in the form of linear combinations of $R^{i}r_{12}^{j}t^{2k}exp(-\zeta R)$ would indeed be sufficient. Our analysis showed (see Appendix) that in this case we can eliminate all the singularities when we also include – similar to Kinoshita $[3]$ – combinations with negative exponents $i$ for $R^{i}$. Restricting the Taylor expansion to even powers of $t$ guarantees the required symmetry with respect to the particle interchange, moreover, odd powers of $t$ mutually cancel in the Taylor expansion of (7) anyhow. However, this Taylor expansion evidently diverges for $t=\pm 1$. The modified expansion

[TABLE]

already correctly reproduces the singular behaviour for $t=\pm 1$, though, the convergence of the series is slower in the vicinity of the singularities. On the other hand, the alternative expansion

[TABLE]

works well anywhere except the vicinity of $t=0$. Based on these potential expansion considerations we suggested the following form for our wave function (para case)

[TABLE]

This form of a function has already a chance to eliminate the potential term $1/\sqrt{1-t^{2}}$. Summations over $l$ in (10) are restricted to [math] and $1$ in order to avoid redundancies caused by combinations of powers of $t^{2k}$.

By linear combinations of the functions of the type (10) a formally simpler form can be obtained

[TABLE]

As (10) can be transformed into (11) it is not surprising that also for a finite truncation of infinite sums in (10), (11) these two functions give for a comparable basis sets almost the same results. In our test calculations the wave functions of the type (10) gave slightly better results, so we decided to work further with $\Psi_{1}$ rather than $\Psi_{2}$ (we will refer to this type of function as Basis A in the text).

For even better cancellation of the nuclear potential we suggested a wave function of a little more complex form,

[TABLE]

Expansion in this basis set (denoted as Basis B) has the advantage that for properly chosen coefficients $A_{ijk},B_{ijk},C_{ijkl}$ acting with (6) on (12) results in the same type of functions on both sides of the Schrödinger equation. Up till now we found nothing that would contradict the statement that (12) is a formal solution of the Schrödinger equation.

In addition, we considered also a modification of (12) with more complex exponential factor $\exp\{-\zeta R(\sqrt{1-t}+\sqrt{1+t})/\sqrt{2}\}$ instead of the simple $\exp(-\zeta R)$. This function we denoted as Basis C.

Notice that each of the suggested forms of wave function (10) - (12) is also automatically antisymmetric with respect to particle interchange.

Substitution of $\Psi$ [either in the form (10), (11) or (12)] into the Schrödinger equation results in a set of coupled equations for the coefficients $C$ (or, $A$, $B$, $C$ in the case of $\Psi_{3}$). Direct solution of these sets of equations is (similarly as was done in excelent work of Pekeris [23]) , however, very difficult, so as a first step the standard variational approach was used.

3 Variational solution of Schrödinger equation

Here we solved the standard matrix equation

[TABLE]

The elements of the matrix $\mathbf{H}$ are $\langle\Phi_{i^{\prime},j^{\prime},k^{\prime},l^{\prime}}\mid H\mid\Phi_{i,j,k,l}\rangle$ and that of the overlap matrix $\mathbf{S}$ are $\langle\Phi_{i^{\prime},j^{\prime},k^{\prime},l^{\prime}}\mid\Phi_{i,j,k,l}\rangle$. Finally the vectors $\mathbf{c}$ consist of the coefficients $C_{i,j,k,l}$ in case of basis A and of $A_{i,j,k}$, $B_{i,j,k}$ and $C_{i,j,k,l}$ in the case of bases B and C.

3.1 Calculation of integrals

When we use coordinates $r_{1}$, $r_{2}$ and $r_{12}$ volume of the integration is $8\pi^{2}r_{1}r_{2}r_{12}\\ .dr_{1}dr_{2}dr_{12}$. After a series of lengthy manipulations we succeeded to show that acting with Hamilton operator on $\Phi_{i,j,k,l}$ results in some linear combination of the type $\Phi_{i^{\prime},j^{\prime},k^{\prime},l^{\prime}}$ again. Consenquently, if we use function $\Psi_{1}$ or $\Psi_{3}$ the only type of integral which we will encounter when solving (10) will be of the type

[TABLE]

The Jacobian of the transformation $r_{1},r_{2}\to R,t$ is $J_{\phi(R,t)}=\frac{R}{2\sqrt{1-t^{2}}}$. Evidently, this transformation is regular except of the point $r_{1}=0,r_{2}=0$. Integration limits of the $r_{12}$ coordinate do not change, thus the region of the integration must be divided into two parts, namely into $t\in\langle-1,0\rangle$, (the $r_{1}>r_{2}$ case) and the region $t\in\langle 0,1\rangle$, (the $r_{2}>r_{1}$ case). The region of integration of coordinate $R$ is independent so when we considered the integrals of kind $\langle\Psi\mid\Psi\rangle$ or $\langle\Psi\mid H\mid\Psi\rangle$ (here $\Psi$ could be $\Psi_{1}$, $\Psi_{2}$ or $\Psi_{3}$) we have these two types of integrals:

[TABLE]

[TABLE]

where

[TABLE]

and in case of basis $C$ $f(R,t,r_{12})$ is the same function the only difference is that it has the exponential function $\exp\{-\zeta R(\sqrt{1-t}+\sqrt{1+t})/\sqrt{2}\}$. Index K is always positive and I,J,L,M and N integers could be also negative. We integrate through coordinate $r_{12}$ first:

[TABLE]

where $J^{\prime}$ is largest odd number from $J+2$. We can change the order of integration and integrate through coordinate $R$ so in case A or B we finally tackle the integrals of type

[TABLE]

and when we are using basis sets of type C, we have integrals of type

[TABLE]

The analytical solutions of integrals of the type (16) and (17) exists.

A software package for evaluation of these integrals in Fortran90 language was written by the author. Calculations were done in real16 (128-bit) arithmetics as the usual real8 (64-bit) arithmetics was found to be unsufficient here. Our experience showed that for more precise calculations even the use of higher precision (real24 or real32) would be advisable. Almost all calculations were done on a PC with 2.8GHz processor with 4GB RAM memory under Linux operating system. The calculation of the integrals was made more effective by pre-calculating of the integrals (16) or (17), so for example for the A-type basis set with $I_{max}=5$ calculation of all integrals of type (14) on this PC is a question of few seconds.

The most time consuming part of the calculation remained the diagonalization of the overlap and transformed Hamiltonian matrices.

3.2 Details and results

The wave functions of the type $\Psi_{1}$, $\Psi_{2}$ and $\Psi_{3}$ (10-12) were tested in our calculations. In practice, the infinite summations over $i,j,k$ in (10-12) have to be obviously truncated to finite values $I_{m},J_{m},K_{m}$, respectively. Unfortunately, for larger powers of $k$ our basis sets are close to linear dependency. We used standard canonical transformation of the basis set with elimination of small eigenvalues below a defined threshold to overcome this problem. In all cases, the choice of the threshold was governed by the criterion to be the minimal threshold where the method still converges. Exact nuclear potential $-Z/r_{i}$ $(i=1,2)$ was used in our variational calculations. For each $\Psi_{i}$ $(i=1,2,3)$ such exponents $\zeta$ were chosen that $E=E(\zeta)=\langle\Psi_{i}(\zeta)|H|\Psi_{i}(\zeta)\rangle$ is minimal with respect to $\zeta$. To find the minimum of $E(\zeta)$ the golden section technique was used. Results of our calculations were compared with reference values for nonrelativistic energy of the ground state of He atom and H*-* ion according to [11].

We used Basis A functions [$\Psi_{1}$ and $\Psi_{2}$, Eqs. (10-11)] first. It turned out that there are almost no differences between using basis sets $\Psi_{1}$ and $\Psi_{2}$ of equivalent sizes. This is not surprising as $\Psi_{1}$ can be in principle generated from $\Psi_{2}$ and vice versa when untruncated infinite summations are considered. Our numerical experience showed that the truncated summations with $\Psi_{1}$ produced slightly better results than with $\Psi_{2}$. Thus, all our results presented with Basis A refer to using $\Psi_{1}$. These are summarised in Tab.1. Each combination of $I_{m},J_{m},K_{m}$ determines the total number of basis functions, $N_{0}$. By $N$ the number of basis set functions after elimination of quasi-linear dependencies is denoted. This represents then the actual basis set size used in the calculation. However, as $N$ is unknown before we start the calculation, $N_{0}$ was used to label our basis sets (preceded by basis type label A,B, or C). The value of optimal exponent $\zeta$ for each basis set is presented as well.

As emphasised in the Appendix, negative powers of $R$ [in combination with positive powers of $r_{12}$ according to (A-15) i.e., $r_{12}^{j}/R^{i^{\prime}}$ up to $i^{\prime}=j$] should be included to eliminate the singularities of the type $R^{-2}$. Unless explicitly specified otherwise, all the basis sets used in our calculations were constructed in this way (both for A and B/C case). To see what happens when we do not include negative powers of $R$ in the wave function a few illustrative calculations with only positive powers of $R$ are included. For basis sets presented in Tab.1 this is the case of the set A438a. Almost no effect on the ground state energy can be seen for this basis set. The effect will be more visible for the basis sets B and C, especially for the cusp parameters.

To test the stability of the method the hydrogen anion energy was considered as well. Evidently, the H*-* anion is a more sensitive system. Here we found to be necessary to use the strict threshold of $10^{-7}$ for elimination of linear dependencies. (For He a threshold of $10^{-12}$ was found to be sufficient for A-type basis sets and even $10^{-16}$ for basis sets B and C.) The difference in minimal threshold for He and H*-* anion explains why our energies for He are closer to the reference value than those for the H*-* anion.

From Tab.1 we see that even with the simple basis set functions of the type (10) the He and H*-* ground state energy with the accuracy of $10^{-8}$ can be achieved with a moderate basis set size.

Evidently, $\Psi_{1}$ and $\Psi_{2}$ can be the exact solutions only for the approximate potential (8) or (9). Therefore, basis sets of the type B or C which satisfy the exact Coulomb potential are expected to produce better results. In Tab.2 results on a hierarchy of these basis sets is presented. When comparing the energies from Tab.1 and Tab.2 we see that the basis sets of the type B and C produce indeed better energies than the A set case. The best result was obtained with basis set C1068. For comparison we present also results with two modified non-standard way generated sets. The set C944a contains only positive powers of $i$ with $I_{m}=7$ ($I_{m}=6$ in all other basis sets with $N_{0}>900$). The basis set C950b contained negative powers of $i$ with $I_{min}=-J_{m}+2$ which corresponds to the constraint that was put on the function in Appendix (A-10).

If we compare the energies corresponding to basis sets C944a and C950b (which are of almost the same size) the energy of C950b is by one order better than the basis set C944a which confirms the importance of inclusion of negative powers of $R$. As the practice shows incorporation of negative powers of $R$ in combination $r_{12}^{j}/R^{i}$ is at least as important as the inclusion of higher positive powers of $R$.

To assess the behaviour of our wave functions in the vicinity of the electron-nuclear and electron-electron coalescence point the cusp parameters $\nu_{1}$ and $\nu_{12}$ were evaluated,

[TABLE]

where $\langle{\bf O}({\bf r}_{1},{\bf r}_{2})\rangle=\int d{\bf r}_{1}d{\bf r}_{2}\,\Psi({\bf r}_{1},{\bf r}_{2})^{*}\,{\bf O}({\bf r}_{1},{\bf r}_{2})\,\Psi({\bf r}_{1},{\bf r}_{2})\,$ and $\delta(\mathbf{r})$ is Dirac delta function. Result are summarised in Tab.3. As expected, the best agreement was obtained for the C-type basis set C1068. The agreement for A-type basis sets was for about two orders worse. Though it is not always justified to compare the A,B and C basis sets directly, these results indicate that better wave functions can be obtained with B and C sets.

4 Discussion

The basis sets could be furthet optimised by omitting the contributions with negligible effect on the result. Therefore, it would probably be possible to find different basis sets with the same total number of functions that could result in better energy. We did not follow this way. Let us notice that almost for each basis for small values of $\zeta$ a ”grey area” exists where the method does not converge due to linear dependency problems. This is more visible for H*-* anion.

It is clear, that functions of type $t^{2k}$ and $t^{2k+2}$ or $t^{2k+1}$ and $t^{2k+3}$, for large $k$ tend to be similar within the interval $\langle-1,1\rangle$. One can easily check that the set $t^{2k},k=0,1,\dots,K_{m}$ is already for $K_{m}\approx 10$ almost linearly dependent on the interval $\langle-1,1\rangle$. Evidently, this limits the basis set size in practice and, consequently, also the precision of our energy (or any other quantity of interest) calculation. From this point of view, basis sets B and C are again superior to A as using the same $K_{m}$ higher precision can be obtained with B or C than with A type basis set (compare Tab.1. and Tab.2.). A remedy how to increase the precision for B or C sets keeping $K_{m}$ the same would be to include more exponents $\zeta_{n}$: $\exp{-\zeta_{n}R}$ or $\exp{-\zeta_{n}.(r_{1}+r_{2})}$. Another possibility would be to use Fock logarithmic terms $R^{i}\log{R}$ to better describe the triple coalescent point of the wave function (i.e., $r_{1}=r_{2}=0$). The large advantage of basis sets of type B is the possibility of exact calculation of integrals of the type

[TABLE]

contrary to the methods that use Hylleraas coordinates. In this case, for the basis of the C type an analytical solution for the integrals over $t$ probably does not exist and the integrals have to be evaluated numerically. An alternative would be a combination of the method with more exponents with Fock logarithmic terms. We expect that each of the aformentioned procedures leads to improvement of basis set linear dependency problems and a higher precision of calculations as a result.

5 Conclusions

In this work a new ansatz for the Helium-like ion wave function was proposed. To our knowledge it is at present the best energy for He atom obtained with the use of hyperradial coordinates. Our preliminary results show that the wave functions of the type B and C are close to the exact solution of the Schrödinger equation. The remarkable property of the functions B and C is that for properly chosen expansion coefficients functions of the same type occur on both sides of the Schrödinger equation. In [26] it was shown that a simple exponential $\exp{-2\zeta R}$ describes the behaviour of the wave functions fairly well in limit cases $r_{1}=0$, $r_{2}=0$ or $r_{12}=0$. We showed that the functions A or B (containing just the aforementioned exponential function) can offer relatively precise ground state energies. However, the calculations of cusp parameters indicate that the use of exponentials $\exp{-2\zeta(r_{1}+r_{2})}$ should be prefered. The important point turnes out to be the inclusion of negative powers of $R$ in combination with positive powers of $r_{12}$ similar to Kinoshita [3] where negative powers of $s=r_{1}+r_{2}$ were included. Incorporation of negative powers of $R$ is at least as important as the use of larger number of positive powers of $R$. For our method the inclusion of the full correction to finite mass of the nuclea (see e.g.[27] p.263) is straightforward, the work is in progress. Our method in its present state could not compete with ICI method of Nakatsuji [16, 17, 18], but our preliminary results urge us to clousure that the wave function of type (12) is close to exact one and one of proposed ways of improvement of the wave function could produce more exact functions. The knowledge of proper functional behavior of three-body wave function could be interesting also for the ICI method.

Finally, there was done some work on generalization of this ansatz to systems with more than two electrons, we believe that we are on the right way to overcome the problems with multielectron integrals that generally depend on $r_{ij}$.

6 Appendix

In this appendix we will show that inclusion of terms $R^{-i}$ (up to $i_{max}=j$) is a necessary condition for a function

[TABLE]

to be (a spinless part of) a correct solution of our problem. We will prove this for a case where we can assume that $V_{nuc}$ (8) can be expanded in a Taylor series around $t=0$,

[TABLE]

Let us act with the Hamilton operator (6) on a function (A-1). We obtain a series of coupled equations

[TABLE]

where the operator $H_{j}$ is

[TABLE]

Setting $\Phi_{j}=0$ for all $j<0$ the first two equations of the set (A-3) read

[TABLE]

Substituting for $\Phi_{0}$ from (A-5) into (A-6) we get

[TABLE]

from which we can get a formal solution

[TABLE]

Generally, a recursive formula for a solution of (A-3) can be found,

[TABLE]

indicating that, using (A-8) and (A-9) all $\Phi_{j}$’s can be, in principle, obtained.

A natural requirement imposed on $\Phi_{j}$’s is that $\Phi_{j}r_{12}^{j}$ [see(1)] are expected to be all finite and decaying to zero with increasing $j$.

As we here consider the situation where the expansion (A-2) can be used we can suggest $\Phi_{j}$’s to be of the form

[TABLE]

Acting with the Hamilton operator on (A-10) can, however, produce singularities of the type $R^{-1}$, $R^{-2}$. Evidently, the suggested form (A-10) must be able to eliminate all these singularities. Let us inspect this more closely.

Substituting (A-10) into (A-8) yields

[TABLE]

For elimination of all singularities of the type $R^{-2}$ from $\Phi_{2}$ it is evidently sufficient to put

[TABLE]

To eliminate the $R^{-1}$ singularities from $\Phi_{2}$ let us set in (A-11) the terms standing at $R^{-1}$ and at all $t^{2k}$ equal to zero. Making use of (A-12) this results in restriction for expansion coefficients (assuming $\Phi_{j}$ is normalized by setting $C_{0,1,0}=1$)

[TABLE]

Thus $\Phi_{1}$ ,$\Phi_{2}$ and of course $\Phi_{0}$ defined by (A-5) can be correctly introduced without any need to incorporate negative powers of $R$. Using the convergence criterion for power series it can be shown that the function $f(t)=\sum_{k=0}^{\infty}C_{1,1,k}t^{2k}$ defined by the sequence (LABEL:App_C11k) converges uniformly to a finite function on its whole interval of definition ($t\in(-1,1)$).

Let us analyse $\Phi_{j}$ for the case $j>2$ now. First, let us notice that if the lowest power of $R$ in $\Phi_{j}$ is $R^{i}$ ($i\neq 0$) then the effect of acting $H_{j}$ on $\Phi_{j}$ (for the case when (A-12) holds) is lowering the lowest power or $R$ to $R^{i-1}$. In accord with this, as $\Phi_{1}$, $\Phi_{2}$ we alredy managed to construct singularity-free, the only singularity in $\Phi_{3}$ (of the type $R^{-1}$) can come from acting $H_{1}$ on $\Phi_{1}$. The only possibility to eliminate this singularity is now to allow the same singularity $R^{-1}$ to occur in $\Phi_{3}$. $\Phi_{4}$ must contain the $R^{-1}$ term again (as both $\Phi_{3}$ and $H_{2}\Phi_{2}$ contain it) and also the $R^{-2}$ due to the fact, that now a non-zero term with $R^{0}$ also occurs.

Generally, notice that while the effect of acting of (A-4) on a term containing no $R$-type singularities ($\Phi_{1}$, $\Phi_{2}$) was the occurence of $R^{-1}$ singularity only, each acting of (A-4) on $R^{-i},i>0$ produces also an $R^{-(i+2)}$ singular term. Consequently, for each $j>2$, if $\Phi_{j}$ contains the $R^{-i}$ term, $\Phi_{j+2}$ must contain $R^{-(i+2)}$ [plus the singularity coming from $\Phi_{j+1}$ which can, however, never be worse than $R^{-(i+2)}$]. Starting from the fact that $\Phi_{3}$ contains $R^{-1}$, $\Phi_{4}$ contains $R^{-1}$ and $R^{-2}$, etc., it is now easy to show recursively that $\Phi_{j+2}$, $j=1,2,\dots$ must contain the $R^{-j}$ term, that is, terms $R^{-i}$ must be incorporated into (A-10) with $i_{max}=j-2$, for $j>2$. In this way, all the $R^{-1}$, $R^{-2}$ singularities have a chance to be eliminated.

Up till now we only showed, that the $R^{-i}$ singularities can be eliminated from the set of equations (A-3). By construction, the singularities are evidently still present in the functions $\Phi_{j}$ and, apparently, also in the wave function (A-1). One has, however, to realise that ${r_{12}}/{R}$ is always finite [see (2)-(5)], as

[TABLE]

For our construction of $\Phi$’s the worst combination of $R$ and $r_{12}$ will evidently be ${r_{12}^{j}}/{R^{j-2}}$, which is, according to (A-15) always finite, as well.

Let us analyze what happens when the infinite series (A-1) is truncated (which is always the case in practice). As an example let us consider the case when in (A-10) $i_{max}=1$ a $j$ is finite In the following we will analyze the most problematic series resulting from (A-9). Let us introduce

[TABLE]

where $J=(j+2)/2$ and $j_{0}=3$ for even $j$ and $J=(j+1)/2$ and $j_{0}=4$ for odd $j$. Acting repeatedly with the $-\frac{\partial^{2}}{\partial R^{2}}$ part of $H_{j}$ on $\Phi_{j}$ yields

[TABLE]

The dominant series resulting from recursive acting of $-\frac{2(1-t^{2})}{R^{2}}\frac{\partial^{2}}{\partial t^{2}}+\frac{6t}{R^{2}}\frac{\partial}{\partial t}$ is

[TABLE]

Finally, acting of $\frac{2jt}{R^{2}}\frac{\partial}{\partial t}$ on $\Phi_{j}$ returns

[TABLE]

All the other parts of $H_{j}\Phi_{j}$ produce series which are minorant to one of series (A-17-A-19). Thus, they need not be analysed. Let us consider the convergence of (LABEL:App_C11k) and the restriction (A-15) now. Evidently, the series (A-17) converges uniformly when $\zeta^{2}/{2}+E<1$. This statement was confirmed also by numerical results. The series (A-18) uniformly converges for $-1<t<1$. Similarly, the series (A-19) converges for finite $j$. As all the other series are minorant to some of the cases (A-17-A-19), also the complete sum of series converges to a finite function for a finite $j$.

7 Acknowledgement

The author is deeply indebted to Pavel Neogrády and Ján Mašek for many inspirating discussions and to Ilja Martišovits for discussions of mathematical aspect of the problems. Special thanks are due to Štefan Varga for help with calculation of some integrals and for help with preparing of the manuscript and also to Štefan Dobiš for technical support. This work is dedicated to the memory of Ladislav Turi-Nagy the former head of Department of Theoretical Chemistry at the Institute of Inorganic Chemistry of Slovak Academy of Sciences, a man who encouraged me very much in my work.