A coupled cluster theory based on quantum electrodynamics: Method for closed shells

Computational and Theoretical Chemistry 1166 (2019) 112574

Contents lists available at ScienceDirect

Computational and Theoretical Chemistry journal homepage: www.elsevier.com/locate/comptc

A coupled cluster theory based on quantum electrodynamics: Method for closed shells

T

Sambhu N. Datta Department of Chemistry, Indian Institute of Technology – Bombay, Powai, Mumbai 400 076, India

ARTICLE INFO

ABSTRACT

Keywords: Coupled cluster QED Relativistic

An electrodynamical coupled cluster (CC) methodology that is based on the standard QED Hamiltonian written using the Dirac-Fock picture of matter fields and Coulomb gauge is discussed here. It employs the combination of a radiative cluster, a pure matter cluster and its pair modification, and relies on the customary CC approach. Averaging over the radiation state is done first, and the cluster operator for radiative effects leads to Lamb, Breit and hyperfine interactions. Relativistic correlation effects are determined next while using the matter cluster in the traditional way of CC. When the matter cluster is extended to include deexcitations to negative-energy levels, vacuum polarization effects are generated from the pair part of Coulomb interaction. The resulting ground state correlation energy includes both relativistic and QED corrections, the latter including Lamb, Breit, hyperfine and pair energy contributions. The many-electron part of the theory is formulated here for closed shell species.

1. Introduction The coupled cluster (CC) treatment has developed into a mature and convenient methodology for the systematic investigation of many-body effects in atoms and molecules. Ever since its inception by Čižek and Paldus [1–3], Mukherjee et al. [4–5] and Lindgren [6] were the first to develop the Multi-reference Coupled Cluster (MRCC) method in Fock space while Jeziorsky and Monkhorst [7] developed the MRCC in Hilbert space. Around the same time Bartlett et al. [8–10] developed molecular applications and computational approaches. Mukherjee and his coworkers were the first to implement a successful form of the Multi-reference Coupled Cluster [11] with a state-specific approach [12]. Bishop [13], Farnell et al. [14] and Kümmel [15] have discussed the coupled-cluster method, its application to and its development in physics. This method was initially developed for nuclear physics by Coester and Kümmel in the 1950 s, while Čižek extended it to atomic and molecular physics in 1966. These are now standard works in manybody theory. As an almost simultaneous event, the concepts and techniques of relativistic quantum chemistry have developed into an interesting and novel subject. Several reviews and monographs have appeared, but the review by Pyykkö [16] and the book by Dyall and Fægri [17] would suffice here. Relativistic effects become pronounced in systems containing heavier atoms, and can alter the electronic structure, thereby causing measurable changes in molecular structure and energetics from the estimated nonrelativistic counterparts. For lighter atoms, intricate spectroscopic features and additional radiative

effects can be observed and compared with theory. A natural outcome of these two achievements has been the development of the relativistic coupled cluster theory and the corresponding method of calculation. The relativistic CC methodology has been prepared by a straight-forward application of the coupled cluster approach to the solution of the wave equation based on the projected Dirac-Coulomb and Dirac-Coulomb-Breit Hamiltonian operators that use phenomenological interactions. These operators are best described as the Hamiltonian operators of matter field, and can be derived from quantum electrodynamics (QED). The presently known formulation of relativistic CC theory has several features. (1) It is normally based on the Dirac-Fock orbitals that can be determined either from the Dirac-Coulomb (DC) Hamiltonian [18–19] or from the Dirac-Coulomb-Breit (DCB) Hamiltonian [20]. In both cases, Breit interaction energy is obtained as expectation value over the Dirac-Fock (DF) ground state wave function [18–25]. In the second case, it also contributes to the determination of the ground state configuration through the SCF process. (2) In practice, variation collapse is avoided at the DF level either by using the matrix representation of operators [26–27] or by employing a so-called kinetically balanced basis set for lighter elements [28]. (3) For a finite basis calculation, spurious spinor solutions of negative energy [29] are not taken into account, apparently to avoid effects akin to continuum dissolution. One must use projected interactions. The use of numerical DF orbitals as basis sets can account for an approximately correct projection [23]. (4) A multi-reference coupled cluster treatment has also been

E-mail address: [email protected]. https://doi.org/10.1016/j.comptc.2019.112574 Received 20 June 2019; Received in revised form 30 August 2019; Accepted 31 August 2019 Available online 11 September 2019 2210-271X/ © 2019 Elsevier B.V. All rights reserved.

Computational and Theoretical Chemistry 1166 (2019) 112574

S.N. Datta

formulated [25]. (5) Some authors base the relativistic CC on DouglasKroll-Hess transformation and use the two-component spinors in order to bypass the two theoretical problems mentioned in (2) and (3) [30]. The ensuing calculation remains approximate through any finite order and requires a large basis set, with the convergence being generally slow and the evaluation of radiative effects tedious. (6) Sometimes authors neglect the spin-orbit splitting of orbitals while selecting the basis spinors. A relativistic treatment is for numerical accuracy that costs computational time, space and effort, and the effects of the spinorbit interaction should be fully retained in the treatment [31]. It is logical to work out a relativistic CC method based on the Hamiltonian of QED rather than starting at the halfway mark. The main objective here is to first obtain the quantum electrodynamical (QED) interactions (Lamb, Breit and hyperfine) from a single procedure based on a radiative cluster and then to implement these in coupled cluster (CC) method. The theoretical background is briefly discussed in Section 2. A summary of the present-day relativistic CC is given in Section 3. QED based clusters are defined in Section 4 while considering the timetested standard form of the Hamiltonian of QED using Coulomb gauge.. It is demonstrated in Section 5 that QED interactions can be obtained from a treatment using the radiative cluster. It becomes possible to get the DF ground state energy with QED corrections to it. The correlation energy and the correlated wave function are obtained in Section 6 from the second step using the matter cluster. The relativistic correlation energy is accompanied by possible QED corrections to it. The matter cluster is allowed to deviate from conventionality so that the Coulombic pair terms give rise to energy corrections due to the creation and annihilation of virtual electron-positron pairs. Concluding remarks are given in Section 7 and the iconic case of noninteracting H2 molecules is inspected in Appendix I. This work provides a justification for the relativistic CC methodologies that have been already developed, and introduces additional interactions to the theoretical treatment.

[

[

DF (r ) =

Km (r )].

p

(r ).

r ),

, t )]+ = 0.

(4)

† (r ) h D, ext

(r ):

(5)

e2 2

d3r1

d3r2

D (r1)

1 r1

D (r2)

r2

(6)

is the field-theoretical density operator,

D (r )

+

D (r ) = † † ( (r ) + r)

† (r )

: +

(r ):

(r )

+ (r )

+:

(r )

†

(r ): .

(7)

Another interesting quantity is the probability current of field theory,

JD (r ) = c:

† (r )

(8)

(r ):

The 4-currents are to be used in writing down interaction between radiation and matter. Sucher [32] recommended utilizing only the part of Coulomb interaction that is projected onto the positive-energy subspace so that the continuum dissolution problem can be avoided. The DC Hamiltonian (in coordinate representation) gives rise to the Hamiltonian operator of quantum field theory (QFT), (9)

HQFT = HD, ext + HC

One may carry out a mean-field treatment with this operator. The electronic ground state configuration in configuration space is represented by the state vector 0N cs and the excited states configurations are written as nN cs (for n ≠ 0). These vectors satisfy the relation

|

n N cs

=|

+ (i ),

+ (i )

= S |um (i ) um (i)|.

+ +

N i=1

=

n N cs ,

(10)

m

The corresponding second quantized state vectors are confined to the N-electron sector of Fock space. The projected interaction is written e2 as + 1 i < j N r r + = HC ++. The corresponding second quantized i

j

operator HC ++ is obtained by using only the first component of density from the second line of (7) in Eq. (6). Similarly, the projected external N Dirac operator is written as + i = 1 hD, ext (i ) + = HD, ext +, the corresponding second quantized operator being written as HD, ext +. The projected Hamiltonian of QFT, HQFT (projected ) = HD, ext + + HC ++ , is the equivalent of the QFT no-pair Hamiltonian no pair no pair no pair HQFT = H D, ext + HC restricted to the N-electron sector of Fock space [32]. The no-pair operators are obtained using only the first and fourth components of the density operator. The corresponding DF ground configuration energy is written as

(1)

(2)

EN0 =

0 N

HQFT (projected)

0 N

=

0 N

HD, ext +

0 N

+

0 N

HC ++

0 N

