Discrete-time construction of nonequilibrium path integrals on the Kostantinov-Perel' time contour

TL;DR

This paper develops a rigorous method for constructing nonequilibrium path integrals on the Kostantinov-Perel' time contour using a discrete temporal mesh, enabling arbitrary initial states in many-body quantum problems.

Contribution

It introduces a discrete-time approach to build nonequilibrium path integrals on the Kostantinov-Perel' contour, overcoming limitations of the Schwinger-Keldysh formalism.

Findings

Enables arbitrary initial state preparation in nonequilibrium path integrals.

Provides a rigorous discrete-time construction method.

Extends the applicability of path integral techniques to more general initial conditions.

Abstract

Rigorous nonequilibrium actions for the many-body problem are usually derived by means of path integrals combined with a discrete temporal mesh on the Schwinger-Keldysh time contour. The latter suffers from a fundamental limitation: the initial state on this contour cannot be arbitrary, but necessarily needs to be described by a non-interacting density matrix, while interactions are switched on adiabatically. The Kostantinov-Perel' contour overcomes these and other limitations, allowing generic initial-state preparations. In this Article, we apply the technique of the discrete temporal mesh to rigorously build the nonequilibrium path integral on the Kostantinov-Perel' time contour.

Equations

$$\displaystyle\hat{\rho}_{0}\equiv\mathrm{e}^{-\beta\hat{\cal H}_{\rm M}}\biggl{/}\mathrm{Tr}(\mathrm{e}^{-\beta\hat{\cal H}_{\rm M}})~{},$$

$${\hat{\mathcal{H}}}_{\mathrm{M}}=\hat{\mathcal{H}}(t_0)-\mu \hat{\mathcal{N}},$$

$$\displaystyle\hat{\rho}_{0}\equiv\hat{\cal U}_{{\gamma_{\rm M}}}\biggl{/}\mathrm{Tr}(\hat{\cal U}_{{\gamma_{\rm M}}})~{}.$$

$$\hat{\mathcal{H}}(t)\equiv {\hat{\mathcal{H}}}_{\mathrm{f}}(t)+{\hat{\mathcal{H}}}_{\mathrm{b}}(t)+\hat{\mathcal{I}}(t).$$

$${\hat{\mathcal{H}}}_{\mathrm{f}}(t)={\hat{d}}^{†}\cdot {\mathbf{T}}_{\mathrm{f}}(t)\hat{d}$$

$${\hat{\mathcal{H}}}_{\mathrm{b}}(t)={\hat{a}}^{†}\cdot {\mathbf{T}}_{\mathrm{b}}(t)\hat{a}$$

$$Z\left[V\right]\equiv \frac{\mathrm{Tr}\left\{{\hat{\mathcal{U}}}_{\gamma }^{(V)}\right\}}{\mathrm{Tr}({\hat{\mathcal{U}}}_{{\gamma }_{\mathrm{M}}})}\equiv \frac{1}{\mathrm{Tr}({\hat{\mathcal{U}}}_{{\gamma }_{\mathrm{M}}})}\int \mathcal{D}\left(\overline{d},d\right)\int \mathcal{D}\left(a^{\ast },a\right){\mathrm{e}}^{i{S}^{(V)}\left[\overline{d},d;a^{\ast },a\right]},$$

$$S^{(V)}[\overline{d},d;a^{\ast },a]\equiv S_{\mathrm{f},Q}^{(V)}[\overline{d},d]+S_{\mathrm{b},Q}^{(V)}[a^{\ast },a]+S_{\mathcal{I}}[\overline{d},d;a^{\ast },a],$$

$$\displaystyle S_{\mathrm{f},Q}^{(V)}[\overline{d},d]\equiv\int_{\gamma}\mathrm{d}z\Big{\{}\overline{d}(z)\cdot i\partial_{z}d(z)-{\cal H}^{(V)}_{\rm f}[\overline{d}(z),d(z);z]\Big{\}}$$

$$\displaystyle S_{\mathrm{b},Q}^{(V)}[a^{*},a]\equiv\!\int_{\gamma}\mathrm{d}z\Big{\{}a^{*}(z)\!\cdot\!i\partial_{z}a(z)-{\cal H}^{(V)}_{\mathrm{b}}[a^{*}(z),a(z);z]\Big{\}}$$

$$S_{\mathcal{I}}\left[\overline{d},d;a^{\ast },a\right]\equiv -{\int }_{\gamma }\mathrm{d}z\mathcal{I}\left[\overline{d}(z),d(z);a^{\ast }(z),a(z);z\right]$$

$${\hat{\mathcal{H}}}_{\mathrm{M}}=\hat{\mathcal{H}}(t_0)-{\mu }_{\mathrm{f}}{\hat{\mathcal{N}}}_{\mathrm{f}}-{\mu }_{\mathrm{b}}{\hat{\mathcal{N}}}_{\mathrm{b}}+\hat{\mathcal{R}},$$

$${\mathcal{H}}_{\mathrm{f}}^{(V)}[\overline{d}(z),d(z);z]\equiv \overline{d}(z)\cdot {\mathbf{T}}_{\mathrm{f}}^{(V)}(z)d(z)$$

$${\mathcal{H}}_{\mathrm{b}}^{(V)}[a^{\ast }(z),a(z);z]\equiv a^{\ast }(z)\cdot {\mathbf{T}}_{\mathrm{b}}^{(V)}(z)a(z),$$

$${\mathbf{T}}_{f/b}^{(V)}(z)\equiv \Theta (t_{0-},z)\left[{\mathbf{T}}_{f/b}(t)+{\mathbf{V}}_{f/b}(z)\right]+\Theta (z,t_{0-}){\mathbf{T}}_{f/b}^{\mathrm{M}}.$$

$$t_j\equiv t_0+(j-1)\delta t,j=1,2,\dots ,N.$$

$${\tau }_j=t_0-i(j-1)\delta t,j=1,2,\dots ,M.$$

$$z_j=t_{j,+},j=1,2,\dots ,N,$$

$$z_{N+j}=t_{(N+1-j),-},j=1,2,\dots ,N,$$

