A second-quantization framework for the unified treatment of relativistic and nonrelativistic molecular perturbations by response theory

A second-quantization framework for the unified treatment of relativistic and nonrelativistic molecular perturbations by response theory

Trygve Helgaker^a兲and Alf Christian Hennum

Department of Chemistry, University of Oslo, P.O. Box 1033 Blindern, N-0315 Oslo, Norway Wim Klopper

Lehrstuhl für Theoretische Chemie, Institut für Physikalische Chemie, Universität Karlsruhe (TH), D-76128 Karlsruhe, Germany

共Received 22 November 2005; accepted 29 March 2006; published online 13 July 2006兲

A formalism is presented for the calculation of relativistic corrections to molecular electronic energies and properties. After a discussion of the Dirac and Breit equations and their first-order Foldy-Wouthuysen关Phys. Rev. 78, 29共1950兲兴transformation, we construct a second-quantization electronic Hamiltonian, valid for all values of the fine-structure constant ␣. The resulting

␣-dependent Hamiltonian is then used to set up a perturbation theory in orders of ␣², using the general framework of time-independent response theory, in the same manner as for geometrical and magnetic perturbations. Explicit expressions are given to second order in␣²for the Hartree-Fock model. However, since all relativistic considerations are contained in the␣-dependent Hamiltonian operator rather than in the wave function, the same approach may be used for other wave-function models, following the general procedure of response theory. In particular, by constructing a variational Lagrangian using the␣-dependent electronic Hamiltonian, relativistic corrections can be calculated for nonvariational methods as well. ©2006 American Institute of Physics.

关DOI:10.1063/1.2198527兴

I. INTRODUCTION

Although, in molecular quantum chemistry, a nonrelativ- istic Schrödinger treatment of the electronic system is often sufficient, there are many situations where relativistic effects cannot be ignored. In such cases, it may sometimes be nec- essary to resort to a fully relativistic treatment of the elec- tronic system—for example, by solving the Dirac-Hartree- Fock self-consistent field equations in place of the simpler nonrelativistic Hartree-Fock equations. This is particularly true for molecular systems containing heavy atoms, where a nonrelativistic treatment may give erroneous results, bearing little relationship to the true system. Often, however, the nonrelativistic treatment may be essentially right but in need of a relativistic correction for agreement with experimental observations or for a reliable prediction of experimental mea- surements. In such cases, the best approach may be to retain the nonrelativistic solution but to improve on it by means of perturbation theory.¹

In the literature, several approaches have been developed for including the effects of relativity in molecular calcula- tions, avoiding a full four-component treatment. Apart from the relativistic direct perturbation theory共DPT兲pioneered by Rutkowski and Kutzelnigg,^2–6 several effective two- component methods have been proposed—in particular, the Douglas-Kroll-Hess共DKH兲method,^7–14the zero-order regu- lar approximation 共ZORA兲 method of Chang et al.¹⁵ and van Lentheet al.,^16,17and Dyall’s method of elimination

of small components.^18–21In the present paper, we present a relativistic perturbation theory within the language of second quantization, using techniques previously developed for the calculation of molecular properties such as forces and force constants.^22,23 By isolating all explicitly relativistic complications in the Hamiltonian operator, our approach provides a uniform treatment of relativistic and nonrelativis- tic properties, applicable to all electronic-structure models for which nonrelativistic response theory has been developed.

In this paper, we discuss in detail the evaluation of the first- and second-order relativistic corrections to the elec- tronic energy in Hartree-Fock theory, for which the calcula- tion of perturbational corrections is particularly simple.

However, for more accurate wave-function models such as those of coupled-cluster theory, the theory of molecular properties 共and therefore relativistic corrections兲 is only slightly more complicated, provided the electronic energy is formally cast in a variational Lagrangian form. The tech- niques presented here in detail for variational wave functions are therefore also applicable to nonvariational models such as coupled-cluster theory.

To calculate relativistic corrections to the molecular electronic energy, we first need to construct an electronic Hamiltonian, valid for all values of the speed of light c in vacuum. We begin in Sec. II by considering the two-electron Breit equation, which we subject to a first-order nonunitary Foldy-Wouthuysen transformation.²⁴ Next, in Sec. III, we use this transformed Breit equation as our basis for

a兲Electronic mail: [email protected]

