# Merging Features from GreenÕs Functions and Time Dependent Density Functional Theory: A Route to the Description of Correlated Materials out of Equilibrium?

###### Abstract

We propose a description of nonequilibrium systems via a simple protocol that combines exchange-correlation potentials from density functional theory with self-energies of many-body perturbation theory. The approach, aimed to avoid double counting of interactions, is tested against exact results in Hubbard-type systems, with respect to interaction strength, perturbation speed and inhomogeneity, and system dimensionality and size. In many regimes, we find significant improvement over adiabatic time dependent density functional theory or second Born nonequilibrium GreenÕs function approximations. We briefly discuss the reasons for the residual discrepancies, and directions for future work.

###### pacs:

71.10.Fd, 71.27.+a, 31.70.Hq, 71.15.MbCurrent address: ] Department of Physics, Nanoscience Center P.O.Box 35 FI-40014 University of Jyväskylä, Finland.

Hybrid methods are a valuable option in physics, to merge concepts and perspectives into a more general and effective level of description. This work adds an item from condensed matter physics to the list; we propose a hybrid method which combines non-perturbative exchange-correlation (XC) potentials from Time Dependent Density Functional Theory (TDDFT) GR84 (); UllrichBook (); Botti () with many-body perturbative self-energy schemes from Non-Equilibrium Green’s Functions (NEGF) KBE (); Keldysh (); StefLeeu (); BalzBon (), to deal with systems with strong electronic correlations and out of equilibrium.

An accurate first-principles description of the real-time dynamics of systems with strong electronic correlations is an important, difficult and basically unsolved problem of condensed matter research. General frameworks like TDDFT and NEGF do indeed allow for an in-principle-exact treatment of strong electronic correlations. However, they both rely on key ingredients [the exchange-correlation (XC) potential for TDDFT and the self-energy for the NEGF] that in general are only approximately known. For TDDFT, a systematic and general way to include non-local, non-adiabatic effects in the XC potential is lacking, while for NEGF a main hindrance is that self-energies based on many-body perturbation theory, already computationally demanding, are usually inadequate for strong electronic correlations. While considerable progress has been made for model system far away from equilibrium (see e.g. CapelleUllrich (); RubioFuksTokatly (); Maitra15 (); BalzerEckstein (); Werner (); Godby (); Romaniello ()) or for the ab initio description of near-equilibrium situations (see e.g. SharmaGross (); Turkovsky ()), a reliable first-principles treatment of the far-from-equilibrium regime is still lacking.

Here, we suggest a step towards the solution of this problem, by a novel combination of TDDFT and NEGF, where perturbative (but systematic) memory-effect corrections augment a non-perturbative local adiabatic treatment of electronic correlations. The approach is fully conserving in the Kadanoff-Baym sense conserv () and, using the so-called generalized Kadanoff-Baym ansatz Lipavsky () (see below), can be made viable for realistic systems.

Putting in practice our proposal at the ab-initio level requires access to continuum non-perturbative XC potentials, and this point is addressed at the end of the paper. However, the scope of our method can already be illustrated here using simple lattice models. This has the advantage of avoiding complex implementations and technicalities that, indispensable to deal with real-world systems, are usually unnecessary (possibly even unwanted) for an explorative assessment of a new methodology. Our results show that in many situations (see also the supplementary material, SM) the hybrid method provides significant progress over adiabatic-TDDFT and perturbative schemes for NEGF, thus holding promise for an improved treatment of the nonequilibrium dynamics of realistic correlated systems.

Models systems.- We consider small Hubbard-type 1D and 3D clusters, isolated or coupled to two 1D semi-infinite non-interacting leads. In the latter case, the cluster consists of 1 site (single impurity). These systems are exposed to time-dependent (TD) local perturbations and/or (where applicable) to electric biases in the leads. The Hamiltonian for the above setups is

(1) |

which has contributions from the cluster, the leads, and the cluster-leads couplings. In standard notation,

(2) |

where labels nearest-neighbour sites in the cluster , is the tunneling amplitude, are time-dependent on-site energies in the cluster, and are contact-interaction strengths. Further, . For the lead Hamiltonian, , where refers to the right (left) lead, and