$$z_{2N+j}={\tau }_j,j=1,2,\dots ,M.$$

$$\mathrm{Tr}\left\{{\hat{\mathcal{U}}}_{\gamma }^{(V)}\right\}=\int \mathrm{d}\left[a_0^{\ast },a_0\right]\int \mathrm{d}\left[{\overline{d}}_0,d_0\right]⟨a_0,d_0|{\hat{\mathcal{U}}}_{\gamma }^{(V)}|a_0,-d_0⟩{\mathrm{e}}^{-{|a_0|}^2}{\mathrm{e}}^{-{\overline{d}}_0\cdot d_0}$$

$$=\int \mathrm{d}\left[a_0^{\ast },a_0\right]\int \mathrm{d}\left[{\overline{d}}_0,d_0\right]⟨a_0,d_0|\prod _{j=1}^{2N+M-1}{\hat{\mathcal{U}}}_{z_j\to z_{j+1}}^{(V)}|a_0,-d_0⟩{\mathrm{e}}^{-{|a_0|}^2}{\mathrm{e}}^{-{\overline{d}}_0\cdot d_0}$$

$$=\int \prod _{j=1}^{2N+M}\mathrm{d}\left[a_j^{\ast },a_j\right]\int \prod _{j=1}^{2N+M}\mathrm{d}\left[{\overline{d}}_j,d_j\right]{\mathrm{e}}^{-{\sum }_{j=1}^{2N+M}{|a_j|}^2}{\mathrm{e}}^{-{\sum }_{j=1}^{2N+M}{\overline{d}}_j\cdot d_j}$$

$$\times \int \mathrm{d}\left[a_0^{\ast },a_0\right]\int \mathrm{d}\left[{\overline{d}}_0,d_0\right]{\mathrm{e}}^{-{|a_0|}^2}{\mathrm{e}}^{-{\overline{d}}_0\cdot d_0}⟨a_0,d_0|a_{2N+M},d_{2N+M}⟩$$

$$\times \left(\prod _{j=1}^{2N+M-1}⟨a_{j+1},d_{j+1}|{\hat{\mathcal{U}}}_{z_j\to z_{j+1}}^{(V)}|a_j,d_j⟩\right)⟨a_1,d_1|a_0,-d_0⟩$$

$$=\int \prod _{j=1}^{2N+M}\mathrm{d}\left[a_j^{\ast },a_j\right]\int \prod _{j=1}^{2N+M}\mathrm{d}\left[{\overline{d}}_j,d_j\right]{\mathrm{e}}^{-{\sum }_{j=1}^{2N+M}{|a_j|}^2}{\mathrm{e}}^{-{\sum }_{j=1}^{2N+M}{\overline{d}}_j\cdot d_j}$$

$$\times {\mathrm{e}}^{a_1^{\ast }\cdot a_{2N+M}}{\mathrm{e}}^{-{\overline{d}}_1\cdot d_{2N+M}}\left(\prod _{j=1}^{2N+M-1}⟨a_{j+1},d_{j+1}|{\hat{\mathcal{U}}}_{z_j\to z_{j+1}}^{(V)}|a_j,d_j⟩\right).$$

$$\int \mathrm{d}\left[a_1^{\ast },a_1\right]{\mathrm{e}}^{-a_1^{\ast }\cdot \left(a_1-a_{2N+M}\right)}\int \mathrm{d}[{\overline{d}}_1,d_1]{\mathrm{e}}^{-{\overline{d}}_1\cdot (d_1+d_{2N+M})}f(a_1,d_1)=f(a_{2N+M},-d_{2N+M}).$$

$$\mathrm{Tr}\left\{{\hat{\mathcal{U}}}_{\gamma }^{(V)}\right\}$$

$${\hat{\mathcal{U}}}_{z_j\to z_{j+1}}^{(V)}\to {\hat{\mathcal{U}}}_{+\delta t}^{(V)}\approx {\mathrm{e}}^{-i\delta t\left[\hat{\mathcal{H}}(t_j)+\hat{V}(t_{j,+})\right]},1\le j\le N-1;$$

$${\hat{\mathcal{U}}}_{z_N\to z_{N+1}}^{(V)}=1;$$

$${\hat{\mathcal{U}}}_{z_j\to z_{j+1}}^{(V)}\to {\hat{\mathcal{U}}}_{-\delta t}^{(V)}\approx {\mathrm{e}}^{i\delta t\left[\hat{\mathcal{H}}(t_{2N+1-j})+\hat{V}(t_{2N+1-j,-})\right]},N+1\le j\le 2N-1;$$

$${\hat{\mathcal{U}}}_{z_{2N}\to z_{2N+1}}^{(V)}=1;$$

$${\hat{\mathcal{U}}}_{z_j\to z_{j+1}}^{(V)}\to {\hat{\mathcal{U}}}_{-i\delta t}^{(V)}={\mathrm{e}}^{-\delta t{\hat{\mathcal{H}}}_{\mathrm{M}}},2N+1\le j\le 2N+M-1.$$

$$\hat{\mathcal{H}}(z_j)=\hat{\mathcal{H}}(t_j),\hat{V}(z_j)=\hat{V}(t_{j,+}),1\le j\le N-1;$$

$$\hat{\mathcal{H}}(z_j)=\hat{\mathcal{H}}(t_{2N+1-j}),\hat{V}(z_j)=\hat{V}(t_{2N+1-j,-}),N+1\le j\le 2N-1;$$

$$\hat{\mathcal{H}}(z_j)={\hat{\mathcal{H}}}_{\mathrm{M}},\hat{V}(z_j)=0,2N+1\le j\le 2N+M-1;$$

Full text

• March 2019

Keywords: Nonequilibrium, Path integrals, Time contours

1 Introduction

The theory of many-body systems brought away from equilibrium strongly relies on the concept of time contours. Introduced in different forms by several authors, including Schwinger [1], Keldysh [2], and Kostantinov and Perel’ [3], time contours provide an elegant way to deal on equal footing with time-ordered and anti-time-ordered products of operators which are both needed in the calculation of time-dependent observables out of equilibrium. For appropriately defined time contours, such cumbersome machinery is replaced by a single time-ordering procedure. The price to pay is that the size of the time domain has to be doubled with respect to the physical real-time axis.