0021-9606/2006/125共2兲/024102/17/$23.00 125, 024102-1 © 2006 American Institute of Physics

constructing ac-dependent second-quantization Hamiltonian that reduces to the standard nonrelativistic Hamiltonian of second quantization asctends to infinity. Finally, in Sec. IV, we set up a relativistic perturbation theory based on this c-dependent Hamiltonian, employing the standard apparatus of molecular response theory, with the Hartree-Fock model as an example. Section V contains some concluding remarks.

II. THE TRANSFORMED BREIT EQUATION

We begin this section by considering the one-electron Dirac equation in Sec. II A, subjecting it to a first-order共non- unitary兲 Foldy-Wouthuysen transformation. Next, we go on to consider, in Sec. II B, the corresponding transformation of the two-electron Breit equation. After a discussion of the associated unitary transformations in Sec. II C and of the nonrelativistic limit in Sec. II D, we conclude by presenting a derivation of the Breit-Pauli Hamiltonian in Sec. II E.

A. The Dirac equation

Consider an electron of massm_eand charge −emoving in the presence of an external electromagnetic field repre- sented by the scalar and vector potentials␾^andA, as well as in the potential generated by point-charge nuclei of charge ZKeat positionsRK. Assuming a static external electromag- netic field, the time-independent Dirac equation is given by

H^D⌿共r兲=E⌿共r兲, 共1兲

where the Dirac Hamiltonian takes the form

H^D=c␣^·␲⁺␷⁺␤^mec². 共2兲 Herecis the speed of light in vacuum; the vector␣contains the three Dirac operators␣x,␣y, and␣z, which together with the operator ␣0=␤ satisfy the anticommutation relations 关␣i,␣j兴+= 2␦ij; the vector ␲ contains the three kinetic- momentum operators

␲^{=p+eA= −}iប⵱ +eA; 共3兲 and␷is the potential-energy term

␷⁼␷nuc+␷ext= − e²

4␲⑀0

兺

whererKis the distance from the electron to nucleusK and

⑀0is the electric constant. In the standard representation, the Dirac operators are represented by the matrices

␣i=

冉

冊

冉

冊

where the two-by-two Pauli spin matrices are given by

␴x=

冉

冊

冉

冊

冉

冊

These matrices are representations of the Pauli spin operators

␴i, obeying the algebra共Einstein summation assumed兲

␴i␴j=␦ij+i⑀ijk␴k, 共7兲

where␦ijis the Kronecker delta and⑀ijkthe Levi-Civita an- tisymmetric permutation symbol.

For the purpose of setting up a second-quantization per- turbation theory, we shall find it convenient to work in a different representation of the Dirac equation. We begin by writing the Dirac equation in terms of spinors ␺^L共r兲 and

␺^S共r兲as

H^D⌿共r兲=E⌿共r兲, 共8兲

where the Hamiltonian H^Dand the bispinor wave function

⌿共r兲are given by

H^D=

冉

冊

冉

冊

Subjecting this equation to the transformation generated by the nonsingular matrix C, we obtain the equivalent general- ized eigenvalue problem

共C^†H^DC兲C⁻¹⌿共r兲=E共C^†C兲C⁻¹⌿共r兲. 共10兲 In particular, using the first-order Foldy-Wouthuysen trans- formation matrix²⁴

C=

冉

冊

H^D⌽共r兲=ES^D⌽共r兲, 共12兲 where the Hamiltonian and overlap matrices are given by

H^D=

冉

冊

S^D=

冉

冊

We have here introduced the kinetic-energy operator t= 2mec²␳^2=共␴^·␲兲²

2me

= ␲² 2me

+ eប 2me

B·␴^, 共14兲 where the magnetic induction Bis related to the vector po- tential asB=⵱⫻A. The transformed bispinor in Eq.共12兲is related to the old bispinor in Eq.共8兲as

⌽共r兲=

冉

冊

冉

冊

In the following, we shall use this first-order Foldy- Wouthuysen-transformed Dirac equation to set up our rela- tivistic perturbation theory.

In the nonrelativistic limit, wherec⁻¹ and␳both tend to zero, the generalized eigenvalue problem Eq.共12兲reduces to a standard Hermitian eigenvalue problem, with a factoriza- tion into positive共electronic兲and negative共positronic兲equa- tions,

