Kadanoff–Baym Equations for Nonequilibrium Dynamics

• Kadanoff–Baym equations are a two-time, non-Markovian framework that defines interacting quantum dynamics via contour-ordered Green’s functions with memory effects and initial correlations.
• They enable the computation of spectral evolution, transport properties, and excitation spectra using approximations like Hartree–Fock, second Born, GW, and T-matrix, while ensuring conservation laws.
• Advanced numerical strategies such as auxiliary-Hamiltonian mappings and memory truncation reduce computational complexity, extending the method’s scope to open and dissipative quantum systems.

Kadanoff–Baym equations are the nonequilibrium Dyson equations for contour-ordered single-particle Green’s functions on the Keldysh contour. They provide a two-time, non-Markovian formulation of interacting quantum dynamics in which the one-body propagator evolves under a one-body Hamiltonian and a time-nonlocal self-energy, while initial correlations are incorporated through the imaginary branch of the contour. In this form, they are a central tool of nonequilibrium Green’s-function theory for correlated dynamics, transport, spectroscopy, and open-system evolution (Stahl et al., 2021, Myöhänen et al., 2010).

1. Contour formulation and real-time structure

The basic object is the contour-ordered Green’s function

$$G(1,2)=-i\langle T_\gamma[\psi(1)\psi^\dagger(2)]\rangle,$$

with $1\equiv(r_1,\sigma_1,t_1)$, $2\equiv(r_2,\sigma_2,t_2)$, and $T_\gamma$ denoting ordering on the Keldysh contour. In compact notation, the Kadanoff–Baym equation can be written as

$$\left(i\partial_{t_1}-h(t_1)\right)G(t_1,t_2)=\delta(t_1,t_2)+\int_\gamma dt\,\Sigma(t_1,t)\,G(t,t_2),$$

together with the adjoint equation differentiated with respect to $t_2$. The same content may be expressed as the contour Dyson equation

$G=G_0+G_0\Sigma G,$

with $\Sigma=\Sigma[G]$ the many-body self-energy functional (Friesen et al., 2010, Balzer et al., 2013).

For nonequilibrium initial-value problems, the contour is typically L-shaped: forward real time, backward real time, and an imaginary branch for the initial thermal state. Projection onto contour sectors yields the standard set of components $G^{<,>}$, $G^{R/A}$, $1\equiv(r_1,\sigma_1,t_1)$0, and $1\equiv(r_1,\sigma_1,t_1)$1. The mixed and Matsubara components encode initial correlations, and KMS boundary conditions enforce consistency with the initial equilibrium density matrix. This decomposition converts the contour equation into coupled integro-differential equations for real-time propagation, but it preserves the causal memory structure: the value at time $1\equiv(r_1,\sigma_1,t_1)$2 depends on the entire earlier history through self-energy convolutions (Myöhänen et al., 2010, Stahl et al., 2021).

This two-time formulation is the defining structural feature of the Kadanoff–Baym framework. It differs from time-local kinetic equations by retaining spectral evolution, memory, and initial-correlation effects at the level of the one-particle propagator itself. That structure is also what makes the equations computationally demanding.

2. Self-energies, conserving structure, and excitation spectra

In practical calculations the Kadanoff–Baym equations are closed by a self-consistent approximation for $1\equiv(r_1,\sigma_1,t_1)$3. The papers considered here repeatedly emphasize the role of $1\equiv(r_1,\sigma_1,t_1)$4-derivable approximations, for which

$1\equiv(r_1,\sigma_1,t_1)$5

When solved self-consistently, Hartree–Fock, second Born, $1\equiv(r_1,\sigma_1,t_1)$6, and $1\equiv(r_1,\sigma_1,t_1)$7-matrix approximations are conserving in the Baym–Kadanoff sense: they enforce exact conservation of total particle number, energy, and momentum during the time evolution (Myöhänen et al., 2010, Friesen et al., 2011).

The diagrammatic content of these approximations matters. Hartree–Fock is time-local and describes only mean-field renormalization. Second Born adds the second-order bubble and exchange diagrams and is the simplest time-nonlocal approximation. $1\equiv(r_1,\sigma_1,t_1)$8 replaces bare interactions by a screened interaction $1\equiv(r_1,\sigma_1,t_1)$9 obtained from bubble resummation. The $2\equiv(r_2,\sigma_2,t_2)$0-matrix approximation resums particle-particle ladder diagrams and is therefore particularly sensitive to repeated two-particle scattering processes (Friesen et al., 2010).

A major consequence is visible in excitation spectra. Time-local Hartree–Fock can give incomplete neutral spectra, while second Born already reproduces many additional excitations characterized as double excitations (Säkkinen et al., 2011). In finite quantum wells, real-time propagation of the full two-time Kadanoff–Baym equations in second Born captures correlation-induced double-excitation features absent at the Hartree–Fock level (Balzer et al., 2012). In the language used there, propagating the Kadanoff–Baym equations with a conserving self-energy is equivalent to solving a Bethe–Salpeter equation with a correspondingly dressed kernel while fulfilling a frequency sum rule (Säkkinen et al., 2011).