Typically, one uses a time contour to evaluate nonequilibrium Green’s functions and associated observables. On both the Schwinger-Keldysh and the Kostantinov-Perel’ contours, equations of motion can be written down straightforwardly and many-body perturbation theory can be developed [2, 4, 5, 6, 7, 8, 9]. In particular, the Kadanoff-Baym equations [4] are the cornerstone of any numerical calculation of nonequilibrium Green’s functions for correlated systems [10, 11, 12, 13, 14].

Complementary to numerical approaches, path integrals constitute, in the context of nonequilibrium physics, a very useful analytical technique to derive effective theories, kinetic equations, and to define transparent physical approximations [15, 16, 17]. These derivations typically rely on nonequilibrium generalizations of procedures which appear also with equilibrium path integrals, such as the Gaussian integration over a subspace of fields (e.g. integration over bosons in a coupled fermion-boson system), the Hubbard-Stratonovich transformation, and the saddle-point approximation.

For the path integral to be meaningful, the action should be defined in a rigorous way, in particular, taking care of the correct operational definition of the inverse single-particle Green’s function operator, whose form depends on the time contour. Since the time contours used for nonequilibrium theories are different than the equilibrium (Matsubara) one, it is necessary to set up the corresponding actions case by case. This has been done for the Schwinger-Keldysh (SK) time contour [15]. Here, we argue that a similarly rigorous derivation for the case of the Kostantinov-Perel’ (KP) time contour is timely.

In what follows, we first summarize the main features of SK and KP time contours, respectively. We denote by $\gamma_{\rm K}=\gamma_{{\rm K},+}\cup\gamma_{{\rm K},-}$ the SK contour. Here $\gamma_{{\rm K},+}$ and $\gamma_{{\rm K},-}$ are the forward and backward branches of the real-time axis $(-\infty,+\infty)$, respectively. These are obtained by doubling the real-time degrees of freedom: $\gamma_{{\rm K},+}$ consists of the real axis traveled forward from $-\infty$ to $+\infty$, while $\gamma_{{\rm K},-}$ consists of the real axis traveled backwards from $+\infty$ to $-\infty$. Nonequilibrium theories on the SK contour require to specify the density matrix, $\hat{\rho}_{-\infty}$, which describes the system in the remote past, i.e. for $t\rightarrow-\infty$. Crucially, this density matrix is required to be non-interacting [15]. Many-body interactions, if present, are switched on adiabatically in such a way that the Hamiltonian coincides with the physical one at a given time $t=t_{0}$.

The KP contour [3, 6] allows a much greater flexibility with respect to the choice of the initial state of the system. In this formalism, one specifies the preparation of the system at an arbitrary time $t_{0}$ ($>-\infty$) via the density matrix $\hat{\rho}_{0}$, whose most general form is [6]

[TABLE]

where $\beta>0$ is a positive constant and $\hat{\cal H}_{\rm M}$ is an arbitrary Hermitian operator. A standard (but not mandatory) choice is

[TABLE]

where $\hat{\cal H}(t_{0})$ is the physical Hamiltonian of the system at time $t=t_{0}$, $\mu$ is the chemical potential, and $\hat{\cal N}$ is the number operator. This describes a system at thermal equilibrium at the initial time, in the grand-canonical ensemble, with inverse temperature $\beta$. However, it should be emphasized that $\hat{\cal H}_{\rm M}$ can be chosen to be a much more general operator [6], which allows to select a specific initial state, avoiding to assign the same weight to degenerate eigenstates of the Hamiltonian $\hat{\cal H}(t_{0})$. This allows to go beyond the grand-canonical ensemble and, for example, to select a broken-symmetry state instead of a thermal mixture, something which cannot be captured with the choice in equation (2). (We will show this by means of a concrete example in section 7.) Importantly, the KP formulation also allows to consider fully interacting systems without resorting to the adiabatic switching on of interactions, and is therefore appropriate for the study of strongly correlated systems. Since the adiabatic switching on is bypassed, one can control the preparation of the system at a finite time $t_{0}$, rather than in the far past. Because of this flexibility, the KP contour is a suitable framework for theories involving coupled fermion-boson systems, such as in the problem of nonequilibrium superconductivity [18] and/or initial broken-symmetry states, such as in nonequilibrium magnetism [19].

In the latter formalism, the real-time domain is $[t_{0},\infty)$, and the KP time contour $\gamma$ is given by the union of three branches: $\gamma=\gamma_{+}\cup\gamma_{-}\cup\gamma_{\rm M}$. The forward ($\gamma_{+}$) and backward ($\gamma_{-}$) branches are analogous to the SK branches, except that their domain of definition is $[t_{0},\infty)$, while the Matsubara branch ($\gamma_{\rm M}$) is the imaginary-time domain needed to describe the initial statistical mixture, $[t_{0},t_{0}-i\beta)$, where $\beta$ is the initial inverse temperature ($\hbar=1$ throughout this Article). One can write equation (1) in terms of an evolution operator along the Matsubara branch, i.e.,

[TABLE]

While the Green’s function problem on the KP contour is thoroughly discussed in the literature [6], a rigorous path-integral formulation seems to be available only in the case of the SK contour [15]. Given the usefulness of path integrals when dealing with problems that can be simplified by field integration [18], it is convenient to set up the same tools for the KP contour. The problem basically consists in deriving the KP nonequilibrium action, which correctly keeps into account the boundary terms arising from the construction of the path integral. This is the goal of the present Article.

Our Article is organized as follows. In section 2 we introduce a generic time-dependent Hamiltonian, including both fermionic and bosonic degrees of freedom. In section 3 we present the nonequilibrium partition function and the associated action that we need to derive from “first principles”. In section 4 we apply the technique of the discrete temporal mesh to derive the discrete form of the action. In section 5 we derive the general form of the free-fermion and free-boson Green’s functions on the KP contour, which allows to treat boundary conditions exactly while moving to the continuous time representation. We then show several important simplifications that occur in relevant particular cases. In section 6 we complete the transition to the continuous-time representation of the action. In section 7 we present the solution of a simple nonequilibrium model using both ordinary quantum mechanics and the KP path integral, with the aim of demonstrating the flexibility of the latter in fixing the initial conditions. Finally, in section 8 we summarize our main results and their applicability.

