Theory Department
Max Planck Institute of Microstructure Physics
Theory Department   >   Research   >   Development of functional theories for static and time-dependent quantum systems

Development of functional theories for static and time-dependent quantum systems

Electronic structure via potential functional approximations

In the original form of density functional theory (DFT), suggested by Thomas and Fermi (TF) and made formally exact by Hohenberg and Kohn, the energy of a many-body quantum system is minimized directly as a functional of the density. Its modern incarnation uses the Kohn-Sham (KS) scheme, which employs the orbitals of a fictitious non-interacting system, such that only a small fraction of the total energy need be approximated. With the goal to simulate ever larger molecular and solid state systems, interest is rapidly reviving in finding an orbital-free approach to DFT. The major bottleneck in modern calculations is the solution of the KS equations, which can be avoided with a pure density functional for the kinetic energy of non-interacting fermions, Ts. Despite decades of effort, no generally applicable approximation for Ts has been found.

We consider the potential as a more natural variable to use in deriving approximations to quantum systems. In exact potential functional theory (PFT) we obtain the ground-state (gs) energy from

(1)
where is the kinetic energy operator, the inter-electronic repulsion, the one-body potential, and Ψ are N-particle wavefunctions. Defining
(2)
with Ψν denoting the gs wavefunction of potential ν(r), the exact gs energy is given by
(3)
n denoting the gs density. In practice, this requires two separate approximations, FA[ν] and nA[ν] yielding a direct approximation

In a recent letter [1] we go beyond those results by considering explicit potential functional approximations to interacting and non-interacting systems of electrons. We show that (i) the universal functional, F[ν], is determined entirely from knowledge of the density as a functional of the potential, such that only one approximation is required, namely nA[ν], (ii) the variational principle imposes a condition relating energy and density approximations, (iii) a simple condition guarantees satisfaction of the variational principle, (iv) with an orbital-free approximation to the non-interacting density as a functional of the potential, Ts is automatically determined, i.e., there is no need for a separate approximation, and (v) satisfaction of the variational principle improves accuracy of approximations.

We deduce an approximation to F from any nA[ν](r) by introducing a coupling constant in the one-body potential, νλ(r) = (1-λ)ν0(r) +λν(r), where ν0(r) is some reference potential. Via the Hellmann-Feynman theorem (and for the choice ν0=0) we obtain

(4)

where

This establishes that the universal functional is determined solely by the knowledge of the density as functional of the potential. Moreover, insertion of nA(r) on the right defines an associated approximate Fcc[nA[ν]], where cc denotes coupling constant.

Furthermore, consider the variational principle. In PFT, this yields[2]

(5)
with a possibly different value from EνA,dir for a given pair of approximations {FA, nA}. We also show that this minimization over trial potentials can be avoided, such that Fcc[nA[ν]] does yield EνA,var, if, and only if,
(6)
where
This condition is satisfied by the exact response function, but not necessarily by an approximation nA[ν]; furthermore, enforcing Eq. (6) on nA[ν] improves upon the accuracy.

Of much more practical use is the application of these results to the non-interacting electrons in a KS potential, mimicking inter-electronic interactions via the Hartree (νH) and exchange-correlation (νXC) contribution:

(7)
For a given which determines this equation can be easily solved by standard iteration techniques, bypassing the need to solve the KS equations via a given where the subscript S denotes single-particle quantities. However, once self-consistency is achieved, we need Ts to extract the total energy of the interacting electronic system. All our derivations apply equally to the non-interacting problem, so we deduce:
(8)
which is the analog of Eq. (4) for non-interacting electrons in the external potential ν(r) (which is called νs(r) when it is the KS potential of some interacting system). This defines a kinetic energy approximation determined solely by the density approximation, Tscc[nsA[v]]. This result eliminates the need for constructing approximations to the non-interacting kinetic energy, Ts.

To illustrate the power of these results, we consider a simple example of non-interacting, spinless fermions in a one-dimensional box with potential ν(x). In Fig. (1) we plot exact and approximate kinetic energy densities (KED) demonstrating that our coupling-constant approach from Eq.(8) evaluated on nsA[ν] from Ref. [3] yields both locally and globally more accurate results, and required no separate approximation for the KED.

