On Fermionic Shadow Wave Functions for strongly-correlated multi-reference systems based on a single Slater determinant

0.62 MB·English

by  F. Calcavecchia

epl draft
On Fermionic Shadow Wave Functions for strongly-correlated
multi-reference systems based on a single Slater determinant
Francesco Calcavecchia1,2,3 and Thomas D. Kühne4,5,6
5 1 Institute of Physics, Johannes Gutenberg-University, Staudingerweg 7, D-55128 Mainz, Germany
1 2 Graduate School of Excellence Materials Science in Mainz, Staudingerweg 9, D-55128 Mainz, Germany
0
3 LPMMC, UMR 5493, Boîte Postale 166, 38042 Grenoble, France
2
4 Institute of Physical Chemistry, Johannes Gutenberg-University, Staudingerweg 7, D-55128 Mainz, Germany

n 5 Center for Computational Sciences, Johannes Gutenberg-University, Staudingerweg 9, D-55128 Mainz, Germany
a
6 Department of Chemistry, University of Paderborn, Warburger Str. 100, D-33098 Paderborn, Germany
J

7

]
h PACS 02.70.Ss–Quantum Monte Carlo methods
p PACS 31.15.vn–Electron correlation calculations for diatomic molecules
- PACS 71.15.-m–Methods of electronic structure calculations
p
m
Abstract – We demonstrate that extending the Shadow Wave Function to fermionic systems
o facilitatestoaccuratelycalculatestrongly-correlatedmulti-referencesystemssuchasthestretched
c H molecule. Thisdevelopmentconsiderablyextendsthescopeofelectronicstructurecalculations
2
s. andenablestoefficientlyrecoverthestaticcorrelationenergyusingjustasingleSlaterdeterminant.
c
i
s
y
h
p
[ Introduction. – One of the most outstanding prob- [7], the formal scaling of quantum Monte Carlo (QMC)

lems of computational physics and quantum chemistry is methods is similar to those of HF and DFT [8–12]. How-
1
theabilitytodeviseaquantitativelyprecise,yetcomputa- ever, since it typically relies on HF or DFT orbitals to
v
8 tionally tractable, method to accurately break a chemical construct the Slater determinant (SD), it only allows to
2 bond across an entire reaction coordinate. A particularly extract the vast majority of dynamic electron correlation,
4 simple example is the H molecule, in particular when the but suffers from exactly the same static correlation error
2
1
covalent bond between the H atoms is stretched. Effec- that is characteristic for single-reference electronic struc-
0
tive single-particle theories, such as the widely employed ture methods [13].
.
1
Hartree-Fock (HF) or Density Functional Theory (DFT) In the present work we demonstrate that extending the
0
methods, describe the covalent bond well, but the energy ShadowWaveFunction(SWF),whichwasfirstintroduced
5
1 is severely overestimated upon dissociation [1]. This well- by Kalos and coworkers [14,15], to fermionic systems al-
: known problem is attributed to the multi-reference char- lows bypass the static correlation problem and permits to
v
acter of the stretched H molecule, or static electron cor- studystrongly-correlatedmulti-referencesystemswithina
i 2
X relation that arises in situations with degeneracy or near- much more efficient single-determinant scheme [16–20].
r degeneracy,asintransitionmetalchemistryandstrongly-
a Method. – Letusbeginbyverybrieflyreviewingthe
correlated systems in general [2]. As a consequence, the
basic principles underlying the VMC method that relies
stretched H molecule and similar problems are typically
2 ontheRayleigh-Ritzvariationalprincipleandimportance
dealt with using multi-determinant wave functions [3].
sampled Monte Carlo (MC) techniques to efficiently eval-
However, for larger systems with many degeneracies, the
uate high-dimensional integrals [21–23]. Due to the fact
number of determinants quickly becomes unfeasible [4,5].
that the exact WF of the electronic ground state is un-
known from the outset, it is approximated by a trial WF
VariationalMonteCarlo(VMC)isanaccuratestochas-
ψ (R,α), where R ≡ (r ,r ,...r ) are the particle co-
tic method that utilizes the full many-body wave func- T 1 2 N
ordinates. The variational parameters α ≡ (α ) ,
tion (WF) and permits to approximately solve the many- i i=1,...n
which corresponds to the lowest variational energy
body Schrödinger equation [6]. In contrast to quantum-
chemical methods, where the computational complexity (cid:82) dRψ∗(R,α)Hψ (R,α)
grows rapidly with the number of considered electrons N E = (cid:82) dRψT∗(R,α)ψ T(R,α) (1)
T T
p-1
F. Calcavecchia and T. D. Kühne
represents the best possible approximation of the elec- 2. Sample from the trial WF, which corresponds to the
tronic ground state within the given trial WF ψ (R,α). initial variational parameters α(0) and estimate the
T
Thus, the accuracy of a VMC simulation critically de- n-dimensional vector δα(0) by means of Eq. 6. Then,
pends on how well the particular trial WF mimics the normalize δα(0):
exact ground state WF.
δα(0)
In order to efficiently evaluate the high-dimensional in- δα(0) :=
tegral of Eq. (1), it is conveniently rewritten as |δα(0)|
(cid:82) dR|ψ (R,α)|2HψT(R,α) (cid:113)
E = (cid:82) dTR|ψ (R,αψ)T|(2R,α) , (2) 3. λ= |α(n0)|2
T
whichenablestocomputeitbymeansofMCmethods. By 4. α(1) =λδα(0)
sampling M points from the probability density function
5. i:=i+1
|ψ (R,α)|2
ρ(R)= (cid:82) T (3) 6. Estimate δα(i) by sampling from the trial WF with
dR|ψ (R,α)|2
T variational parameters α(i) and normalize it:
usingtheM(RT)2algorithm(alsoknownastheMetropolis
δα(i)
algorithm)[24],thevariationalenergycanbeestimatedas δα(i) :=
|δα(i)|
M
1 (cid:88)
E (cid:39) M Eloc(Ri), (4) 7. If i>2 then λ:=λ(cid:0)1+0.1×δα(i)·δα(i−1)(cid:1)
i=1
where 8. If i>3 then
Hψ (R,α)
E (R)≡ T (5) (cid:32) (cid:33)
loc ψ (R,α) |α(i)−α(i−2)|
T λ:=λ 0.85+0.3×
is the so-called local energy. (cid:80)2 |α(i+1−j)−α(i−j)|
j=1
Stochastic Reconfiguration. Having shown how the
high-dimensionalintegralofEq.(1)canbeefficientlycom- 9. If i>5 then
putedusingtheM(RT)2 algorithm,itisstillnecessaryde-
(cid:32) (cid:33)
terminetheoptimalvariationalparametersα, whichmin- |α(i)−α(i−4)|
λ:=λ 0.75+0.5×
imizesthevariationalenergy. Herethisissolvedbydevis- (cid:80)4 |α(i+1−j)−α(i−j)|
j=1
ing a modified version of the Stochastic Reconfiguration
(SR)methodfromSorella[25]. TheSRschemeprescribes
10. α(i+1) =λδα(i)
that the variational parameters are varied according to
n 11. Repeat steps 5-10 until convergence.
δα =λ(cid:88)f (cid:0)s−1(cid:1) , (6)
l k kl
ThemostcommonlyemployedtrialWFtodescribethe
k=1
electronic structure within VMC is the so-called Jastrow-
where
  Slater (JS) WF and consists of a single SD that is multi-
s = (cid:104)O O (cid:105)−(cid:104)O (cid:105)(cid:104)O (cid:105)
 lk k l l k  pliedbyasimplecorrelationfactoroftheJastrowformto
f = (cid:104)H(cid:105)(cid:104)O (cid:105)−(cid:104)O H(cid:105) (7)
k k k
recover most of the dynamic correlation effects [26]:
 O = ∂ ln(ψ ) 
k ∂αk T
and(cid:104)·(cid:105)≡(cid:104)ψ |·|ψ (cid:105). Oncethedirectioninthevariational ψ (R)≡det(φ (r↑))det(φ (r↓))J (R)J (R,Q), (8)
T T JS α β α β ee ep
parameters space that minimizes the variational energy
has been computed, the step length λ along this direction where α and β are the row and column indexes of the
needs to be determined. Due to the fact that determin- SDs for the spin-up and spin-down electrons, while φ
α
ing the new direction δα is computational approximately are the single-particle orbitals that are typically deter-
equally expensive than calculating the the variational en- mined by mean-field theories, such as HF or DFT. The
ergy, we found it convenient to start with a rather small Jastrow correlation factor J =e−(cid:80)i,ju(rij), where u(rij)
valueforλandcontinuouslyadjustingitontheflyduring is a two-body pseudopotential, for the electron-electron
theoptimization. Specifically,λisdecreasedwheneverthe and electron-proton interactions are denoted as J and
ee
search direction changes, which means that we are "mov- J , respectively. Here we have implemented the Yukawa-
ep
ing too fast", and increased otherwise. The modified SR Jastrow pseudopotential for both J and J , which is
ee ep
algorithm then reads as follows: defined as
1−e−Fr
1. i=0 uYUK(r)≡A r (9)
p-2
Fermionic Shadow Wave Functions
and where A and F are both variational parameters. The where a is a normalization factor, the eventual expression
Yukawa-Jastrow pseudopotential is able to satisfy Kato’s for the improved trial WF reads as
cusp condition, since
(cid:90)
AF2 e−τHψT(R)=e−τ2V(R) dSe−τ2V(S)e−(R−4τS)2(cid:104)S|ψT(cid:105).
u (r)−r−→→0 AF − r+O(r2). (10)
YUK 2 (17)
However, throughout our derivation we have assumed
However, we have not exploited the cusp condition to fix
that τ (cid:28) 1. Hence, the imaginary-time propagation is
oneofthetwoparameters,butinsteadwehavedetermined
rather short and the trial WF only slightly improved. For
both of them by means of the modified SR algorithm.
thepurposetoelongatethepropagationinimaginary-time
Shadow Wave Function. Any arbitrary trial WF can and to solve the Schrödinger equation exactly, the de-
be systematically improved using the SWF formalism of scribed procedure needs to be applied repeatedly, which
Kalos and coworkers [14], which can be derived by ap- eventually results in a formalism similar to the path-
plying the imaginary-time propagator e−τH that projects integral approach [28,29]. Yet, there is no explicit impor-
ψT (cid:54)⊥ ψGS onto the ground state WF ψGS. In order to tance sampling in path-integral MC methods [30]. There-
demonstrate this, let us decompose the trial WF into fore, following our original intention to find an improved
and computational efficient trial WF, we rather truncate
+∞
(cid:88) theprojectionafteronestepandrefinetheobtainedfunc-
ψ = c φ , (11)
T n n
tionalformatthevariationallevel. Inotherwords,instead
n=0
ofapproachingthelimitτ →0,weinsteadsubstituteτ by
where φn are the eigenfunctions of the Schrödinger equa- avariationalparameterC inthegaussianterm. Moreover,
tion, i.e. Hφn = Enφn for all n ∈ N, with En indi- we interpret the exponential e−V(R) as the Jastrow corre-
cating the associated energy eigenenvalues. Applying the lation factor J (R) for the protons and likewise e−V(S)
p
imaginary-time propagator onto ψ, we obtain as the corresponding two-body correlation term J (S) for
s
the shadows. The identity (cid:104)S|ψ (cid:105)=ψ (S) implies that
+∞ T T
e−τHψ = (cid:88)c e−τEnφ . (12) the original trial WF has to be evaluated on the shadow
T n n
coordinates S ≡ (s ,s ,...r ). The latter is particularly
n=0 1 2 N
important for the term that dictates the symmetry of the
Theprojectore−τH causesthatallexcitedcomponentsare SWF,whichisaproductoforbitalsforabosonicandaSD
exponentially decaying 1, so that eventually the ground for a fermionic system, respectively. From this it follows
state energy E0 is projected out, i.e. that any trial WF ψT can be systematically improved by
shadow formalism. The resulting SWF then reads as
+∞
(cid:88)
lim e−τHψ = lim c e−τEnφ ∝φ . (13) (cid:90)
τ→∞ T τ→∞n=0 n n 0 ψSWF(R)=Jp(R) dSe−C(cid:80)Ni=1(ri−si)2Js(S)ψT(S).
As a consequence, any arbitrary trial WF ψ can be sys- (18)
T
tematically improved by From the discussion above it is apparent that the SWF
can also be thought of as an one-step Variational Path
e−τHψ (R) = (cid:104)R|e−τH|ψ (cid:105) (14a) Integral [31].
T T
(cid:90) Although the implementation of the SWF for bosons is
= dS(cid:104)R|e−τH|S(cid:105)(cid:104)S|ψ (cid:105), (14b)
T relatively straightforward, the extension to fermionic sys-
temsisnontrivialduetotheantisymmetryrequirementof
where we have introduced an integral over a complete set
theWFtoobeythePauliexclusionprinciple. Thenatural
of Dirac deltas |S(cid:105) and omitted the inessential normaliza-
way to devise an antisymmetric version of the SWF is to
tionfactor. Assumingthatτ (cid:28)1,wenowusetheTrotter
introduceaSDforeachofthespinsasafunctionofS,i.e.
formula to approximate det(φ (s↑)) and det(φ (s↓)), respectively. This results in
α β α β
e−τ(K+V) ∼e−τ2Ve−τKe−τ2V, (15) the so-called Fermionic Shadow Wave Function (FSWF)
(cid:90)
where K represents the kinetic term of the Hamiltonian ψ (R) = J (R)J (R,Q) dSe−C(R−S)2J (S,R)
FSWF ee ep se
[27]. Using the equality
× J (S,Q) det(φ (s↑))det(φ (s↓)), (19)
sp α β α β
(cid:104)x|e−τK|y(cid:105)= e−(x−4τy)2 , (16) where J (S,R) is the electron-shadow and J (S,Q) the
a se sp
shadow-proton Jastrow correlation factor [16,18–20]. The
1If some energy eigenvalues En are negative, the corresponding FSWF,however,suffersfromasignproblem[19,20],which
term is exponentially increasing instead of decaying. Nevertheless, differs from the infamous fermion sign problem of projec-
itisalwayspossibletoaddanappropriatelychosenconstantenergy-
tion QMC methods such as Green’s function or diffusion
shifttoH,suchthatallexcitedcomponentsareagainexponentially
decaying. MC [32,33], but limits its applicability to relatively small
p-3
F. Calcavecchia and T. D. Kühne
systems. A simple ansatz to bypass the sign problem is 0.1
the Antisymmetric Shadow Wave Function (ASWF)
0.05
ψASWF(R) = Jee(R)Jep(R,Q)det(φα(r↑β))det(φα(r↓β)) y]
R
× (cid:90) dSe−C(R−S)2Jse(S,R)Jsp(S,Q), (20) ergy [ 0
n
E-0.05
where det(φα(r↑β)) and det(φα(r↓β)) are SDs as a function nal
of the electronic coordinates only [17]. Even though the o
ati -0.1 Full CI
AanSyWorFdearl,retahdeyFiSnWcluFdiesssumpaenryio-rbosidnyceciotrraeclcaotuionntsenffoetctosnolyf Vari JJASSS--WpDwFFT-pw/FSWF-pw
for symmetric, but moreover also for asymmetric, three- -0.15 ASWF-DFT
FSWF-DFT
body and backflow correlation effects [34,35].
-0.2
Application to the H molecule. – The effective- 0 2 4 6 8
2
r [Bohr]
nessofthevariousSWFsbymeansoftheVMCmethodis
demonstrated on the H molecule, whose Hamiltonian (in
2
atomic units) reads as
Figure1: VariationalbindingenergycurvesoftheH molecule,
2
1 (cid:88) 1 1 obtained by subtracting the variational energy with r =
H = −2 ∇2i + ||r −r || + ||q −q || − 10 Bohr (dissociated energy) to the total variational energy.
1 2 1 2
i=1,2 Weremarkthat,asaconsequenceofsuchoffset,alowerbind-
(cid:88) 1 ingenergydoesnotimplyalowervariationalenergy,sincesuch
, (21)
||r −q || offset is different for each trial wave function.
i j
i=1,2
j=1,2
verysmall,butmoreimportantlyapproximatelyconstant.
where r represents the electronic coordinates and q the
Thisistosaythatdespitethestrongmulti-referencechar-
protonic coordinates. Since H includes the bare Coulomb
acter of the stretched H molecule, the FSWF is not only
potential, spin interactions are neglected. Due to the fact 2
capable to recover the dynamic but also the static corre-
that the electrons of the H molecule possess antiparal-
2
lation energy, using a single SD only. As such, even with
lel spins, the SDs can be replaced by the orbitals them-
a nearly minimal basis set, the FSWF is very competitive
selves. To that extend we have considered two possibili-
with exact or highly accurate, but computational much
ties. The first one is to use simple translational invariant
more demanding, electronic structure techniques such as
plane waves (pw) orbitals, by setting φ=1 (since k =0).
1
full CI QMC [4,5], projection QMC [32,33] and the Cou-
In this way, only the Jastrow correlation factor accounts
pled Cluster method, which for the H molecule is exact
for all the relevant physics. Alternatively, more accurate 2
and equivalent to the full CI method [40].
orbitals can be computed at the DFT level. Here we have
employed the PWscf program of the Quantum Espresso
Conclusions. – To summarize, in the present work
package together with a pw cutoff of just 2 Ry, the bare
we have extended the SWF to fermionic systems and de-
CoulombpseudopotentialandthePBEexchangeandcor-
rivedtheunderlyingconnectionbetweentheSWFformal-
relation functional [36,37]. Since for the H molecule the
2 ism and projection QMC methods via the imaginary-time
signproblemisirrelevant, itispossibletodirectlyemploy
propagator. Moreover, we have demonstrated that the
both, the ASWF and the more accurate FSWF. The spe-
commonly used JS-DFT trial WF is able to accurately
cific trial WFs we have considered in the present work are
describe the covalent H bond, but fails to recover most
2
summarized in Tab. 1.
of the static correlation energy of the stretched dimer,
Results and Discussion. – In Fig. 1 the H binding whichpossessasizablemulti-referencecharacterandtypi-
2
energy curves as obtained using the various trial WFs are cally would require the usage of a multi-determinant WF.
shown together with the exact full configuration interac- However, the ASWF and especially the FSWF permits
tion (CI) reference data for comparison [38,39]. The dif- to study strongly-correlated multi-reference systems and
ferences between the considered trial WFs and the exact is able to quantitatively reproduce the exact H2binding
full CI reference is plotted in Fig. 2. It is evident that energy curve within an very efficient single-determinant
although the DFT orbitals are generally superior to sim- QMC method.
ple pw orbitals, in either case the commonly employed
JS trial WF substantially overestimates the potential en- ∗∗∗
ergyoftheH moleculeupondissociation. Incontrast,the
2
ASWFandFSWFareabletoquantitativelyreproducethe T.D.K acknowledges financial support from the IDEE
exact full CI reference data. As can be seen in Fig. 2, in project of the Carl Zeiss Foundation, as well as the Gauss
particularfortheFSWF,theenergydifferenceisnotonly Center for Supercomputing (GCS) for providing comput-
p-4
Fermionic Shadow Wave Functions
JS-pw J (R)J (R,Q)
ee ep
JS-DFT J (R)J (R,Q)φDFT(r↑)φDFT(r↓)
ee ep
(cid:90)
ASWF-pw/FSWF-pw J (R)J (R,Q) dSJ (S,R)J (S,Q)
ee ep se sp
(cid:90)
ASWF-DFT J (R)J (R,Q)φDFT(r↑)φDFT(r↓) dSJ (S,R)J (S,Q)
ee ep se sp
(cid:90)
FSWF-DFT J (R)J (R,Q) dSJ (S,R)J (S,Q)φDFT(s↑)φDFT(s↓))
ee ep se sp
Table 1: Employed trial WFs for the unpolarized H molecule.
2
[10] KolorencJ.andMitasL.,Rep.Prog.Phys.,74(2011)
0.06 JS-pw
JS-DFT 026502.
ASWF-pw/FSWF-pw [11] Austin B. M., Zubarev D. Y. and Lester W. A.,
0.05 ASWF-DFT
y] FSWF-DFT Chem. Rev., 112 (2012) 263.
R [12] Morales M. A., Clay R., Pierleoni C. and Ceper-
e [0.04 ley D. M., Entropy, 16 (2014) 287.
c
en [13] Lyakh D. I., Musial M., Lotrich V. F. and
er0.03 Bartlett R. J., Chem. Rev., 112 (2011) 182.
Diff [14] Vitiello S., Runge K. and Kalos M. H., Phys. Rev.
y 0.02
g Lett., 60 (1988) 1970.
Ener0.01 [15] Reatto L. and Masserini G. L., Phys. Rev. B, 38
(1988) 4516.
[16] Kalos M. H. and Reatto L., Monte carlo variational
0
theory for fermionsinProgress in Computational Physics
ofMatter,editedbyReatto L.andManghi F.,(World
2 4 6 8
Scientific, Singapore) 1995.
r [Bohr]
[17] Pederiva F., Vitiello S. A., Gernoth K., Fantoni
S. and Reatto L., Phys. Rev. B, 53 (1996) 15129.
[18] Pederiva F. and Chester G. V., J. Low Temp. Phys.,
Figure 2: Difference between the binding energy curves of the
113 (1998) 741.
H2molecules as obtained using the various trial WFs and the [19] Calcavecchia F., Pederiva F. and Kühne T. D.,
exact full CI reference [38,39].
Journal of Unsolved Questions, 1 (2011) 13.
[20] Calcavecchia F., Pederiva F., Kalos M. H. and
Kühne T. D., Phys. Rev. E, 90 (2014) 053304.
ing time through the John von Neumann Institute for [21] Binder K., Rep. Prog. Phys., 60 (1997) 487.
Computing (NIC) on the GCS share of the supercom- [22] KalosM.H.andWhitlockP.A.,MonteCarloMeth-
puter JUQUEEN at the Jülich Supercomputing Centre ods (Wiley-VCH, Weinheim) 2008.
(JSC). F.C. would like to acknowledge Markus Holzmann [23] Landau D. P.andBinder K.,AGuidetoMonteCarlo
for useful comments and the Nanosciences Foundation of Simulations in Statistical Physics (Cambridge University
Grenoble for financial support. Press, Cambridge) 2013.
[24] Metropolis N., Rosenbluth A. W., Rosenbluth
M. N., Teller A. H. and Teller E., J. Chem. Phys.,
21 (1953) 1087.
References [25] Sorella S., Phys. Rev. B, 71 (2005) 241103.
[26] Jastrow R., Phys. Rev., 98 (1955) 1479.
[1] KochW.andHolthausenM.C.,AChemist’sGuideto [27] Trotter H. F., Proc. Am. Math. Soc., 10 (1959) 545.
DensityFunctionalTheory(Wiley-VCH,Weinheim)2001. [28] Feynman R. P. and Hibbs A. R., Quantum Mechanics
[2] Cohen A. J., Mori-Sánchez P. and Yang W., Sci- and Path Integrals (McGraw-Hill, New York) 1965.
ence, 321 (2008) 792. [29] Kleinert H., Path Integrals in Quantum Mechan-
[3] Köhn A., Hanauer M., Mück L. A., Jagau T.-C. ics, Statistics, Polymer Physics, and Financial Markets
and Gauss J., WIREs Comput. Mol. Sci., 3 (2013) 176. (World Scientific, Singapore) 2009.
[4] Booth G. H., Thom A. J. W.andAlavi A.,J.Chem. [30] Ceperley D. and Alder B., Science, 231 (1986) 555.
Phys., 131 (2009) 054106. [31] Ceperley D. M., Rev. Mod. Phys., 67 (1995) 279.
[5] Cleland D., Booth G. H. and Alavi A., J. Chem. [32] KalosM.H.,LevesqueD.andVerletL.,Phys.Rev.
Phys., 132 (2010) 041103. A, 9 (1974) 2178.
[6] McMillan W. L., Phys. Rev., (1965) . [33] Ceperley D. M.andAlder B. J.,Phys.Rev.Lett.,45
[7] Pople J. A., Rev. Mod. Phys., 71 (1999) 1267. (1980) 566.
[8] Lüchow A., WIREs Comput. Mol. Sci., 1 (2011) 388. [34] Feynman R. P., Phys. Rev., 94 (1954) 262.
[9] Needs R. J., Towler M. D., Drummond N. D. and [35] Feynman R. P.andCohen M.,Phys. Rev.,102(1956)
RiosP.L.,J.Phys.: Condens.Matter,22(2010)023201. 1189.
p-5
F. Calcavecchia and T. D. Kühne
[36] Giannozzi P., J. Phys.: Condens. Matter, 21 (2009)
395502.
[37] Perdew J. P., Burke K. and Ernzerhof M., Phys.
Rev. Lett., 77 (1996) 3865.
[38] Kolos W. and Wolniewicz L., J. Chem. Phys., 49
(1968) 404.
[39] Kolos W., Szalewicz K. and Monkhorst H. J., J.
Chem. Phys., 84 (1986) 3278.
[40] Bartlett R. J. and Musiał M., Rev. Mod. Phys., 79
(2007) 291.
p-6