2 Hamiltonian

We consider a general system of fermions and bosons described by the following time-dependent Hamiltonian:

[TABLE]

Here,

[TABLE]

is the single-fermion Hamiltonian,

[TABLE]

is the single-boson Hamiltonian, and $\hat{\cal I}(t)$ includes all the remaining terms. The creation/annihilation operators, $\hat{d}^{\dagger}$/$\hat{d}$ for fermions and $\hat{a}^{\dagger}$/$\hat{a}$ for bosons, respectively, are grouped into vectors whose components are distinguished by single-particle quantum numbers (e.g. wavevector and spin projection for electrons on a lattice). We do not make any assumption on the physical nature of the fields that were just introduced. Fermions can be e.g. electrons, holes, or Nambu fermions. Bosons can be e.g. photons or phonons. The term $\hat{\cal I}(t)$ may include any form of interaction between these fields, that is, both between particles of the same species (e.g. electron-electron Coulomb interaction) and between particles of different species (e.g. electron-phonon coupling). The form of the interaction terms does not affect the derivation of the single-particle inverse Green’s function operator, which is the main mathematical issue that we solve here. The single-particle Hamiltonian matrix $\mathbf{T}_{\mathrm{f/b}}(t)$ can include a time-independent hopping (or, if the hopping is diagonalized, the single-particle energy spectrum), as well as the effect of any external time-dependent fields. Setting either $\mathbf{T}_{\mathrm{f}}(t)=0$ in (5) or $\mathbf{T}_{\mathrm{b}}(t)=0$ in (6) allows one to restrict the theory to the particular cases of purely bosonic or fermionic systems, respectively. In the following, we will always use bold fonts to denote matrices in the basis of single-particle quantum numbers.

3 Nonequilibrium path integrals

We denote the contour variable by $z\in\gamma$. For every real time coordinate $t\in[t_{0},\infty)$, there are two values of $z$: one belonging to $\gamma_{+}$, which we denote by $t_{+}$, and one belonging to $\gamma_{-}$, which we denote by $t_{-}$. For every imaginary-time coordinate on the Matsubara branch $[t_{0},t_{0}-i\beta)$ there is a single value of $z$.

The nonequilibrium partition function is expressed as a coherent-state path integral over Grassmann variables $(\overline{d},d)$ (representing the fermion fields) and complex variables $(a^{*},a)$ (representing the boson fields) [15]. As in the previous section, with a single Latin letter, such as “$d$”, we intend an array of variables, one for each different value of the set of single-particle quantum numbers. In the coherent-state representation, these variables acquire a further dependence on the contour variable $z$. The general expression for the partition function is

[TABLE]

where $\hat{\mathcal{U}}^{(V)}_{\gamma}$ is the evolution operator along the contour $\gamma$ , namely $\hat{\mathcal{U}}^{(V)}_{\gamma}=\hat{\mathcal{U}}_{\gamma_{\rm M}}\hat{\mathcal{U}}^{(V)}_{\gamma_{-}}\hat{\mathcal{U}}^{(V)}_{\gamma_{+}}$ and we have included some source potential $V(z)$ depending on the contour variable $z$. One has $\hat{\mathcal{U}}^{(0)}_{\gamma}=\hat{\mathcal{U}}_{\gamma_{\rm M}}$ and therefore $Z[V=0]=1$. The same identity, i.e. $Z=1$, applies if $V(t_{+})=V(t_{-})$, $\forall t$. For simplicity, we assume that the sources are quadratic in the particle fields. The nonequilibrium action is

[TABLE]

where

[TABLE]

and

[TABLE]

contain all the terms which are quadratic in the fermion or boson fields, respectively, while

[TABLE]

includes all the other terms, such as interactions between particles of the same or different kind. It is important to notice that, differently from the case of the SK contour, in the present case we must include the Matsubara branch and account for the fact that the Matsubara Hamiltonian $\hat{\cal{H}}_{\rm M}$, introduced in (1), is, in general, different from the physical Hamiltonian at the initial time, $\hat{\cal H}(t_{0})$. In general, one has

[TABLE]

where $\mu_{\rm f(\rm b)}$ is the chemical potential for the fermion (boson) subsystem, $\hat{\cal N}_{\rm f(\rm b)}$ is the fermion (boson) number operator, and $\hat{\cal R}$ includes all the remaining fermionic and bosonic terms. If the equilibrium state is a grand-canonical thermal distribution, then $\hat{\cal R}=0$. If, instead, one wants to implement a non-thermal initial state [6], a non-vanishing $\hat{\cal R}$ is needed.

As in the previous discussion about the physical Hamiltonian, we separate the terms of the Matsubara Hamiltonian which are quadratic in the particle fields from those which are not. The quadratic terms have the form $\hat{d}^{\dagger}\cdot\mathbf{T}^{\mathrm{M}}_{\mathrm{f}}\,\hat{d}$ and $\hat{a}^{\dagger}\cdot\mathbf{T}^{\mathrm{M}}_{\mathrm{b}}\,\hat{a}$, where we have introduced the Matsubara hopping matrices $\mathbf{T}^{\mathrm{M}}_{\mathrm{f/b}}$. After moving to the path integral formulation, all the terms which are not quadratic are collected into the quantity ${\cal I}\left[\overline{d}(z),d(z);a^{*}(z),a(z);z\right]$, introduced in equation (11), while the quadratic terms are included into equations (9) and (10) as

[TABLE]

and

[TABLE]

where

[TABLE]

