Landauer’s formula with finite-time relaxation: Kramers’ crossover in electronic transport

Landauer’s formula is the standard theoretical tool to examine ballistic transport in nano- and meso-scale junctions, but it necessitates that any variation of the junction with time must be slow compared to characteristic times of the system, e.g., the relaxation time of local excitations. Transport through structurally dynamic junctions is, however, increasingly of interest for sensing, harnessing fluctuations, and real-time control. Here, we calculate the steady-state current when relaxation of electrons in the reservoirs is present and demonstrate that it gives rise to three regimes of behavior: weak relaxation gives a contact-limited current; strong relaxation localizes electrons, distorting their natural dynamics and reducing the current; and in an intermediate regime the Landauer view of the system only is recovered. We also demonstrate that a simple equation of motion emerges, which is suitable for efficiently simulating time-dependent transport.

Let us a consider a single electronic level of energy ω k connected to a manifold of non-interacting states that comprise the implicit reservoir E k . The index k denotes the extended reservoir state that these states are coupled to. The E L(R) from the main text are composed of all E k for k ∈ L(R). The Hamiltonian of this partial system is (A1) This Hamiltonian describes the systems shown in Fig.  1(b) in the main text. For an implicit reservoir state α, the isolated Green's functions are and g a α (t, t ) = −θ(t − t) g > (t, t ) − g < (t, t ) ; (A5) or, in terms of their Fourier transforms, g < α (ω) = 2πıf (ω α )δ(ω − ω α ), (A7) g r α (ω) = 1/(ω − ω α + ıη), and g a α (ω) = 1/(ω − ω α − ıη), where η is the infinitesimal positive number and f (ω) is the Fermi-Dirac distribution. The subscript α on g α is used to distinguish it from an extended reservoir state Green's functions (g k ) or the full system Green's functions (G ij ) as used for the L − S − R system. For the interaction of one level with all the other levels, Eq. (A1), we have (symbolically, on the Keldysh contour) g k = g 0k + g 0k Σ k g k , with g 0k being the single isolated extended reservoir site. Using the relation between on-contour and real-time non-equilibrium Green's functions, 1 for the retarded Green's function we have g r k (ω) = g r 0k (ω) + g r 0k (ω)Σ r k (ω)g r k (ω), where Here, γ(ω) = 2π α t 2 α δ(ω−ω α ). Evaluating the integral, one obtains where E k (ω) is a frequency-dependent energy shift and γ(ω) is a frequency-dependent relaxation rate. Similarly, the lesser self-energy is evaluated as The retarded Green's function of the extended reservoir site then becomes As is seen, "friction" is in general non-Markovian, since the dephasing rate γ(ω) is frequency-dependent. The frequency shift, E k (ω), is also non-Markovian in the same sense. Now let us find the the full lesser Green's function.
To this end we will use the Keldysh equation 1,2 resulting in The important fact here is that the Fermi-Dirac distribution here is not evaluated at ω k or at any other fixed frequency. Instead, it is evaluated at ω and, therefore, this factor is the same (at fixed ω) for any Green's function of any site within the same reservoir. In particular, it guarantees that the current vanishes when the Fermi-Dirac distribution is the same for the two extended reservoirs.
The often used wide band approximation would result in ω-independent γ and vanishing E k . In this approximation, the result for the lesser Green's function is In a similar manner, the retarded Green's function is Accordingly, one has, using the identity g r 2. Landauer-like Formula We are following a notation similar to the original Meir-Wingreen paper. 3 The Hamiltonian of the system is with k = ω k . Unlike the Meir-Wingreen scenario where such a Hamiltonian fully describes the steady-state since the sizes of left and right reservoirs are assumed to be infinite, here we take a finite number of extended reservoir states, but each such site k is connected to an implicit reservoir according to the Hamiltonian (A1).
The time-dependent current from the left reservoir to the system can be written as 3 When the steady state is established, we can take the Fourier transform, Since the reservoir sites are non-interacting (in a twoor more-electron sense) we have the following Dyson equation 2 (A24) Using these identities, Eq. (A22) can be rewritten as We emphasize that even though g k (ω) are single-particle non-interacting Green's functions, the corresponding quasiparticles do have finite lifetime because of their coupling to the implicit reservoirs, which is different from the original Meir-Wingreen formulation. For example, g ≷ k (ω) are not delta functions with respect to ω, but rather Lorentzians (when the wide band limit is taken, Eq. (A17)).
To proceed further we use (i) the Keldysh equation, 1,2 G ≷ = G r Σ ≷ G a , and (ii) that the self-energy due to the interaction with reservoir sites is Σ , and (iii) the linear relation between real-time Green's functions G > −G < = G r −G a . 1 Using these identities the current can be rewritten as Considering the equilibrium property of the isolated state k, Eq. (A19), one gets It is clearly seen that once the Fermi-Dirac distribution becomes identical on the left and on the right, the current vanishes. Actually, any partial current also vanishes, i.e., the current for a specific choice of indices m, n, α, β, k, l and frequency ω, as of course is expected due to the necessary detailed balance at equilibrium. A concise expression for the current can be written by introducing the spectral density and the expression for the current becomes Using the matrix notation, This expression is Eq. (3) from the main text.
3. Small γ, Nr → ∞ As an example, we examine a system that is a single noninteracting state. In this case, the continuous form for the extended reservoir self-energy from integrating the single site Green's functions, Eq. (12) from the main text, is Σ r(a) (ω) = ω ∓ ıγ/2 − ı 4J 2 − (ω ∓ ıγ/2) 2 /2, (A31) when W = 4J, and can be broken into the real and imaginary components Re Σ r (ω) = 1 4 8 ω 2 + 4J 2 γ 2 + 16 ω 2 − 4J 2 2 + ω 4 The term in the integral of Eq. (A30) can then be written in terms of these components, and then expanded in γ and ω. Since both γ and ω appear as the same order in the expansion and the self-energy is small outside the band edge, the error in the ω integral is small when γ is small. In practice, the bias is typically also taken to be small. For a single central site, The first term is the transmission coefficient for a 1D lattice in the Landauer formula, so the correction is of order γ 2 .