(11) using the state vector 0N of Fock space. Because of the restriction (10) on the configuration space vectors, the pair part of Coulomb inPair no pair teraction HC = HC HC makes zero contribution to the energy of the DF ground configuration. The Dirac-Fock Hamiltonian is the sum of the Fock operators for all N electrons,

m

p

d 3r :

† + (r ) + (r )

=

t ) = S am (t ) um (r ),

( r , t ) = S bp ( t )

± (r

where

0 where Aext is the zeroth component of the external field 4-vector, and 0 eAext usually represents the sum of the Coulombic energy of interaction of the electron with all nuclei. The corresponding matter field is written as (r , t ) = + (r , t ) + † (r , t ) where ψ+ is the field operator for the bound and scattered states of positive energy for the attractive interaction between the particle and nuclear centers, and ψ– is the operator for the scattered states of the positive-energy antiparticle, (charge conjugated to the eigenstates of the negative-energy electron). In diagonal representation, + (r ,

t ),

HC =

The external-field Dirac operator is given by 0 hD, ext = mc 2 + c ·p + eAext (r )

± (r ,

3 (r

and the interparticle Coulomb interaction as

1 (r ), 2 DF

[Jm (r )

, t )]+ =

HD, ext =

A field theoretical formulation needs to start from the choice of a specific physical picture. The mean field picture is adopted here, with umσ (υnσ) being the positive-energy (negative-mass) eigenspinors of the N-electron relativistic Fock operator

occ um

† ± (r

The field operators are used to write down different components of Hamiltonian. To begin, the external-field electronic Hamiltonian operator is written as

2. Theoretical background

F = hD, ext +

t ),

± (r ,

(3)

where a, a† (and b, b†) are particle (and antiparticle) destruction and creation operators, and the summation holds over both discrete and continuous spectra. The fermionic field operators obey the equal time anticommutation rule,

HDF =

N F i=1

VDF =

(i )

d 3r :

HDF = HD, ext + VDF , † (r )

DF (r )

(r ): .

(12)

The eigenvalues of the Dirac-Fock Hamiltonian operator are written 2

Computational and Theoretical Chemistry 1166 (2019) 112574

S.N. Datta n as EDF , n N

HDF

3. Relativistic coupled cluster theory n =EDF

(13)

n N

Current relativistic CC formalisms are based on HQFT . There are two slightly different methods. In both cases the correlation corrections are calculated using the purely electronic cluster operators. In the first form, the mean field is derived only from the (projected) Coulomb term, and a nonrelativistic-like CC is carried out. One may also consider Breit interaction along with the Hamiltonian of QFT. Breit interaction energy can be obtained either as an expectation value over the HF ground state configuration, or as an expectation value over the wave function that results from the CC treatment. The difference between the two would represent Breit interaction correction to correlation energy. In the second and refined version, the mean field is determined by both (projected) Coulomb and (projected) Breit interactions, and subsequently the correlation corrections are calculated. In this case the orbitals (Dirac-Fock-Breit or DFB spinor eigenvectors) as well as the correlated wave function (written using the coefficients of the cluster operators) are influenced by Breit operator. The total correlation energy here would be slightly different from the sum of the Coulomb correlation and the Breit term induced correction to it. These two approaches, makeshift and refined, have been adopted in earlier work as found in references 18–25 and 30–31. Both represent a straight-forward application of the nonrelativistic CC theory for many-electron systems with Dirac operator replacing the traditional one-electron Schrödinger Hamiltonian, and with the possibility of adding Breit operator to Coulomb interaction. Here the intermediately normalized ground state wave function is related to the DF ground state configuration by an exponential operator containing the cluster operatorT ,

One may notice that the DF ground state energy differs from the DF ground configuration energy EN0 by an additional average of the projected interaction energy. The quantized radiation field is prescribed now. Operators for the creation and destruction of a photon are written as Ak† and Ak where k is the wave vector and λ is the unit vector in one of the two directions of transverse polarization. These operators follow Bose particle commutation rules [Ak , A k ] = 0 and [Ak , A k† ] = k, k , . Furthermore, ωk = ck, Nk = Ak† Ak is the number operator and Nk is the number state such that 0 EDF

Nk Nk

=

k,k

,

k, k

,

(14)

while

Ak Nk Ak†

=

Nk

=

k, k

c2

(2

,

c2

(2

1 2 (N k) k 1 2 (N k) k

)1 2 Nk

+

1)1 2

1 ,

Nk + 1 ,

(15)

and Nk = (2 c k ) Nk . In transverse gauge, the electromagnetic 4potential operators are

d3r '

A0 (r , t ) = 1

A (r , t ) =

k

0 ' D (r , t ) , '

r

r

ei (k ·r

[Ak

kt)

+ h. c. ].

(16)

where Ω is the volume in which the photons are counted. The Hamiltonian for the quantized radiation field per unit volume is 0

Hrad =

1

(k 2

c ) Nk

k

(17)

k

and the state vector is 0 Erad =

1

Nk

0 rad

0 N

= {Nk } such that the energy density is

k

such that

(18)

k

0 N

The covariant interaction of the matter 4-current with the radiation 4-potential is (1)

(2)

d3r Jµ Aµ = Hint + Hint

Hint = e

(1)

d 3r

D (r ,

(2)

(20)

T = T2 =

e c

d3rJD (r , t )· A (r , t )

respectively. In Coulomb gauge, eA0 (r , t ) =

(21)

d3r

the pair interactions have been double counted retain only one interaction between each electron pair, that is, only (1) HC = Hint 2 , thereafter causing a departure from the approximate covariance of equal time representation. The vector interactions em(2) bodied in Hint are responsible for known QED effects such as Breit interaction of order mc2α4Z2 (sum of the electron-electron magnetic interaction or Gaunt term and the retarded interaction) which is in general frequency-dependent though the time-independent form is normally utilized for a bound state calculation, Lamb shift of order mc2α5Z4 (mainly a part of the electron self-energy [33]), and hyperfine interaction of order (m2/M)c2α4Z3 (magnetic interaction between the electron and nucleus). The quantity α is fine structure constant. One writes the Hamiltonian operator of QED as 0

(2)

eT 0 N

0 N

(24)

=1

is assumed to be the true ground state wave function,

EN0 )

0 N

=Ecorrel

0 N

(25)

m< n r
rs Cmn ar† as† an am

(26)

that is a linear combination of the double excitations (while m, n, etc. are occupied orbitals and r, s, etc. are virtuals), one may write (23) as

e D (r , t ) . However, all r r (1) in Hint . One needs to

HQED = Hrad + HD, ext + HC + Hint .

0 N

Moreover, the exponential operator may be expanded in a series containing powers of T that leads to a coupling among the terms in the cluster. A large number of terms in the second and higher orders make a nonzero contribution. Using

and

Hint =

=

(HD, ext + HC ++

(19)

t ) A0 (r , t )

0 N

Because

where the interaction involving a longitudinal photon and a transverse virtual photon are

Hint = e 2

(23)

0 N

=eT

0 N

=(1 + T2 + T4 + ...)

(27)

0 N

where

T4 =

rstu m < m < n < n Cmm nn r< s< t
ar† as† at† au† an an am am , rs Cmm

tu Cnn ,

(28)

etc. The doubles cluster T2 makes the most prominent contribution to correlation energy, m
0 N

rs HD, ext + HC ++ rs mn Cmn = Ecorrel,

(29)

while in principle T4 improves the calculation by supplying quadruples and determines a better set of coefficients for the doubles from the relation

(22)

This Hamiltonian is used in the present work. 3

Computational and Theoretical Chemistry 1166 (2019) 112574

S.N. Datta rs 0 N

HD, ext + HC ++

p < q rs mn t
+

mn 0 N

p
HD, ext + HC ++

rs HD, ext + HC ++ tu pq Cmn

EN0

tu Cpq =0

tu pq

(2)

T1,int

tu Cpq

e

= (30)

These are the CCD equations involving the ground state wave function in (27). The choice T = T2 may be modified by adding a cluster of single excitations

Cmr ar† am

T1 =

=(1 + T1 + T2 + T3 + T4 + ...)

0 N

The operator

0

=eTmat

0

0 Erad =

0

0

0

Hrad

0

0

+ Ak†

e

i (k·r

k t )].

is of order αZ while the interaction (2) n1

T1,int

(2) n2

and

0 N

=(1 + T1 + T2 + T3 + T4 + ...)

T1,int

(2) Hint

0 rad

has the

are non-

(37)

0

(2)

eT1,int

=

0

(2)

such that

eT1,int

=

0

eTmat

(2)

0

eT1,int

0

(38)

0 rad

= 1.

5. Effects of radiative cluster A few observations can be made at this juncture. These are as follows: (1) Averaging over the reference state of photons yields (2)

(2)

0

0 = [H D + HC + Erad ] eT1,int

HQED eT1,int

(34)

0 rad

(35)

0

k t)

0 rad

(33)

HD, ext + HC ++

(2) T1,int

eT1,int

and at least two clusters for a CC treatment, one for the radiative effect and the other for the matter correlation. The following relations are observed:

EN0 =

ei (k·r

(2)

When one starts from the Hamiltonian operator of QED, one needs to consider a product state vector 0 rad

1

where T1 , T2 , T3, T4 , etc., are given in Eqs. (26)–(32). As the matter clusters consist of excitations from the Dirac-Fock ground state configuration, 0 0 = 1, an obvious variant of Eq. (24). The final intermediately normalized state can be written as

(32)

4. QED based coupled cluster

0 N

Nk ) ck ]