In the above equations, $\Theta(z,t_{0-})$ equals $1$ if $z$ lies on the Matsubara branch, and [math] otherwise ($t_{0-}$ is the initial time taken on the backward branch). On the other hand, $\Theta(t_{0-},z)$ equals $1$ if $z$ lies on either $\gamma_{+}$ or $\gamma_{-}$, and [math] otherwise. We have also included the quadratic sources into (15). Non-quadratic sources, if present, must be included in (11). To better grasp the structure of $\mathbf{T}^{(V)}_{\mathrm{f/b}}(z)$ in (15), we emphasize that: 1) $\mathbf{T}_{\mathrm{f/b}}(t)$ is the physical time-dependent hopping matrix, so it depends on the real time coordinate $t$ and, if seen as a $z$-dependent quantity, it has the same value for $z=t_{+}$ and $z=t-$; 2) $\mathbf{V}_{\mathrm{f/b}}(z)$ is the source matrix, so it contributes only when $z$ is on the real-time branches and, for a given physical time $t$, it must satisfy $\mathbf{V}_{\mathrm{f/b}}(t_{+})\neq\mathbf{V}_{\mathrm{f/b}}(t_{-})$, so that $Z[V]\neq 1$; 3) $\mathbf{T}^{\mathrm{M}}_{\mathrm{f/b}}$ is constant on the Matsubara branch and does not contribute when $z$ is on the real-time branches, being replaced in this case by $\mathbf{T}_{\mathrm{f/b}}(t)$.

Once this issue about the Matsubara Hamiltonian is taken into account, writing down the explicit form of equation (11) is formally analogous to what is done in the case of the SK contour [15]. This term therefore presents no difficulties. In this Article, instead, we deal with the derivation of the explicit and unambiguous forms of (9) and (10). Indeed, the operator $i\partial_{z}$ appearing in these equations is a shorthand for an object that needs to be defined with great care within a discrete-time formulation. The difference between the KP $\gamma$ and SK $\gamma_{\rm K}$ time contours requires a generalization of the procedure given in [15], which employs a discrete temporal mesh to properly take into account boundary conditions. Besides this, it should be kept in mind that, starting from the discrete temporal mesh, the path integral is well defined only in the limit of a vanishingly small time step, which is taken at the end of the derivation (see section 6) and is not optional. With respect to [15], in our derivation we also include an arbitrary matrix structure of the hopping matrix.

4 The KP action

We now proceed to derive the path integral defining the nonequilibrium partition function on the KP time contour. We initially use the finite interval $[t_{0},t_{\infty}]$ on the real-time axis (with $t_{\infty}>t_{0}$) and take the limit $t_{\infty}\rightarrow+\infty$ at the end of the derivation. This interval is divided into $N-1$ arbitrarily small sub-intervals of width $\delta t$, so that $t_{\infty}-t_{0}=(N-1)\delta t$. We take the continuum $N\rightarrow\infty$ limit only at the end. We are therefore led to introduce a discrete set of time values:

[TABLE]

Note that $t_{1}=t_{0}$ and $t_{N}=t_{\infty}$. To each of the physical time values in (16) we assign two distinct contour coordinates, $t_{j,+}\in\gamma_{+}$ and $t_{j,-}\in\gamma_{-}$.

We then consider the Matsubara branch $[t_{0},t_{0}-i\beta]$. We split this interval into $M-1$ infinitesimally small parts, such that $\beta=(M-1)\delta t$, and we introduce

[TABLE]

We give a common name to the full set of discrete contour variables by introducing $2N+M$ contour coordinates $z_{j}$ defined as

[TABLE]

To determine the action $S^{(V)}[\overline{d},d;a^{*},a]$ in equations (7) and (8), which is a function of the fermionic (Grassmann) fields $(\bar{d},d)$ and the bosonic (complex) fields $(a^{*},a)$, we use a standard procedure that involves the decomposition of the identity operator of the full Hilbert space over fermionic and bosonic coherent states [15, 20, 21]. Denoting by $a_{j}$, $a^{*}_{j}$, $\overline{d}_{j}$, $d_{j}$ the vectors collecting the fields needed to specify the Hamiltonian in the $j$-th decomposition of the identity, we write

[TABLE]

Since $\left|a_{1},d_{1}\right>$ depends only on $a_{1}$ and $d_{1}$, and not on $a^{*}_{1}$ and $\overline{d}_{1}$, we notice that $\overline{d}_{1}$ enters the integral in the combination $\mathrm{e}^{-\overline{d}_{1}\cdot\left(d_{1}+d_{2N+M}\right)}$ only, while $a^{*}_{1}$ enters the integral in the combination $\mathrm{e}^{-a^{*}_{1}\cdot\left(a_{1}-a_{2N+M}\right)}$ only. Using the representations of ordinary and Grassmann $\delta$ functions [20], we then observe that, for any function $f(a_{1},d_{1})$,

[TABLE]

We then obtain

[TABLE]

Taking $\delta t$ to be infinitesimally small, we have

[TABLE]

We now introduce $\hat{\cal H}(z_{j})$, $\hat{V}(z_{j})$ and $\delta z_{j}$ in such a way that

[TABLE]

and

[TABLE]

We then write the contour evolution operator compactly as

[TABLE]

and we obtain

[TABLE]

where $a_{1}=a_{2N+M}$, and $d_{1}=-d_{2N+M}$ in the arguments of the exponentials. In the continuous-time representation, one formally puts

[TABLE]

obtaining

[TABLE]

where the action is given by equation (8). The boundary conditions $a_{1}=a_{2N+M}$ and $d_{1}=-d_{2N+M}$ make the definition of the operator $\partial_{z}$ on the contour non-trivial, something which is fully captured by using the discrete temporal mesh.

In passing, we mention that discrete-temporal-mesh methods have also been used to tackle other problems, related, for example, to the correct definition and use of time-derivative-dependent variable transformations in path integrals [22].

5 Green’s functions for the quadratic part of the action

5.1 Derivation of the general formula

In the discrete-time representation, the quadratic ($Q$) actions, (9) and (10), read as

[TABLE]

and

[TABLE]