Large γ
As γ increases, the Green's function for a single reservoir site, Eq. (A18), approaches g r k (ω) = −2/γ in the relevant region of integration (i.e., where the difference in Fermi distributions is non-negligible). The transmission coefficient in that region, then, is Using the explicit form for the total Green's function, with S consisting of a single non-interacting state with frequency ω s , this simplifies to where we used that k v 2 k = J 2 since the couplings v k come from a unitary transformation times the total coupling to the system. When the bias window includes this peak and at zero temperature, the current in the large γ regime becomes This expression is independent of the bias, as long as it is nonzero. As we note in the main text, the transition to this value of the current, though, does depend on the bias.

Appendix B: Markovian Master Equation
In the main text and above, we derived a Landauer formula for the L − S − R system, which results in a Markovian master equation for the real-time dynamics in the small-to-intermediate γ regime. As we show below and mention in the main text, outside of this regime this formalism gives physically invalid results. In the appropriate region, though, it allows for the direct calculation of the full time dynamics and can be readily expanded to include many-body interactions or time-dependent terms. Here we will derive the full solution-in all regimes-to the Markovian master equation. Given the Hamiltonian the starting point is the Markovian master equation, This form for γ k + and γ k − ensures that in the absence of the interaction H I , the extended reservoirs will relax into an independent equilibrium of their own Hamiltonians, H L and H R . For simplicity, we have assumed that γ is the same in both the left and right regions, and for all k.
This master equation describes the evolution of a system in the presence of explicit reservoir states, with a different mechanism than used in the main text. The following will examine the behavior of Eq. (B2) over the full range of the relaxation γ.