Figure 1: Top panel: Exact (ex) KED of the coupling-constant (cc) and the direct (dir) approach (black and blue) with the corresponding semiclassical (sc) approximations (red and green) for one particle in ν(x) = -5 sin2(πx), 0 < x < 1 . Bottom panel: Absolute errors.

Potential functional approximations to the density, nsA[ν], are presently being developed via a systematic asymptotic expansion in terms of the potential, which has already been found in simple cases [3], [4]. The leading terms are local approximations of the TF type, and the leading corrections, which are relatively simple, nonlocal functionals of the potential, greatly improve over the accuracy of local approximations in a systematic and understandable way. However, a practical realization of the presented approach awaits general-purpose approximations to nsA[v](r) for an arbitrary three-dimensional case. We have shown here that it no longer awaits the corresponding kinetic energy approximations.

References

[1]
A. Cangi, D. Lee, P. Elliott, K. Burke and E. K. U. Gross, Phys. Rev. Lett. 106, 236404 (2011).
[2]
E.K.U. Gross and C.R. Proetto, J. Chem. Theory Comput. 5, 844 (2009).
[3]
P. Elliott, D. Lee, A. Cangi, and K. Burke, Phys. Rev. Lett. 100, 256406 (2008).
[4]
A. Cangi, D. Lee, P. Elliott, and K. Burke, Phys. Rev. B 81, 235128 (2010).

To top

Functionals for non-colinear magnetism

Most formulations of spin density functional theory (SDFT) restrict the magnetization vector field to have global collinearity. Nevertheless, there exists a wealth of strong non-collinearity in nature, for example molecular magnets, spin-spirals, spin-glasses and all magnets at finite temperatures.

The local spin density approximation (LSDA) can be extended to these non-collinear cases [1] but this extension has the undesirable property of having the exchange-correlation (XC) field parallel to the magnetization density at each point in space. When used in conjunction with the equation of motion for the spin magnetization in the absence of spin currents and external fields [2], [3] , this local collinearity eliminates the torsional term, resulting in no time evolution. This severe shortcoming of LSDA, where the physical prediction is qualitatively wrong, opens up an important new direction for the development of functionals where this time evolution is correctly described.

Towards this goal we have taken two approaches. The first extends the Kohn-Sham optimized effective potential (OEP) method to the non-collinear case [3], [4] , while the second is to develop a new XC functional[5] based on spin-spiral waves which goes beyond the LSDA by using gradients of the magnetization density. In both cases, we require the Kohn-Sham (KS) equation for two-component spinors Φi , which has the form of a Pauli equation. For non-interacting electrons moving in an effective scalar potential vs and a magnetic vector field Bs it reads as (atomic units are used throughout):

(1)

This equation can be derived by minimizing the total energy which, in SDFT, is given as a functional of the density

and the magnetization density

Fig. 1: Fully non-collinear magnetization density and B field obtained using the LSDA and exchange-only EXX functionals for an unsupported Cr-monolayer in Néel state. Arrows indicate the direction and information about the magnitude (in atomic units) is given in the colour bar.

For a given external scalar potential vext and magnetic field Bext this total energy reads

(2)
where is the Hartree energy. The XC potential and XC magnetic field are given by
     and      (3)
respectively. The exact functional form of Exc[ρ,m] is unknown and has to be approximated in practice.

Assuming that the densities [ρ,m] are non-interacting (v,B)-representable one may, equivalently, minimize the total-energy functional (2) over the effective scalar potential and magnetic field. Thus the conditions

and
(4)
must be satisfied.

If we were to use a functional that depends explicitly on spinor valued wavefunctions we can stay within a single global reference frame, in contrast to the case where functionals originally designed for collinear magnetism are used in a non-collinear context by introducing a local reference frame at each point in space. The most commonly used orbital functional is the EXX energy given by

(5)