共mec²+t±␷兲␸±共r兲= ±E␸±共r兲, 共16兲 where we have retained the共infinite兲additive rest-mass term.

In the following, we shall refer to bispinors with ␸−共r兲⬅0 and ␸+共r兲⬅0 in Eq. 共15兲 as electronic and positronic bis- pinors, respectively.

B. The Breit equation

For many-electron systems, we must also consider two- electron interactions. The two-electron analog of the Dirac equation is the Breit equation

H^B⌿共r₁,r₂兲=E⌿共r₁,r₂兲. 共17兲 Here⌿共r1,r₂兲is a two-electron wave function and the two- electron Breit Hamiltonian is given by

H^B=H₁^D+H₂^D+G12+B12, 共18兲 whereHi

Dis the Dirac Hamiltonian Eq.共2兲of electroniand G12= e²

4␲⑀0

1

r₁₂, 共19a兲

B12= − e² 8␲⑀0

r₁₂² 共␣1·␣2兲+共␣1·r₁₂兲共r₁₂·␣2兲

r₁₂³ , 共19b兲

are the Coulomb and Breit interaction operators, respec- tively. To transform the Breit Hamiltonian in the same man- ner as the Dirac Hamiltonian in Sec. II A, we write the Breit equation in the form

H^B⌿共r1,r₂兲=E⌿共r1,r₂兲, 共20a兲

H^B=H₁^D丣H₂^D+G12+B12, 共20b兲

where the one-electron operatorH₁^D丣H₂^Dis the direct sum of the Dirac Hamiltonians Eq. 共9兲 of the two electrons and the two-electron Coulomb and Breit operators are given by

G12=g₁₂I₄, g₁₂= e² 4␲⑀0

1

r₁₂, 共21a兲

B12=b₁₂J₄,

共21b兲 b₁₂= − e²

8␲⑀0

r₁₂²共␴1·␴2兲+共␴1·r₁₂兲共r12·␴2兲

r₁₂³ .

The four-dimensional exchange matrixJ₄that appears in the Breit operator is obtained from the corresponding unit matrix I₄ by reversing the order of the columns. The Breit operator Eq.共21b兲is often decomposed into two distinct contributions

b₁₂=b₁₂^G +b₁₂^g , 共22a兲

b₁₂^G = − e² 4␲⑀0

␴1·␴2

r₁₂ , 共22b兲

b₁₂^g = − e² 8␲⑀0

共␴1·⵱1兲共␴2·⵱2兲r₁₂, 共22c兲 whereb₁₂^G is the Gaunt term and b₁₂^g is the gauge term.共We include no positive-energy projection operators in our Hamiltonian, restricting ourselves to relativistic perturbation theory for calculations without such projectors applied.兲The two-electron wave function in the Breit equation Eq. 共20兲 may be expanded in antisymmetrized direct products of bis- pinors

⌿共r₁,r₂兲=_␮⬎␯

兺

共23兲 where the bispinors⌿_␮共ri兲 may be taken as the eigenfunc- tions of the Dirac HamiltonianHiD.

We now subject the Breit equation to the same first-order Foldy-Wouthuysen transformation as the Dirac equation in Eq. 共10兲. Using the transformation matrix

C₁₂=C₁丢C₂, 共24兲

where Ci is obtained from Eq. 共11兲 by replacing ␳ ^with␳i

=␴i·␲i/ 2mec, we obtain

H^B⌽共r₁,r₂兲=ES^B⌽共r₁,r₂兲, 共25兲 where the transformed two-electron wave function is given by⌽共r₁,r₂兲=C₁₂⁻¹⌿共r₁,r₂兲and the transformed Breit Hamil- tonian and overlap matrices by

H^B=C₁₂^†H^BC₁₂=H₁^D丣H₂^D+G₁₂+B₁₂, 共26a兲

S^B=C₁₂^† C₁₂=S₁^D丢S₂^D. 共26b兲

The Dirac matrices H₁^DandH₂^Dare those given in Eq.共13兲 and the two-electron Coulomb and Breit operators are given by

G₁₂=C₁₂^†G12C₁₂=共C₁^†丢C₂^†兲g12共C1丢C₂兲, 共27a兲