(3) |

Here, is the (site-independent) bias in lead , the tunnelling amplitude and . The coupling between the leads and the cluster (impurity) are given by

(4) |

All energies units are expressed in terms of the hopping parameter (for the 1-site impurity cluster we use instead), and time is measured in the units of the inverse hopping parameter (assuming atomic units). We now switch to continuum variables for generality and notational convenience, and provide some elements of TDDFT and NEGF relevant to our approach.

NEGF.- The nonequilibrium propagator satisfies the equation of motion (and a similar one for ). Here, is the single-particle Hamiltonian, with kinetic energy , Hartree potential ,
and external potential . is the self-energy, which introduces a memory dependence. We integrate over the Keldysh contour KBE (); Keldysh (). is an embedding self-energy which accounts for the leads (if present),
while accounts for XC effects Myohanen2008 ().
Standard approximations for are Second-Born (2BA), T-matrix (TMA) and screened interaction (GW) StefLeeu (); BalzBon (). For real time, the lesser part of (denoted ) gives the density and the current.

TDDFT.- The time-dependent density is obtained in terms of the Kohn-Sham (KS) orbitals . These obey the KS equation , where
, and accounts for XC effects. Then, . Within a NEGF treatment, the KS density can be obtained from , with
. In practical implementations, the functional
dependence of on is often replaced by an Adiabatic Local Density Approximation (ALDA), i.e. .

A hybrid TDDFT-NEGF approach.-
Our proposal is to augment a perturbative self-energy from
a conserving many-body scheme with a non-perturbative XC potential , local in space/time.
Alternatively, this prescription
can be be seen as recasting an ALDA-TDDFT based on in a NEGF approach, but augmenting it with a non-local, non-adiabatic
perturbative self-energy .
To avoid double counting we subtract an ALDA potential
obtained from the same approximation as was used for .
The basic equation of our approach is

(5) |

To actually proceed with Eq.(Merging Features from GreenÕs Functions and Time Dependent Density Functional Theory: A Route to the Description of Correlated Materials out of Equilibrium?), at we solve for , i.e. we find , and self-consistently on the imaginary-time track; then we propagate self-consistently on the Keldysh contour, thus fulfilling the conservations laws of Kadanoff and Baym. The hybrid scheme involves no additional computational costs compared to standard NEGF time propagation. Since the augmentation is of the form of a time-local potential, our scheme can similarly be implemented in a density matrix formalism. This means that a Generalized Kadanoff-Baym Ansatz (GKBA) Lipavsky (); Bonitz (); Stefanucci () can be employed to reduce computational costs allowing for first-principles calculations of realistic systems.

The non-perturbative XC potentials- For lattice systems, depends on the system’s dimensionality. In 1D, we describe the non-perturbative, adiabatic local correlations in terms of Verdozzi08 (), and in 3D in terms of Karlsson (). is computed with the Bethe-ansatz from the 1D Hubbard model GunSchonNoack (); LimCape (), and with DMFT DMFT (); DMFT-Hubb.ref () from the 3D homogeneous Hubbard model Karlsson ().

The correction.- For concreteness, in this paper and are computed in the 2BA (some results in the TMA are also shown). The calculation and use of for Hubbard-type interactions in a NEGF time evolution has been discussed before (see e.g. Puig ()) and is not repeated here. Rather, we provide additional details of the perturbative correction . For the homogeneous (Hubbard) reference system, we use , where , and the three terms on the RHS respectively are the total energy in the second Born approximation, the non-interacting kinetic energy and the Hartree energy for the 1D (or 3D) homogeneous Hubbard model. We compute in -space:

with the retarded propagator, the statistical Fermi factor (we consider zero temperature), the single-particle energies, and . In Fig. 1 we plot for the 1D and 3D Hubbard model, for different interaction values. We also show the non-perturbative potentials used in Eq. (2). They exhibit a discontinuity at half-filling, which is always present in 1D but only for large values in 3D, reflecting the Mott-Hubbard metal-insulator transition Karlsson (). The discontinuity is absent in the 2BA. Note that, at exactly half-filling, and are both zero. Finally, from the TMA is shown. The discontinuity is absent also in this case, and at low/high filling approaches .