Single Extended Reservoir State Green's Function
In the absence of S, the extended reservoir states decay into the equilibrium state, exp(−βH L(R) )/Z, i.e., the occupations decay to a Fermi-Dirac distribution as e −γt/2 . This can be shown from the exact time-dependent lesser Green's function in a reservoir site uncoupled from S, which is given by where Lρ = dρ/dt can be found from Eq. (B2) and the superoperator L is the Lindbladian. The equilibrium state of the site has filling f k , so tr c † k c k ρ eq = f k . In the Fock basis of a single state, the equilibrium state is A Jordan-Wigner transformation maps the electron creation and annihilation operators onto spin operators, which allows us to write c k ρ eq = (f k σ x + ıf k σ y )/2. For the equation of motion, Eq. (B2), the Lindblad operator is block diagonal in the the σ x , σ y subspace and σ I , σ z subspace. This means we can separately solve for the dynamics using L for these two subspaces. If we wish to calculate the action on a generic operator O = a 0 σ I + a x σ x + a y σ y + a z σ z , for the σ x , σ y subspace: Acting on both sides with (1/2)σ x tr gives da Solving these equations of motion we obtain (B6) Then acting with c k and taking the trace gives for t ≥ 0. This can be readily employed to find the retarded and advanced Green's functions for the single state or its Fourier transform, g r(a) k (ω) = (ω − ω k ± ıγ/2) −1 . Physically, the reservoir sites are exchanging electrons with a larger external reservoir with an infinite number of electrons and states without memory. The lesser Green's function is also found to be This expression recovers the original lesser Green's function, Eq. (A17), except that the distribution is evaluated at the ω k of the state, rather than being a continuum over ω. As we discuss below, this has the effect of broadening the density of states after they are occupied, rather than occupying after broadening.

Steady-State Current
The general solution to the steady-states of the master equation (Eq. (B2)) can also be found in an analogous way to Jauho, Meir, and Wingreen and follows the same process as the derivation for Eq. (A30). In this case we use the equilibrium relation from the previous section, Eq. (B9), rather than using the one that has an ω-dependent distribution f (ω). In practice, this derivation is the same but with a filling, f k , dependent on the reservoir energy, k .
That is, after applying the Markovian equilibrium property, Eq. (B9), the expression for the current, Eq. (A27), is instead Again, a concise expression for the current can be written by introducing the spectral density and the population-weighted spectral densitỹ where the Fermi-Dirac distributions are included. Then, the current can be written as Note that, in this case, the integrand in Eq. (B13) does not include the Fermi-Dirac distribution as a separate prefactor as in the Meir-Wingreen / Landauer formula, but rather appears as a convolution with a Lorentzian due to the inclusion of a finite relaxation time from the Markovian master equation.
When the reservoir states are symmetrically coupled to the system, a simplified expression results: For this symmetric case, when the two extended reservoirs have the same chemical potential, the calculated current I is always equal to zero. However, in the asymmetric case, Eq. (B10) can yield a non-zero current even with no applied bias. We will examine this fact in Appendix C in order to develop a bound for when the Markovian master equation is consistent with physical expectations.

Small and Large γ Limits
In similar fashion to the main text, we can derive the limiting expressions for the steady-state current in the low and high relaxation rate regimes. The spectral function for a single reservoir site connected to an implicit bath is given by , so long as f (ω) changes little over the width of A k (ω), and the expression for the current can be rewritten as This recovers Eq. (3) from the main text, where the spectral density and Green's functions include the relaxation rate γ. What this means is that the master equation can be used to simulate transport, provided that γ is in a suitable range. Even the intermediate regime can be accurately simulated, so long as N r is sufficiently large. Similarly for large γ, the integral within Eq. (B14) can be calculated by transforming out from the frequency domain: which applies to both interacting and non-interacting Green's functions. In the non-interacting case, this can be thought of an integration over the density of states for each site in S. This yields an expression for the current in the large γ regime to be where again the sum is over just a single set of k (either in the left or right extended reservoir, which are identical). When S consists of a single site, then the large-γ current is found to be The sum of the v 2 k terms is the transformation that diagonalizes the extended reservoir's single particle Hamiltonian. 4 Thus, where n L(R) is the occupation of the extended reservoir state in real space at the site immediately adjacent to the system on the left (right). This makes a correspondence with a setup with just a single extended reservoir site on each end (N r = 1), I 3 ≈ 2eJ 2 (n L − n R )/γ. Essentially, the sites further away in the N r > 1 case are effectively decoupled from the system as any flow of electrons away from those sites is suppressed due to the strong relaxation. Therefore, the current in the large γ regime is also independent of N r with the exception of discretization effects.
In contrast to the Landauer-like formalism, the Markovian master equation allows for an explicit derivation showing only the occupation at the boundaries matter, i.e., all other electrons are prevented from flowing to the system by the strong relaxation. This is possible due to the well-defined occupation of each reservoir site. In the full model, the states are broadened and then occupied, so there is no equivalent transformation to a single edge site occupation like there is in the Markovian case.
We can, however, go further with the Markovian procedure. Within a small bias window, the reservoir coupling is approximately constant v 2 k ≈ J 2 /N r , and the sum of the occupation terms is approximately the total number of states within the bias window Substituting this in, we find I 3 ≈ eV [2πJ 2 /(γW )](2π ). Comparing with Eq. (10) in the main text, the large γ current here is additionally inversely related to the bandwidth and grows linearly with V . Thus, while there is similar physical behavior to the full model, the large γ current is quantitatively very different for the Markovian master equation, which reflects its lack of validity in this regime.

