Chemical Physics Letters 610–611 (2014) 167–172
Contents lists available at ScienceDirect
Chemical Physics Letters journal homepage: www.elsevier.com/locate/cplett
Orbitals of the dipositronium A.J.C. Varandas a,∗ , M. Brajczewska b , J. da Providência b , J.P. da Providência c a
Departamento de Química, and Centro de Química, Universidade de Coimbra, 3004-535 Coimbra, Portugal CFC, Departamento de Física, Universidade de Coimbra, 3004-516 Coimbra, Portugal c Departamento de Física, Universidade da Beira Interior, P-6201-001 Covilhã, Portugal b
a r t i c l e
i n f o
a b s t r a c t
Article history: Received 5 June 2014 In final form 4 July 2014 Available online 14 July 2014
The stability of dipositronium is discussed. Our purpose is not to compete with highly accurate calculations that are available but point out that its dissociation energy can be understood on the basis of simplified models if satisfying the important requirement of consistency in the description of the ‘atoms’ and molecules of positronium. The importance of correlations in positronium and dipositronium is analyzed. Also studied is the replacement of positrons by heavier Fermions in positronium and dipositronium. It is found that the Born–Oppenheimer approximation is valid provided the mass of the heavier Fermions exceeds 200 times the mass of the electron. © 2014 Elsevier B.V. All rights reserved.
1. Introduction The existence of diatomic positronium was predicted by Wheeler in 1946 [1]. In the following year, the stability of this system was confirmed by a variational calculation [2]. A very accurate upperbound for the dipositronium binding energy (−0.5159767 Eh ) is already available [3], with similarly accurate calculations being available for other related species such as positronium hydride [4]. It is not our aim to compete with such very precise although highly technical numerical simulations. Instead, we point out that the stability of the dipositronium molecule can be understood on the basis of simple models satisfying the important requirement of consistently describing ‘atoms’ and molecules of positronium. By consistent description, which is an important ingredient in our argument, we mean that the model wave function of the dipositronium should contain as factors the model wave functions of the positronium atoms, at least for some particular value of parameters specifying correlations between electrons and positrons belonging to distinct positronium atoms. 2. One-center model The dipositronium is a system of two electrons and two positrons. In atomic units, its Hamiltonian is given by H = T + V,
T=
4 p2 i
i=1
2
,
V=
4 ei ej i
|ri − rj |
∗ Corresponding author. E-mail address:
[email protected] (A.J.C. Varandas). http://dx.doi.org/10.1016/j.cplett.2014.07.015 0009-2614/© 2014 Elsevier B.V. All rights reserved.
where e1 = e3 = −1, e2 = e4 = 1. Let us assume, in the framework of a variational approach, that positrons and electrons have a common orbital that is described by the un-normalized wave function 2
(r) = e−(˛/2)r , the center of the orbital being the center of mass of the four particles. According to the Pauli Principle, the spin states must be antisymmetric for both species. The space part of the wave function is of the form (r1 , r2 , r3 , r4 ) = K(r1 )(r2 )(r3 )(r4 ), where K is a normalization constant. The expectation value of the Hamiltonian reads
E = |H| = 3˛ − 2
2˛ ,
with atomic units being used throughout, unless specified otherwise. This is twice the expectation value of the positronium Hamiltonian with a wave function of the form (r1 , r2 ) = K(r1 )(r2 ). It then suggests that the dipositronium is better described in terms of a two-center model. 3. Two-center model
,
(1)
Consider now a two-center model. Let us assume that the wave function of dipositronium is the product of two positronium wave functions with separate centers. The un-normalized wave function of one of the positrons and electrons is 2
1 (r) = e−(˛/2)r .
(2)
168
A.J.C. Varandas et al. / Chemical Physics Letters 610–611 (2014) 167–172
The un-normalized wave function of the other positron and electron is 2
2 (r) = e−(˛/2)(r−q) ,
(3)
where q is the vector between the centers of the two positronium atoms, and q = |q| is its modulus. The Hamiltonian is given by Eq. (1). We denote by ri , and i (i = 1, 2, 3, 4), respectively, the space coordinates and spin quantum numbers of the particles. The wave function must be antisymmetric under the separate interchange of the electron coordinates (r1 , 1 ) ↔ (r3 , 3 ) and positron coordinates (r2 , 2 ) ↔ (r4 , 4 ). If it is antisymmetric under the interchange 1 ↔ 3 , it is symmetric under the interchange r1 ↔ r3 , and viceversa. The same holds for the positron coordinates (r2 , 2 ), (r4 , 4 ). The total spin of electrons and positrons is 0 (singlet state) or 1 (triplet state) according to whether the wave function is antisymmetric or symmetric with respect to the interchange of the respective spin quantum numbers. We assume that electrons and positrons are in singlet spin states, which is energetically the most favorable situation. The space factor of the positronium wave function is therefore obtained by symmetrizing the wave function DD (r1 , r2 , r3 , r4 ) = 1 (r1 )1 (r2 ) × 2 (r3 )2 (r4 )
(4)
with respect to the electron coordinates r1 , r3 and with respect to the positron coordinates r2 , r4 . In turn, we denote by ED the wave function obtained from DD by interchanging r1 , r3 , and by DE the wave function obtained from DD by interchanging r2 , r4 . Similarly, we denote by EE the wave function obtained from DD by interchanging r1 , r3 and r2 , r4 . Explicitly, we write: ED (r1 , r2 , r3 , r4 ) = DD (r3 , r2 , r1 , r4 ), DE (r1 , r2 , r3 , r4 ) = DD (r1 , r4 , r3 , r2 ), EE (r1 , r2 , r3 , r4 ) = DD (r3 , r4 , r1 , r2 ). Denoting now by i ( j ), i = ↑ , ↓ , j = 1, 2, 3, 4 the Fermion spin wave function, the spin factor of the positronium wave function assumes the form (↑ (1 )↓ (3 ) − ↑ (3 )↓ (1 )) × (↑ (2 )↓ (4 ) − ↑ (4 )↓ (2 )). We find, for ˛ = 1, NDD =
4
3
d rj 2DD = 6 ,
4 j=1
NDE =
2 /2)
d rj DD DE = 6 e−(q 3
,
j=1 4
NEE =
2
d rj DD EE = 6 e−q , 3
4 j=1
KDD =
4
3
d rj
4 j=1
KDE =
i=1
3
4
d rj j=1 4
KEE =
1 6 2
(∇ i DD ) · (∇ i EE ) =
1 6 2 (6 − q2 )e−q , 2
6−
1 2 q 2
(∇ i DD ) · (∇ i DE ) =
2 /2)
e−(q
,
i=1 4
3
d rj
4 j=1
VDD =
(∇ i DD ) · (∇ i DD ) = 36 ,
i=1
√ d rj 2DD V = −2 211/2 ,
4
√ 2 3 d rj DD EE V = −2 211/2 e−q ,
4 j=1
VDE =
K=
√
˛
KDD + 2KDE + KEE . NDD + 2NDE + NEE
The potential energy is given by V=
VDD + 2VDE + VEE . NDD + 2NDE + NEE
The corresponding expressions for the case in which electrons and positrons are not simultaneously in singlet spin states are presented in Appendix A. The dipositronium energy is obtained by minimizing E = K + V with respect to ˛ and q . It is convenient √ to replace the parameters ˛, q by ˛ and ˛ q. The minimization √ with respect to ˛, for fixed ˛ q, is immediate, and leads to −(V2 /(4K))˛=1 . The minimization of −(V2 /(4K))˛=1 with respect √ to ˛ q is straightforward. We find, Emin = −0.227952 Eh For the binding energy of the dipositronium molecule with respect to two positronium atoms we obtain −0.0157454 Eh , which is a quite remarkable result. Naively, we might expect that the wave function (4) describes the ground-state of the dipositronium molecule. This is not true, because the ground-state of the dipositronium molecule has spherical symmetry and the wave function in Eq. (4) does not. Although there is no reason why we should restrict the description of the dipositronium molecule to its ground-state wave function, it is certainly instructive to look into that question, which is considered in the next section. The wave functions 1 (r1 )1 (r2 ) and 2 (r3 )2 (r4 ), as defined in Eqs. (2) and (3), describe atoms of positronium with localized centers of mass. Localizing the centers of mass is unavoidable if we wish to bring the atoms together to produce a molecule. As a consequence, the molecule which is produced will not be in its true ground-state, but in a non-stationary state, which is described by Eq. (4), and there is nothing wrong with that. As a consequence, it is quite probable that the molecule is produced, not exactly in its true groundstate, but in a non-stationary state, such as the one which is described by Eq. (4). This argument suggests that the number −0.0157454 Eh , although remarkable, is not completely accidental. It should be emphasized that the above does not reflect in any way a limitation of the variational principle. What we just wish to express is that the wave function associated with a specific molecular formula, e.g., the hydrogen molecule, conveys information that goes beyond the ground-state properties, covering in the case of the given example, the rotational spectrum. Indeed, chemistry is not restricted to ground-state properties. For example, the groundstate of levorotary and dextrorotary compounds is neither one nor the other. Such properties are an outcome of the production process. Translated to the case under consideration, it is remarkable that the wave functions which are used to approximate the groundstate favor the tendency for asymmetric structures related to a specific separation between the positronium atoms, although the exact ground-state is, of course, symmetric.
3
j=1
VEE =
potential and kinetic energies are multiplied by extra factors, and ˛, respectively. The kinetic energy is given by
3
d rj DD DE V
j=1 2 /2)
= 11/2 e−(q
√
√ √ q q 8 Erf √ Erf 2+ − , √ q q 2 2 2
where V has been given in Eq. (1). For ˛ = / 1, the correspond√ ing results are easily obtained, if q is replaced by ˛ q and the
4. Correlations Correlations play a crucial role in the description of the groundstate dipositronium. In order to take them into account, we assume that the space part of the un-symmetrized dipositronium wave function is of the form
DD (r1 , r2 , r3 , r4 ) = exp
ˇ [(r1 − r2 )2 + (r3 − r4 )2 + ((r1 + r2 )2 2
+(r3 + r4 − 2q)2 )]
,
(5)
A.J.C. Varandas et al. / Chemical Physics Letters 610–611 (2014) 167–172
where ˇ, , q are variational parameters. This wave function reduces to Eq. (4) for = 1, ˇ = ˛/2 (keeping in mind Eq. (2), and Eq. (3)). Define now:
169
3.0
2.5
ED (r1 , r2 , r3 , r4 ) = DD (r3 , r2 , r1 , r4 ), DE (r1 , r2 , r3 , r4 ) = DD (r1 , r4 , r3 , r2 ),
2.0
EE (r1 , r2 , r3 , r4 ) = DD (r3 , r4 , r1 , r2 ). The following matrix elements are then required: NDD = DD |DD ,
NED = DD |ED ,
NDE = DD |DE ,
NEE = DD |EE ,
1.5
1.0
TDD = DD |T |DD ,
TED = DD |T |ED ,
TDE = DD |T |DE ,
TEE = DD |T |EE ,
VDD = DD |V |DD ,
VED = DD |V |ED ,
VDE = DD |V |DE ,
VEE = DD |V |EE .
0.5
0.0 0.0
Obviously, NDE = NED , TDE = TED , VDE = VED . We obtain, for ˇ = 1, NDD
2
e−2q 6 , 64 3
TDD
2q2 e 1 + 6
8 3/2 (1 + ) 3(1 + )6 = , 64 3
TEE =
e
−2q2
VDD
−
, 3
6 (6 + 6 − 4 2 q2 ) , 128 3 2 2 q2 e 1 + 6
−
TED = TDE =
0.6
0.8
ˇq. Horizontal axis, . Vertical axis,
16 3/2 (1 + )3
3(1 + 6 + 2 ) 8 2 q2 − 1+ (1 + )2
,
11/2 =− , 16 3
1.0
ˇq.
Thus, the minimum of the energy E(ˇ, , q), for fixed , with rspect to ˇ, is
−
NED = NDE =
0.4
Figure 1. Contour plot of the minimum of E(ˇ, , q), with respect to ˇ, for fixed
6 = , 64 3
NEE =
0.2
(V(,
4K(,
ˇq))
ˇq and
2
ˇq)
.
The numerical value of this quantity is displayed in Figure 1. It shows that the absolute minimum, −0.424413 Eh is obtained for = 0 and is independent of the value of q. We observe that the absolute minimum, −0.424413 Eh is twice the energy of the positronium, namely −0.212207 Eh . If, consistently with the dipositronium wave-function in Eq. (5), the positronium is described by the function
(r1 , r2 ) = exp
ˇ [(r1 − r2 )2 + (r1 + r2 )2 ] 2
.
(6)
2
VEE =
e−2q 11/2 , 16 3 2q2 1 + e −
VED = VDE =
3+
16 5/2 (1 + )5/2
3 + q
√ 2 211/2 q + 1 + 6 Erf
√ √ 2q 2q 6 − 8 1 + Erf 1+
1 + 4 + 3 2
.
The corresponding expressions for ˇ = / 1 are easily obtained from the previous ones, by an obvious choice of variables. The energy may be written as E(ˇ, , q) = ˇK(,
ˇq) +
where K(, q) =
TDD + 2TED + TEE , NDD + 2NED + NEE
V(, q) =
VDD + 2VED + VEE , NDD + 2NED + NEE
ˇV(,
ˇq),
The wave function in Eq. (5) improves the description of the dipositronium ground-state since −0.424413 Eh is much closer to the theoretical value of −0.5159767 Eh (cf. [3]) than the value −0.227952 Eh which is associated with the wave function in Eq. (4). However, this does not necessarily mean that the wave function in Eq. (5) is more appropriate to describe the formation of the dipositronium molecule starting from the positronium atoms. In fact, the wave function in Eq. (5), which describes a nonstationary state by representing a superposition of a huge number of eigenstates, may under such a circumstance be advantageous.
170
A.J.C. Varandas et al. / Chemical Physics Letters 610–611 (2014) 167–172
3.0
exchanged wave function ED (r1 , r2 , r3 , r4 ) = DD (r3 , r2 , r1 , r4 ) . We find,
2.5
ND =
4
V(13)D =
2.0
V(12)D =
1.5 NE =
0.0
0.2
0.4
0.6
0.8
fixed
ˇq. Horizontal axis, . Vertical axis,
2
3
d ri
i=1 4
3
d ri i=1 4
KE =
(DD ) 11/2 = 3/2 |r1 − r2 | 2 (4 − ˇ )(−ˇ2 + 4 2 ) 6
i=1 4
1 2
3
d ri
3
1 2
i=1
2(3 + ) (ˇ + 2)11/2 DD ED = 5/2 |r1 − r2 | (2 + ˇ)((1 − ˇ + )(3 + )(ˇ + 2))
4
3(1 + )6
2
(∇ j DD ) =
((−4 + ˇ2 )(ˇ2 − 4 2 ))
(∇ j DD )(∇ j ED )
3
2
4(1 − ˇ + ) ((2 + ˇ)(1 − ˇ + ) (ˇ + 2))
The value of ˇ which minimizes E(ˇ, , q) for fixed
ˇ=
(V(,
4(K(,
ˇq))
ˇq is given by
2
ˇq))
2
3/2
.
The norm is N = ND + NE . The expectation value of the potential energy is V = VD + VE , where VD = −2V(12),D and VE = 2V(13),E − 4V(12),E . The expectation value of the kinetic energy is K = KD + KE . The minimum of the binding energy for fixed , is
. W=−
and is displayed in Figure 2. For = 0, we find ˇ = 4/(9) = 0.141471 . The drawback of the wave function in Eq. (5) is that it predicts no binding of the dipositronium relative to two positronium atoms, not only with respect to energy but also to spatial separation. In the next section a new variational ansatz is presented which overcomes this deficiency.
5. A new variational ansatz The wave function in Eq. (5) does not account for the binding energy of the dipositronium with respect to the positronium atoms, and is not sensitive to the distance between these atoms. It is therefore necessary to replace the wave function in Eq. (5) by a more flexible one which is still compatible with the positronium description. In this section, we investigate a trial wave function of the form DD (r1 , r2 , r3 , r4 ) = exp(−TDD ),
(7)
where TDD
3/2
j=1 2
3(1 − ˇ + ) (1 + 6 + 2 − 2ˇ(1 + ))6
=
3/2
√ 2(2 + ˇ)(ˇ + 2)11/2 DD ED = 5/2 |r1 − r3 | ((2 + ˇ)(1 − ˇ + )(ˇ + 2))
j=1 4
d ri
ˇq.
2
((2 + ˇ)(1 − ˇ + ) (ˇ + 2))
2
3
d ri
i=1 4
1.0
Figure 2. Contour plot of ˇ, such that E(ˇ, , q) is minimized with respect to ˇ, for
i=1 4
3/2
√ 2 2(1 − ˇ + )11/2 (DD ) = 3/2 |r1 − r3 | (−2 + ˇ)(ˇ − 2)((2 + ˇ)(1 − ˇ + )(ˇ + 2))
3
V(12)E =
0.0
3
d ri
((−4 + ˇ2 )(ˇ2 − 4 2 ))
d ri D ED =
V(13)E =
KD =
i=1 4
i=1 4
1.0
0.5
6
2
3
d ri (DD ) =
ˇ = {(r1 − r2 )2 + (r3 − r4 )2 + [(r1 + r2 )2 + (r3 + r4 )2 ] 2 + 2(r1 · r3 + r2 · r4 )},
with 2 ≥ ≥ −2. For = 0, this wave function reduces to the one in Eq. (5) for q = 0 . The wave function in Eq. (7) is still consistent with the one in Eq. (6). To compute the energy, we also need the
V2 . 4NK
The absolute minimum −0.447608 Eh occurs for = 0.0886, = −2. The energy difference between two positronium atoms and a dipositronium molecule is now −0.023195 Eh . As shown, the wave function in Eq. (7) slightly overestimates the dipositronium binding energy. Conversely, the wave function in Eq. (5) underestimates it. The wave function in Eq. (4) is the one which works better. It is also the wave-function in Eq. (4) which better satisfies the criterion of consistency between the positronium and dipositronium descriptions. According to the model in this section, the dipositronium radius is 1.69329 a0 which is too small a value [5]. Notice that the positronium radius is 1.67927 ˛0 on the basis of a consistent wave function. The choice made in Eq. (7) was only done for simplicity reasons, although it is clearly preferable to consider TDD =
ˇ 2
(r1 − r2 )2 + (r3 − r4 )2 + [(r1 + r2 )2 + (r3 + r4 − q)2 ]
+2
r1 −
q 2
·
r3 −
q 2
+ r2 −
q 2
·
r4 −
q 2
,
keeping q as an extra variational parameter. Since the parameters , must satisfy 2 ≥ ≥ −2, an important consequence of includ/ 0, ing the parameter is that the minimum energy occurs for = which renders the wavefunction sensitive to the value of |q|. On the basis of the model in Section 3, we expect the separation between the centers of the positronium atoms to be q ∼ 3.631523 a0 . It is encouraging that this may lead to a size of the dipositronium molecule in agreement with the results in [5] (Figure 3).
A.J.C. Varandas et al. / Chemical Physics Letters 610–611 (2014) 167–172
171
-0.1 -0.415 -0.15 -0.42 -0.2 -0.425 -0.43
-0.25
-0.435
-0.3 0.05
0.1
0.15
-0.35
0.2
-0.445 1 Figure 3. Dipositronium energy as a function of for = −2 (Section 5).
√ 1 = √ 2
− [r21 + r22 + (r3 − q)2 + (r4 − q)2 ] 2
q √ 2
2 /2
+ (0)e−q
− 2 (0) +
5
6
7
7. Born–Oppenheimer approximation
It is challenging to investigate the replacement the positrons by Fermions with the same charge but different mass. For instance, the positrons might be replaced by anti-muons, anti-tauons or protons. It is very important to use orbitals with different widths for Fermions with different masses. We have considered exp(− ˛r2 ) for electrons and exp(− ˛ r2 ) for positively charged Fermions, treating both ˛ and as variational parameters. Center of mass corrections (CMC) were not taken into account, so that the descriptions of positronium and dipositronium are consistent. Explicitly, the dipositronium wave function under consideration is written as
˛
4
3
Figure 4. Ground-state energy of the system electron-anti-Fermion as a function of log10 (M), where M is the anti-Fermion mass in atomic units. With CM corrections, full line. Without CM corrections, dashed line.
6. Replacement of the positrons by Fermions with the same charge but different mass
DD (r1 , r2 , r3 , r4 ) = exp
2
Consider the protons fixed, one at the origin and the other at the point q, and assume that the electrons are in a singlet state. In coordinate space the wave function is symmetric, 12 (r1 , r2 ) = 1 (r1 )2 (r2 ) + 1 (r2 )2 (r1 ), with 1 (r) = exp(− ˛r2 ), 2 (r) = exp [− ˛(r − q)2 ]. We obtain, for ˛ = 1, 2 /2
12 |12 = 3 (1 + e−q 12 | −
. 12 |
q 2
√ The optimizing satisfies = M. The optimizing values of the parameters ˛, , q are such that, for different M, the product √ ˛q ranges between 1.3 and 1.5 a0 , showing that in the molecular systems the atoms retain a remarkable individuality. Correlation effects, which are expected to become less important as the mass of positively charged Fermions increase, were not considered. The importance of correlation effects for the atomic system (which are accounted for by CMC) are illustrated in Table 1.
2 /2
e−q
),
1 2 1 2 ∇ − ∇ |12 = 3 2 1 2 2
3 2
+
3 2
−
1 2 2 q e−q /2 , 4
1 1 1 1 1 1 − − − − + |12 |r1 − r2 | |r1 | |r2 | |r1 − q| |r2 − q| |q|
+ (q) +
q 2
2 /2
e−q
+
1 √ 2 1 + e−q /2 q
.
The corresponding expressions for ˛ = / 1 are easily obtained from the previous ones, by an appropriate choice of variables. We find that, for M = ∞ , the BO minimum energy is −0.986901 Eh , √ √ and occurs for ˛q = 1.31037, ˛ = 0.848054 a−1 . The value 0 −0.986901 Eh is the asymptotic limit of the model described in Section 6 for M, → ∞, as may be seen from the third column of Table 1. For M finite, the BO minimum energy is the ground-state eigen-energy of the Schrödinger equation which determines the wave-function of the protons, under the interaction predicted by the BO approach. The BO results are presented in the last column of Table 1 and compare favorably with the results in the third column. Acknowledgements
Table 1 energies of the atomic and molecular systems, Numerical values of the ground-state √ without CM corrections ( = M). For the atomic system with CM corrections, = M (cf. Figure 4). log10 M
Atomic no CMC
Atomic with CMC
Molecular no CMC
BO approx.
0 1 2 3 4 5 6 7 8 ∞
−0.106103 −0.244978 −0.350755 −0.398793 −0.416051 −0.421742 −0.423566 −0.424145
−0.212207 −0.38583 −0.420211 −0.423989 −0.424371 −0.424409 −0.424413 −0.424413
−0.227952 −0.531224 −0.794730 −0.919499 −0.964838 −0.979846 −0.984663 −0.986192 −0.986677 −0.986901
−0.278702 −0.762948 −0.916081 −0.964505 −0.979818 −0.984661 −0.986192 −0.986677 −0.986901
This work is supported by Fundac¸ão para a Ciência e a Tecnologia, Portugal, under contracts PTDC/CEQ-COM/3249/2012 and PTDC/AAG-MAA/4657/2012. The support to the Coimbra Chemistry Centre through the project PEst-OE/QUI/UI0313/2014 is also acknowledged. Appendix A. Let us assume that the spin factor of the dipositronium wave function is symmetric under the exchange of positron spin quantum numbers and antisymmetric under exchange of electron spin quantum numbers, in which case it is given by (↑ (1 )↓ (3 ) − ↑ (3 )↓ (1 )) × (↑ (2 )↓ (4 ) + ↑ (4 )↓ (2 )),
172
A.J.C. Varandas et al. / Chemical Physics Letters 610–611 (2014) 167–172
or (↑ (1 )↓ (3 ) − ↑ (3 )↓ (1 )) × ↑ (2 )↑ (4 ), or
K=
(↑ (1 )↓ (3 ) − ↑ (3 )↓ (1 )) × ↓ (2 )↓ (4 ). Then the space factor of the dipositronium wave function is given by DD + ED − DE − EE . In this case (or the opposite one), the kinetic energy is given by K=
and of electron spin quantum numbers the kinetic energy is given by
KDD − KEE NDD − NEE
KDD − 2KDE + KEE NDD − 2NDE + NEE
and the potential energy is given by V=
VDD − 2VDE + VEE . NDD − 2NDE + NEE
References
and the potential energy is given by V=
VDD − VEE . NDD − NEE
Similarly, if the spin part of the dipositronium wave function is symmetric under the exchange of positron spin quantum numbers
[1] [2] [3] [4] [5]
A. Wheeler, Ann. N. Y. Acad. Sci. 48 (1946) 219. E.A. Hylleraas, A. Ore, Phys. Rev. 71 (1947) 493. D.B. Kinghorn, R.D. Poshusta, Phys. Rev. A 47 (1993) 3671. S. Bubin, K. Varga, Phys. Rev. A 84 (2011) 012509. L.J. Dunne, J.N. Murrell, Int. J. Quantum Chem. 96 (2004) 512.