VU Research Portal

(1)

Quantum fluctuations and kinetic correlation in the strongly Interacting Limit of

Density Functional Theory

Grossi, J.

2020

document version

Publisher's PDF, also known as Version of record

Link to publication in VU Research Portal

citation for published version (APA)

Grossi, J. (2020). Quantum fluctuations and kinetic correlation in the strongly Interacting Limit of Density Functional Theory.

General rights

Copyright and moral rights for the publications made accessible in the public portal are retained by the authors and/or other copyright owners and it is a condition of accessing publications that users recognise and abide by the legal requirements associated with these rights. • Users may download and print one copy of any publication from the public portal for the purpose of private study or research. • You may not further distribute the material or use it for any profit-making activity or commercial gain

• You may freely distribute the URL identifying the publication in the public portal ?

Take down policy

If you believe that this document breaches copyright please contact us providing details, and we will remove access to the work immediately and investigate your claim.

E-mail address:

vuresearchportal.ub@vu.nl

(2)

8

Z E R O P O I N T O S C I L L AT I O N S I N H A R T R E E - F O C K T H E O R Y

I suppose it is tempting, if the only tool you have is a hammer, to treat everything as if it were a nail. — Abraham H. Maslow[181] 8.1 introduction

The reason for dealing with Hartree-Fock theory at this stage of a thesis about Density Functional Theory is at least twofold. On one hand Density Functional Theory has a long tradition, which dates back at least to the first hybrid functionals[182–187], of incursions in the realm of Hartree-Fock based methods in search of ingredients for improving approximations. On the other hand, an intriguing parallelism between Hartree-Fock theory and Density Functional Theory has been probed very recently [188]. It turns out in fact that , via a coupling l, one can connect the real physical system, eq. (2.1), to the Hartree-Fock system via a coupled hamiltonian ˆHHF

l (see eq. (8.8) below). Its l!• limit can be proven to be a functional of the Hartree-Fock density $HFalone:

it corresponds in fact to the minimum total electrostatic energy of N equal negative point charges in a positive background with a smooth density $HF, and it is a rigorous lower bound for W•DFT[$]1. The proof

of this relies on specific assumptions concerning the scaling of the Fock expectation value as a function of l. Since Zero Point oscillations are the common thread of this thesis, we would like to export some of the related ideas, described in the previous Chapters, to this setting.

In particular, in the first part of this Chapter we will argue that the subleading nature of the Fock exchange expectation value which makes the theory possible breaks down at lower dimensionality (d < 3). In the second part, we explore the possibility of using this system as an approximation to the ground state energy of the Jellium at low densities. Finally, we conclude the Chapter outlining the (open) problem of minimizing variationally the subleading term of the strong coupling limit of Hartree-Fock theory.

1 Throughout this chapter, an extra label (DFT\HF) will be attached to the quantities of interest, when needed, to highlight the distinction.

(3)

8.2 bird’s eye survey on adiabatic connection in hartree-fock theory

We begin by assuming that we have solved the Hartree-Fock equations for a given system (we consider here closed-shell systems for simplicity and restricted Fock), so that we are given a set of N/2 Hartree-Fock orbitals fj(r), from which we can define the Hartree-Fock density:

$HF(r) =2 N/2

Â

j=1|

fj(r)|2 (8.1)

Next, we fix the Hartree operator ˆJ ˆJ= N

Â

i=1 Z Rddr 0$HF(r0) |ri r0| | {z } :=vJ₍_r_i₎ , (8.2)

and the Exchange operator ˆK, which acts on a trial function x(r) according to ˆK= N

Â

i=1 vF(ri), vF(ri):= ˆvFx(ri):= N/2

Â

j=1 fj(ri) Z Rd f?_j(r0₎_x₍_r0₎ |ri r0| dr 0 (8.3) The Exchange operator can be written as a kernel:

ˆv(r, r0_{) =} N/2

_Â

j=1

f?_j(r)fj(r0)

|r r0_| . (8.4)

In fact, writing the one body reduced density matrix (1RDM)2

g(r, r0) = N

Â

i=1 gi(r, r0) = N

Â

i=1 Z Vdr . . . dri 1dri+1. . . Y ?₍_{r, r} i 1, ri, ri+1. . .)Y?(r, ri 1, r0i, ri+1. . .), (8.5) the Exchange expectation value can be written as

hY|ˆK|Yi =

Â

N

i=1

Z

Rd_⇥Rddrdr

0_ˆv₍_{r, r}0₎_g_i₍_{r, r}0₎_. _(8.6)

Clearly, in terms of the 1RDM the Hartree expectation reads instead hY|ˆJ|Yi = N