where the label occ indicates that the summation runs only over occupied states. From the conditions of Eq. (4) and using the functional form given in Eq. (5) (although generalization to other orbital functionals is straightforward), the OEP equations for the non-collinear systems were derived[3, 4].

In order to explore the impact of treating non-collinear magnetism in the way outlined above and at the same time to ensure that our numerical analysis be as accurate as possible, we implemented the OEP equations for the fully non-collinear case within the full-potential linearized augmented plane wave (FP-LAPW) method implemented within the ELK code[6] This method is then applied to study the non-collinear spin magnetism in an unsupported Cr (111) monolayer. In Fig.1 are shown the magnetization density and B field calculated using both the LSDA and the OEP method.

A striking feature of the OEP B field is that, unlike its LSDA counterpart, it is not locally parallel to the magnetization density and this will produce manifestly different spin-dynamics. This is because the equation of motion for the spin magnetization reads

(6)

where Js is the spin current and γ the gyromagetic ratio. In the time-independent LSDA and conventional GGA, m ( r ) and Bxc ( r) are locally collinear, as is clear from Fig. 1, and therefore m(r) × Bxc(r) vanishes. This also holds true in the adiabatic approximation of time dependent SDFT which, by Eq. (6), implies that these functionals cannot properly describe the dynamics of the spin magnetization. In contrast, already at the static level, for the EXX functional m(r) × Bx(r) does not vanish. In fact, in the ground state of a non-collinear ferromagnet without external magnetic field, m(r) × Bxc (r) exactly cancels the divergence of the spin current,∇·Js, i.e. these terms are equally important, and it is essential to have a proper description of m(r) × Bxc (r). These results indicate that a time-dependent generalization of our method could open the way to an ab-initio description of spin dynamics. How well this functional really performs in describing the spin dynamics remains a question for future investigations.

To top

 

 

Fig. 2: Illustration of spin spiral waves along one spatial coordinate for two different choices of wavevector q=k 1,2. The angle Θ is determined by s while the angle φ=q · r. Image reproduced from the thesis of F. Eich.

Unfortunately OEP calculations are computationally demanding, leading us to develop[5] a new semi-local XC functional for non-collinear magnetism. This is analogous to the case in DFT, where going from the LDA to the generalized gradient approximations (GGAs) led to additional accuracy. To derive a new functional, our starting point was the spin-spiral wave (SSW) phase of the electron gas, which has magnetization:

(7)

where s is the spin projection on the z-axis and q is the SSW wavevector, it is illustrated in Fig. 2. Using quantum monte-carlo or the random phase approximation (RPA), the XC energy for such a state can be calculated and parameterized. To then construct a functional, we will use gradients of the magnetization to create effective s and q variables which can then be inserted into our parameterization. One particular choice is to use quantities:

Fig. 3: The magnitude (heat map) and direction (arrows) of the magnetization density (left) and XC magnetic field (right) of the Cr monolayer.
(8)
giving effective s and q:
(9)

The functional can then be written as:

(10)
where is the non-collinear LSDA functional discussed earlier. Sxc is the enhancement factor based on the parameterization of the SSW energies, and is designed so that the total functional reduces to LSDA in the appropriate limit. It can also be shown that this functional respects the zero-torque theorem.

To test this functional, it was implemented in ELK [6] and applied to the same case studied by OEP EXX above. In Fig. 3, we plot the self-consistent magnetization density m(r) and the XC magnetic field BXC(r) for a 2d slice in the plane of the Cr monolayer. To demonstrate that this functional does lead to locally non-collinear terms, in Fig.4 we plot the z-component of m(r) × BXC(r) in the region around an atom of Chromium. Unlike the LSDA case and similar to the OEP case, this quantity is non-zero implying non-zero spin currents in the ground-state. For TDDFT this functional shows much promise as it contains non-zero local torques on the magnetization density (missing in LSDA) while still a semi-local functional and hence computationally less demanding that OEP.

Fig. 4: Heat map of the z-component of

References