B₁₂=C₁₂^† B12C₁₂=共C₁^†丢C₂^†兲b12共J2C₁丢J₂C₂兲, 共27b兲 where J₂=␴x is the two-by-two exchange matrix. We now introduce the notation

A=

冢

冣

where superscripts eand p denote electronic and positronic components, respectively: the first two for particle one共row and column兲, the last two for particle two共row and column兲. After some algebra, we find that the elements of the Cou- lomb operator are given by

G₁₂^eeee=G₁₂^ppee=G₁₂^eepp=G₁₂^pppp=g₁₂+␳1g₁₂␳1+␳2g₁₂␳2

+␳1␳2g₁₂␳1␳2, 共29a兲 G₁₂^epee= −G₁₂^peee=G₁₂^eppp= −G₁₂^pepp=␳1g₁₂−g₁₂␳1

+␳1␳2g₁₂␳2−␳2g₁₂␳1␳2, 共29b兲 G₁₂^eeep= −G₁₂^eepe=G₁₂^ppep= −G₁₂^pppe=␳2g₁₂−g₁₂␳2

+␳1␳2g₁₂␳1−␳1g₁₂␳1␳2, 共29c兲 G₁₂^pepe= −G₁₂^eppe= −G₁₂^peep=G₁₂^epep=␳1␳2g₁₂−␳1g₁₂␳2

−␳2g₁₂␳1+g₁₂␳1␳2, 共29d兲

and the elements of the Breit operator by

B₁₂^eeee= −B₁₂^ppee= −B₁₂^eepp=B₁₂^pppp=␳1␳2b₁₂+␳1b₁₂␳2

+␳2b₁₂␳1+b₁₂␳1␳2, 共30a兲

B₁₂^epee=B₁₂^peee= −B₁₂^eppp= −B₁₂^pepp=␳2b₁₂+b₁₂␳2

−␳1b₁₂␳1␳2−␳1␳2b₁₂␳1, 共30b兲 B₁₂^eeep=B₁₂^eepe= −B₁₂^ppep= −B₁₂^pppe=␳1b₁₂+b₁₂␳1

−␳2b₁₂␳1␳2−␳1␳2b₁₂␳2, 共30c兲 B₁₂^pepe=B₁₂^eppe=B₁₂^peep=B₁₂^epep=b₁₂−␳1b₁₂␳1−␳2b₁₂␳2

+␳1␳2b₁₂␳1␳2. 共30d兲

Whereas the Coulomb operator is unaffected by changing subscriptseetoppbut changes sign whenpe is changed to ep, the Breit operator behaves in the opposite manner, chang- ing sign wheneeis changed toppbut not whenpechanges toep. We also note the following conjugation properties:

共G₁₂^eeee兲^†=G₁₂^eeee, 共G₁₂^peee兲^†=G₁₂^epee,

共31a兲 共G₁₂^eepe兲^†=G₁₂^eeep, 共G₁₂^pepe兲^†=G₁₂^epep,

共B₁₂^eeee兲^†=B₁₂^eeee, 共B₁₂^peee兲^†=B₁₂^epee,

共31b兲 共B₁₂^eepe兲^†=B₁₂^eeep, 共B₁₂^pepe兲^†=B₁₂^epep,

which are identical forG₁₂andB₁₂.

C. Unitary transformations of the Dirac and Breit equations

In Secs. II A and II B, we introduced a nonunitary first- order Foldy-Wouthuysen transformation of the Dirac and Breit equations, based on the nonsingular matrix C of Eq.

共11兲. For some purposes, it is more convenient to consider the related standard eigenvalue problems generated by the the unitary first-order Foldy-Wouthuysen matrix

U=C共C^†C兲^−1/2=C关S^D兴^−1/2=共1 +␳²兲^−1/2

冉

冊

共32兲 Carrying out the transformation Eq.共10兲withCreplaced by U, we obtain the Dirac eigenvalue problem

H˜^D⌽共r兲=E⌽共r兲, 共33兲

where the transformed Dirac Hamiltonian is given by

H˜^D=

冉

冊

in the following notation for renormalized operators:

˜a=共1 +␳²兲^−1/2a共1 +␳²兲^−1/2. 共35兲 Transforming the Breit equation Eq. 共20兲 with U₁₂

=C₁₂关S^B兴^−1/2in place ofC₁₂in Eq.共26a兲, we likewise obtain