where the free-fermion and free-boson inverse Green’s functions (GFs) are matrices in the discrete-time coordinates ($j,j^{\prime}$), as well as in the free-particle indexes. Note that we are using the term “free” in the sense of non-interacting. External fields coupling with the particles are included in the definition of $\mathbf{G}^{-1}_{Q}$, but interactions are left into the action term $S_{\cal I}$, equation (11). Our goal here is to find the direct GFs by inverting the matrices $\mathbf{G}^{-1}_{Q}$.

In order to lighten the notation, in the following steps we will omit subscripts and superscripts of the GFs, except for those referring to the discrete-time coordinates. To distinguish between bosons and fermions, we introduce the index $\xi=\pm 1$, where the $+$ sign is for bosons and the $-$ sign is for fermions.

The free-particle inverse GF matrices are then written compactly as

[TABLE]

where the dependence of $\mathbf{T}_{j}$ on the contour coordinate is detailed in (15). The only formal difference between fermions and bosons is in the $\xi$-dependent term appearing in the second line of (31), which reflects the different boundary conditions arising from the construction of the path integral.

The matrix in equation (31) is invertible, i.e. the left inverse is the same as the right inverse. We have

[TABLE]

where

[TABLE]

To solve equation (32), we first introduce the Ansatz

[TABLE]

in the second line of (32); separating the case $j=2$ from $j\neq 2$, we obtain

[TABLE]

The second equation ($j\neq 2$) is solved by iteration: for some $n<j-1$,

[TABLE]

where $\theta_{p}$ is the discrete step function ($\theta_{p}=1$ if $p$ is true, $\theta_{p}=0$ if $p$ is false), and we have introduced the contour-anti-ordered product

[TABLE]

which is accompanied by the contour-ordered product

[TABLE]

From equation (36), we take $n=j-2$ and, using the first among equations (35), we find

[TABLE]

We take $j=2N+M$ and solve for every $j^{\prime}$:

[TABLE]

where we have introduced the quantities

[TABLE]

Inserting (40) into (39), we solve the latter for $3\leq j\leq 2N+M-1$, while we solve for $j=2$ using (35). The solutions for these cases can be combined into a single one:

[TABLE]

It can be seen that (42) coincides with (40) for $j=2N+M$. Therefore, (42) represents the full solution for all values of $j,j^{\prime}$. From (34), we find the free-particle GF:

[TABLE]

We now take the continuum limit. We first introduce

[TABLE]

These identities give the real (if $1\leq j\leq 2N$) or complex (if $2N+1\leq j\leq 2N+M$) time coordinates on the contour corresponding to the discrete index $j$.

Taking $\delta t$ as infinitesimally small, we find

[TABLE]

and

[TABLE]

where $\mathcal{T}_{\gamma}$ is the contour-ordering operator along $\gamma$, and $\overline{\mathcal{T}_{\gamma}}$ is the analogous contour-anti-ordering operator. In particular,

[TABLE]

where the simplification is possible due to the contour structure of the hopping detailed in equation (15). In the second step of (5.1), the contour-anti-ordered exponential is different from $\mathbf{1}$ only in the presence of sources, while the term $\exp\left(\beta\mathbf{T}^{\mathrm{M}}\right)$ is a matrix generalization of the inverse Boltzmann factor. This suggests to introduce a generalized occupation-number matrix,

[TABLE]

We will see that, in the simplest cases, this quantity reduces to the standard occupation number.

We now switch to the full continuum notation for the GF, obtaining

[TABLE]

where $\Theta(z,z^{\prime})$ is the step function on $\gamma$, with $\Theta(z,z)=1$.

Equation (49) represents the explicit nonequilibrium free-fermion and free-boson GFs on the KP time contour $\gamma$, for the most general case of a single-particle Hamiltonian with a non-trivial matrix structure, including time-dependent fields and sources. We note that, in the absence of sources ($V=0$), equation (49) reduces to the noninteracting GF derived in [6] with the method of equations of motion. We will now consider several relevant cases in which (49) can be significantly simplified.

5.2 GFs in the absence of sources

Most nonequilibrium problems require functional differentiation with respect to the sources, followed by the evaluation of the result at $\mathbf{V}(z)=\mathbf{0}$. Therefore, although nonzero sources are needed to make the path integral meaningful, at some point in the calculation one typically needs to compute the GFs in the absence of sources. In this case, $\mathbf{T}(z)=\mathbf{T}(t)$ if $z\in\gamma_{+}\cup\gamma_{-}$, while $\mathbf{T}(z)=\mathbf{T}^{\mathrm{M}}$ if $z\in\gamma_{\mathrm{M}}$. Equations (45) and (46) simplify to

[TABLE]

[TABLE]

[TABLE]

and

[TABLE]

where $\mathcal{T}$ and $\overline{\mathcal{T}}$ are the standard time-ordering and anti-time-ordering operators on the real axis, respectively, $t_{j}$ are the time coordinates defined as in (5.1), and $0<\tau<\beta$ parameterizes the Matsubara branch as $t_{0}-i\tau$. In particular,

[TABLE]

The occupation number matrix introduced in (5.1) becomes

[TABLE]

which is the matrix generalization of the occupation number resulting from either the Bose-Einstein or the Fermi-Dirac statistics. Once $z$ and $z^{\prime}$ are specified, equation (49) is simplified using (50)-(53) and (55).

5.3 GFs in the equilibrium case for a diagonal single-particle Hamiltonian

Let us consider the case of an equilibrium single-particle GF and assume that the single-particle Hamiltonian (and, therefore, the GF as well) is diagonal in the single-particle indexes. It should be noted that treating an equilibrium single-particle GF does not mean that the system must be at equilibrium. It only means that one chooses to exclude the nonequilibrium features from the definition of the single-particle Hamiltonian, i.e. the time-dependent fields are included into (11), despite being quadratic contributions to the action. Let $k$ be the set of single-particle quantum numbers. For a system of electrons on a lattice, such set would consist of the wavevector and a spin projection.

Including a possibly $k$-dependent chemical potential term (e.g. a spin-dependent chemical potential for electrons), which enters the definition of the Hamiltonian on the Matsubara branch, we find $T_{k,k^{\prime}}(z)=\delta_{k,k^{\prime}}\varepsilon_{k}(z)$, with