[1]
J. Kuebler, K.-H. Hoeck, J. Sticht and A. R. Williams, J. Phys. F18, 469 (1993).
[2]
K. Capelle, G. Vignale and B. L. Gyoerffy, Phys. Rev. Lett.87, 206403 (2001).
[3]
S. Sharma, J. K. Dewhurst, C. Ambrosch-Draxl, S. Kurth, N. Helbig, S. Pittalis, S. Shallcross, L. Nordstroem and E.K.U. Gross Phys. Rev. Lett.98, 196405 (2007)
[4]
S. Sharma, S. Pittalis, S. Kurth, S. Shallcross, J. K. Dewhurst and E.K.U. Gross Phys. Rev. B76, 100401 (Rapid Comm.) (2007)
[5]
F. G. Eich and E. K. U. Gross, arXiv:1212.3658 [cond-mat.str-el] (2013).
[6]
http://elk.sourceforge.net

To top

Basic rules for functional construction at finite temperature

Because of the small mass ratio between electrons and nuclei, standard electronic structure calculations treat the former as being in their ground state, but routinely account for the finite temperature of the latter, as in ab initio molecular dynamics[1]. But as electronic structure methods are applied in ever more esoteric areas, the need to account for the finite temperature of electrons increases. Phenomena where such effects play a role include rapid heating of solids via strong laser fields [2], dynamo effects in giant planets [3] , magnetic [4] and superconducting phase transitions [5], , shock waves [7] , warm dense matter [8] , and hot plasmas [9] .

Within density functional theory, the natural framework for treating such effects was created by Mermin [10] . Application of that work to the Kohn-Sham (KS) scheme at finite temperature also yields a natural approximation: treat KS electrons at finite temperature, but use ground-state exchange-correlation (XC) functionals. This works well in recent calculations [7] , [8] , where inclusion of such effects is crucial for accurate prediction. This assumes that finite-temperature effects on exchange-correlation are negligible relative to the KS contributions, which may not always be true.

In the present work we establish the basic rules that will allow the building of finite-temperature exchange-correlation functionals neyond these standard approximations [11] .

Central to the thermodynamic description of many-electron systems is the grand-canonical potential, defined as the statistical average of the grand-canonical operator

(1)
where Ĥ, Ŝ, ^N, τ and μ are the Hamiltonian, entropy, and particle-number operators, temperature and chemical potential, respectively. In detail,
where ^T and ^Vee are the kinetic energy and the Coulomb electron-electron interaction operators, and ^V represents an external scalar potential v(r). The entropy operator is given by
where k is the Boltzmann constant and
is a statistical operator, with | Ψ N,i⟩ and wN,i being orthonormal N-particle states and statistical weights, respectively, with the latter satisfying the (normalization) condition ΣN,iwN,i = 1. The statistical average of an operator ^A is obtained as
(2)
To create a DFT at finite temperature, Mermin [10] rewrites this as (in modern parlance)
(3)
where the minimizing n(r) is the equilibrium density n0(r), and
(4)
is the finite-temperature analog of the universal Hohenberg-Kohn functional, defined through a constrained search [15]. This depends only on τ and not on μ. We denote as the minimizing statistical operator in Eq. (4), and define the density functionals:

i.e., each density functional is the trace of its operator over the minimizing ^Γ for the given τ and density. Because it arises so often in this work, we have defined the "kentropy" Kτ[n]. Next consider a system of non-interacting electrons at the same temperature τ and denote its one-body potential as vs(r).
All the previous considerations apply, and we choose vS(r) to make its density match that of the interacting problem. This defines the KS system at finite temperature. The corresponding functionals will be labeled with the subscritp s.

We have found that the exchange (x) and correlation (c) parts of the grand-canonical potential, as well as the kentropy and the potential energy (U) are bound to fulfil the exact inequalities [11] :

(5)
and no approximation should violate these basic rules.

Some of the most important results in ground-state DFT come from uniform scaling of the coordinates [17] . In the following considerations, when we refer explicitly to wavefunctions, we shall restrict to wavefunctions having finite norm on their entire domain of definition. Under norm-preserving homogeneous scaling of the coordinate r→γr, with γ > 0, to the scaled wave function [17] .