H˜^B⌽共r1,r₂兲=E⌽共r1,r₂兲, 共36a兲

H˜^B=H˜

1 D丣H˜

2 D+G˜

12+B˜

12, 共36b兲

where the elements of G˜

12=关S^B兴^−1/2G₁₂关S^B兴^−1/2 and B˜

12

=关S^B兴^−1/2B₁₂关S^B兴^−1/2 are obtained by replacing g₁₂ with˜g₁₂ in Eq.共29兲andb₁₂ by˜b

12in Eq.共30兲.

Comparing the standard Dirac and Breit eigenvalue problems Eqs.共33兲and共36兲with the corresponding general- ized eigenvalue problems Eqs. 共12兲 and 共25兲, we note that the Hamiltonians differ only in that the spinor operators in the standard eigenvalue problems have been subjected to the normalization Eq.共35兲. It is possible to construct a relativis- tic perturbation theory based on either representation. How- ever, we shall find it convenient to use the representation based on un-normalized operators. For a general discussion of the relationship between the normalized and un- normalized Foldy-Wouthuysen transformations, see Ref.1.

D. The nonrelativistic limit of the Breit equation Let us now consider the nonrelativistic limit of the Breit equation Eq. 共36兲. From the structure of the transformed Dirac Hamiltonian Eq. 共34兲 and of the Coulomb and Breit operators Eqs. 共29兲 and 共30兲, we find that the transformed Breit equation Eq. 共36兲, expressed in terms of bispinors 共omitting arguments for brevity兲, takes the form

冢

冣

⫻

冢

冣

冢

冣

The zero-order diagonal terms arise from the diagonal part of the Dirac Hamiltonian Eq. 共34兲 and from the Coulomb op- eratorg₁₂in Eq.共29a兲, whereas the zero-order skew-diagonal terms arise from the zero-order Breit operator b₁₂ in Eq.

共30d兲. To extract the nonrelativistic limit for two electrons, we follow DPT, aligning the energy scale with the nonrela- tivistic limit

Es=E− 2mec² 共38兲

and introducing the scaled wave function,

冢

冣

^冢

^冣

In terms of the scaled energy and wave function, the Breit equation Eq.共37兲becomes

冢

冣 ^冢

^冣

=Es

冢

冣 冢

冣

where the purely electronic diagonal block of the Hamil- tonian is given by

H˜

eeee

B − 2mec²=t₁+t₂+␷1+␷2+g₁₂+O共c⁻²兲. 共41兲 As c tends to infinity, the scaled Breit equation Eq. 共40兲 factorizes, reducing to the nonrelativistic two-electron Schrödinger equation

共t1+t₂+␷1+␷2+g₁₂兲␸s ee=Es␸s

ee. 共42兲

We also note that the lowest-order relativistic correction to the nonrelativistic energy is given, in the usual manner of perturbation theory, by the expectation value of those terms inH˜

eeee

B − 2mec² that are proportional toc⁻², with respect to the zero-order wave function ␸s

ee. Finally, we note that, in Sec. IV, we shall introduce a similar scaling of the Hartree- Fock wave function in second quantization, in setting up relativistic perturbation theory.

E. The Breit-Pauli Hamiltonian

A number of derivations of the Breit-Pauli Hamiltonian have been published over the years—see Refs. 25–27. An alternative derivation is given here, based on the transformed Breit equation Eq.共36兲.

To determine the Breit-Pauli Hamiltonian, we consider, according to the discussion in Sec. II D, the purely electronic part of the unitarily transformed Breit Hamiltonian Eq.共36兲. From Eqs. 共29兲, 共30兲, and 共34兲 and from 共1 +␳i

2兲^−1/2= 1

−¹₂␳i2+O共c⁻⁴兲, we obtain H˜

eeee

B =

兺

i=1

2

冉

冊

−1

2关␳1,关␳1,g₁₂兴兴−1

2关␳2,关␳2,g₁₂兴兴

+关␳1,关␳2,b₁₂兴+兴++O共c⁻⁴兲, 共43兲 where the operators t, ␷, g₁₂, andb₁₂ are found in Eqs.共4兲, 共14兲, and 共21兲. Ignoring the rest-mass term 2mec² and the terms proportional toc⁻⁴, we arrive at the Breit-Pauli Hamil- tonian,