vanishing only when n1 and n2 are even positive numbers including zero. For the matter cluster Tmat , the double excitation operator T2 is a staple in both nonrelativistic and relativistic CC treatments. It can be fortified by adding the singles, triples, quadruples, etc. The clusters T1 , T2 , T3 and T4 are neither hermitean nor anti-hermitean. The net cluster in (2) the exponential is written as the sum (T1,int + Tmat ) . The first intermediately normalized state is

etc. The coefficients in the cluster of any given order are to be partly determined by the lower order coefficients.

=

n EDF + (Nk

×

order mc2α3Z3. The moments

m< m < n r< s< t

0

0 [EDF

Nk

(36)

rst † † † Cmm n ar as at an am am

T3 =

n N

Nk

d3r JD (r )·[Ak

to the argument of the exponential operator. This gives T = T1 + T2 in the exponent and is known in literature as the CCSD. Of course the coefficients involved in the operator T2 become modified from those in CCA. CCSD is computationally affordable. It works better than MP2 and CISD. Nevertheless, it is not accurate, and accuracy requires some amount of triples in some forms to be included in a calculation near the equilibrium geometry. Triples must be considered for a calculation near the breaking of a single bond. Diradicals are treated by MRCC, and again some triples are needed. Mixing of quadruples is felt for double bond breaking. In general one considers CCSDTQ with selective slices of triples and quadruples: 0 N

k

n N

Nk

(31)

m, r

n

c

0 rad

+

(2) operatorHint

linearly varies with the photon creaThe interaction tion and destruction operators. Its expectation value over the photon 0 ground state rad vanishes. It is bilinear in matter field, and because of the presence of Dirac α matrix operator, it accommodates single-par(2) (2) ticle excitations. The radiative cluster T1,int is to work withHint at least as a linear factor and the product should give nonzero average values over the states of radiation as well as matter states. This condition re(2) quires the simplest T1,int to be formed from single-particle excitations and to linearly vary with photon operators. Thus it needs to differ from (2) Hint only by a multiplicative factor for each intermediate state. Rayleigh-Schrödinger perturbation theory (RSPT) gives the second order energy correction V Q V that indicates the cluster operator to be 1 Q V whereQ = S I I (E H0 ) , the summation is over the inter-

+ 0 rad

(2)

(

(2) Hint T1,int 1 +

(2)2 1 T 3 ! 1,int

)

+ ...

0 rad

.

(39)

(2) Also, (2)

eT1,int

=1+ 0 rad

1 (2)2 T1,int 2!

0 rad

+ O ( 4Z 4).

(40)

(3) This leads to (2)

HQED eT1,int 0 rad (2) eT1,int

0 = Erad + HD, ext + HC + [

0 rad

(2)

(2)

Hint T1,int

0 rad

0 rad

+ O (mc 2 6Z 6)].

I

mediate states (I), H0 being the zeroth-order Hamiltonian and V being the perturbation. Cluster Q V creates a difference at higher orders. For example, the third order energy correction in RSPT is VQVQV V V Q Q V instead of V Q V Q V 2 that is obtained from an exponential operator using Q V . When the higher order terms become zero or negligibly small, the cluster operator Q V can be chosen. (2) (2) The operatorHint T1,int (=V Q V ) is to be manifestly hermitean. Following this argument, the radiative cluster is written as

(41)

(4) Eq. (41) represents the second order energy correction due to the (2) perturbationHint , the result being of order mc2α4Z4. In the perturbation theory based picture. The next non-vanishing contribution comes at the (2) fourth order in Hint , that is, at order mc2α6Z6, which can be safely ne4 4 glected for α Z ≪ 1. This justifies the choice of the radiative cluster for use in most of the atomic and molecular cases. (5) When one puts all the pieces together, one obtains the effect of 4

Computational and Theoretical Chemistry 1166 (2019) 112574

S.N. Datta

the interaction of one or two electron(s) absorbing (emitting) and subsequently emitting (absorbing) a transverse virtual photon, (2)

0 rad

k

n N

0 rad

0 [EDF

e2 c2

=

e i (k·r

Nk [Ak

Nk n

(2)

Hint T1,int

k t)

d3r

ei (k·r

Nk [Ak

k t)

d3r JD (r )·

+ Ak†

n EDF + (Nk

e

Nk ) ck ]

+ Ak†

e

that linearly diverges. When it is translated into the nonrelativistic limit, a renormalization of mass is deemed necessary and the same is done by subtracting a similar correction for the free electron. Thus the visible part of this self-energy is given by the difference

i (k·r 1

i (k·r

k t) ]

n N

k t)

ext

JD (r )·

1 8

k

(42)

] Nk .

d 3k

3

2 3 c

=

HLamb =

(43)

Furthermore, when the polarization vectors are considered in real forms, one obtains (2)

(2) Hint T1,int

0 rad e2 4 2c

=

d3k k

[ d3r Nk + 1 n N E0 n DF EDF

J (r )· n D

d3r

+

d3r e

n N

ck

ik·(r r )

Nk n + ck EDF

n N E0 DF

J (r )· n D

d3r e

×

JD (r )·

HLamb =

(2)

(2)

Hint T1,int

JD (r )·

.

(44)

ext

e2 4 2c

(45)

HSE , vac + HBreit

d3k k

d3r

JD (r )· n

0 EDF

n N

n EDF

n N

ck

e2 4 2c 2

d3k k2

Nk

[(Nk + 1) d3r

d3r eik·(r

d3r

ik·(r r )

d3r e

r)

JD (r )· J D (r )· ].

=

=

J (r ) n D

2 3 c n N

0

·

dk n N

d3r

0 EDF may be kept in the denominator

JD (r ) En

DF

n N

J (r ) n ( 0) D

0 EDF ) ln

d3r

n EDF

(49)

·

n N

JD (r )×

0 + ck EDF CO n EDF

(50)

0 EDF

·

n N

JD (r )×

mc 2 n EDF

(51)

n EDF

e2 4 2c 2

d 3k k2

d3r eik·(r

d3r

r)

JD (r )·kJ D (r )·k k2

(52)

a sum of magnetic and retarded interactions. Singer transformation gives the first part as

HMagnetic =

e2 2c 2

d3r

d 3r

JD (r )·JD (r ) r r

(53)

while the equality

1 2

0 EDF

0 + ck EDF

n N

J (r ) n ( n) D n EDF ) ln

JD (r )· JD (r )

[JD (r )·JD (r ) n EDF

0 EDF

0 + ck EDF

HBreit

while retaining the one-particle self-energy accompanied by in scattering, that is, one scattering of the electron by the external potential preceded by the emission (absorption) of a virtual photon and followed by the absorption (emission) of the same photon. The energy ext contribution that evolves from HSE , vac is already of order greater than αZ, and it eventually leads to the renormalization of mass and gives rise to the effect historically known as Lamb shift. Averaging over the transverse polarization in Eq. (46) yields ext HSE , vac

n EDF

DF

For any specific k vector, the two space integrals over the exponential functions in Eq. (47) are equal by symmetry. A sum over the polarization vector is carried out to get

JD (r )· J D (r )·

5.1. Lamb shift (2) T1,int

JD (r ) En

5.2. Breit interaction

(47)

n The energy difference EDF

n N

JD (r )· (46)

HBreit =

·

n n EDF where use has been made of ħckco = mc2 ≫ EDF . The k-integral is in reality a principal value integral, and this leads to the absolute value in (51). The Lamb operator is of order mc2α5Z4, and it is easily amenable to the calculation of Lamb shift as average over the DiracFock nth bound state configuration. As a classic example (albeit in the one-electron case), the 2S1/2–2P1/2 shift in hydrogen atom is about 1057.8 MHz. In this work we have used the self-energy interpretation due to Bethe and Salpeter [33] which accounts for more than 98% of 2 s1/ 2–2p1/2 shift in hydrogen atom (1052 MHz). Lamb shift can also arise from vacuum polarization (−27 MHz), vertex corrections and higher order corrections [33–34]. The calculated net shifts for H, D and He+ differ from the experimental ones only in the sixth significant digit.

where

HSE , vac =

d3r

n (EDF

ext

0 rad

2 3 c2

×

n N

n N

J (r ) n ( 0) D

A general expression valid for an arbitrary DF state, (say, for the nth state), can be written as

Most of the electron self-energy in presence of the external potential can be discarded, keeping only the part that exists even in radiation ext vacuum (HSE , vac ). Considering that for the remaining term the virtual photon energy is much greater than the excitation energy, n 0 ck > > EDF EDF , the denominators can be approximated. Furthermore, the sum over N-electron states can be replaced by unit operator, leading to the operator HBreit . Thus, after correcting for the electron self-energy, 0 rad

2 3 c2 n (EDF

0 rad ik·(r r )

d3r

dk

0

HSE , vac

which is only logarithmically divergent. The logarithmic divergence is removed by using the cut-off kco = mc/ħ as the upper boundary of the kintegral such that for a real-life calculation, one is left with an effective Hamiltonian operator for Lamb shift,

The photon matrix elements can be easily calculated by using Eqs. (14) and (15). The sum over the discrete variable k for a finite Ω can be replaced by an integral over the continuous variable k in the limit of infinite volume as shown below, 1

free

HSE , vac = HSE, vac

Nk ×

(48)

2

d3k ik·r a · k b·k 1 e = a ·b k2 k2 2r

a ·r b · r r2

(54)

reduces the second (retarded) part of Breit interaction into the form 5

Computational and Theoretical Chemistry 1166 (2019) 112574

S.N. Datta

many-electron treatment. Four main many-electron methods are CI, MBPT, Green’s function (GF) and CCA. In this work the CC theory is pursued, as it represents the most dominant and successful many-body approach of the present times.

HRetarded e2

=

d 3r

4c 2

1

d3r

JD (r )·(r

r

JD (r )· JD (r )

r

r ) JD (r )·(r r r 2

r)

.

(55)

6. Matter clusters

Therefore, Breit interaction can be written as the sum

6.1. Negative energy solutions

HBreit e2 4c 2

= +

d3r

JD (r )·(r

d3r

1 r

It is well known that the free particle 4-component wave functions need both positive-energy and negative-energy eigenspinors for completeness. Similarly, in Furry or Dirac-Fock picture, both positive-energy spinors (bound state and continuum solutions) and negative-energy spinors (continuum solutions) together form an orthonormal complete set. A trial square-integrable spinor can be written as a linear combination of the eigenvectors of both positive energy (usually representing bound state components only) and negative energy (spurious solutions). This leads to the possibility of variation collapse [39] and a min-max principle for solving the involved wave equation [40–43]. In this work it is assumed that the Dirac-Fock orbitals have been obtained from the min-max principle discussed in ref. 41. It is often taken for granted that the positive-energy Dirac or DiracFock eigenspinors representing bound states form a complete space for normalizable solutions. Though this assumption is wrong, it has been deeply entrenched in quantum chemical calculations. In the free particle picture, positive and negative energy solutions are distinctly known and the completeness relations hold separately for them. Sucher showed that the non-perturbative use of the interaction associated with the Feynman gauge photon propagator in place of interaction associated with the Coulomb gauge propagator leads to energy levels that are incorrect at the level of atomic fine structure [44]. Basically the

JD (r )·JD (r )

r

r ) JD (r )·(r r r 2

r)

.

(56)

After self-energy corrections, it reduces to its usual form in coordinate representation

HBreit =

B (i , j )

(57)

1 i
B (i , j ) =

e2 2rij

i· j

+

i· rij j·rij rij2

(58)

The Breit interaction energy in the ground state of a light atom such as helium is of the order of 105 MHz, and for neon it is about 108 MHz. To compare, the 2p1/2-2p3/2 fine structure in hydrogen atom is 1.095 × 104 MHz. Breit interaction cannot be directly observed whereas the fine structure is obtained from spectroscopy. 5.3. The hyperfine interaction The hyperfine correction is another QED effect, an additional magnetic interaction to be accommodated within the Dirac-CoulombBreit Hamiltonian. Following the two-fermion formulation of Chraplyvy [35–36] and a subsequent development by Barker and Glover [37], this interaction between the electrons and the fermion nuclei in a molecule is written as 2 2 fermion Zn e ge gn n 4mMn c 2 (nucleus )

Hhf = n

(

Dn

r

Rn 3

3

(r

Rn) Dn·(ri ri Rn 5

× Rn)

d3r 8 3

negative energy projector

† + (r ) el D, el + (r )·

Dn

3 (r i

)

Rn) .

Pair

HC

eff

(2)

= (HD, ext + HC +

0 Erad

0 rad

+ HLamb + HBreit , ++ + Hhf )

np

HC

(61)

Additional QED correction terms are known to arise from the polarization of vacuum due to the creation and annihilation of virtual Pair electron-positron pairs using the operatorHC . Such corrections are mostly blocked on the ground of the exclusion principle. After the Pauli blocking, the 1-pair and 2-pair contributions to energy appear as tiny positive corrections of orders mc2α6Z6 and mc2α8Z8, respectively. The pair terms do not appear in a relativistic configuration interaction (RCI) calculation that is based on the configurations prepared from only the DF positive-energy eigenvectors (PERCI). The 1-pair (and 2-pair) term (s) appear(s) when the all-energy eigenvectors are considered (AERCI), that is, the spurious solutions of negative energy from the DF calculation are included to obtain de-excitations from the ground state configuration in the RCI. The AERCI corresponds to a many-electron

eT1,int 0 rad

= HC

6.2. Earlier work on pair terms

(2)

HQED eT1,int

p

For a finite basis calculation in Fock space, an approximation to the no-pair Coulomb operator can be obtained by defining ψ+ in terms of positive energy eigenvectors and ψ– using the spurious eigenvectors of negative energy. The pair operators would be accordingly defined. Correlation effects are determined from the matter clusters. However, the cluster now includes not only the excitations from the DF ground state configuration to the conventional virtual orbitals but also a mixture of excitations to the virtuals and de-excitations to the spurious levels (indicated by primes) in addition to double deexcitations.

(59)

In the above σD stands for the Dirac spin matrix vector. The first two terms within the bracket in (59) give the dipolar interaction between the electron spin and the nuclear spin at a finite distance. The third term represents the Fermi contact interaction. The hyperfine splitting is of order m2c2α4Z1-3/M ~ mα2Z1-3/M hartree, M being the nuclear mass. To give an estimate of the order of magnitude, the hyperfine splitting of the hydrogen atom is 1420 MHz in its ground state, 177 MHz in 2S1/2 state and 59 MHz in 2p1/2 and 2p3/2 states. Hyperfine structure of Cd+ has been calculated by Li et al. using the relativistic CC [38]. These results allow the definition of an effective Hamiltonian from equation (41) that stands for a formulation on combining the standard Hamiltonian of QED with the radiative cluster,

HQED =

(i ) = S vp (i ) vp (i) contributes to a trial

spinor at order (p/mc) ~ αZ so that the energy levels can be incorrect at order mc2α4Z4. It is easy to realize that the negative energy solutions (in practice, the eigenvectors) must be included in the calculation of correlation energy. In field theory, the Coulomb pair operator is obtained as

(60)

where the hyperfine interaction has been added to complement the electronic Breit operator [37], and following Sucher’s suggestion, Breit operator has been considered in the projected form. The import is that by using an exponential cluster operator for radiative effects, an effective Hamiltonian can be obtained. After separating the radiation energy, the operator can be used for any QED based 6

Computational and Theoretical Chemistry 1166 (2019) 112574

S.N. Datta

min–max procedure, and the (AERCI – PERCI) energy difference was shown in 1992 and 1994 [45–46] to be in excellent agreement with an analytical estimate made by Sucher in 1987 of the 1-pair energy for helium-like species [47] when αZ ≤ 0.2. The vacuum polarization effect on energy is fundamentally a correlation effect, and it can be realized from the cluster operator technique if one considers the more complete Coulomb interaction while exploring the influence of a more detailed matter cluster. However, the pair effects would be prominent for large αZ.

rs mn

p < q rs mn t< u

+

0 N

p
0 rad

=

0

Hrad

0 N

+

(2)

+ Hint

0 rad

+

0 N

eff

0 N

HQED eff

= (HQED

(63)

0 N

]

0 0 0 [(EN0 + Erad ) + (EBreit + ELamb + Ehf0 )]) eTmat

0 0 N

0 N

(64)

where Tmat is extended to cater for the pair operator. A hierarchy of matrix equations can be derived from (64) and then solved. 6.4. CCD as example

E1 1 pair

1 pair

T2

=

2 pair

+ T2

,

rs T2 = m < n Cmn ar† as† an am , r
E2

The primed indices (p' and q') represent the negative energy eigenvectors. These de-excitations were included in AERCI back in 1992 and 1994 [45–46]. Recent authors fail to cite either this reference or even the earlier one where the 1-pair term was analytically pursued by Sucher [47]. Eq. (29) is translated in the present treatment as

m
m
0 N

rs HD, ext + HC ++ + HBreit + HLamb rs mn Cmn = Ecorrel

(66) 0 N

HD, ext +

Pair HC

+ HBreit + HLamb

Pair

HD, ext + HC

+ HBreit + HLamb

Pair

HD,ext + HC ++ + HC rs mn

+

+ HBreit + HLamb Pair

HD, ext + HC ++ + HC

0 EBreit

+

0 ELamb

HD, ext + HC ++ +

+

Ehf0 ) tpqu

Pair HC

t u pq

t u Cpq

rs Cmn

+ HBreit + HLamb + Hhf t u Cpq

+ HBreit + HLamb tpqu

rs Cmn

t u Cpq = 0.

rp mn

rp Cmn = E1

pq mn

pair

1 2mc 2

= =

1 4mc 2

m
| rp mn |2 ~O (mc 2 6Z 6), | p q mn |2 ~O (mc 2 8Z 8).

(70)

The CCD treatment can be improved by including singles clusters, 1 pair 2 pair Tmat = T1 + T2 + T2 + T2 , and the resulting method is called CCSD. The energy Eq. (64) is written as

pair

pq Cmn = E2

pair

6.5. Ccsd

(67) 0 N

+ HBreit + HLamb + Hhf

See Fig. 1 for the evaluation of pair effects through second order. Fig. 1(a) shows the outcome of MP2 and Delta RCI (AERCI–PERCI), using the negative-energy eigenvectors. The second order energy correction is positive here because the energy denominator is positive. Fig. 1(b) is directly from the QED perturbation theory where the electron-positron pair term is treated as the perturbation operator. In this case, the second order correction to energy is negative as usual, but these diagrams are blocked by Pauli’s exclusion principle. The net contribution is positive. This was first discussed in references 45 and 46, and later in 42 and 43. Fig. 1

pq † † Cmn a p aq an am .

(65)

m
Pair

HD, ext + HC ++ + HC

tu Cpq

Similar equations can be written down for other correlation coeffirs cients, one each for every coefficient of type Cmn . Two more similar sets of equations can be written down for the determination of coefficients t u t u like Cmn for single pairs and Cmn for double pairs. These three sets of rs }, equations can be solved together to obtain the sets of coefficients {Cmn t u t u {Cpq } and {Cpq }. The first three lines of Eq. (69) show the influence of rs other correlation coefficients on Cmn , and constitute the dominant effect. The last four lines reveal the minor influence of pair coefficients. As the 1-pair coefficient is at least of order α3Z3 and the 2-pair coefficient starts at order α4Z4, the last four lines in (69) can be neglected for low αZ, that is, correlation, 1-pair and 2-pair effects can be practically independent of one another. When the pair creation operators are placed in an exponential these will also create higher order pair terms. However, the higher order pair energy is non-negligible only for high αZ. The CCD pair energies can be determined by using the MBPT expressions for the involved coefficients. Very good estimates are obtained from the second order perturbation rp pq = rp mn [ m + n = p q mn theory, Cmn and Cmn r p] [m+ n ] [45–46]. For low αZ (< 0.1), these are well app q rp pq = rp mn 2mc 2 and Cmn = p q mn 4mc 2 proximated by Cmn [45–47]. It is easy to find

(62)

0 [(EN0 + Erad ) + (HBreit + HLamb + Hhf )

Tmat = T2 + T2

rs + HBreit + HLamb tu pq Cmn

(69)

0 N

0 0 0 = (EN0 + Erad ) + (ELamb + EBreit + Ehf0 )

= Ecorrel eTmat

+ HBreit + HLamb + Hhf

Pair

0 0 ELamb , EBreit and Ehf0 being the expectation values ofHLamb , HBreit andHhf respectively over the Dirac-Fock ground state configuration 0N . Combining Eqs. (38), (41) and (60), one obtains

0 rad

0 N

p
whereas (34) and (60) yield

HQED

p
(EN0

HD, ext + HC ++ + HC

0 = Erad + EN0

0 N

0 N

t u 0 0 (EN0 + EBreit + ELamb + Ehf0 ) tpqu Cpq

Correlation effects are determined from the matter clusters. Eqs. (22), (34) and (35) together give 0

Pair

HD, ext + HC ++ + HC

p < q rs mn t, u

+

6.3. Correlation treatment

HQED

HD, ext + HC ++ +

+ HBreit + HLamb

Pair HC

0 0 tu (EN0 + EBreit + ELamb + Ehf0 ) tu pq Cpq

p
0

Pair

HD,ext + HC ++ + HC

0 N

pair

(68)

eff

|(HQED

0 0 0 [(EN0 + Erad ) + (EBreit + ELamb + Ehf0 )]) eTmat |

= Ecorrel + E1 in CCSD, and it gives

whereas the equation corresponding to (30) appears as 7

pair

+ E2

pair

0 N

(71)

Computational and Theoretical Chemistry 1166 (2019) 112574

S.N. Datta

Fig. 1. (a) Second order energy diagrams representing the admixing of deexcited configurations in MP2 or during AERCI, equivalent to the min-max procedure. A negative-energy electron line (u or v) proceeds upward, and a hole line (K or K′) runs downward. (b) Second order propagator diagrams with an electron line (K or K′) proceeding upward and a positron line (u) moving downward. (Reproduced and modified from references 45 and 46).

(

eff

0 N

2 pair

+ T2

natural recourse. Nonzero matrix elements of the charge current can be obtained from two states differing by a single excitation. Hence singly excited intermediate electronic states in HLamb can contribute to both (63) and (64). The electron-nucleus hyperfine interactions are oneelectron effects and they can contribute to the correlation energy and the correlated wave function through the second term in (69), thereby modifying the coefficients and subsequently updating the correlation energy in (66). These contributions would be of course more extensive in CCSD and its derivative procedures that include some of the higher order excitations.

1 pair

0 0 0 [(EN0 + Erad ) + (EBreit + ELamb + Ehf0 )]) T1 + T2 + T2

(HQED

1

2

+ 2 T1

)

0 N

= Ecorrel + E1

+ E2

pair

pair .

(72) The operator in Eq. (60) can be partitioned as eff

0 (Erad + HLamb + Hhf )

HQED =

no pair (H D, ext

no pair

+ HC

Pair

Pair

+ HBreit , ++) + (H D, ext + HC

(73)

).

7. Conclusions

Neglecting pair energy contributions at greater than second order, one obtains the CCSD correlation energy as m, r

1 4

(hD, ext )rm Cmr +

rs ( rs mn + rs mn B ) Cmn + m, n, r , s

The hyperfine splitting can be estimated with relative ease, and it is detectable from magnetic resonance spectroscopies even for heavier atoms and molecules. For lighter atoms, intricate spectroscopic features (such as energy ordering of electronic states and the spin-orbit splitting) can be observed. Lamb shift can be observed for atoms. Breit interaction energies can be estimated for both atoms and molecules. However, most of the radiative effects remain concealed in a molecule because of extensive rotational and vibrational energy contributions as well as rovibronic and vibronic interactions at a finite temperature. A few observations can be made now: (1) The effective cluster considered in the present work has been

1 2

rs mn Cmr Cns = Ecorrel.

(74)

m, n, r , s

The amplitude contributions are to be determined from equations of the types r

(

eff

(HQED

1 pair

2 pair

ECCSD) T1 + T2 + T2

m

+ T2

1

2

+ 2 T1

)

0 N

= 0,

Teff

0 0 0 ECCSD = (EN0 + Erad ) + (EBreit + ELamb + Ehf0 )

+ (Ecorrel + E1

+ E2

pair

(2)

(75)

eff

(

HQED 1 + T1 + T2 + T2

1 pair

+ T2

(

2 pair

+ 2 T1

1 pair

E T1 + T2 + T2

+ T2

2 pair

1

2

)

1

2

+ 2 T1 0 N

=(1 + Teff ) such that . (2) The tactic employed has been to first calculate an average over the radiation state so that the radiative, matter and pair effects become separated. Formulations in the earlier sections involve

)

= 0.

+ T3 (77)

0 N

rs

2 pair

+ T2

+ T4 + ...)

and

mn

1 pair

= T1,int (1 + T1 + T2 + T3 + T4 + ...) + (T1 + T2 + T2

pair ),

(76)

HQED Teff

The main application of CCSD is to CCSD(T) where expansion up to 3 third order is considered using triples T1 T2 and T1 in the wave function. CCSD(T) is considered as a gold standard methodology for quantum chemical calculations on the ground state when a complete basis set extrapolation is utilized. It would also be possible to construct a state specific MRCC (Mk-MRCC [12]) that includes the QED interactions. In summary, while Breit interaction was added to QFT-based CC as an afterthought, here it is directly involved. Its action is at par with that of Coulomb interaction, though smaller in the absolute magnitude by an order of α2Z2. The coefficients of the exponential cluster T2 are determined by it. The pair clusters do not greatly affect Breit interaction n 0 EDF < < 1, a situation genthat relies on the inequality ( ck ) 1 EDF erally not conducive to pair creation and annihilation. The Lamb shift, generally of order mc2α5Z4, is updated to the many-body level as a

0 rad

0 N

0 = Erad + (HLamb + HBreit + Hhf )(1 + [T2 + T4 ])

+ (HD, ext + HC ++)(1 + T2 + T4 ) Pair

+ HC

1 pair

[(T2

2 pair

+ T2

pair

)(1 + T2 + T4 ) + T4

].

(78)

(3) The factor 1 gives the mean field energy plus lowest order radiative and pair corrections. (4) The matter clusters basically contribute to correlation energy. What is new in the present work? One newness is to get the three QED interactions (Lamb, Breit and hyperfine) from a single procedure based on the radiative cluster. Another is to get the pair energy from extended matter cluster formalism. The third novelty is that additional correlation energy contributions arise from radiative effects and pair terms, in, for example, Eq. (72). Fourth, higher order pair energy can appear here and these would be evidenced in high-Z atoms. Consider the example of N' noninteracting minimal-basis H2 8

Computational and Theoretical Chemistry 1166 (2019) 112574

S.N. Datta

molecules with N (=2N') electrons. This system has been treated at the nonrelativistic level (in the limit c → ∞) in the text by Szabo and Ostlund [48]. A relativistic and electrodynamical version is discussed here. Each molecule has two sets of doubly degenerate Dirac-Fock 4component spinor orbitals (ψ1↑, ψ1↓) and (ψ2↑, ψ2↓) corresponding to the bonding and antibonding sigma molecular orbitals of the nonrelativistic theory. The nature of these spinors is shown, all the associated terms are defined, and integrals are given in Appendix I. The bonding spinors are fully occupied in the DF ground state configuration. Because the molecules do not interact with one another, there are only 2 i 2i N' doubles , (i = 1,2,…, N'), with equal coefficients C for each 1i 1i double in the expanded matter cluster. It is transparent that there is no intermediate state to connect with DF ground state configuration through the 3-current operator and the contribution of HLamb is zero in Eqs. (63) and (66). Also, the hyperfine corrections for two different electron spins cancel each other in the ground state. Therefore, the correlation energy is determined only from Coulomb and Breit interactions: B Ecorrel = N C (K12 + K12 )=N[

2 1 2

B + (K12 + K12 )]

2

[

given in Appendix I. The nonrelativistic energies (energies in the limit c → ∞) are shown in equation (I.10). Familiar relativistic corrections such as the kinetic energy correction and the Darwin term are given in equation (I.11). Because the orbital angular momentum is zero in each orbital of the minimal basis calculation, the spin-orbit interaction is absent in this case. The QED corrections to energy values appear only in the form of Breit integrals as shown in Eq. (I.12). As mentioned above, the Lamb corrections do not materialize because of the want of nonzero orbital angular momentum states, while the total of hyperfine interaction energy given in Eqs. (I.13) through (I.15) becomes zero for a closed shell. Relativistic correction to correlation energy is given in Eqs. (I.16) and (I.17). To order mc2α4Z4, the only QED correction to correlation energy appears from Breit interaction as shown in Eqs. (I.18) and (I.19). The lowest order vacuum polarization effects in this example are shown in (I.20). This work has been strictly limited to the basic theory. Detailed treatments necessary for the open-shell CC (multireference CC) or a state-specific CC are still to be worked out. Also, application has been limited to the simplest exemplary system of the minimal basis hydrogen molecules. As mentioned earlier, methodologies have been established for relativistic extension of CCA, and numerical results have been generated by different researchers in this field [18–25,30–31,38]. It would be interesting to evaluate the QED contribution to correlation effects and to compare the net QED effects with the molecular energetics at a sufficiently low temperature where the rotational, vibrational, vibronic and ro-vibronic activities mostly remain frozen. Declarations A preprint of this manuscript has been published in arXiv:1904.11936 [physics.atom-ph]. The Author declares no conflict of interest.

(79)

]

It is easy to determine the coefficient C from (69). For pair correc2i 1i tions, there are N' 1-pair doubles of type , another N' 1-pair 1i 1i 22 1 2i doubles of type i , and N' 2-pair doubles of type i i , where the 1i 1i 1i 1i primed orbitals are the corresponding negative energy spinors. This example is very familiar from the nonrelativistic theory, and the only importance of it here is to get estimates of relativistic and QED corrections to various energy values and wave functions. The latter are Appendix I

For N' noninteracting minimal basis H2 molecules, the upper and lower components of the 4-component spinors ψ1↑i, ψ1↓i, ψ2↑i and ψ2↓i (for each molecule numbered as i where i = 1, … , N') are written as

u1 = N1

(10 ), l = c [

u1 = N1

(10 ), l = c [

u2 = N2

(10 ), l = c [

u2 = N2

(1 0 ) , l

1

1

2

2

1

+ mc 2

eA ext ]

1

1

+ mc 2

eAext ]

1

·pu1

pz (1 ) N1 2mc (px + ipy )(1 )

· pu1

N1 (px 2mc

ipy )(1 ) pz (1 )

+ mc 2

eAext ]

1

· pu2

pz (1 ) N2 2mc (px + ipy )(1

= c [ 2 + mc 2

eAext ]

1

· pu 2

N2 (px 2mc

2

(I.1)

)

ipy )(1 pz (1

) )

where 1σ and 1σ* are nonrelativistic-type orthonormal molecular orbitals of appropriate symmetries, and σ is the Pauli spin matrix vector. The normalization constants are N1,2 = (1+ < p2 > 1σ,1σ* /4m2c2)–1/2 ~ (1– < p2 > 1σ,1σ* /8m2c2). 0 0 =0 The electronic configuration {ψ1↑i ψ1↓i} with N = 2 N' gives the Dirac-Coulomb-Breit (DCB) mean field energy ENDFB = EN0 + EBreit while ELamb and Ehf0 = 0 . Here, 0 B EN0 = N E20, EBreit = N J11 1 2

= (hD, ext )22 + 2J12 E20

B 1

= (hD, ext )11 + J11, K12,

= 2(hD, ext )11 + J11,

B 2

E20B

B = J11 ,

B = 2J12

=

B J11 ,

B 1

=

1

B K12 ,

E20DCB

=

B 1

+ DCB 2

E20

=

+

2

+

B 2

E20B

B E2 = 2(hD, ext )22 + J22, E2 B = J22 , E2 DCB = E2 + E2 B J11= < 11 11 > =K11, J12= < 12 12 > , K12 = < 12 21 > , J22 = < 22 22> B B B B J11 = < 11 B (1, 2) 11 > =K11 , J22 = < 22 B (1, 2) 22 > = K22 B B J12 = < 12 B (1, 2) 12 > , K12 = < 12 B (1, 2) 21>

There are N' doubles of type

2 i 2i 1i 1i

(I.2)

for i = 1,2,…, N', with equal coefficients C in the expanded cluster. The operator HLamb has zero contribution 9

Computational and Theoretical Chemistry 1166 (2019) 112574

S.N. Datta

in the minimal basis case, and the net contributions of operator Hhf to Eqs. (63), (66) and (69) also vanish. As αZ ≪ 1, the influence of pairs on correlation coefficients is neglected here. Therefore, the correlation energy is given by DC Ecorrel = N CDC K12, DCB Ecorrel

B = N CDCB (K12 + K12 ).

(I.3)

The coefficients can be determined from (69) to obtain

CDC = K121 [ 1 (E 2 2

=

DC

E20 ) = ( 2 1

B CDCB = (K12 + K12 ) [ DCB

1

E20DCB ) = (

= 2 (E2 DCB

2 DC + 1 ) + [J 1 2 11

[

DC

B 2

+ 2 [(J11 +

=

2J12 + K12,

2 B 2 1 2 DCB + (K12 + K12) ] ], B B J11 ) + (J22 + J22 )] 2(J12

1 B (J 2 22

+

DC

+ J22]

[

DCB 1

B 1 )

2 1 2 K12 ] ],

B B + J12 ) + (K12 + K12 )

B J11 ).

(I.4)

The correlation energy calculation is manifestly size-extensive, DC Ecorrel /N = DCB Ecorrel /N

=[

2 DCB

[

DCB

2 1 2 + K12 ) ~

2 DC

(

DC

+ (K12 +

2 4 K12 2 DC + K12 8 3DC ..., B 2 1 2 B 2 K12) ] ]~ (K12 + K12 ) 2 DCB

Thus the contribution of Breit interaction to DF energy is DCB Ecorrel

DC Ecorrel

1

=

N

2

1

B B K12 (2K12 + K12 )

DC

2

DC

2 B K12 (J22

B J11 )

B J11

1

(I.5)

+ ...

per diatomic unit while the contribution to correlation energy is even smaller, B J22 2

B J11

+ ...

(I.6)

DC

+N Also, the DCI energy is not size-extensive, the Dirac-Coulomb energy per unit being N [ DC ( MP2, though size-extensive, has 2 2( 2 energy K12 1) , the higher order cluster contribution being absent in both DCI and MP2. The molecule is strictly in the nonrelativistic limit as Z = 1. The overall effective nuclear charge is Zeff|e| where Zeff somewhat varies from Z. Henceforth while showing the orders, Z will be written in place of Zeff. The normalized orbitals in the nonrelativistic limit are 1

(1 or 2),( or )

= (1 or 1

)

() ()

2 1 2 K12 ) ].

,

or

1 = or 0 . 0 1

or

2 DC

(I.7)

One obtains the expansion 2 1

p2

(hD, ext )11 = hnonrel

1

8m3c 2

2Z 4 [1 (R1)]2 + O (mc 2 5Z 5) m2e 2

+

(I.8)

the spin-orbit interaction being absent as the orbital angular momentum is zero. The operator hnonrel is the Schrödinger Hamiltonian operator for the H2 molecule. A similar relation with 1σ* in place of 1σ holds for (hD)22. The two-electron integrals are found as

(

p2 u1

J11 = K11 = 1

2m2c 2

)J

u1, u1

+ 2Ju1, l1 + O (mc 2 6Z 6);

J12 = Ju1, u2 + J12, p2 1 + p2 1 4m2c 2

J12 =

J1

+ Ju1, l2 + Jl1, u2 + O (mc 2 6Z 6 );

,1

K12 = Ku1, u2 + K12, p2 1 + p2 1 4m2c 2

K12 =

K1

+ Ku1, l2 + Kl1, u2 + O (mc 2 6Z 6).

,1

B B B J11 , J12 , J22

and The Breit integrals The nonrelativistic energies

are of order mc α Z and higher orders in αZ.

EN0 , nonrel = N E02,nonrel , E02,nonrel = 2 hnonrel 1, nonrel

= hnonrel

1

(I.9)

2 4 4

B K12

+ J1 ,1 , 2, nonrel = hnonrel E2, nonrel = 2 hnonrel 1 + J1

1 ,1

1

+ J1

+ 2J1 ,

,1

,1

, K1

,1

, (I.10)

are supplemented by the familiar relativistic corrections 0 EN , rel

N

2

p2 1

= E2,0 rel = 1, rel

E2, rel = 4

p2 12 8m3c 2

=

2, rel

4m3c 2

=

2Z 4

+

m2e 2

p2 12

+

8m3c 2

p2 12 4m3c 2

+

2

+

2

2Z 4 [1 m2e 2

[1 (R1)]2

2Z 4 [1 m2e 2

2Z 4

m2e 2

p2 1 J 2m2c 2 1 ,1

(R1)]2 p2

1

+ 2Ju1, l1 + O (mc 2 5Z 5),

J 2m2c 2 1 ,1

*(R1)]2 + 2 J12

[1 *(R1)]2

p2

1

J 2m2c 2 1

+ 2Ju1, l1 + O (mc 2 5Z 5),

K12 + O (mc 2 5Z 5), ,1

+ 2Ju2, l2 + O (mc 2 5Z 5).

3 2

(I.11)

The p /m c term is the kinetic correction, the contact term is the Darwin interaction, and additional corrections are obtained from the twoelectron interaction. These are accompanied by the QED corrections

10

Computational and Theoretical Chemistry 1166 (2019) 112574

S.N. Datta 0 EN , QED

B = E2,0 QED = J11 ,

N

B = J11 ,

1, QED

2, QED

B E2, QED = J22 , B = 2J12

B K12 .

(I.12)

(10), l

()

= 0 and uP = 0 , lP = 0 such that the 1 Dirac matrix βP may be replaced by the unit matrix of rank 4. The total nuclear spin states are written in this notation as Furthermore, the hyperfine splitting of orbitals is calculated using the proton spinors uP =

Singlet : I = 0,Iz = 0 > = 2 Triplet : 1, 1 > = I1 I2 ; 1, 0> = 2

1 2(

1 2(

I1 I2

P

I1 I2 );

I1 I2 + I1 I2 ); 1,

(I.13)

1 > = I1 I2 . 2

The hyperfine interaction contributes to the orbital energies when the nuclear spin states are |1,1 > and |1,–1 > . To order α Zme/MP, contributions to the 1σ spinor energies are hf 1,(11)

hf 1,(00)

=

hf 1,(00)

=

hf 1,(11)

hf 1,(10)

=

hf 1,(11)

=

=

hf 1,(10) hf 1,(11)

=

= 0,

=

hf

(1 ),

(I.14)

where |Δhf(1σ)| is the hyperfine splitting of each electronic spin-orbital due to coupling with nuclear spins, hf

e 2 2ge gP

(1 ) =

1 1 r R1 3

1

2mMP c 2

31

(z Z1)2 1 r R1 5

8 [1 (r = R1)]2 . 3

(I.15)

Similar expressions are obtained for hyperfine corrections to ε2. However, the hyperfine corrections to orbitals do not contribute to total energy in the ground state of a closed shell molecule, as the hyperfine energies of two different electron spins cancel each other. Finally, the ground state correlation energy exhibits the trends (through order mc2α4Z4): rel

=

DC

nonrel

Ecorrel, rel = N QED

=

rel

DCB

Ecorrel, QED = N

(

1 ( E2, rel 2

=

2 DC

DC

QED

[

2 1 + K12 )

=

1 B (J22 2

2 DCB

E2,0 rel ) 2

+(

(I.16) nr 2 1 2 + K12 ) ,

2 nonrel

(I.17)

B J11 )

(I.18) 1 2

B 2 + (K12 + K12 )]

+(

2 DC

2 1 2 + K12 ) .

(I.19)

After some calculations, the pair energy values are found: E1 pair N E2 pair N

=

=

1 [ mc 2

1 [ 4mc 2

< 11 21 > 2 + < 11 22 > 2 ] K11 2 + K12

2

+ 2 < 11 1 2 > 2 ]

O (mc 2 6Z 6), O (mc 2 8Z 8),

(I.20)

the primed orbitals being the negative energy eigenvectors. The 2-pair correction is negligibly smaller than the 1-pair term. Even the 1-pair term is smaller than the relativistic and other QED corrections in absolute magnitude by an amount of order α2Z2. The size extensive character is obvious at every step of calculation – not only at the known nonrelativistic level but also in relativistic corrections, radiative effects, relativistic correlation energy and pair energies.

References [13]

[1] J. Čižek, J. Paldus, L. Sroubkova, Cluster expansion analysis for delocalized systems, Int. J. Quant. Chem. 3 (1969) 149–167. [2] J. Paldus, J. Čižek, Stability Conditions for the Solutions of the Hartree-Fock Equations for the Simple Open-Shell Case, Chem. Phys. Lett. 3 (1969) 1–3. [3] J. Čižek, J. Paldus, Coupled cluster approach, Phys. Scr. 21 (1980) 251–254. [4] D. Mukherjee, R.K. Moitra, A. Mukhopadhyay, Correlation problem in open-shell atoms and molecules: A non-perturbative linked cluster formulation, Mol. Phys. 30 (1975) 1861–1888. [5] D. Mukherjee, R.K. Moitra, A. Mukhopadhyay, Applications of a non-perturbative many-body formalism to general open-shell atomic and molecular problems: calculation of the ground and the lowest π-π* singlet and triplet energies and the first ionization potential of trans-butadiene, Mol. Phys. 33 (1977) 955–969. [6] I. Lindgren, A coupled‐cluster approach to the many‐body perturbation theory for open‐shell systems, Int. J. Quant. Chem. 14 (S12) (1978) 33–58. [7] B. Jeziorski, H.J. Monkhorst, Coupled-cluster method for multideterminantal reference states, Phys. Rev. A 24 (1981) 1668–1681. [8] R.J. Bartlett, G.D. Purvis III, Molecular applications of coupled cluster and manybody perturbation methods, Phys. Scr. 21 (1980) 255–265. [9] R.J. Bartlett, Many-body perturbation theory and coupled cluster theory for electron correlation in molecules, Ann. Rev. Phys. Chem. 32 (1981) 359–401. [10] R.J. Bartlett, C.E. Dykstra, J. Paldus, Coupled cluster methods for molecular calculations, Reidel, Dordrecht, 1984, pp. 127–159. [11] D. Mukherjee, S. Pal, Use of cluster expansion methods in the open-shell correlation problem, Adv. Quant. Chem. 20 (1989) 291–373. [12] U.S. Mahapatra, B. Datta, D. Mukherjee, A size-consistent state-specific

[14] [15]

[16] [17] [18] [19] [20] [21] [22] [23]

11

multireference coupled cluster theory: formal developments and molecular applications, J. Chem. Phys. 110 (1999) 6171–6188. R.F. Bishop, An overview of coupled cluster theory and its applications in physics, Theor. Chim. Acta 80 (1991) 95–148. D.J.J. Farnell, R.F. Bishop, K.A. Gernoth, Coupled cluster method calculations of quantum magnets with spins of general spin quantum number, J. Stat. Phys. 108 (2002) 401–428. H.G. Kümmel, A biography of the coupled cluster method. In: Recent progress in many-body theories, in: R.F. Bishop, T. Brandes, K.A. Gernoth, N.R. Walet, Y. Xian (Eds.), Proceedings of the 11th international conference, World Scientific, Singapore, 2002, pp. 334–348. P. Pyykkö, Relativistic effects in structural chemistry, Chem. Rev. 88 (1988) 563–594. K.C. Dyall, K. Fægri Jr., Introduction to Relativistic Quantum Chemistry, Oxford University Press, Oxford, 2007. E. Ilyabaev, U. Kaldor, Relativistic coupled-cluster calculations for closed-shell atoms, Chem. Phys. Lett. 194 (1992) 95–98. E. Ilyabaev, U. Kaldor, Relativistic coupled-cluster calculations for open-shell atoms, Phys. Rev. A 47 (1993) 137–143. E. Eliav, U. Kaldor, Y. Ishikawa, Open-shell relativistic coupled-cluster method with Dirac-Fock Breit wave functions: energies of the gold atom and its cation, Phys. Rev. A 49 (1994) 1724–1729. R. Wesendrup, T. Hunt, P. Schwerdtfeger, Relativistic coupled cluster calculations for neutral and singly charged Au3 clusters, J. Chem. Phys. 112 (2000) 9356–9362. B.K. Sahoo, R. Chaudhuri, B.P. Das, D. Mukherjee, Relativistic coupled-cluster theory of atomic parity nonconservation: Application to 137Ba+, Phys. Rev. Lett. 96 (2006) 163003. B.K. Sahoo, B.P. Das, Relativistic coupled-cluster analysis of parity nonconserving

Computational and Theoretical Chemistry 1166 (2019) 112574

S.N. Datta

[24]

[25] [26] [27] [28] [29] [30]

[31] [32] [33]

amplitudes and related properties of the 6s 2S1/2−5d 2D3/2.5/2 transitions in 133Cs, Mol. Phys. 115 (21-22) (2017) 2765–2774, https://doi.org/10.1080/00268976. 2017.1317859. L.F. Pašteka, E. Eliav, A. Borschevsky, U. Kaldor, P. Schwerdtfeger, Relativistic coupled cluster calculations with variational quantum electrodynamics resolve the discrepancy between experiment and theory concerning the electron affinity and ionization potential of gold, Phys. Rev. Lett. 118 (2017) 023002. A. Ghosh, R.K. Chaudhuri, S. Chattopadhyay, Relativistic state-specific multireference coupled cluster theory description for bond-breaking energy surfaces, J. Chem. Phys. 145 124303 (2016), https://doi.org/10.1063/1.4962911. J. Wood, I.P. Grant, S. Wilson, The dirac equation in the algebraic approximation. IV. application of the partitioning technique, J. Phys. B – At. Mol. Opt. Phys. 18 (1985) 3027–3042. I.P. Grant, Variational methods for Dirac wave equations, J. Phys. B – At. Mol. Opt. Phys. 19 (1986) 3187–3206. R.E. Stanton, S. Havriliak, Kinetic balance: A partial solution to the problem of variational safety in Dirac calculations, J Chem. Phys. 81 (1984) 1910, https://doi. org/10.1063/1.447865. G.W.F. Drake, S.P. Goldman, Application of discrete-basis-set methods to the Dirac equation, Phys. Rev. A 23 (5) (1981) 2093–2098, https://doi.org/10.1103/ PhysRevA.23.2093. S. Hirata, T. Yanai, W.A. De Jong, Third-order Douglas-Kroll relativistic coupledcluster theory through connected single, double, triple, and quadruple substitutions: applications to diatomic and triatomic hydrides, J. Chem. Phys. 120 (2004) 3297–3310, https://doi.org/10.1063/1.1639361. A. Shee, L. Visscher, T. Saue, Analytic one-electron properties at the 4-component relativistic coupled cluster level with inclusion of spin-orbit coupling, J. Chem. Phys. 145 (2016) 184107, https://doi.org/10.1063/1.4966643. J. Sucher, Foundations of the relativistic theory of many-electron atoms, Phys. Rev. A 22 (1980) 348–361. H.A. Bethe, E.E. Salpeter, Quantum mechanics of one- and two-electron atoms,

Springer, Berlin, 1957, pp. 90–100. [34] W. Greiner, J. Reinhardt, Quantum electrodynamics, Second Ed., Springer, Berlin, 1996, pp. 233–246. [35] Z.V. Chraplyvy, Reduction of relativistic two-particle wave equations to approximate forms, I. Phys. Rev. 91 (1953) 388–391. [36] Z.V. Chraplyvy, Reduction of relativistic two-particle wave equations to approximate forms, II. Phys. Rev. 92 (1953) 1310–1315. [37] W.A. Barker, F.N. Glover, Reduction of relativistic two-particle wave equations to approximate forms, III. Phys. Rev. 99 (1955) 317–324. [38] C.B. Li, Y.-M. Yu, B.K. Sahoo, Relativistic coupled-cluster-theory analysis of energies, hyperfine-structure constants, and dipole polarizabilities of Cd+, arXiv:1892.02740v1 [physics.atom-ph] (2018); Phys. Rev. A 97 (2018) 022512. [39] S.N. Datta, Variational stability in Dirac-Hartree-Fock theory, Chem. Phys. Lett. 74 (1980) 568–572. [40] S.N. Datta, S. Jagannathan, Analysis and merit of the constrained-component variation in Dirac theory, Pramana – J. Phys. 23 (1984) 467–473. [41] S.N. Datta, G. Devaiah, The minimax technique in relativistic Hartree-fock calculations, Pramana – J. Phys. 30 (1988) 387–405. [42] S.N. Datta, A review of the min-max approach to the solution of relativistic electron wave equation, arXiv:1709.07061 [quant-ph] (2017). [43] S.N. Datta, Min-max and max-min approaches to the solution of wave equation in relativistic quantum chemistry, Ind. J. Phys. 93 (2019) 285–299. [44] J. Sucher, On the choice of the electron-electron potential in relativistic atomic physics, J. Phys. B – At. Mol. Opt. Phys. 21 (1988) L585–L592. [45] S.N. Datta, Bound state solutions of the two-electron Dirac-Coulomb equation, Pramana – J. Phys. 38 (1992) 51–75. [46] S.N. Datta, Relativistic quantum chemistry and rigorous variational analysis, Proc. Ind. Acad. Sci., Chem. Sci. 106 (1994) 445–466. [47] J. Sucher, Relativistic many-electron Hamiltonians, Phys. Scr. 36 (1987) 271–281. [48] A. Szabo, N.S. Ostlund, Modern Quantum Chemistry, Dover, New York, 1996, pp. 287–319.

12