(6)
corresponds the scaled density nγ(r) = γ3n(γr). A simple 1D sketch of the scaling procedure is given in this figure:
Writing Ψγ(r1,....rN) = ⟨r,...,rN|Ψγ⟩ in terms of the (representation-free) element |Ψγ⟩ of Hilbert space, the scaled statistical operator is defined as [11]
(7)
where the statistical weights are hold fixed, i.e., the scaling only acts on the states.

With the above definition, the statistical average of an operator whose pure-state expectation value scales homogeneously , scales homogeneously [17] as well. In particular, we have:

The scaling behavior of the density functionals is, however, more subtle. First consider the non-interacting functionals in some detail. We observe that
(8)
which implies
(9)
For non-interacting electrons, the statistical operator at a given temperature that is the minimizer for a given scaled density is simply the scaled statistical operator, but at a quadratically scaled temperature, an effect that is obviously absent in the ground-state theory.
 Next, we consider the interacting case. The analysis of the scaling relations allows to derive a set of inequalities that provide tight constraints on the functionals, and are used in non-empirical functional construction in the ground state. We obtain the following:
(10)
(11)
(12)
which are reversed if γ < 1.

Lastly, we consider the adiabatic coupling constant for finite temperature, its relationship to scaling, and derive the adiabatic connection formula. Define

(13)
with being the corresponding minimizing .

Of course, non-interacting functionals are not affected by a coupling constant modification. Eq. (10) implies that the exchange and Hartree density functionals have a linear dependence on λ. Employing the minimization property of Eq. (13) and the Hellmann-Feynman theorem, we find [11]

(14)
where
(15)

Eq.(14) is the finite-temperature adiabatic connection formula, whose zero-temperature limit played a central role in ground-state DFT.

(Eq. 5), and the scaling inequalities can be combined, analogously to Ref. [17] to show that Uxcτ[n](λ) is monotonically decreasing in λ , as shown in the schematic drawing below.
Fig. 1. Geometrical interpretation of the adiabatic connection formula. The shaded area between the curve and the λ axis represents Ωxcτ[n]. All magnitudes are negative, in agreement with the inequalities in Eq. (5).

So far, all results presented have been exact. To see them in practice, consider the finite-temperature local density approximation (LDA) to Ωxcτ[n]

(16)
where ωxcunifτ(n) is the XC grand canonical potential density of a uniform electron gas of density n. Because a uniform electron gas is a quantum mechanical system, its energies satisfy all our conditions, guaranteeing by construction that LDA satisfies all the exact conditions listed here. In the Jacob's ladder of functional construction [13] , more sophisticated approximations should also satisfy these conditions. To give one simple example of the usefullness of our results, Eq. (10) implies
(17)

where

is a dimensionless measure of the local temperature. Thus the largest fractional deviations from ground-state results should occur (in LDA) in regions of lowest density, but these contribute less in absolute terms.

In summary, there is a present lack of approximate density functionals for finite temperature. We have derived many basic relations needed to construct such approximations, and expect future approximations to either build these in, or be tested against them. In principle, such approximations should already be implemented in high-temperature DFT calculations, at least at the LDA level, as a check that XC corrections due to finite temperature do not alter calculated results. If they do, then more accurate approximations than LDA will be needed to account for them.

To top

References