This establishes a central methodological point: the Kadanoff–Baym equations are not only transport equations for one-particle observables. Through the choice of $2\equiv(r_2,\sigma_2,t_2)$1, they define a self-consistent approximation to neutral and charged excitation spectra, damping, and correlation satellites.

3. Higher-order observables and positivity constraints

Although the Kadanoff–Baym equations evolve a one-particle Green’s function, equations of motion permit extraction of certain two-particle observables. For a local Hubbard interaction, the local double occupancy

$2\equiv(r_2,\sigma_2,t_2)$2

can be written as

$2\equiv(r_2,\sigma_2,t_2)$3

This relation is exact at the formal level; the approximation enters through the chosen self-energy $2\equiv(r_2,\sigma_2,t_2)$4 (Friesen et al., 2010, Friesen et al., 2011).

That identity exposes an important limitation of conserving approximations. In a Hubbard dimer driven out of equilibrium, second Born and $2\equiv(r_2,\sigma_2,t_2)$5 can yield negative double occupancy at finite times, even though both approximations are $2\equiv(r_2,\sigma_2,t_2)$6-derivable and conserve particle number and energy. Since $2\equiv(r_2,\sigma_2,t_2)$7 is the expectation value of a positive operator, negative values are unphysical. Among the tested schemes, only the $2\equiv(r_2,\sigma_2,t_2)$8-matrix approximation yields double occupancies that remain non-negative in all studied cases (Friesen et al., 2010).

The same issue appears in entanglement calculations. For a non-magnetic spin-$2\equiv(r_2,\sigma_2,t_2)$9 system, the local single-site entanglement entropy can be reconstructed from the density $T_\gamma$0 and double occupancy $T_\gamma$1. In a transport setup with an interacting impurity, $T_\gamma$2 obtained from Kadanoff–Baym propagation within second Born or $T_\gamma$3 can again become negative, which makes the entropy formula ill-defined. In the $T_\gamma$4-matrix approximation, by contrast, the pair correlation function is proven to remain non-negative; this provides a controlled route from Kadanoff–Baym dynamics to entanglement entropy in the low-density regime (Friesen et al., 2011).

These examples sharply delimit what conservation does and does not guarantee. Conserving approximations preserve macroscopic balance laws, but they do not automatically preserve microscopic positivity constraints on pair densities or related local observables. The Kadanoff–Baym framework therefore requires approximation diagnostics beyond conservation alone.

4. Finite, open, and dissipative systems

The formalism applies to both isolated and embedded systems. In quantum transport, the central interacting region is coupled to noninteracting leads that are integrated out into an embedding self-energy,

$T_\gamma$5

The resulting embedded Kadanoff–Baym equation contains $T_\gamma$6 and consistently incorporates initial correlations, initial embedding, and arbitrary AC or DC biases. In this form, the framework yields time-dependent densities, currents, dipole moments, and nonequilibrium spectral functions for open correlated conductors (Myöhänen et al., 2010).

The same structure extends to explicitly open quantum systems coupled to thermal environments. For fermions in an attractive one-dimensional potential interacting with a bosonic heat bath through elastic $T_\gamma$7 collisions, the nonequilibrium Kadanoff–Baym equations describe equilibration, thermalization with the bath, and decoherence of off-diagonal density-matrix elements. The retarded propagator simultaneously determines the evolving spectral properties of bound and scattering states (Neidig et al., 2023). A later extension generalizes this construction to arbitrary fermionic or bosonic systems and baths, including three spatial dimensions, different traps, and different interaction potentials, and uses it to study bound-state formation, decay, regeneration, and thermodynamics in open systems (Neidig et al., 16 Jun 2025).

A different extension reformulates Keldysh theory for Lindbladian dynamics. In that setting, a contour-ordered NEGF can still be defined, but the contour Hamiltonian acquires branch-dependent non-Hermitian and cross-branch terms. The result is a generalized Kadanoff–Baym framework for interacting systems with dissipative Lindblad evolution, preserving a diagrammatic many-body structure while accommodating loss and gain processes (Stefanucci, 2024).

Together these formulations show that the Kadanoff–Baym equations are not restricted to isolated Hamiltonian dynamics. They also organize embedding, reservoirs, heat baths, and even Lindbladian dissipation within closely related two-time structures.

5. Memory kernels, computational complexity, and algorithmic reformulations

The main numerical obstacle is the memory kernel. On a time grid with $T_\gamma$8 points, the direct two-time formulation stores $T_\gamma$9 on an $$\left(i\partial_{t_1}-h(t_1)\right)G(t_1,t_2)=\delta(t_1,t_2)+\int_\gamma dt\,\Sigma(t_1,t)\,G(t,t_2),$$0 mesh, with memory scaling as $$\left(i\partial_{t_1}-h(t_1)\right)G(t_1,t_2)=\delta(t_1,t_2)+\int_\gamma dt\,\Sigma(t_1,t)\,G(t,t_2),$$1, while evaluation of the history integrals leads to a computational cost that typically scales as $$\left(i\partial_{t_1}-h(t_1)\right)G(t_1,t_2)=\delta(t_1,t_2)+\int_\gamma dt\,\Sigma(t_1,t)\,G(t,t_2),$$2 (Stahl et al., 2021, Balzer et al., 2013). This is the canonical cost barrier of full Kadanoff–Baym propagation.