H^BP=H^NR+H^MV+H^ext+H^C1+H^C2+H^G+H^g. 共44兲 We have here introduced the nonrelativistic Hamiltonian

H^NR=t₁+t₂+␷1+␷2+g₁₂, 共45兲 the mass-velocity operator

H^MV= − 2mec²共␳14+␳24兲= − 1

2mec²共t₁²+t₂²兲, 共46兲 the external-field and one- and two-electron electrostatic Coulomb relativistic operators

H^ext= −¹₂关␳1,关␳1,␷ext兴兴−¹₂关␳2,关␳2,␷ext兴兴, 共47a兲 H^C1= −¹₂关␳1,关␳1,␷nuc兴兴−¹₂关␳2,关␳2,␷nuc兴兴, 共47b兲 H^C2= −¹₂关␳1,关␳1,g₁₂兴兴− ¹₂关␳2,关␳2,g₁₂兴兴, 共47c兲 and finally the Gaunt and gauge-term operators

H^G=关␳1,关␳2,b₁₂^G兴+兴+, 共48a兲

H^g=关␳1,关␳2,b₁₂^g兴+兴+, 共48b兲 obtained by decomposing the Breit operatorb₁₂according to Eq. 共22兲.

To evaluate the Coulomb corrections, we note that the double commutator of the potential f共r1兲with␳1is given by 共in the Einstein summation convention, with indicesi,j, and kfor the Cartesian directions兲

关␴1i␲1i,关␴1j␲1j,f共r₁兲兴兴

=␴1i␴1j关␲1i,关␲1j,f共r1兲兴兴+关␴1i,␴1j兴关␲1j,f共r1兲兴␲1i

=共␦ij+i⑀ijk␴1k兲共p1ip_1jf共r1兲兲+ 2i⑀ijk␴1k共p1jf共r1兲兲␲1i

= −ប²共ⵜ₁²f共r1兲兲− 4s₁共⵱1f共r1兲兲⫻␲1, 共49兲 where we have introduced the spin operator as

si=ប

2␴i, 共50兲

and made use of the relation Eq. 共7兲. Substituting ␷=␷nuc

+␷ext andg₁₂ for f共r₁兲, we find that the leading relativistic electric-field and Coulomb corrections are given by

H^ext= eប² 8m_e²c²

兺

i=1 2

共⵱i·Ei兲

+ e

4m_e²c²

兺

i=1 2

si·共Ei⫻␲i−␲i⫻Ei兲, 共51a兲

H^C1= e²ប² 8⑀0m_e²c²

兺

K

ZK

兺

i=1 2

␦共riK兲 + e²

8␲⑀0m_e²c²

兺

兺

i=1

2 s_i·r_iK⫻␲i

r_iK³ , 共51b兲

H^C2= − e²ប²

4⑀0m_e²c²␦共r₁₂兲− e² 8␲⑀0m_e²c²

⫻

冉

r₁₂³

冊

where the first and second terms in the operators are the Darwin and spin-orbit corrections, respectively, andEiis the external electric field evaluated at the position of electroni.

Although these expressions have been obtained for static fields E= −⵱␾ ^andB=⵱⫻A, they are, in fact, correct also in the time-dependent case, where E= −⵱␾⁻⳵^A/⳵^{t. How-} ever, for a consistent derivation of the Breit-Pauli Hamil- tonian in time-dependent situations, we would need to con- sider the time-dependent Breit equation rather than the time- independent eigenvalue problem studied here.

Next, to evaluate the Gaunt contributionH^Gto the Breit- Pauli operator Eq.共44兲, we obtain by substitutingb₁₂^G of Eq.

共22b兲into Eq.共48a兲the expression

H^G= − e²

16␲⑀0m_e²c²共␴1i␴1k␴2j␴2k␲1i␲2jr₁₂⁻¹

+␴1i␴1k␴2k␴2j␲1ir₁₂⁻¹␲2j+␴1k␴1i␴2j␴2k␲2jr₁₂⁻¹␲1i

+␴1k␴1i␴2k␴2jr₁₂⁻¹␲1i␲2j兲. 共52兲 Inserting everywhere Eq. 共7兲, we find that there are three distinct Gaunt contributions: the orbit-orbit Gaunt operator H_OO^G with no spin operators, the spin-other-orbit operator H_SoO^G with single spin operators, and the spin-spin operator H_SS^G with pairs of spin operators,