Â

i=1 Z Rdv J₍_r₎_g i(r, r)dr. (8.7) 2 We keep the summation over i, because it will turn useful for our calculations on classical distributions that arise in the l !• limit of the Hartree-Fock adiabatic connection.

(4)

8.2 bird’s eye survey on adiabatic connection in hartree-fock theory 107

8.2.1 The strong interacting limit in HF

Consider the following l-dependent operator : ˆ

HHF

l = ˆT+l(Vee ˆJ ˆK) + ˆVext+ ˆJ+ ˆK. (8.8) For each fixed l, its ground state YHF

l allows to compute the different contributions to the total ground state energy, EHF

l = hHˆlHFi. In fact, for l=0 we have the Hartree-Fock hamiltonian, while for l=1 we have the physical one.

We are now ready to introduce the Adiabatic connection integrand in Hartree-Fock theory. First, recall

8 < : EHF 0 :=minYhY|ˆT+ ˆVext+ ˆJ+ ˆK|Yi = EHF+UH[$HF] +Ex EHF 1 :=minYhY|ˆT+ ˆVee+ ˆVext|Yi (8.9) In Quantum Chemistry, the quantityEHFis referred to as Hartree-Fock energy:

EHF=_hYHF₀ |ˆT+ ˆVee+ ˆVext|YHF0 i. (8.10)

The Hellman-Feynman theorem allows to write

Z ₁ 0 hY HF l |ˆVee ˆJ ˆK|YHFl idl=EHF1 EHF0 =E_|HF1 _{zEHF_} :=EHF c UH[$HF] EHFx . (8.11) Rearrange the terms to write finally

EHF_c +EHF_x | {z } :=EHF xc = Z ₁ 0 ⇣ hYHF_l |ˆVee ˆJ ˆK|YlHFi +UH[$HF] +2EHFx ⌘ | {z } :=WHF l dl, (8.12) where we added on both sides EHF

x for a more direct comparison with

eq. (2.45) The quantity WHF

l has several interesting features. First, as pointed out in ref. [188], for small arguments it recovers the Møller-Plesset expansion: W_lHF⇠EHF_x + •

Â

m=2 mE_cMPmlm 1, l⌧1. (8.13) When l_!• one might expect, in analogy with the DFT case, that YHF

l minimizes the operator l ˆVee ˆJ ˆK alone.

To proceed further, we need to make a conjecture, based on physi-cal assumptions. Both ˆJ and ˆK are clearly one-body operators (their expectation value can be written in terms of the 1RDM alone), and the main difference lies in the fact that while ˆJ is local, ˆK is not (compare eqs. (8.6),(8.7)). If the electron-electron repulsion becomes infinite, we expect the off-diagonal elements of the 1RDM to go to 0 as the particles

(5)

become more and more separated. Therefore, we make the assumption that

hYHF_l | ˆK|YHF_l i =o⇣hYHF_l |ˆVee ˆJ|YHFl i ⌘

, l 1. (8.14)

eq. (8.14) implies that YHF • := lim l!•Y HF l =arg inf_YhY|_{| {z }}ˆVee ˆJ :=Hˆ_•HF |Yi. (8.15) Analogously to the DFT case, the infimum indicates that the lower bound is reached outside of the space of allowed wavefunctions. In-deed, being ˆHHF

• a classical operator, furthermore without any density

constraint, we expect YHF

l to collapse on a sum of delta peaks, centered on the minima ui of ˆH•HF3 |YHF_• (r)_|2= 1 N N

Â

i=1 d(ri ui) = _N1d(r u). (8.16)

To justify eq. (8.14), we are now going to illustrate a variational argu-ment proposed in ref. [188]. As in the DFT case, at large but finite l we expect the true ground state wavefunction YHF

l to be spread around the equilibrium positions u_i. We therefore make a variational ansatz YHF_l (r)⇠Ya(l)(r) = a(l)dNdet(A) (2p)dN !1 4 e a(2l) (r u)T·A·(r u), l 1 (8.17) where a(l) is a monotonically increasing function of l in order to have liml!•Ya(l) = YHF• . The relevant contributions to the energy

with t, h, k>0 (see Appendix Ffor a full derivation of (8.18c)). Retain-ing the leadRetain-ing term for a(l), i.e. a(l)_⇠ lq, we can write

EHF

l ⇠tlq+lE•HF+hl1 q+klq(

1 d

2 )+1, l 1. (8.19) We can minimize this expression with respect to q. The minimizing exponents q are summed up in table8.1as a function of the dimension. Notice that for d=1, one should use an effective interaction to avoid the divergence of the Coulomb integral. This compromises the simple scaling of the coulomb interaction. It is interesting to notice that the main assumption of this treatment, eq. (8.14), breaks down in d= 1