Closed systems: the 3D case. - We start our analysis with a 3D cubic cluster with sites, open boundary conditions, and a single interacting (and perturbed) site at the cluster centre (Fig. 2c). We compare time-dependent densities from the hybrid-approach, 2BA and ALDA, against exact results. The system is highly inhomogeneous, and despite the local character of the interaction and external perturbation, non-local effects are important: the exact (not shown) can have large nonzero components at all sites Karlsson (). Using symmetry, we map the cluster to a 10-site one (Fig. 2d), as in Karlsson (). We consider both weak (, panels a,b) and strong correlations (, panels c,d). The temporal shape of the external fields we use is Gaussian (bottom-row panels, red curves), with a slower or faster onset/offset (in the following, referred to as fast or slow perturbations). For additional time profiles we refer to the SM.

For the weakly correlated, slowly perturbed case (panel a), all approximations follow the exact solution. For the fast perturbation (panel b), non-adiabatic effects emerge, and this leads to the failure of the ALDA; the remaining approximations perform well, with the hybrid method marginally better than 2BA. In contrast, for the slow perturbation and stronger correlations (panel c), the agreement of the 2BA is poor, while the other treatments still follow the exact solution. For the most unfavourable and extreme regime of strong correlations and fast perturbations (panel d), ALDA and 2BA are largely out of phase, and only the hybrid approximation reproduces the main structures of the exact solution with the correct phase. Overall, the hybrid approximation exhibits a fairly good agreement in all regimes, and is superior to the others in the most extreme regime.

Closed systems: the 1D case.- We next consider when all sites are interacting and exposed to a space- and time-dependent perturbation. A 3D system for this situation which is also an exactly solvable benchmark is not easily accessible, due to the unfavourable scaling of the configuration space. We thus turn to a numerically more convenient 1D test-case (this also makes possible to assess the hybrid approach at low-dimensionality), choosing a 1D ring with 8 interacting sites (Fig. 3). To explore the role of space inhomogeneity, we resort to a (rather artificial) perturbation sinusoidally modulated in space: , where () and is temporal profile. The phase guarantees that the sine nodes are between sites and the amplitude at site has always the same sign. For the time profile, (step, s), (ramp, r) or (gaussian, g). Results are shown in Fig. 3 (for a more systematic study see the SM).

With highly inhomogeneous fields () no approximation reproduces the exact dynamics. Moreover for the hybrid method shows artificial density oscillations. The latter, also present in the the TDDFT-ALDA based on , are induced by the sharp discontinuity in and are not removed by the 2BA self-energy (thus, non-local, non-adiabatic effects beyond the 2BA should be also taken into account). For more homogeneous fields () the different approximations compare more favourably to the exact dynamics with superiority of the hybrid method. Looking at , the hybrid approximation is in phase with the exact curve but, for densities changing across half-filling, it still exhibits the artificial oscillations (see the SM). Further, ALDA does not perform well, and 2BA tends to be out-of-phase with the exact solution. Finally, for a slowly varying-in-space perturbation () the hybrid approach (in contrast to the other approximations) is in excellent agreement with exact results. This applies for all time profiles .

Open systems -
Finally, we test the hybrid method in open systems (Fig. 4). Specifically, using
a single-orbital Anderson impurity coupled to two 1D semi-infinite leads Thygesen () (system shown in Fig. 4e),
we consider
i) the conductance in the wide-band limit (WBL), Fig. 4a; ii) the finite-bias, finite-lead-width regime, Fig. 4b-e.
Starting with i), we find the exact density (and thereby the exact linear conductance via the Friedel sum rule) in the WBL Wiegmann (); Burke (); StefanucciKurth (); Evers (). Fig. 4a displays for the absolute deviation from the exact as function of and for different approximate treatments. We consider stronger correlations (, see the plateau in the conductance); here except for , the hybrid method performs as the best compared to 2BA or ALDA, and it is significantly better in the range (symmetrical considerations apply above half-filling).

