Table Of Content

Time-dependent localized Hartree-Fock potential V. U. Nazarov Research Center for Applied Sciences, Academia Sinica, Taipei 11529, Taiwan (Dated: January 22, 2013) By minimizing the difference between the left- and the right-hand sides of the many-body time-dependent Schrödinger equation with the Slater-determinant wave-function, we derive a non- adiabatic and self-interaction free time-dependent single-particle effective potential which is the generalization to the time-dependent case of the localized Hartree-Fock potential. The new potential can be efficiently used within the framework of the time-dependent density-functional theory as we demonstrate by the evaluation of the wave-vector and frequency dependent exchange kernel fx(q,ω) and theexchange shear modulus µx of thehomogeneous electron gas. 3 1 0 Time-dependent (TD) density-functional theory exchangekernelf (q,ω)andtheexchangeshearmodulus x 2 (DFT) [1] is in the perpetual search of effective single- µ [14] of the homogeneous electron gas (HEG). x n particle potentials which accurately (in ideal - exactly) Let us consider a N-electrons system with the Hamil- a map the propagation of an interacting many-body tonian J system onto that of the non-interacting one. While 0 the exchange-correlation (xc) functionals based on Hˆ(t)= −1∆ +v (r ,t) + 1 . (1) 2 (cid:20) 2 i ext i (cid:21) |r −r | the Local-Density Approximation (LDA) [2, 3] and Xi Xi<j i j ] its semi-local refinement of the Generalized Gradient h Approximation(GGA) [4–6] have provenvery successful The many-body wave-function Ψ(t) satisfies the p in the ground-state DFT, in TDDFT the usefulness of Schrödinger equation - t (semi-) localapproachesis limited due to the fundamen- n iΨ˙(t)=Hˆ(t)Ψ(t), tal spatial non-locality of the exact time-dependent xc a u functional [7]. wherethedotdenotesthetimederivative. Weareasking q Beyond LDA and GGA, the concept of the opti- the question: What is the potential v (r,t) such that [ eff mized effective potential (OEP) [8, 9] plays one of the the functional 1 key roles in the systematic non-heuristic construction v of DFT. In the ground-state case, OEP is defined as a iΨ˙ (t)−Hˆ(t)Ψ (t) 2dr ...dr (2) 5 single-particlepotentialwhichminimizesthe many-body Z (cid:12) s s (cid:12) 1 N 4 (cid:12) (cid:12) 6 HamiltonianexpectationvalueontheSlater-determinant is minimal at(cid:12)an arbitrary time t(cid:12)for the wave-function 4 wave-function. From the DFT perspective, OEP is the Ψ (t) being the Slater determinantbuilt with the single- . firsttermintheadiabaticconnectionseriesinthepowers s 1 particle orbitals ψ (r,t) which satisfy the single-particle of the interaction constant (exact exchange) [10]. The α 0 Schrödinger equation 3 generalization of the OEP to the time-dependent case 1 has been proposed [11], which, however, is impractical 1 v: for applications due to the formidable complexity of the iψ˙α(r,t)=(cid:20)−2∆+veff(r,t)(cid:21)ψα(r,t). (3) i integral equation involved. X In this work we propose and implement an alterna- Let us assume that the TD part of v (r,t) was absent ext r tive approachto the development of the time-dependent at t < 0, and we have already solved the static problem a single-particle effective potential for many-body prob- of determining v (r) at t < 0 and the corresponding eff lems. This is based on the variational principle of the orbitals ψ (r,t) at t≤0. At t=0, the time-dependence α minimizationofthedifferencebetweentheleft-andright- of the potential is switched on. Knowing the orbitals hand sides of the time-dependent Schrödinger equation, ψ (r,0), we find v (r,0) which, determining ψ˙ (r,0) α eff α which principle we had introduced almost three decades by Eq. (3), minimizes the functional (2) at t = 0. We ago [12]. By this and with no further approximations or then find ψ (r,∆t)≈ψ (r,0)+ψ˙ (r,0)∆t, where ∆t is α α α ad hocassumptions,wederiveatime-dependenteffective a small time increase. The procedure is repeated up to potentialwiththefollowingusefulproperties: (i)Itsatis- an arbitrary time t. fies the exact large-separationasymptotic condition and Inthe∆t→0limit,thisschemereads: Withfixed(but is self-interaction free; (ii) It is expressed in terms of an yetunknown)orbitalsψ (r,t),wearelookingforthepo- α equation easily solvable for both a finite and an infinite tential v (r,t) which, determining ψ˙ (r,t) by Eq. (3), eff α periodic (or homogeneous) problems; (iii) In the static minimizes the functional (2). This gives v (r,t) as a eff case, our effective potential reduces to the earlier known functionalof the orbitals,and, finally, the orbitalsthem- localized Hartree-Fock(HF) potential [13]. As an imme- selvesarefoundbytheself-consistentsolutionofEqs.(3). diateapplicationofthisapproach,wederivethedynamic Wenotethataprocedureoftheminimizationofthesame 2 functional(2)withrespecttoψ˙ asindependentlyvaried where α functions retrieves the TD HF equations [12]. N Since Ψs(t) obviously satisfies the equation n(r,t)= |ψα(r,t)|2, (9) αX=1 N 1 ρ(r,r ,t)= ψ (r,t)ψ∗(r ,t) (10) iΨ˙ (t)= − ∆ +v (r ,t) Ψ (t), (4) 1 α α 1 s (cid:20) 2 i eff i (cid:21) s αX=1 Xi are the particle density and the single-particle density- matrix, respectively [16]. Equation (8) is our main re- the functional (2) can be rewritten as sult. We note that the only difference of Eq. (8) from the earlier known equation for the static localized HF potential (see Refs. 13 and 17 for spin-neutral and spin- 2 polarized cases, respectively) is the time-dependence of 1 Z  v˜(ri,t)− |r −r | |Ψs(t)|2dr1...drN, (5) all the quantities involved. It must, however, be em- Xi Xi<j i j phasized that without the derivation from the time-   dependentvariationalprinciple,the generalizationofthe static localized HF potential to the time-dependent case wherev˜=v −v . Equatingtozerothefirstvariation eff ext by just inserting the time variable into the static equa- of Eq. (5) with respect to δv˜=δv, we find tion would have been ungrounded. Trivially, Eq. (8) reduces to the equation for the localized HF potential in the time-independent case. 1 Solution of Eq. (8) in the case of a few-body system Z  v˜(ri,t)− |r −r ||Ψs(t)|2 doesnotpresentadifficultproblem,ashavealreadybeen Xi Xi<j i j  (6) pointed out in Ref. 13 in conjunction with the time- × δv˜(ri,t)dr1...drN =0, independentcase: DuetoEq.(10),thekernel|ρ(r,r1,t)|2 Xi of the integral equation (8) is a separable function with respect to r and r′ variables. Neither presents it a problem in the case of infinite periodic systems, when the which can be rewritten using the permutational symme- equation reduces to the matrix one. We now apply the try of the wave-function as TDlocalizedHFpotentialtoobtainthewave-vectorand frequency dependent exchange kernel f (q,ω) of HEG, x which, being a fundamental quantity by itself, is also an 1 important input in the theory of optical response of a Z  v˜(ri,t)− |r −r ||Ψs(t)|2δv˜(r1,t)dr1...drN=0, weakly ingomogeneous interacting electron gas [18]. Xi Xi<j i j  DynamicexchangekernelofHEG–InthecaseofHEG and a weak externally applied potential δv (r,t) = ext δv (q,ω)ei(q·r−ωt), we linearize Eq. (8) with respect to and, due to the arbitrariness of δv˜, ext the latter potential. The zero-order orbitals are plane- waves[19], andto the zerothandfirstorderswe havefor the density-matrix 1 Z Xi v˜(ri,t)−Xi<j |ri−rj||Ψs(t)|2dr2...drN =0. ρ0(r,r1)= 2 f(ǫk)eik·(r−r1),   (7) V Xk Z Straightforward but rather lengthy transformations car- δρ(r,r ,ω)=δv (q,ω)eiq·r 2 f(ǫk)−f(ǫk+q) eik·(r−r1), ried out in Ref. [15] lead from Eq. (7) to the following 1 s V ω−ǫ +ǫ +iη Xk k+q k equation for the exchange potential v = v˜−v , where x H (11) v (r) is the Hartree potential, H whereǫ =k2/2arefree-particleeigenenergies,f(ǫ )are k k their occupation numbers, δv (q,ω) is the perturbation s 1 n(r,t)v (r,t)= v (r ,t)− |ρ(r,r ,t)|2dr of the Kohn-Sham (KS) potential, V is the normaliza- x Z (cid:20) x 1 |r−r1|(cid:21) 1 1 tion volume, and η is an infinitesimal positive. After + ρ(r,r1,t)ρ(r1,r2,t)ρ(r2,r,t)dr dr , linearization, Eq. (8) yields Z |r −r | 1 2 1 2 (8) 3 δv (q,ω) 4π 1 1 f(ǫ )B(q,k)−f(ǫ )[f(ǫ )−1]C(k) x = − k k+q k dk, (12) δvs(q,ω) (2π)3 Z (cid:20)ω−ǫk+q+ǫk+iη ω+ǫk+q−ǫk+iη(cid:21) n0−A(q) where 0 2 A(q)= f(ǫ )f(ǫ )dk, (13) (2π)3 Z k k+q ) rs = 2 u. B(q,k)= (2π2)3 Z f(ǫ|kk1)−f(kǫ1k|12+q)dk1, (14) )(a.-0.05 rs = 5 ω C(k)= 2 f(ǫk1) dk . (15) 2, (2π)3 Z |k−k |2 1 / 1 F k -0.1 ( With the use of Eq. (12), the exchange kernel is now hx Re f f x found as Im f x fxh(q,ω)≡ δδvnx((qq,,ωω)) = -0.15 0 0.25 0.5 0.75 1 (16) ω/ωpl δvx(q,ω)δvs(q,ω) = δvx(q,ω) χh −1(q,ω), δvs(q,ω) δn(q,ω) δvs(q,ω) (cid:0) s(cid:1) FIG.1. (coloronline)Exchangekernelfxh(q,ω),q=0.5×kF, ofHEGofrs =2(blackcurvesonline)andrs=5(redcurves whereχhs(q,ω)isLindharddensity-responsefunction[20] online). Solidanddashedcurvesarerealandimaginaryparts of fh. x 2 f(ǫ )−f(ǫ ) χh(q,ω)= k k+q . (17) s V ω−ǫ +ǫ +iη Xk k+q k 80 Re ǫ In Ref. [15], we evaluate integrals (13) and (15) analyti- Im ǫ cally and reduce the integral (14) to a single-fold one. 60 It can be seen from Eqs. (12) and (16) that Imfx is ω) nonzero inside the single particle-hole excitation contin- 2,40 ǫh uum only. This deficiency of f derived from TD local- / x F ǫh ized HF potential is not limited to HEG, but, as can be k s ( easily seen, persists for any extended (periodic) system. hǫ20 Therefore, such subtle effect as the high-frequency tail of Imf of HEG [21] cannotbe accountedfor within the x 0 present approach. Instead, we will now show that f de- x -10 rived from the TD localized HF potential significantly 0 0.25 0.5 0.75 1 corrects the Lindhard dielectric function of HEG within ω/ω pl the single particle-hole continuum. InFig.1,theexchangekernelfxh(q,ω)obtainedbythe FIG. 2. (color online) Dielectric function ǫh(q,ω), q = 0.5× use of Eqs. (12) and (16) is plotted at q = 0.5×kF for kF, of HEG of rs =5 evaluated by the use of EQ. (18) with r = 2 and 5. With the inclusion of f , the dielectric the exchange kernel fh included (red curves online) and its s x x function of HEG can be written as Lindhardcounterpartǫhs(q,ω)(blackcurvesonline). Solidand dashed curves are real and imaginary parts of the dielectric 4π χh(q,ω) function, respectively. ǫh(q,ω)=1− s . (18) q2 1+χh(q,ω)fh(q,ω) s x Thefollowinglimiting casescanbe furtherworkedout In Fig. 2, the dielectric function of HEG of r = 5 s from Eqs. (12)-(15) [23] obtained through Eq. (18) is plotted together with the Lindhard dielectric function. From this we judge that 3π dynamic exchange plays a significantrole at this density limfh(q,ω 6=0)=− , (19) q→0 x 4k2 and can hardly be considered as a weak perturbation. F π In Fig. 3, the static exchange kernel fh(q) is plotted lim lim fh(q,ω)=− , (20) x q→0ω→0 x k2 for HEG of r = 5 as a function of the wave-vector. We F s 2π finda qualitativeagreementwiththe MonteCarlo(MC) lim lim fh(q,ω)=− , (21) simulations of Ref. 22. q→∞ω→0 x q2 4 iscomparativelyeasyforevaluation,whichisin contrast 0 toearlierknownTDoptimizedeffectivepotential. These present properties open a way to efficiently use this potential MC −2π/q2 within the context of time-dependent density-functional ) -10 theory, as we demonstrate by the derivation of the ex- u. a. change kernel fx(q,ω) of the homogeneous electron gas. ( ) I acknowledgesupport from NationalScience Council, q ( hfx−π-/2k0F2 Taiwan, Grant No. 100-2112-M-001-025-MY3. -30 0 1 2 3 4 5 [1] A.ZangwillandP.Soven,Phys. Rev.A 21, 1561 (1980); q/k E. Runge and E. K. U. Gross, Phys. Rev. Lett. 52, 997 F (1984); E. K. U. Gross and W. Kohn, ibid. 55, 2850 (1985). FIG.3. (coloronline)Staticexchangekernelfh(q)ofHEGof x [2] W.KohnandL.J.Sham,Phys.Rev.140,A1133(1965). rs =5. Solid line (red online) is the present result. Symbols [3] J. P. Perdew and Y. Wang, with error bars are fh(q) by MC simulations from Ref. [22]. xc Phys. Rev.B 45, 13244 (1992). Thedashed line(blueonline) shows theasymptoticbehavior [4] J.P.PerdewandW.Yue,Phys. Rev.B 33, 8800 (1986). at large q, as stipulated byEq. (21). [5] J. P. Perdew, K. Burke, and Y. Wang, Phys. Rev.B 54, 16533 (1996). [6] J. P. Perdew, K. Burke, and M. Ernz- and since [14] erhof, Phys. Rev.Lett. 77, 3865 (1996); Phys. Rev.Lett. 78, 1396 (1997). 4µ lim lim−lim lim fh(q,ω)= xc, (22) [7] G. Vignale, Physics Letters A 209, 206 (1995). (cid:20)ω→0q→0 q→0ω→0(cid:21) x 3 n2 [8] R. T. Sharp and G. K. Horton, Phys. Rev.90, 317 (1953). where µ is the exchangecontributionto the shear mod- x [9] J. D.Talman and W. F. Shadwick,Phys. Rev.A 14, 36 ulus of HEG, we can write by virtue of Eqs. (19), (20), (1976). and (22) [10] A. Görling and M. Levy,Phys. Rev.A 50, 196 (1994). [11] C. A. Ullrich, U. J. Gossmann, and E. K. U. Gross, 3 3 2/3 0.0091 Phys. Rev.Lett. 74, 872 (1995). µ = 2 ≈ , (23) x 64(cid:0)π5(cid:1)/3r4 r4 [12] V.U.Nazarov,Mathematical Proceedings of the Cambridge Philosophical s s [13] F. Della Sala and A. Gorling, The Journal of Chemical which is smaller than the high-density result of Ref. 24: Physics 115, 5718 (2001). µ = 3(32)2/3 ≈ 0,0146. In Table I, we compare µ [14] Z.QianandG.Vignale,Phys.Rev.B65,235121(2002). xc 40π5/3r4 r4 x [15] See EPAPS Document No ... s s obtained via Eq. (23) with µxc of Refs. 14, 25, and 24. [16] InEq.(8),theintegrationoverthespacecoordinatesalso implies thesummation overspin indices. [17] Z. Zhou and S.-I.Chu, Phys.Rev.A 71, 022513 (2005). TABLE I. Exchange shear modulus µx in units of 2ωpln. [18] V.U.Nazarov,G.Vignale, andY.-C.Chang,Phys.Rev. Presentresultsareshowntogetherwiththoseforµxc ofRefs. Lett. 102, 113001 (2009). 14 (QV),25 (NCT), and 24 (CV). [19] Cautionmustbeexercisedwhensolvingtheground-state problem for HEG with Eq. (8): Although the ground- rs 1 2 3 4 5 state effective potential is constant, it occurs to be infi- present 0.01102 0.01559 0.01909 0.02204 0.02464 nitefortheinfinitesystem.Aproperlimitingprocedure, QV 0.00738 0.00770 0.00801 0.00837 0.00851 however, of starting from a finite volume with periodic NCT 0.0064 0.052 0.0037 0.0020 0.0002 boundary conditions, resolves this difficulty unambigu- ously. CV 0.01763 0.02494 0.03054 0.03527 0.03943 [20] J. Lindhard, K. Dan. Vidensk. Selsk. Mat.-Fys. Medd. 28, 1 (1954). Inconclusion,within the well-definedprocedureofthe [21] K. Sturm and A. Gusarov, minimizationofthedifferencebetweentheleft-andright- Phys. Rev.B 62, 16474 (2000). hand sides of the time-dependent Shrödinger equation, [22] S. Moroni, D. M. Ceperley, and G. Senatore, wehavederivedatime-dependentsingle-particleeffective Phys. Rev.Lett. 75, 689 (1995). [23] ItmustbenotedthatinEq.(20)wehavenotbeenableto potential for a system of arbitrary number of electrons evaluate the numerator analytically, but rather, having undertheactionofatime-dependentexternalfield. This evaluated it numerically to 3.1415, surmised it to be π. potentialisnon-local,non-adiabatic,andself-interaction This,however,doesnotaffectthefollowingdiscussionat free, and it satisfies the exact large-separation asymp- all. totic condition. Atthe sametime,oureffectivepotential [24] S. Conti and G. Vignale, Phys.Rev.B 60, 7966 (1999). 5 [25] H. M. Böhm, S. Conti, and M. P. Tosi, (1997); R. Nifos`ı, S. Conti, and M. Tosi, Journal of Physics: Condensed Matter 8, 781 (1996); Physica E: Low-dimensional Systemsand Nanostructures1, 188 (1997) S. Conti, R. Nifos`ı, and M. Tosi, Jour- R. Nifos`ı, S. Conti, and M. P. Tosi, Phys. Rev. B 58, nal of Physics: Condensed Matter 9, L475 12758 (1998). SUPPLEMENTAL MATERIAL TO THE PAPER BY V. U. NAZAROV “TIME-DEPENDENT LOCALIZED HARTREE-FOCK POTENTIAL” Derivation of Eq. (8). Using the permutational symmetry, we can rewrite Eq. (7) as (N −1) (N −1)(N −2) v˜(r ,t) |Ψ (t)|2dr ...dr + (N −1)v˜(r ,t)− − |Ψ (t)|2dr ...dr =0, (S.1) 1 Z s 2 N Z (cid:20) 2 |r −r | 2|r −r | (cid:21) s 2 N 1 2 2 3 or v˜(r ,t)n(r ,t) 1 N −2 1 1 +(N −1) v˜(r ,t)− − |Ψ (t)|2dr ...dr =0, (S.2) N Z (cid:20) 2 |r −r | 2|r −r |(cid:21) s 2 N 1 2 2 3 where n(r,t) is the particle-density. Further simplifications give v˜(r ,t) 1 ρ (r ,r ;r ,r ;t) (N−2) ρ (r ,r ,r ;r ,r ,r ;t) 1 + v˜(r ,t)− 2 1 2 1 2 dr − 3 1 2 3 1 2 3 dr dr =0, (S.3) N(N −1) Z (cid:20) 2 |r1−r2|(cid:21) n(r1,t) 2 2 Z n(r1,t)|r2−r3| 2 3 where the k-particle density-matrix is ρ (r ,...r ;r′,...r′;t)= Ψ (r ...r ,r′′ ...r′′ ,t)Ψ∗(r′...r′,r′′ ...r′′,t)dr′′ ...dr′′ . (S.4) k 1 k 1 k Z s 1 k k+1 N s 1 k k+1 N k+1 N For the Slater determinant wave-function Ψ (t) the equalities hold s 1 ρ (r ,r ;r ,r ;t)= n(r ,t)n(r ,t)−|ρ(r ;r ,t)|2 , (S.5) 2 1 2 1 2 1 2 1 2 N(N−1) (cid:2) (cid:3) 1 ρ (r ,r ,r ;r ,r ,r ;t)= n(r ,t)n(r ,t)n(r ,t)−n(r ,t)|ρ(r ;r ,t)|2 −n(r ,t)|ρ(r ;r ,t)|2 3 1 2 3 1 2 3 1 2 3 2 1 3 3 1 2 N(N−1)(N−2) (cid:2) −n(r ,t)|ρ(r ;r ,t)|2+ρ(r ;r ,t)ρ(r ;r ,t)ρ(r ;r ,t)+ρ(r ;r ,t)ρ(r ;r ,t)ρ(r ;r ,t)], 1 2 3 1 2 2 3 3 1 1 3 3 2 2 1 (S.6) where ρ(r;r′,t)=Nρ (r;r′,t)= ψ (r,t)ψ∗(r′,t). (S.7) 1 α α Xα Therefore 1 |ρ(r ;r ;t)|2 n(r ,t)|ρ(r ;r ,t)|2 v˜(r ,t)−V (r ,t)− v˜(r ,t)− 1 2 dr + 2 1 3 dr dr 1 H 1 Z (cid:20) 2 |r −r |(cid:21) n(r ,t) 2 Z n(r ,t)|r −r | 2 3 1 2 1 1 2 3 (S.8) ρ(r ;r ,t)ρ(r ;r ,t)ρ(r ;r ,t) − 1 2 2 3 3 1 dr dr +C(t)=0, Z n(r ,t)|r −r | 2 3 1 2 3 where n(r′,t) V (r,t)= dr′, (S.9) H Z |r−r′| 1 n(r ,t)n(r ,t) 1 |ρ(r ;r ,t)|2 C(t)= v˜(r ,t)n(r ,t)dr − 2 3 dr dr + 2 3 dr dr =0, (S.10) Z 2 2 2 2Z |r −r | 2 3 2Z |r −r | 2 3 2 3 2 3 and the second equality in Eq. (S.10) follows from Eq. (S.8). Eqs. (S.8) and (S.10) give immediately Eq. (8). 2 Integrals (13)-(15). For integrals (13) and (15) we have straightforwardly 1 A(q)= Θ(2k −q)(2k −q)2(4k +q), (S.11) 48π2 F F F 2 C(k)= H(k,k ), (S.12) (2π)3 F where Θ(x) is the Heaviside’s step function and π p+k H(k,p)= (p2−k2)log +2kp . (S.13) k (cid:20) (cid:12)p−k(cid:12) (cid:21) (cid:12) (cid:12) (cid:12) (cid:12) (cid:12) (cid:12) For integral (14) we have 2 2 kF k2 −k2−q2 k·q B(q,k)= Θ(k −q)H(k,k −q)+ Θ(2k −q) P k,k , F 1 , dk , (S.14) (2π)3 F F (2π)2 F Z (cid:18) 1 2k q kq (cid:19) 1 1 |kF−q| where k P(k,k ,x,y)= 1 log k4+4k2k2x2+2kk (kk (2y2−1)−2xy(k2+k2))+k4+2kk x−y k2+k2 1 2k (cid:26) (cid:20)q 1 1 1 1 1 1 1 (cid:21) (S.15) (cid:0) (cid:1) −log (1−y)(k−k )2 1 (cid:0) (cid:1)(cid:9)