(6)

8.3 the jellium model from the adiabatic connection in hartree-fock theory 109

Table 8.1: Optimal exponent

d 1 2 3

q ⇡1 2/3 1/2

q(1 d₂ ) +1 ⇡1 2/3 1/2

and, at d = 2, hˆKil becomes already dominant with respect to the kinetic energy, complicating the analysis of the subleading terms. Therefore, from now on we shall assume d = 3, throughout this Chapter.

Finally, the high coupling limit of WHF

• is explicitly known, since we

have W_lHF⇠min Y hY|ˆVee ˆJ|Yi +UH[$ HF_] | {z } :=Eel[$HF] +2EHF_x , l 1 (8.20)

as a consequence of the assumption (8.14). The following chain of inequalities holds[188]:

W_•HFEel[$HF]W•DFT[$HF]. (8.21)

It is interesting to notice that the first two terms of eq. (8.20) define Eel[$HF], a functional of the Hartree-Fock density only. For it, there is a

straightforward interpretation: it is the minimum electrostatic energy of N point charges in a positive background with density $HF_{, plus}

the background-background repulsion, UH[$HF].

8.3 the jellium model from the adiabatic connection in hartree-fock theory

In the Jellium model[189], the detailed structure of the positive back-ground binding the electrons (typically, the atomic lattice) is replaced by a homogeneous continuum of positive charge, $+. Correspondingly,

the Hamiltonian (2.1) is replaced by ˆ Hjel= 1₂ N

Â

i=1r 2 ri+ 1 2 N

Â

i6=j 1 |ri rj| N

Â

i=1 Z V $+ |ri r0|dr 0₊1 2 Z V⇥V $2+ |r r0_|dr0dr. (8.22) The Jellium model plays an essential role in Chemistry and Quantum Physics[166, 190]. In particular, in Density Functional Theory it’s ubiquitous as main ingredient of the Local Density Approximation[6,

31,171,191]. In this approximation the density is modeled as if locally behaves as a Uniform Electron Gas, i.e. a system of electrons with a fixed, constant density. Despite their apparent difference (the Uniform Electron Gas has a constant density, while the Jellium has a constant background), very recently the equivalence between the Jellium and the

(7)

Uniform Electron Gas ground state has been proved independently in[36] and[192].

In the thermodinamic limit, the ground state energy per electron

ejel associated with ˆHjelis a function of the Bohr radius rs, defined as

the radius of a sphere containing on average one electron4

rs=✓ 3_4pn ◆1 3 , n := lim N,V!• N V. (8.23)

Upon a scaling r 7! rrs, it can be shown that the kinetic energy operator scales as ⇠r 2

s , whereas the electrostatic operator (last three terms in

ˆ

Hjel, eq.(8.22)) scales as ⇠rs1 [166]. This different scaling determines

two very different regimes, depending on whether we have rs ⌧ 1

(high density) or rs 1 (low density). When rs ! 0, we have the

particles behave as a non-interacting electron gas. For small finite rs,

the energy per electron is well approximated by[166]

e(rs)⇠ 1.11_r₂ s

0.458

rs +0.031 log rs rs⌧1 (8.24)

On the contrary, by virtue of the predominant electrostatic energy, in the low density regime rs!•, the electrons are believed to localize

in lattice points to form a Wigner crystal[189, 193]. A perturbation treatment, first proposed by Carr[18], suggests the expansion

e(rs)⇠ 0.896_r s + 1.33 r32 s rs 1 (8.25)

The coefficient of the leading term is the Madelung constant of a bcc lattice, and is obtained by minimising the electrostatic energy alone. The subleading term is obtained by assuming that the electrons perform zero-point oscillations around their equilibrium positions and performing a normal modes analysis. We now observe that

min Y ⌧ Y 1 2 N

Â

i6=j 1 |ri rj| N

Â

i=1 Z R3 $+ |ri r0|dr 0₊1 2 Z R3_⇥R3 $2+ |r r0_|dr0dr y =Eel[$+]. (8.26)

In other words, the strong interacting limit of the adiabatic connection integrand in Hartree-Fock theory, for a uniform Hartree-Fock density

$HF=$+, yields the leading term of the low density expansion of the

Jellium.

It would be interesting to see whether the subleading orders of ˆHHF

l can provide some estimate to the subleading term in eq. (8.25) too. Similar to the DFT case, this would imply to take into account the

4 In general dimension d we have rs =

✓_G (d 2+1) npd2 ◆1 d