[TABLE]

and the number matrix reduces to the familiar occupation number,

[TABLE]

In this case, the free-particle GF is diagonal, with diagonal components given by

[TABLE]

where $t$ and $t^{\prime}$ are the (complex) time coordinates corresponding to the contour coordinates $z$ and $z^{\prime}$, respectively. Recall that $\xi=\pm 1$ distinguishes between the bosonic and fermionic cases, respectively. Equation (58) represents the generalization of the results given in [15] for the SK time contour, to the case of the KP time contour.

To better illustrate the main features of equation (58), we explicitly consider the various combinations obtained when the positions of $z$ and $z^{\prime}$ on the contour are specified. First, consider the case in which both $z$ and $z^{\prime}$ belong to the real-time branches (either to $\gamma_{+}$ or to $\gamma_{-}$). In this case, one obtains

[TABLE]

if $z>z^{\prime}$ (on the contour) and

[TABLE]

if $z<z^{\prime}$. These coincide with the standard greater and lesser nonequilibrium GFs which are also found in the SK formalism [15].

On the KP contour, one also obtains the other Langreth components [6]. If $z$ belongs to one of the two real-time branches, while $z^{\prime}$ is on the Matsubara branch,

[TABLE]

In the opposite case, when $z$ belongs to the Matsubara branch and $z^{\prime}$ to one of the real-time branches,

[TABLE]

Finally, if both $z$ and $z^{\prime}$ belong to the Matsubara branch,

[TABLE]

where $\theta(\tau-\tau^{\prime})$ is the ordinary step function.

5.4 Determinant

In several applications which involve field integration over either the fermionic/bosonic degrees of freedom, it is necessary to know the determinant of $(-i\mathbf{G}^{-1}_{j,j^{\prime}})$. In general, explicit calculations of this quantity are difficult.

However, in the case when the hopping is diagonal, the calculation can be done easily and directly from equation (31). The determinant can be evaluated using the Laplace theorem combined with the fact that the determinant of a triangular matrix is the product of its diagonal elements. The result is independent of $M$ and $N$ being even or odd, implying that it is well defined in the limit $M,N\rightarrow\infty$. We obtain (omitting the single-particle quantum numbers)

[TABLE]

6 From the discrete to the continuum representation

We now apply the continuum representation for the time domain, in order to make the path-integral construction meaningful. We define the operator

[TABLE]

which satisfies

[TABLE]

where $\mathbf{G}_{\xi}(z^{\prime},z^{\prime\prime})$ is given by (49) and $\delta(z,z^{\prime\prime})$ is the Dirac delta on $\gamma$, namely

[TABLE]

In the continuum, equation (29) reads as

[TABLE]

while equation (30) reads as

[TABLE]

The dependence on the sources $V$ is implicit. The inverse GFs appearing in (68) and (69) are to be considered as merely symbolic objects. Any calculation (involving, for example, field integration or extremization of the action) must, at some stage, rely on a transformation from the inverse GF to the direct GF in equation (49). This encodes the boundary conditions arising from the procedure of construction of the path integral, while being fully well-defined in the continuum limit of the time domain.

7 An example: Spin evolution in a time-dependent Zeeman field

To show how the KP time contour can be used to specify the initial state in practice, we consider a simple and well-known problem that can be solved analytically with elementary means, allowing for a straightforward comparison of the result with that obtained in the realm of path integrals.

The model involves a single electronic orbital subjected to a time-dependent magnetic field along a certain direction $x$, which couples to the spin through the Zeeman coupling. We neglect interactions, which would affect the doubly-occupied state and choose a basis for the electronic states where spin is quantized along $z$, with $z\perp x$. The Hamiltonian is

[TABLE]

where $\hat{s}_{x}=(\hat{d}^{\dagger}_{\uparrow}\hat{d}_{\downarrow}+\hat{d}^{\dagger}_{\downarrow}\hat{d}_{\uparrow})/2$ and $B_{x}(t)$ is the external magnetic field, in appropriate units. We assume $B_{x}(t_{0})=0$. Our goal is to compute the time-dependent density matrix

[TABLE]

from which single-particle observables (particle number and spin components) can be read off. The quantity (71) is obviously very sensitive to the initial preparation of the system, demonstrating spin precession if the spin in the initial state has a non-vanishing component perpendicular to $x$, or no evolution at all otherwise. In section 7.1 we compute (71) by using tools of basic quantum mechanics. Then, in section 7.2 we compute (71) by using the KP partition function that we have derived before, to show that an appropriate choice of $\hat{\cal H}_{\mathrm{M}}$ allows to choose the initial state and therefore to reproduce all possible dynamical evolutions of the system.

7.1 Solution of the problem with ordinary quantum mechanics

We first solve the problem using ordinary quantum mechanics. The Fock space consists of four states: the vacuum $\left|0\right>$, the singly-occupied states $\hat{d}^{\dagger}_{\sigma}\left|0\right>\equiv\left|\sigma\right>$, and the doubly-occupied state $\hat{d}^{\dagger}_{\uparrow}\hat{d}^{\dagger}_{\downarrow}\left|0\right>\equiv\left|\uparrow\downarrow\right>$. The evolution operator, satisfying $i\,\partial_{t}\,\hat{\cal U}(t,t_{0})=\hat{\cal H}(t)\,\hat{\cal U}(t,t_{0})$ and $\hat{\cal U}(t_{0},t_{0})=\hat{1}$, is obtained directly as

[TABLE]

where

[TABLE]

The expectation value in (71) can be calculated from

[TABLE]

where $W_{n}$ is the statistical weight of the Fock state $\left|n\right>$. The final result (after a number of simple trigonometric manipulations) is

[TABLE]