ii) Next, we consider 1D tight-binding leads (of bandwidth ). We fix a static to to be away from
the particle-hole symmetric ground state (where ). As benchmark, we use open-ended, Anderson-impurity finite chains with up to sites treated with tDMRG tDMRGref (); montangero1 (). When (panels b,d), the agreement between hybrid and tDMRG densities/currents is fairly good, especially in the transients ( and from tDMRG never fully reach a steady state within the simulation time, in contrast to hybrid, 2BA, and ALDA. ones However, for stronger correlations and lower transparency (panels c,e), the impurity density from the hybrid scheme is closest to the tDMRG one than other schemes, whilst for the currents ALDA performs best.
The unconvincing performance of the hybrid approximation for comes probably from multiple-scattering processes, neglected by 2BA.
To corroborate this conjecture we have tested the hybrid method also using the TMA, which includes such processes. In Fig. 4c and e) the TMA hybrid method shows an improvement over the ALDA and the pure TMA calculation and thus supports the conjecture (for an expanded discussion and additional results, see the SM).

Conclusions and outlook.- By merging elements of TDDFT and NEGF, we proposed a simple, easy to implement, nonequilibrium scheme aimed to improve the treatment of local non-perturbative correlation effects and, at the same time, to incorporate non-local, non-adiabatic effects. Results from Hubbard-type systems are quite encouraging. Taking a mildly optimistic stand, we can argue that our approach extends the applicability of ALDA-TDDFT and NEGF based on perturbation theory, thus providing a way forward to merge (strong) correlations and memory effects in general. On the other hand, one can certainly envisage situations where non-perturbative and non-local correlations are very important, and this is where perhaps corrections beyond the 2BA (e.g., GWA or TMA or mixed, or other) could be employed. We note that Hubbard-type systems usually are challenging benchmarks to perturbative approximations such as 2BA, TMA or GWA. The latter generally perform much better for continuum systems with long-range interactions. Thus, we speculatively suggest that our hybrid method could perform even better for realistic systems. This is where the real merits of our proposal could possibly be: Using continuum XC potentials tailored for strong correlations (obtained from e.g. the strictly correlated approachGoriGiorgi1 (); GoriGiorgi2 (); GoriGiorgi3 (), where the discontinuities in manifest in a different way) and simplifications for perturbative self-energies (such as the GKBA Lipavsky (); Bonitz (); Stefanucci ()), our approach would be a leeway to an improved first principle treatments of realistic systems in nonequilibrium when strong local electronic correlations and memory effects play a role.

###### Acknowledgements.

We wish to acknowledge M. Puig von Friesen for discussions in the early stages of this work.## References