, G being the Euler’s Gamma function.

(8)

8.3 the jellium model from the adiabatic connection in hartree-fock theory 111

kinetic energy as a perturbation. However, recall from table8.1and the subsequent considerations that in d = 3 the kinetic energy operator and the exchange operator enter to the same order, ⇠pl. Therefore

we are lead to consider a l-dependent Hamiltonian ˆ˜HHF

l , modified to take into account the background-background interaction of the Jellium5 ˆ˜HHF l = 1 2 N

Â

i r 2 ri+l 1 2 N

Â

i6=j 1 |ri rj| N

Â

i=1 Z V $+ |ri r0|dr 0₊ 1 2 Z Vdrdr 0 $2+ |r r0_| ˆK ! (8.27)

where the remaining missing terms from eq. (8.8) have been dropped as they contribute to subleading orders in the high coupling limit. What we would like to compute is

minYhY| ˆ˜HlHF|Yi N +l 0.896 rs ⇠ p l a HF ZP r32 s (8.28) and compare the value aHF

ZP with the value in expression (8.25).

8.3.1 Hartree-Fock approximation for the rs 1 regime

We use as ansatz the Hartree product Yl(r1, . . . rN) = N

’

i=1 fW_l_,i(ri), flW,i(r) = ✓ wl p ◆3 4 e wl2 (|r ui|)2_{, (8.29)}

ui being the lattice sites. This very simple model does not include

spin effects, which for arguments analogous to the ones provided in Chapter4should fall exponentially (indeed, the electrons are distin-guishable as they are labelled according to the lattice site they stick to), and neglects the effect of the anisotropy of the lattice by using spherically symmetric Gaussians, analogously to the Einsten model for phonons. The frequency wl needs to be computed variationally by minimising the energy per electron associated with eq. (8.27). The kinetic energy reads simply

hYl|ˆT|Yli = N 3

4wl. (8.30)

5 It is a well-known result (see e.g. [166]) that the divergent terms in (8.22) arising from the electron interaction cancel exactly with the ones coming from the electron-background interaction and the electron-background-electron-background in the thermodinamic limit, as a result of the charge screening operated by the background.

(9)

Focusing on the electrostatic terms, we take the result from Drum-mond and co-workers[19], which is re-derived in appendix G for completeness: ⌧ Yl 1 2 N

Â

i6=j 1 |ri rj| N

Â

i=1 Z V $+ |ri r0|dr 0₊ 1 2 Z Vdrdr 0 $2+ |r r0_| Yl = N2k 3 f 3p _G

Â

6=~0 e 2wlG2 |G|2 N r wl 2p (8.31) where we defined k_f := 3 2 2/3p₃ p rs , (8.32)

and G are the vectors of the reciprocal lattice.

Next, we move to the computation of the Exhange expectation value. Consider the HF orbitals

fHF_j (r) = p1 Ve

ikj·r_, _(8.33)

which provide the density $+via equation (8.1). We have

hYl|ˆK|Yli = N

Â

i=1 Z Rd_⇥Rddrdr 0_ˆv₍_{r, r}0₎_g_i₍_{r, r}0₎ = 1 2V N/2

Â

j=1 N

Â