[1]
R. Car and M. Parrinello, Phys. Rev. Lett. 55, 2471 (1985).
[2]
J. Gavnholt, A. Rubio, T. Olsen, K. S. Thygesen, and J. Schiøtz, Phys. Rev. B 79, 195405 (2009).
[3]
R. Redmer, T. R. Mattsson, N. Nettelmann, and M. French, Icarus 21, 798 (2011).
[4]
N. M. Rosengaard and B. Johansson, Phys. Rev. B 55, 14975 (1997).
[5]
G. Profeta, C. Franchini, N. N. Lathiotakis, A. Floris, A. Sanna, M. A. L. Marques, M. Lüders, S. Massidda, E. K. U. Gross, and A. Continenza, Phys. Rev. Lett. 96, 047003 (2006).
[6]
P. Cudazzo, G. Profeta, A. Sanna, A. Floris, A. Continenza, S. Massidda, and E. K. U. Gross, Phys. Rev. Lett. 100, 257001 (2008).
[7]
S. Root, R. J. Magyar, J. H. Carpenter, D. L. Hanson, and T. R. Mattsson, Phys. Rev. Lett. 105, 085501 (2010).
[8]
A. Kietzmann, R. Redmer, M. P. Desjarlais, and T. R. Mattsson, Phys. Rev. Lett. 101, 070401 (2008).
[9]
M. W. C. Dharma-wardana and M. S. Murillo, Phys. Rev. E 77, 026401 (2008).
[10]
N. D. Mermin, Phys. Rev. 137, A1441 (1965).
[11]
S. Pittalis, C. R. Proetto, A. Floris, A. Sanna, C. Bersier, K. Burke, and E. K. U. Gross, Phys. Rev. Lett. **** , ****** (2011).
[12]
R. G. Dandrea, N. W. Ashcroft, and A. E. Carlsson, Phys. Rev. B 34, 2097 (1986).
[13]
C. Fiolhais, F. Nogueira, and M. A. Marques, eds., A Primer in Density Functional Theory (Springer, 2003).
[14]
D. C. Langreth and J. P. Perdew, Solid State Commun. 31, 567 (1979).
[15]
R. Parr and X. Yang, Density-functional Theory of Atoms and Molecules (Oxford University Press, 1989).
[16]
M. Greiner, P. Carrier, and A. Görling, Phys. Rev. B 81, 155119 (2010).
[17]
M. Levy and J. P. Perdew, Phys. Rev. A 32, 2010 (1985).

To top

Correlation potentials for molecular bond dissociation within the self-consistent random phase approximation

The adiabatic connection fluctuation dissipation (ACFD) formula within time-dependent density functional theory (TDDFT) is an exact expression for the correlation energy within DFT. The first approximation, obtained by setting the exchange-correlation (XC) kernel to zero, yields the so-called random phase approximation (RPA) which incorporates long-range correlation effects such as the van der Waals forces. It also removes static correlation errors or fractional spin errors when breaking chemical bonds.

For a unique solution the RPA should be solved self-consistently, which involves the construction of the corresponding RPA correlation potential. In this work we have calculated self-consistent RPA potentials for diatomic molecules and investigated the RPA correlation potentials at different inter-atomic distances. By increasing the inter-atomic distance correlation effects become more important and the exact correlation potential is known to exhibit strong correlation features such as a sharp step and a peak at the bond-midpoint.

Within the ACFD framework the exact correlation energy is written as