- (1) E. Runge and E. K. U. Gross, Phys. Rev. Lett. 52, 997 (1984).
- (2) C. A. Ullrich, Time-Dependent Density-Functional Theory: Concepts and Applications (Oxford University Press, Oxford, 2012).
- (3) S. Botti, A. Schindlmayr, R. Del Sole, and L. Reining, Rep. Prog. Phys. 70 357 (2007).
- (4) L.P. Kadanoff and G. Baym, Quantum Statistical Mechanics (Benjamin, New York, 1962).
- (5) L.V. Keldysh, Sov. Phys. JETP 20, 1018 (1965).
- (6) G. Stefanucci and R. van Leeuwen, Nonequilibrium Many-Body Theory of Quantum Systems: A Modern Introduction (Cambridge University Press, Cambridge, 2013).
- (7) K. Balzer and M. Bonitz, Nonequilibrium Green’s Functions Approach to Inhomogeneous Systems, Lecture Notes in Physics Vol. 867 (Springer, Berlin, Heidelberg, 2013).
- (8) E. Orestes, K. Capelle, A. B. F. da Silva, and C. A. Ullrich, J. Chem. Phys. 127, 124101 (2007).
- (9) J. I. Fuks, M. Farzanehpour, I. V. Tokatly, H. Appel, S. Kurth, and A. Rubio, Phys. Rev. A 88, 062512 (2013).
- (10) J. I. Fuks, K. Luo, E. D. Sandoval, and N. T. Maitra, Phys. Rev. Lett. 114, 183002, (2015).
- (11) K. Balzer and M. Eckstein, Phys. Rev. B 89, 035148 (2014).
- (12) H. Aoki, N. Tsuji, M. Eckstein, M. Kollar, T. Oka, and P. Werner, Rev. Mod. Phys. 86, 779 (2014).
- (13) M. J. P. Hodgson, J. D. Ramsden, T. R. Durrant, and R. W. Godby Phys. Rev. B 90, 241107(R) (2014).
- (14) A. Stan, P. Romaniello, S. Rigamonti, L. Reining, and J.A. Berger, New J. Phys. 17, 093045 (2015).
- (15) S. Sharma, J. K. Dewhurst, A. Sanna, and E. K. U. Gross, Phys. Rev. Lett. 107, 186401 (2011).
- (16) V. Turkowski and T. S. Rahman, J. Phys.: Condens. Matter 26, 022201 (2014).
- (17) G. Baym and L. P. Kadanoff, Phys. Rev. 124, 287 (1961).
- (18) P. Lipavsky, V. Spicka, and B. Velicky, Phys. Rev. B 34, 6933 (1986).
- (19) See Supplemental Material at http://link.aps.org/ supplemental/10.1103/PhysRevLett.116.236402, for details and additional results of our systematic studies.
- (20) P. Myöhänen, A. Stan, G. Stefanucci, and R. van Leeuwen, EPL 84, 67001 (2008)
- (21) L. Hedin, Phys. Rev. 139, A796 (1965).
- (22) S. Hermanns, N. Schlünzen, and M. Bonitz, Phys. Rev. B 90, 125111 (2014).
- (23) S. Latini, E. Perfetto, A.-M. Uimonen, R. van Leeuwen, and G. Stefanucci, Phys. Rev. B 89, 075306 (2014).
- (24) C. Verdozzi, Phys. Rev. Lett. 101, 166401 (2008).
- (25) D. Karlsson, A. Privitera, and C. Verdozzi, Phys. Rev. Lett. 106, 116401 (2011).
- (26) N. A. Lima, M. F. Silva, L. N. Oliveira, and K. Capelle, Phys. Rev. Lett. 90, 146402 (2003).
- (27) K. Schönhammer, O. Gunnarsson, and R. M. Noack, Phys. Rev. B 52, 2504 (1995).
- (28) W. Metzner and D. Vollhardt, Phys. Rev. Lett. 62, 324 (1989).
- (29) A. Georges, G. Kotliar, W. Krauth, and M. J. Rozenberg, Rev. Mod. Phys. 68, 13 (1996).
- (30) M. Puig von Friesen, C. Verdozzi, and C.-O. Almbladh, Phys. Rev. Lett. 103, 176404 (2009).
- (31) This model is a frequent choice for testing perturbative NEGF approaches to quantum transport, see e.g. K. S. Thygesen and A. Rubio, Phys. Rev. B. 77, 115333 (2008).
- (32) P. B. Wiegmann and A.M. Tsvelick, J. Phys. C 16, 2281 (1983).
- (33) J. P. Bergfield, Z.-F. Liu, K. Burke, and Ch. A. Stafford, Phys. Rev. Lett. 108, 066801 (2012).
- (34) P. Tröster, P. Schmitteckert, and F. Evers Phys. Rev. B 85, 115409 (2012).
- (35) G. Stefanucci and S. Kurth, Phys. Rev. Lett. 107, 216401 (2011).
- (36) G. De Chiara, M. Rizzi, D. Rossini, and S. Montangero, J. Comput. Theor. Nanosci. 5, 1277 (2008), open source code at: http://qti.sns.it/dmrg/phome.html
- (37) F. Heidrich-Meisner, A. E. Feiguin, and E. Dagotto, Phys. Rev. B 79, 235336 (2009).
- (38) P. Gori-Giorgi, M. Seidl, and G. Vignale, Phys. Rev. Lett. 103, 166402 (2009).
- (39) F. Malet and P. Gori-Giorgi, Phys. Rev. Lett. 109, 246402 (2012).
- (40) A. Mirtschink, M. Seidl and P- Gori-Giorgi Phys. Rev. Lett. 111, 126402 (2013)