i=1 Z Rd_⇥Rddrdr 0eikj·(r r 0₎ |r r0_| gi(r, r0) (8.34)

The Exchange expectation value then reads hYl|ˆK|Yli = _2V1 N

Â

i=1 ⇣ wl p ⌘3 2

Â

|k|<kf Z Vdrdr 0e wl 2 (|r ui|)2e wl2 (|r0 ui|)2eik·(r r0) |r r0_| (8.35) = N 2V ⇣ wl p ⌘3 2

Â

|k|<kf Z Vdtdt 0e wl 2 (|t|)2e wl2 (|t0|)2eik·(t t0) |t t0_| (8.36) Rotate the coordinates:

8 < : t t0 _{= ~}_t t+t0 ₌₂~_T, (8.37) hYl|ˆK|Yli = N 2V ⇣ wl p ⌘3 2

Â

|k|<kf Z Vd ~_Td_~_t e wl 2 (|~T+~t2|)2e wl2 (|~T ~t2|)2eik·~t |~t| (8.38)

(10)

8.3 the jellium model from the adiabatic connection in hartree-fock theory 113 Using e wl2 (|~T+~2t|)2e wl2 (|~T ~t2|)2 ₌e wl|~T|2e wl4 |~t|2, write hYl|ˆK|Yli = N 2V _|_k

Â

_|<_k f Z Vd~t e wl4 |~t|2eik·~t |~t| = N 2V _|_k

Â

_|<_k f Z • 0 2pdt Z 1 1dy e wl 4 t2eikty_t (8.39) = N 2V _|_k

Â

_|<_k f 4FD ⇣ k p_w l ⌘ kpwl = N 2 Z WF dk (2p)3 4FD ⇣ k p_w l ⌘ kpwl (8.40) = Nkf p_w lFD ⇣ _k_f p_w l ⌘ 2p2 (8.41)

where FD(x) denotes the Dawson integral. So the total energy per

electron reads el = h Yl|Hˆl|Yli N = 3 4wl+l 2k3 f 3p _G

Â

6=~0 e |2wlG|2 |G|2 l r wl 2p +l k_f pwlFD ⇣ _k f p_w l ⌘ 2p2 (8.42) Approximating the sum with an integral and deriving w.r.t. wl we find del dwl ⇡ 3 4 l 0 @ k 3 f 3pw2 l p_w lkf ⇣ 2k2 f +wl ⌘ FD ⇣ _k f p_w l ⌘ 4p2_w3/2 l 1 A_⇡ 3 4 l (1+p)k3_f 3p2_w2 l (8.43) Define wl:= ˜w p lto write del d ˜w ⇡ p l✓ 3 4 (1+p)k3_f 3p2_˜w2 ◆ . (8.44)

To find the minimum value of the energy, we set del

d ˜w =0, yielding ˜w= 2 p 1+p 3p k3/2f = r 1+ 1 p ✓ 1 rs ◆3/2 = 1.15 r3/2s (8.45) Then we integrate the energy derivative:

el = Z d ˜wdel d ˜w _˜w₌1.15 r3/2s =pl✓ 1+2p 4p3_r3 s ˜w+ 3 ˜w 4 ◆ ˜w=1.15 r3/2s +c=pl 1.72 r3/2s +c (8.46) We conclude that aHF

ZP ⇡1.72. This is qualitative comparable with the

result in eq. (8.25), i.e. 1.33. It is clearly higher than the value obtained in[19], i.e. 1.5, where no Exchange (which enters as a positive definite quantity in our treatment) was taken into account. In the future, we shall explore more the possibility of improving this estimate, e.g. by using the normal modes of the lattice in the ansatz (8.29). Moreover, the rigor of the justification of the computation outlined in this section is still unsatisfying in the eye of the author, and needs further pondering.

(11)

8.4 variational computation of the zero point hartree-fock energy

The variational argument used to support eq. (8.14) does not determine per se the matrix A in eq. (8.17), i.e. it is not trivial which are the frequencies that minimise the expectation value of ˆHHF

l at large finite

l. The purpose of this last section is to formulate the problem in a

more formal way. We start expanding HHF • :

HHF

• ⇡lE•+ l₂(r u)T·H· (r u), (8.47)

Hbeing the Hessian of HHF_• computed in its minimum, u.

The next term in the expansion of the energy and the corresponding minimizing wavefunction is obtained by minimizing the operator

ˆ

HZP:= ˆT+ l₂(r u)T·H· (r u) l ˆK. (8.48)

If ˆK was not present in eq (8.48), the minimizer could be written in terms of harmonic oscillators with frequencies provided by the eigenvalues of H, similarly to the DFT case. However, given the presence of the Exchange operator, a different set of oscillators could be occupied. We reconsider therefore the wavefunction (8.17):

YZPE_l (r) = 0 @l N3 2 _det₍_A₎ (2p)dN 1 A 1 4 e p2l (r u)T·A·(r u), (8.49)

A being a symmetric real 3N_⇥3N matrix. We then evaluate the

energy, making use of the expression derived in Appendix, eq.F.14: EZPE_l =hYZPE_l |HˆZP|YZPEl i =

p l 4 Tr(A)+ p l 4 Tr(A 1H) + p l N

Â

i=1 $HF_p(ui) 2 Z R3d 3_t e 1 2 ~tTsi~t |~t| , (8.50) where si are N 3⇥3 matrices which constitutes N distinct diagonal

blocks in A: A = 0 B B B @ s1 . . . . ... s2 ... ... 1 C C C A, (8.51)

The optimal matrix A can be determined, in principle, by setting[194] dEZPE

l

dA =0 (8.52)

(notice that neglecting ˆK one would get A = pH). Unfortunately,

(12)

8.4 variational computation of the zero point hartree-fock energy 115

problem. Moreover, a preliminary study on the Hydrogen atom seems to suggest that the minimizations of ˆT+l ˆVee ˆJ ˆK could lead

to an alternative set of occupied states, namely the excited states corresponding to the minimization of ˆT+l ˆVee ˆJ alone.

(13)