Appendix C: Bounds on the Validity of the Markovian Master Equation
This section quantifies the regimes where the master equation is physically valid. The Markovian master equation is always mathematically valid in that it gives proper quantum evolution. However, as we will see, it does not accurately represent the equilibrium state at larger values of γ, which, e.g., leads to spurious currents and a break down of detailed balance.

Broadening and the Fermi Level
The Markovian master equation broadens extended reservoir states across a wide range of ω when the relaxation rate γ is large. That is, even with the Fermi level fixed in the isolated extended reservoir, there is excessive electron occupation beyond this level in the open system (succinctly, the Markovian equation occupies the states then broadens them, rather than broadening then occupying). As an example, Fig. 1 shows the density of states times the occupation for both the full and Markovian approaches. For small γ, Fig. 1(a), the relaxation is weak enough that the states are still relatively localized in energy. However, for large γ, Fig. 1(b), the electronic occupation is smeared too much, and this allows current to flow even without a drop in chemical potential. This difference is most apparent for states at the Fermi level, ω k ≈ ω F , as this is where the distribution is most rapidly changing.
In the two approaches, only the lesser Green's functions are different, and these only differ by the distribution function, f (ω) compared to f (ω k ). We quantify the error, ∆, between the two by the integrated absolute difference, To upper bound this error, we will use two features of ∆. One is that the error is maximal when the state ω k is at the Fermi level, ω k = ω F . Two is that the Fermi distribution can be replaced by the piecewise continuous in order to obtain a bound of the error at the Fermi level. That is, this replacement has a greater absolute difference to f (ω F ) = 1/2 than the original distribution function for all ω. Using these two features, the error for any ω k is bounded by The condition for this error to be small is then γβ 1. This can be interpreted as requiring that the broadening due to the relaxation must be smaller than the broadening caused by thermal processes. As well, it has a simple, intuitive mathematical meaning: The maximal slope of the Fermi-Dirac distribution should be much smaller than γ, so that the γ-induced smearing has no significant effect on the occupation.
Expanding on the above brief account, we can show that the error is maximal when ω k is at the Fermi level by extremizing ∆, which can be rewritten aŝ When ω k is at the Fermi level, the term containing f (ω k ) integrates to zero, as it is a symmetric function multiplied by an antisymmetric function around ω k . Also when ω k is at the Fermi level, f (ω) is symmetric around ω F (which comes from the relationship f Thus, the integrand is antisymmetric and it evaluates to zero. This therefore gives an extremum in the error. Moreover, when ω k < ω F , the slope is positive and when ω k > ω F , the slope is negative, therefore the error is maximal at the Fermi level. In the limiting cases, ω k → ±∞, the total error ∆ is zero, meaning that there is agreement when the state is far away from the Fermi level, as excess smearing has no consequences when the distribution is flat. Additionally, bounds on the relaxation rate and the number of states required can be found for a given system by numerically integrating the spectral density above the Fermi level for the zero temperature case. This yields a direct measure of the improperly occupied high-energy states, even when they do not directly contribute to the electronic current.

Asymmetric Reservoirs
Next, we will examine a system in between asymmetric reservoirs. When the asymmetry of the L and R extended reservoir is due to a shift in their relative energies by δ, then we have ω k = ω k + δ/2 for k ∈ L and ω k = ω k −δ/2 for k ∈ R. Solving the equation of motion, Eq. (B2), shows that there is a non-zero current for the zero applied bias (V = 0). Note that δ is a parameter that quantifies the asymmetry of the system.
Physically, the Fermi level on both reservoirs are identical, so the steady-state current between them should be zero. When the reservoir density of states is broadened by γ, the replacement of f (ω) by f (ω k ), as the Markovian model does, results in some electronic occupation above the Fermi level. This can give rise to an electronic current. Figure 2 shows the steady-state current I as a function of the relaxation γ for increasing values of the applied bias V .
To simplify the analytic forms for the steady-state current, we will introduce a system where the extended reservoirs are "Markovian" rather than 1D. The density of states for the 1D case is approximately constant near ω = 0 as long as the bandwidth is large enough. For the following, the states in the reservoirs are equally spaced between frequencies ω min and ω max and then shifted by the asymmetry parameter δ, given k ∈ {1, . . . , N r }: ω k = (ω min − ∆/2) + ∆k + δ/2 where the state spacing is ∆ = (ω max − ω min )/N r and the couplings are constant v k = 8J 2 /(2πN r ). Taking the limit as N r → ∞, the integral of the single particle Green's functions yields for L and similarly for R. Here the bandwidth is W = ω max − ω min . The calculation of the steady-state current then continues as previously. At zero temperature (β → ∞), expanding the current for large bandwidth (and small V , δ), the leading terms are The first term gives the linear response current. For small γ, it is eV /2π . The second term gives a residual, unphysical current. So long as the asymmetry is small compared to the bandwidth or γ is sufficiently small, this term will be negligible compared to the linear response contribution. Moreover, to accurately calculate results in the Landauer regime for steady states, γ needs to be of the order of W/N r (this condition depends on the details of the setup, as discussed in the main text regarding the small γ regime). Thus, so long as then the simulation will accurately predict steady state behavior. In other words, one needs N r J 2 δ /W 2 V (or, for our particular example in the main text, N r δ/V ), which is a condition that is basically always fulfilled, to ensure the simulation gives accurate results in the Landauer regime. This covers both small and intermediate regimes discussed throughout the text.
The value of σ is related to the coupling between the reservoir and system and describes the particle flow rate. In the general case, this would be a function of the system-reservoir coupling and the total Green's function. For the example system, however, σ scales as 1/N r . In linear response, the conductance is independent of J and W , so the sum of the total current from all N r states must be a constant.
In the steady state, all four of these currents-I L k , I R k , I LS k , I SR k -must be equal, giving the solution The quantity γ L γ R /(γ L +γ R ) is the "reduced γ" between the two reservoirs and is simply γ/2 when the relaxation rates are equal. The small γ (γ σ) approximation for Eq. (E3) with equal relaxation rates is I k ≈ e(γ/2)(f L k − f R k ), as used in the main text. This rate argument can be extended to the full N r = 1 system by solving a system of 4N r equations for each of the incoming and outgoing currents and then equating the total current into n S , k I LS k , with the total current leaving k I SR k . When γ L = γ R = γ, this yields the current through the system as For small γ, this recovers the expression in the main text. Additionally, the occupation of the central site, n S , is found to be equal to the mean target filling, k (f L k + f R k )/(2N r ). With a symmetrical distribution of reservoir state energies and bias window (as in the case of the main text), then the central site onsite density is 1/2.
Lastly, a full solution involving a different σ k for each reservoir mode is obtainable. Assuming that the σ k are symmetric between the left and right sides, the total current flow recovers Eq. (E4) to the lowest order of γ. The relaxation strength limits the total current through the system and, so long as γ is sufficiently small, then the current is independent of the system-reservoir coupling. * mpz@nist.gov 1 H. Haug and A.-P. Jauho, Quantum kinetics in transport and optics of semiconductors (Springer, 1996). 2 A. Jauho, N. Wingreen, and Y. Meir, Phys. Rev. B 50, 5528 (1994).