It is important to notice that the freedom on the choice of the weights $W_{n}$ allows to prepare the system in a non-trivial initial state and to implement symmetry breaking. For example, if we make the rotationally symmetric choice $W_{\uparrow}=W_{\downarrow}$, we get $\rho_{\sigma,\sigma^{\prime}}(t)=\delta_{\sigma,\sigma^{\prime}}\left(W_{\uparrow\downarrow}+W_{\uparrow}\right)$: the density matrix becomes time-independent, encoding only information about the (constant) population of the orbital, as all the spin components are zero. This choice, in more sophisticated cases, is not uncommon at all: it is exactly what one would obtain from a grand canonical thermal mixture, where $W_{n}=\mathrm{e}^{-\beta(E_{n}-\mu N_{n})}/\mathcal{Z}$ and $\mathcal{Z}=\sum_{n}\mathrm{e}^{-\beta(E_{n}-\mu N_{n})}$. In such case, the equality between $W_{\uparrow}$ and $W_{\downarrow}$ comes from the degeneracy of the states $\left|\uparrow\right>$ and $\left|\downarrow\right>$ with respect to the Hamiltonian $\hat{\cal H}(t_{0})$. It is by allowing $W_{\uparrow}\neq W_{\downarrow}$ (hence, deviating from a thermal mixture) that rotational symmetry is broken, the $z$ direction is selected as the one along which the spin is aligned in the initial state, and a non-trivial spin evolution occurs in response to the application of a magnetic field perpendicular to $z$.

7.2 Solution of the problem by means of the KP path integral approach

Can we obtain the same result (75) within the KP path integral formulation? In particular, do we have the same freedom on the preparation of the initial state? To see this, let us follow the procedure outlined above. Recalling that $B_{x}(t_{0})=0$, we take the following Hamiltonian on the Matsubara branch:

[TABLE]

where $\hat{n}_{\sigma}=\hat{d}^{\dagger}_{\sigma}\hat{d}_{\sigma}$. The inclusion of the second term ($\propto\Delta$) in equation (76), which is equivalent to assuming a spin-dependent chemical potential, is the key tool that allows to break rotational symmetry in the initial-state preparation. We emphasize that (76), for $\Delta\neq 0$, is crucially different from the grand-canonical Matsubara Hamiltonian, which would be equal to $-\mu\left(\hat{n}_{\uparrow}+\hat{n}_{\downarrow}\right)$. Therefore, we are not choosing a thermal mixture.

The nonequilibrium partition function, from equation (7), is

[TABLE]

where

[TABLE]

is the equilibrium partition function. Because $Z[0]=1$, we have

[TABLE]

The action is quadratic, as in (68). The free-electron GF has the general form given by (49), with $\mathbf{T}$ being a matrix in spin space only:

[TABLE]

where $\mathbf{1}$ is the $2\times 2$ identity and $\bm{\sigma}_{n}$ are ordinary $2\times 2$ Pauli matrices. The source term in the action is

[TABLE]

so that the density matrix is obtained from functional differentiation with respect to the source fields as

[TABLE]

Carrying out the path integral in (77), we obtain

[TABLE]

The explicit calculation of (82) requires the following identity:

[TABLE]

where $V=0$ is intended in the last line. We finally find

[TABLE]

To calculate the GF, let us consider equation (49) specialized to our case. Since the GF needs to be calculated in the absence of sources, we can use the simplifications discussed in section 5.2. We find

[TABLE]

where we have exploited the fact that $\mathbf{T}(t)$ commutes with itself at different times. The exponentials of the Pauli matrix give the usual result:

[TABLE]

Equation (86) includes information about the preparation of the initial state via the occupation number matrix:

[TABLE]

The components of (86), after simple trigonometric manipulations, can be written as

[TABLE]

where we have identified the components of the density matrix $\rho_{\sigma,\sigma^{\prime}}(t)$ according to (85) (minding the matrix transposition).

We immediately see that the resulting density matrix, calculated from the path integral on the KP contour, has the same form as that computed with ordinary quantum mechanics—see equations (75)—once the following identifications are made:

[TABLE]

where the weights $W$ in the right-hand sides are those appearing in (75). One solution of equations (90) is given by the grand-canonical statistical weights

[TABLE]

where $\mathcal{Z}$ is the equilibrium partition function given in (78). Note that $\Delta\neq 0$ (i.e. a non-thermal choice for $\hat{\cal H}_{\mathrm{M}}$) is essential to realize spin symmetry breaking in our example, as it allows for $W_{\downarrow}\neq W_{\uparrow}$ and to tune the relative weights at will. Moreover, taking the limit $\Delta\rightarrow+\infty$ allows to prepare the system in the $\left|\downarrow\right\rangle$ state ($W_{\downarrow}=1$, $W_{\uparrow}=0$) and, viceversa, taking the limit $\Delta\rightarrow-\infty$ allows to prepare the system in the $\left|\uparrow\right\rangle$ state ($W_{\downarrow}=0,W_{\uparrow}=1$).

This simple example illustrates the power of the path integral approach on the KP contour, which gives us the full freedom of defining the Hamiltonian on the Matsubara branch in an arbitrary way.

8 Summary and conclusions

In summary, we have presented a discrete temporal procedure for the rigorous construction of the nonequilibrium path integral on the KP time contour. Our main result is that we have rigorously converted the general (contour-independent) expressions given by equations (9) and (10) into the forms of equations (68) and (69), via the expression of the GF, equation (49), that is specific of the KP contour.

Our procedure generalizes the one used for the SK contour [15], allowing us to include the imaginary-time (Matsubara) branch in addition to the real-time branches. Consequently, path-integral theories on the KP time contour can account, in principle, for arbitrarily correlated initial states and/or statistical mixtures, or particular chosen pure states.

Of course, many of the results that can be usually achieved with path integrals can also be obtained by applying other methods, such as the equations of motion for nonequilibrium GFs. However, path integrals are extremely useful when part of the fields entering the action can be integrated away, and/or when a suitable Hubbard-Stratonovich transformation on the interaction term yields a meaningful physical theory which can then be studied, e.g., via a saddle-point approximation. Interested readers can find the application of all these procedures in the context of a theory of laser-induced nonequilibrium superconductivity [18].

Acknowledgements

This work was supported by the European Union’s Horizon 2020 research and innovation programme under grant agreement No. 785219 - GrapheneCore2. We thank Alex Kamenev for very stimulating comments on our manuscript.