H^G=H_OO^G +H_SoO^G +H_SS^G. 共53兲 共We here use the termorbit-orbitin a loose sense, noting that it should strictly speaking be reserved for the corresponding contribution fromH^g兲. First, retaining only the spin-free part of Eq.共7兲, we obtain for the Gaunt orbit-orbit operator

H_OO^G = − e²

16␲⑀0m_e²c²␦ik␦jk共␲1i␲2jr₁₂⁻¹+␲1ir₁₂⁻¹␲2j

+␲2jr₁₂⁻¹␲1i+r₁₂⁻¹␲1i␲2j兲

= − e²

16␲⑀0m_e²c²共2␲1ir₁₂⁻¹␲2i+ 2␲2ir₁₂⁻¹␲1i

+关␲1i,关␲2i,r₁₂⁻¹兴兴兲

= − e²

8␲⑀0m_e²c²共␲1·r₁₂⁻¹␲2+␲2·r₁₂⁻¹␲1兲 + e²ប²

4⑀0m_e²c²␦共r12兲. 共54兲 Next, retaining upon substitution of Eq.共7兲only those terms in Eq.共52兲that are linear in the spin operators, we generate the spin-other-orbit operator,

H_SoO^G = − e²

16␲⑀0m_e²c²关␦iki⑀jkl␴2l共␲1i␲2jr₁₂⁻¹−␲1ir₁₂⁻¹␲2j

+␲2jr₁₂⁻¹␲1i−r₁₂⁻¹␲1i␲2j兲+P₁₂兴

= − e²

16␲⑀0m_e²c²兵i⑀jil␴2l共关␲1i,关␲2j,r₁₂⁻¹兴兴 + 2关␲2j,r₁₂⁻¹兴␲1i兲+P₁₂其

= − e²

4␲⑀0m_e²c²

冉

r₁₂³

冊

whereP₁₂indicates the presence of a second term, obtained from the first by permutation of the particle indices. Finally, retaining in Eq.共52兲all terms quadratic in the spin operators, we arrive at the spin-spin operator,

H_SS^G = e²

16␲⑀0m_e²c²⑀ikm⑀jkn␴1m␴2n共␲1i␲2jr₁₂⁻¹−␲1ir₁₂⁻¹␲2j

−␲2jr₁₂⁻¹␲1i+r₁₂⁻¹␲1i␲2j兲

= e²

16␲⑀0m_e²c²共␦ij␦mn−␦in␦jm兲␴1m␴2n关␲1i,关␲2j,r₁₂⁻¹兴兴

= e²

4␲⑀0m_e²c²

冋

−8␲

3 ␦共r12兲共s1·s₂兲

册

where we have used 共ⵜ1iⵜ1jr₁₂⁻¹兲= −共4␲/ 3兲␦ij␦共r12兲−共␦ijr₁₂²

− 3r_12ir_12j兲r₁₂⁻⁵in the final step.

It remains to consider the gauge-term contributionH^gto the Breit-Pauli Hamiltonian Eq.共44兲. Substituting Eq.共22c兲 into Eq.共48b兲, we obtain

H^g= − e²

32␲⑀0m_e²c²␴1i␴1k␴2j␴2l关␲1i␲2j共ⵜ1kⵜ2lr₁₂兲 +␲1i共ⵜ1kⵜ2jr₁₂兲␲2l+␲2j共ⵜ1iⵜ2lr₁₂兲␲1k

+共ⵜ1iⵜ2jr₁₂兲␲1k␲2l兴. 共57兲 The operator in the brackets is symmetric with respect to

permutations ofiandkand ofjandl, as seen by expanding the␲1iand␲2jdifferentiations. Unlike for the Gaunt opera- tor Eq. 共52兲, only the spin-free part of Eq.共7兲 thus contrib- utes upon substitution into Eq. 共57兲, giving the orbit-orbit operator

H_OO^g = − e²

32␲⑀0m_e²c²兵2␲1i共ⵜ1iⵜ2jr₁₂兲␲2j

+ 2␲2j共ⵜ1iⵜ2jr₁₂兲␲1i+关␲1i,关␲2j,共ⵜ1iⵜ2jr₁₂兲兴兴其.