Several algorithmic responses have been developed. One is an auxiliary-Hamiltonian representation for spatially local self-energies. There the nonequilibrium Dyson equation is mapped to a noninteracting auxiliary Hamiltonian with explicit bath orbitals chosen so that the bath hybridization reproduces the original self-energy kernel. The original non-Markovian dynamics is then recast as a time-local evolution in an enlarged Hilbert space, avoiding explicit memory-integral evaluation and substantially reducing storage demands in nonequilibrium DMFT applications (Balzer et al., 2013).

A second strategy is direct memory truncation. If the self-energy decays beyond a finite correlation time $$\left(i\partial_{t_1}-h(t_1)\right)G(t_1,t_2)=\delta(t_1,t_2)+\int_\gamma dt\,\Sigma(t_1,t)\,G(t,t_2),$$3, one may truncate the memory integrals to $$\left(i\partial_{t_1}-h(t_1)\right)G(t_1,t_2)=\delta(t_1,t_2)+\int_\gamma dt\,\Sigma(t_1,t)\,G(t,t_2),$$4. In the truncation scheme analyzed for the Hubbard model, the cost drops from cubic to linear in $$\left(i\partial_{t_1}-h(t_1)\right)G(t_1,t_2)=\delta(t_1,t_2)+\int_\gamma dt\,\Sigma(t_1,t)\,G(t,t_2),$$5 for fixed $$\left(i\partial_{t_1}-h(t_1)\right)G(t_1,t_2)=\delta(t_1,t_2)+\int_\gamma dt\,\Sigma(t_1,t)\,G(t,t_2),$$6, and simulation times up to two orders of magnitude longer become accessible in DMFT calculations. The truncation is especially well controlled when the self-energy is local or momentum independent (Stahl et al., 2021).

A third direction explores global-in-time iterative solvers. Fixed-point, Jacobian-free, and Newton–Krylov variants have been tested in nonequilibrium DMFT. Several global methods remain stable at large propagation times, but a standard forward fixed-point iteration does not. Even for the stable methods, the number of iterations required for fixed time step and target accuracy scales roughly linearly with the number of time steps, a behavior associated with a propagating front in the residual error. This identifies a concrete obstacle to making global solvers competitive with causal time-stepping schemes (Gašperlin et al., 12 Dec 2025).

These developments do not alter the formal content of the equations. They address how the same two-time physics can be represented or approximated more efficiently.

6. Reduced descriptions, practical scope, and known artifacts

The generalized Kadanoff–Baym ansatz replaces the full two-time evolution of $$\left(i\partial_{t_1}-h(t_1)\right)G(t_1,t_2)=\delta(t_1,t_2)+\int_\gamma dt\,\Sigma(t_1,t)\,G(t,t_2),$$7 by a reconstruction from equal-time densities and retarded/advanced propagators. In the scenarios tested by Reeves and collaborators, the terms neglected when deriving the GKBA are explicitly computed and shown to be orders of magnitude smaller than the terms retained, so that they provide only a small correction in the full Kadanoff–Baym equations. In the same study, both GKBA and Kadanoff–Baym calculations capture the dynamics of interacting systems with moderate and even strong interactions well (Reeves et al., 2023).

More drastic reductions lead toward kinetic theory. In first-principles carrier dynamics for bulk silicon, the generalized Baym–Kadanoff ansatz combined with the complete collision approximation yields a scattering term similar to that of the semiclassical Boltzmann equation (Sangalli et al., 2015). At a more formal level, analytical WKB solutions of the Kadanoff–Baym equations show that a generalized Boltzmann equation emerges whenever the WKB approximation holds, including far from equilibrium and in time-dependent backgrounds (Drewes et al., 2012).

At the same time, the full two-time framework has known pathologies when combined with approximate self-energies. In finite Hubbard clusters, self-consistent many-body approximations can generate correlation-induced damping and steady states under strong driving, even for isolated finite systems where such damping is purely artificial. The same study also finds that, for isolated clusters, the steady state reached is not unique but depends on how the external field is switched on. This behavior is tied to approximate self-energies based on infinite partial summations (Friesen et al., 2010).

A plausible implication is that the practical value of the Kadanoff–Baym equations depends on two distinct choices: whether the two-time structure must be retained, and whether the chosen self-energy respects the microscopic constraints relevant to the observable of interest. The formalism is exact at the contour level; its successes and failures in applications are governed by the approximation used for $$\left(i\partial_{t_1}-h(t_1)\right)G(t_1,t_2)=\delta(t_1,t_2)+\int_\gamma dt\,\Sigma(t_1,t)\,G(t,t_2),$$8 and by how the memory structure is represented numerically.