(1)
where χs is the non-interacting KS density response function and χλ is the scaled interacting density response function. The scaling refers to a system with a linearly scaled Coulomb interaction λv(r,r') plus a fictitious potential which keeps the density fixed as λ is changed. The parameter λ runs between 0 (non-interacting KS case) and 1 (fully interacting case). We have used the short hand notation
for any two-point functions f and g. Within TDDFT the function χλ reads
(2)
The scaled XC kernel fxcλ is defined as the functional derivative of the scaled XC potential vxcλ with respect to the density n, evaluated at the ground state density.

In the RPA, fxcλ = 0 and thus corresponds to the simplest approximation within the ACFD formalism. Within RPA the λ-integral in Eq. (1) can be evaluated analytically with the result

(3)
Diagrammatically Eq. (3) is equal to an infinite summation of ring-diagrams.

The RPA correlation potential vc can be obtained as the functional derivative of Eq. 3 with respect to the density. If we let Vs signify the total KS potential, Gs the non-interacting KS Green function and χs= -iGsGs, the functional derivative can be obtained via the chain rule

(4)
The result is the well-known linearized Sham-Schlüter (LSS) equation
(5)
Here, we have used the notation (r1,t1) = 1 etc. and introduced Λ(3,2;1)=-iGs(3,1)Gs(1,2). The correlation part of the self-energy Σc in the RPA is given by
(6)

where
(7)

In Fig. 1 we show the RPA correlation potential of 1D LiH (ZLi=3.6, ZH =1.2). The solid, dotted and dashed curves represent the correlation potentials for 2, 3, and 8 Bohr interatomic distances, respectively. The build up of the peak at the bond midpoint is apparent. However, a step due to the derivative discontinuity, as is present in the exact KS correlation potential, is not observed.

Fig. 1: The correlation potential along the bond with the bond midpoint at zero. The system is 1D LiH with a soft coulomb potential. The Li atom is located at -1, -1.5, and -4 Bohr, respectively. The H atom is located at 1, 1.5, and 4 Bohr, respectively.
Fig. 2 shows the correlation potential of the three dimensional LiH for bond distances 3 (solid), 6 (dotted), and 8 Bohr (dashed). The same features as in the 1D case are also found in the results for the 3D correlation potentials. In the region of the Li atom (-10 to 0) the potential qualitatively resembles that of the 1D system in Fig. 1. We see a well with a peak close to the nucleus. Also in the region of the H atom (0 to 10) figures 1 and 2 show the same structure. A well emerges with increasing bond distance. The difference between 1D and 3D is found only at the bond midpoint. A peak emerges only for the 1D system (Fig. 1). In contrast, the 3D system (Fig. 2) exhibits a peak at zero only for small and intermediate bond distances (3 and 5 Bohr). The orbital and potential basis functions are simply not diffuse enough to extend to the bond midpoint.
Fig. 2: The correlation potential of LiH along the bond axis. The bond midpoint is at zero. The Li atom is located at -1.5, -2.5, and -4 Bohr, respectively. The H atom is located at 1.5, 2.5, and 4 Bohr, respectively.

We have further analyzed the RPA functional in the context of fractional charge and spin by calculating the total energy at different particle numbers. Fig. 3 shows the total energy as a function of the number of electrons at the H atom. The number of electrons at the Li atom will then be four minus the x value. The exact functional (black) has a minimum at 1.0, because it dissociates LiH into a neutral H atom and a neutral Li atom.

Fig. 3: The total energy of dissociated LiH as a function of number of electrons at the H atom for EXX (blue), RPA (red), and exact (black).
In Fig. 3 we see that EXX (blue) and RPA (red) do not dissociate LiH into the neutral atoms. In both cases there is a surplus of electronic charge at the H atom. However, in the case of RPA the surplus (0.16 electrons) is much smaller than in the case of EXX (0.4 electrons). The large improvement may be related to the peak that is present in the RPA correlation potential. The total energy in the dissociation limit is shown to be largely improved in the RPA showing only a small fractional spin error.

References

[1]
M. Hellgren, D.R. Rohr, E.K.U. Gross, Journal of Chemical Physics 136, 034106 (2012).

To top

Discontinuities of the exchange-correlation kernel and charge-transfer excitations in time-dependent density functional theory

Density functional theory (DFT) owes its great success to the fact that the highly intricate many-body wave function is replaced by the much simpler electronic density. In addition, the density is calculated from a system of non-interacting electrons known as the Kohn-Sham (KS) system. The effective KS potential is composed of the external, the Hartree and the so-called exchange-correlation (XC) potential vxc, where the latter is given by the functional derivative of the XC energy Exc with respect to the density, and incorporates all the effects of the Coulomb interaction beyond the Hartree level. The difficulties in using DFT lies in finding good approximations to Exc. Many works have therefore been devoted to study its exact properties. One of the most intriguing features is the existence of derivative discontinuities at integral particle numbers [1], a property which turns out to have several important physical implications.

It is not hard to see that a derivative discontinuity in Exc gives rise to a singularity in vxc in the form of a jump by a constant Δxc, i.e., vxc+(r)=vxc -(r)+Δxc, where ± refers to N→ N0±, and N0 is an integer. A classic example where the jump in vxc becomes important is in the dissociation of closed-shell molecules composed of open-shell atoms. In order to get the correct dissociation limit a step in vxc has to develop over the atom with the lower ionization energy (a similar situation but with closed-shell fragments is depicted in Fig. 1).

Fig. 1: Left: fxc and vxc of a dissociating He-Be2+ system. Right: Discontinuity of the Be2+ atom.

With the advent of time-dependent DFT (TDDFT) the variation of vxc with respect to the density, i.e., the XC kernel fxc(r,r',ω), has become another central quantity of interest. The reason is the possibility of extracting excitation energies from the linear density response function, which in TDDFT can be expressed[2] in terms of the KS response function χs, the Coulomb interaction v, and fxc

(1)

In this work we have evaluated the discontinuities of fxc in order to increase the understanding of its exact behavior [3] . The discontinuities of fxc turned out to be of a more complex nature than those of vxc , where the general structure is an r -and ω-dependent function gxc(r,ω), which may diverge. In order to obtain an exact expression for the discontinuity we study the limiting behavior of Exc as a function of particle number N around an integer N0. The constant jump in vxc can then formally be written as

, (2)

where f(r)=∂n(r)/ ∂N is the Fukui function.The same idea can then be used for the static XC kernel and we find an expression for gxc
(3)

After performing a common denominator approximation to the response functions it is easy to see that the second term exhibits a diverging behavior
which depends on the difference between the ionization energy I and the KS affinity As. This divergency can be seen in f xc of the dissociating system in Fig. 1 where the calculations have been made in the exact-exchange approximation, which is known to have a discontinuity. The right panel of Fig. 1 displays g xc, calculated in the same approximation but defined on an ensemble which allows for fractional charges. A striking similarity can be seen and we can conclude that the peaked structure in the XC kernel of a dissociating system is just the discontinuity of the kernel.

A diverging behavior of the kernel is exactly what is needed in order to describe charge-transfer (CT) excitations. The reason is that Eq. (1) involves only matrix elements of fxc between so-called excitation functions, i.e., products of occupied and unoccupied KS orbitals. As the distance between the fragments increases these products vanish exponentially and thus there is no correction to the KS eigenvalue differences unless the kernel diverges. The formulation of TDDFT for fractional charges has to allow the particle number to change in time and we have therefore proposed an ensemble of the form


where Ψk(t) is the k-particle many-body state at time t and αk(t) are given TD probabilities. A Runge-Gross theorem can be proved also in this case allowing for an analysis of the same kind as in the static case. The function gxc acquires an ω-dependence and Fig. 2 illustrates the difference compared to the static case showing mainly an increase of the exponential growth. This increase is, however, crucial for the description of CT excitations, as seen in Fig. 3, where AEXX uses fx(ω=0) and TDEXX uses fx(ω KS).
Fig. 2: vxc and fxc at ω=0 and ωKS of a dissociating He-Be system.
Fig. 3: CT excitation energies of the system in Fig. 2 in different approximations.

Apart from playing an important role for the excitation energies also the orbitals are affected by the discontinuity. Within chemical reactivity the Fukui function f(r) as defined above is an important reactivity indicator. The Fukui function can be expressed as [4]


and we see that a discontinuity of f xc does not cancel when determining f+. This contribution can be very large for certain molecules as shown in Fig. 4.
Fig. 4: The Fukui function f+ of a 1D Be2+ atom.

References

[1]
J. P. Perdew et al., PPhys. Rev. Lett. 49, 1691 (1982)
[2]
E. K. U. Gross and W. Kohn, Phys. Rev. Lett. 55, 2850 (1985)
[3]
M. Hellgren, E.K.U. Gross, Phys. Rev. A 85, 022514 (2012).
[4]
M. Hellgren, E.K.U. Gross, J. Chem. Phys. 136, 114102 (2012)

To top

Copyright © 2010 Max-Planck-Gesellschaft zur Förderung der Wissenschaften e.V., München, All Rights Reserved
created by MPI für Mikrostrukturphysik | Jobs | Research School | University of Halle | Halle (Saale) | MPG | Impressum