共58兲 Next, invoking the identities ⵜ1iⵜ2jr₁₂= −␦ijr₁₂⁻¹+r_12ir_12jr₁₂⁻³ andⵜ₁²ⵜ₂²r₁₂= −8␲␦共r₁₂兲, we arrive at the following expres- sion for the orbit-orbit operator:

H_OO^g = e²

16␲⑀0m_e²c²关共␲1·r₁₂⁻¹␲2−␲1·r₁₂r₁₂⁻³r₁₂·␲2兲 +P₁₂兴− e²ប²

4⑀0m_e²c²␦共r12兲. 共59兲 Adding this operator to its Gaunt counterpart in Eq.共54兲, we obtain the full Breit orbit-orbit operator

H_OO^B =H_OO^G +H_OO^g = − e²

16␲⑀0m_e²c²关共␲1·r₁₂⁻¹␲2

+␲1·r₁₂r₁₂⁻³r₁₂·␲2兲+P₁₂兴. 共60兲 In the final Breit-Pauli Hamiltonian Eq.共44兲, the total Breit correction is given by

H^B=H^G+H^g=H_OO^B +H_SoO^G +H_SS^G, 共61兲 where the total orbit-orbit operator is given by Eq.共60兲, the spin-other-orbit operator by Eq. 共55兲, and the spin-spin op- erator by Eq.共56兲.

III. THE SECOND-QUANTIZATION HAMILTONIAN Having presented the Breit Hamiltonian for two elec- trons in Sec. II, we shall now construct a corresponding many-electron Hamiltonian, using the formalism of second quantization. Since we shall use this Hamiltonian to set up a relativistic perturbation theory in Sec. IV, we shall pay spe- cial attention to its dependence onc. For a uniform treatment of relativistic and other perturbations, we shall follow the same procedure as done for geometrical and magnetic per- turbations in Refs. 22 and23. Adopting henceforth atomic units, we shall use the dimensionless fine-structure constant

␣=e²/ 4␲⑀0បc⬇1 / 137 as our perturbation parameter.

A. The spinor and bispinor bases

In setting up the second-quantization many-electron mo- lecular Hamiltonian for relativistic perturbation theory, we shall assume that, in a prior nonrelativistic molecular calcu- lation, we have obtained a set ofnone-electron functions or molecular orbitals共MOs兲␸p共r兲, which are taken to be ortho- normal

具␸p兩␸q典=␦pq. 共62兲

From each such MO, we may generate two independent spinors

␸P_␣共r兲=

冉

冊

冉

冊

We here denote MOs by lowercase subscripts ␸p共r兲 and spinors by uppercase subscripts␸P共r兲, sometimes with a spin label attached as in␸P_␣共r兲. Likewise, we may from each MO generate four independent bispinors

␸_Pe␣共r兲=

冢

冣

冢

冣

共64兲

␸P_p␣共r兲=

冢

冣

冢

冣

where the first two are electronic bispinors and the last two positronic bispinors. For bispinors, we use uppercase calli- graphic subscripts ␸_P共r兲, to which we sometimes attach a further label as in␸_Pe␣共r兲or␸_Pe共r兲. In the same manner, we may from each spinor␸P共r兲 generate two bispinors

␸_Pe共r兲=

冉

冊

冉

冊

representing electronic and positronic states, respectively. In passing, we note that, if the electronic bispinors ␸P_e共r兲 are taken to be the electronic eigenfunctions of the nonrelativis- tic eigenvalue problem Eq.共16兲共suitably generalized to bis- pinor form兲, then the corresponding positronic spinors␸P_p共r兲 are not the positronic eigenfunctions of Eq. 共16兲.

It is often convenient to represent spinors in terms of spin orbitals. A given spinor

␸P共r兲=

冉

冊

is then written as a linear combination of alpha and beta spin orbitals

␸P共x兲=␸P_␣共r兲␣共ms兲+␸P_␤共r兲␤共ms兲, 共67兲 wherex is a composite set of the spatial coordinatesr and the dichotomous spin coordinatems, which may take on the discrete valuesms= −¹₂ andms=¹₂. The spin space is spanned by the two functions ␣共ms兲 and␤共ms兲. The functional form of the spin functions is given by the equations

␣

共

兲

共

兲

␤

共

兲

共

兲

and they are therefore orthogonal