Multi-terminal far-from-equilibrium thermoelectric nano-devices in the Kondo regime

TL;DR

This paper investigates the thermoelectric properties of quantum dot devices under far-from-equilibrium conditions, revealing how strong electron interactions influence efficiency and power output in the Kondo regime.

Contribution

It introduces generalized Seebeck coefficients for non-linear regimes and analyzes the impact of Coulomb interactions on thermoelectric performance.

Findings

Maximum efficiency occurs far from equilibrium.

Strong Coulomb interactions reduce maximum power and efficiency.

Kondo regime effects are significant in thermoelectric performance.

Abstract

The quest for good thermoelectric materials and/or high-efficiency thermoelectric devices is of primary importance from theoretical and practical points of view. Low-dimensional structures with quantum dots or molecules are promising candidates to achieve the goal. Interactions between electrons, far-from-equilibrium conditions and strongly non-linear transport are important factors affecting the usefulness of the devices. This paper analyses the thermoelectric power of a two-terminal quantum dot under large thermal $\Delta T$ and voltage $V$ biases as well as the performance of the three-terminal system as a heat engine. To properly characterise the non-linear effects under these conditions, two different Seebeck coefficients are introduced, generalizing the linear response expression. The direct calculations of thermally induced electric and heat currents show, in agreement with recent work, that the efficiency of the thermoelectric heat engine as measured by the delivered power is maximal far from equilibrium. Moreover, the strong Coulomb interactions between electrons on the quantum dot are found to diminish the efficiency at maximum power and the maximal value of the delivered power, both in the Kondo regime and outside of it.

Figures

Equations

$$S=-{\left(\frac{V}{\Delta T}\right)}_{I=0},$$

$$H=\sum _{\lambda ,k,\sigma }{\epsilon }_{\lambda k}{n}_{\lambda k\sigma }+\sum _{\sigma }{\epsilon }_{\sigma }{n}_{\sigma }+U{n}_{\uparrow }{n}_{\downarrow }+\sum _{\lambda k\sigma }(V_{\lambda k\sigma }{c}_{\lambda k\sigma }^{†}{d}_{\sigma }+V_{\lambda k\sigma }^{\ast }{d}_{\sigma }^{†}{c}_{\lambda k\sigma }),$$

$$I_{\lambda }=-e⟨\frac{d{N}_{\lambda }}{dt}⟩=\frac{ie}{\hslash }⟨[N_{\lambda },\hat{H}]⟩$$

$$J_{\lambda }=\frac{i}{\hslash }⟨[H_{\lambda },\hat{H}]⟩-{\mu }_{\lambda }\frac{i}{\hslash }⟨[N_{\lambda },\hat{H}]⟩,$$

$$\displaystyle I_{\lambda}(t)=\frac{2e}{\hbar}\sum_{{k}\sigma}{\rm Re}\bigg{[}{V}_{\lambda{k}\sigma}G^{<}_{\sigma,\lambda{k}\sigma}(t,t)\bigg{]},$$

$$\displaystyle J_{\lambda}(t)=\frac{2e}{\hbar}\sum_{{k}\sigma}(\varepsilon_{\lambda{k}}-\mu_{\lambda}){\rm Re}\bigg{[}{V}_{\lambda{k}\sigma}G^{<}_{\sigma,\lambda{k}\sigma}(t,t)\bigg{]}.$$

$$I_{\lambda }$$

$$J_{\lambda }$$

$$J_{\lambda }=J_{\lambda }^E-{\mu }_{\lambda }{I}_{\lambda }.$$

$$⟨n_{\sigma }⟩=-i\int \frac{dE}{2\pi }{G}_{\sigma }^{<}(E).$$

$$0=\frac{d⟨n_{\sigma }⟩}{dt}=⟨\frac{d{n}_{\sigma }}{dt}⟩=-⟨\frac{i}{\hslash }[n_{\sigma },H]⟩.$$

$$⟨n_{\sigma }⟩=i\int \frac{dE}{2\pi }\frac{{\Gamma }_{\sigma }^L{f}_L(E)+{\Gamma }_{\sigma }^R{f}_R(E)}{{\Gamma }_{\sigma }^L+{\Gamma }_{\sigma }^R}[G_{\sigma }^r(E)-G_{\sigma }^a(E)].$$

$$I=\frac{2e}{\hslash }\sum _{\sigma }{\tilde{\Gamma }}_{\sigma }\int \frac{dE}{2\pi }[f_L(E)-f_R(E)]\mathrm{Im}{G}_{\sigma }^r(E),$$

$$J_L=\frac{2e}{\hslash }\sum _{\sigma }{\tilde{\Gamma }}_{\sigma }\int \frac{dE}{2\pi }(E-{\mu }_L)[f_L(E)-f_R(E)]\mathrm{Im}{G}_{\sigma }^r(E),$$

$$J_R=\frac{2e}{\hslash }\sum _{\sigma }{\tilde{\Gamma }}_{\sigma }\int \frac{dE}{2\pi }(E-{\mu }_R)[f_R(E)-f_L(E)]\mathrm{Im}{G}_{\sigma }^r(E).$$

$$J_L+J_R+({\mu }_R-{\mu }_L)I=0$$

$$J=J_H+J_L+J_R,$$

$$I=I_L+I_R.$$

$$P=IV,$$

$$\eta =\frac{IV}{J_H}.$$

$$⟨⟨d_{\sigma }|d_{\sigma }^{†}⟩{⟩}_E=\frac{1+I_d(E)[⟨n_{\bar{\sigma }}⟩+b_{1\bar{\sigma }}-b_{2\bar{\sigma }}]}{E-{\epsilon }_{\sigma }-{\Sigma }_{0\sigma }+I_d(E)[{\Sigma }_1^T+{\Sigma }_2^T-(b_{1\bar{\sigma }}-b_{2\bar{\sigma }}){\Sigma }_{0\sigma }]},$$

$$I_d(E)=\frac{U}{E-{\epsilon }_{\sigma }-U-{\Sigma }_{0\sigma }-{\Sigma }_{\sigma }^{(1)}-{\Sigma }_{\sigma }^{(2)}}.$$

$${\tilde{\gamma }}_{\alpha }=2\pi \sum _{|f⟩}|⟨T(E_{\alpha })⟩{|}^2\delta (E_{\alpha }-E_f),$$

$${\tilde{\epsilon }}_d={\epsilon }_d+{\Sigma }_1^T({\tilde{\epsilon }}_d)+{\Sigma }_2^T({\tilde{\epsilon }}_d).$$

$$S_n=-{\left(\frac{V}{\Delta T}\right)}_{I(\Delta T,V)=0},$$

$$I(\Delta T,V)=L_{11}V+L_{12}\Delta T,$$

$$S=\frac{L_{12}}{L_{11}}.$$

$$S_{d}=-\left(\frac{\partial V}{\partial\Delta T}\right)_{I=\mathrm{const}}=-\left(\frac{\partial I}{\partial\Delta T}\right)\bigg{/}\left(\frac{\partial I}{\partial V}\right).$$

$$T_K=\frac{\sqrt{U{\Gamma }_N}}{2}\exp \frac{\pi {\epsilon }_d({\epsilon }_d+U)}{{\Gamma }_{N}U},$$

$$\left(\frac{\partial I}{\partial \Delta T}\right)$$

$$D(V)=[f_L(E,T)-f_R(E,T)]=-D(-V)$$

Full text

1 Introduction

Thermoelectricity, the invention of the 19th century, is still at the forefront of research due to its importance for space exploration and automotive industry [1], and many more branches of modern technology both at large [2] and small [3] scales. Attempts to find high performance thermoelectric bulk materials [4], including those with topologically non-trivial [5] band structure, have seen limited progress. In the last few decades, a lot of attention has been put forward towards nano-devices [6] and molecular systems [7, 8], utilizing quantum effects to boost their thermoelectric performance towards the thermodynamic limit [9].

When a temperature gradient $\nabla T$ (or a temperature difference $\Delta T$) is established across a bulk material, the voltage $V$ is generated. The response of the isotropic system is quantitatively characterized by the Seebeck coefficient [10, 11]

[TABLE]

defined under the condition of zero charge current. The same formal definition is valid for a nano-structure with two external electrodes (see Fig. 1(a)). However, Eq. (1) has to be generalized for more complicated geometries with several electrodes. In fact the thermoelectric characterization of nano-structures with three electrodes (as, e.g., the one shown in Fig. 1(b)) requires the definition of the whole matrix of Seebeck coefficients. Such geometries also allow studies of non-local effects [12, 13].

In bulk systems the linear approximation is generally a valid simplification [14]. Under such a condition the thermoelectric figure of merit, $ZT=GS^{2}T/\kappa$, where $G$ is the conductance, $\kappa$ thermal conductance, and $S$ the Seebeck coefficient, is viewed as the most important factor deciding on the usefulness of the material as heat to electricity converter: namely, the efficiency of the thermoelectric heat to electricity converter is given by $\eta=\eta_{C}(\sqrt{ZT+1}-1)/(\sqrt{ZT+1}+1)$, where $\eta_{C}$ is the Carnot efficiency. In nano-structures containing quantum dots, however, we are virtually always dealing with a non-linear situation, as mentioned earlier [15, 16, 17] and carefully discussed recently [6]. The small ratio of the thermalization length to the sample length in bulk systems, and the opposite limit in nano-structures is responsible for their different behaviour. This has a profound effect on the analysis of small heat engines, and implies that the thermoelectric figure of merit $ZT$ is not a useful parameter to judge the efficiency of a device. This also means that nano-systems with a large figure of merit [18] $ZT$ may in fact feature a small efficiency [12, 15, 16, 19, 20].

From a basic physics point of view, the Seebeck coefficient provides additional and novel information about the investigated system compared to that obtained from the electrical conductivity. In the simplest case, the latter depends on the value of the density of states $N(E_{F})$ at the Fermi level, while the former ‘measures’ its slope. The thermopower has been shown to be directly related to the entropy flowing between different parts of the system [11]. In fact, the entropy of nanosystems has been recently measured [21] by means of thermoelectric transport. A strong increase of $S$ in nano-devices and nano-structured materials has been observed [22], in agreement with the earlier proposal [23]. Recent studies [24] have shown that doping or nano-structuring bulk thermoelectric materials may lead to the required modifications of the energy spectrum close to the Fermi energy but also to localization of states which deteriorates the systems performance. The details of the nano-structures may be important, $e.g.$, for systems based on molecules the actual value of the thermopower depends on the length of the molecule [26, 25].

It has been found [9, 16] that non-linear working conditions can favourably affect the performance of heat engines. Indeed recent years have witnessed increased activities in the theoretical studies of transport beyond the linear approximation [27, 28, 29, 30, 31, 32, 33]. The main motivation of that kind of work is related to the desire to find devices with improved thermoelectric performance, which often can be achieved in the non-linear regime only [34]. The early developments in non-linear quantum transport driven by thermal gradients and/or voltage biases have been reviewed recently [36]. Even more recent work includes Refs. [37, 38, 39, 40, 41, 42, 43].

The aim of the present work is to study thermoelectric transport of quantum dot based nano-devices by means of the non-equilibrium Green function approach, taking Coulomb interactions and non-linear effects into account. In the non-linear regime more general definitions of the Seebeck coefficient than that given in Eq. (1) are needed. These are especially important if an externally applied voltage $V$ is present while the system is thermally biased. In the presence of interaction and at low temperature, the Kondo effect is expected to dominate the transport of the system at hand. The existence of the Kondo effect in quantum dots has been predicted a long time ago [44, 45], and later observed experimentally [46, 47, 48].

As the theoretical treatment of an interacting quantum dot is of importance in itself, we also present in some detail the semi-analytical technique proposed recently by Lavagna [49]. It is based on the equation of motion method [50] for the non-equilibrium (or Keldysh) Green functions [51]. She has proposed a few important additions, which allow to properly describe the Kondo effect [52] in linear and non-linear regimes, i.e., under large voltage and temperature differences between the electrodes and for the particle-hole symmetric model.

In the present work we extend our previous studies of non-equilibrium screening effects [53], which is important for a better understanding of the interaction effects in far-from-equilibrium situations as encountered in three-terminal two-quantum-dot high-efficiency heat engines [34, 35].

The use of the equation of motion (EOM) technique to study the single impurity Anderson model has a long history. It started with the work of Anderson himself [54], and has been pursued by others [55, 56, 57, 58], mainly in the context of single impurities in metals. Some of the attempts at generalizing the original EOM and the decoupling procedures have been summarized in [59]. The technique was later applied to the Kondo effect in quantum dots coupled to external electrodes [44, 45]. More recent work includes studies of more complicated systems and geometries like those with one normal and one ferromagnetic or superconducting electrode [60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76].

It has to be mentioned that there is a vast amount of experimental data measuring the thermopower in various nano-systems [77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103]; these provide additional motivation for the theoretical work. Very recently the three terminal two quantum dot system identical to one of the systems studied here has been found [104] to generate a power of $0.13$ fW with efficiency of $\approx 0.15$ of the Carnot value.

The organization of the paper is as follows. In Sec. 2 we present the model and our approach of calculating the charge and heat currents by means of the Keldysh non-equilibrium Green function (GF) technique. The resulting spectral function of the quantum dot is discussed in Sec. 3 at low temperatures and for various conditions including particle-hole symmetry and non-equilibrium. The thermally induced currents are the subject of Sec. 4 where also rectification properties of the non-symmetrical device are mentioned. Three possible definitions of the Seebeck coefficients, including two valid in the non-linear regime, are proposed in Sec. 5. The Coulomb interaction between electrons on the three-terminal quantum dot heat engine is found (Sec. 6) to diminish the performance of the device in question. The appearance of the Kondo effect in the heat engine shows up as a two-leaf structure of the performance diagram on the efficiency-power plane. We end with summary and conclusions in Sec. 7.

2 Modeling the device and calculating currents

Here we discuss the simple geometry where the system consists of a quantum dot tunnel-coupled to a few normal electrodes as shown in Fig. 1. The Hamiltonian of the system is written as

[TABLE]

where $n_{\lambda{k}\sigma}=c^{\dagger}_{\lambda{k}\sigma}c_{\lambda{k}\sigma}$ and $n_{\sigma}=d^{\dagger}_{\sigma}d_{\sigma}$ denote particle number operators for the leads and the dot, respectively. The operators $c^{\dagger}_{\lambda k\sigma}(d^{\dagger}_{\sigma})$ create electrons in respective states $\lambda{k}\sigma$ $(\sigma)$ in the leads $\lambda=L,R,H$ (on the dot). The wave vector $k$ denotes the Bloch state in the electrodes, the spin is $\sigma=\pm 1\;(\uparrow,\downarrow)$, and $\varepsilon_{\sigma}=\varepsilon_{d}+\sigma\mu_{B}B$, where $B$ is the magnetic field introducing Zeeman splitting, $\mu_{B}$ denotes the Bohr magneton, and $\varepsilon_{d}$ is the dot electron energy level. The Hubbard parameter $U$ describes the repulsion between two electrons on the dot.

The charge current in the electrode $\lambda$ is calculated as the time derivative of the average charge in that electrode $\langle N_{\lambda}\rangle=\sum_{k\sigma}\langle n_{\lambda{k}\sigma}\rangle$, and reads

[TABLE]

where $\langle...\rangle$ denotes the statistical average. The calculation of the heat fluxes, $J_{\lambda}$, follows that of the charge:

[TABLE]

where $H_{\lambda}=\sum_{{k},\sigma}\varepsilon_{\lambda{k}\sigma}n_{\lambda{k}\sigma}$ is the energy operator for the electrode $\lambda$. Calculating the commutators and introducing appropriate GFs, one obtains[51]

[TABLE]

The final expressions for the (stationary) currents can easily be written in the following general form: [51]

[TABLE]

The parameters $\Gamma_{\sigma}^{\lambda}(E)=2\pi\sum_{{k}}|V_{\lambda{k}\sigma}|^{2}\delta(E-\varepsilon_{\lambda{k}})$ describe the coupling between the dot and the respective electrodes. The Green functions $G^{i}_{\sigma}(E)=\langle\langle d_{\sigma}|d^{\dagger}_{\sigma}\rangle\rangle^{i}_{E}$ with $i=r,a,<$ determine the spectral properties of the quantum dot as well as the transport properties of the whole system.

Having in mind non-equilibrium transport induced by a voltage bias and/or a temperature difference, we keep the dependence of the Fermi distribution functions $f_{\lambda}(E)$ on the electrode $\lambda$ via its chemical potential $\mu_{\lambda}$ and its temperature $T_{\lambda}$. The heat current (8) can be written as the difference between the energy current $J^{E}_{\lambda}$ and the charge current $I_{\lambda}$:

[TABLE]

Importantly, in the steady state and in the wide band approximation ($i.e.$, for energy independent couplings: $\Gamma_{\sigma}^{L,R}(E)=\Gamma_{\sigma}^{L,R}$), one can derive [49] exact expressions for the currents. First, from $\langle n_{\sigma}\rangle=\langle c^{\dagger}_{\sigma}(t)c_{\sigma}(t)\rangle$ and the definition of the lesser Green function, one has [51]

[TABLE]

The derivation then makes use of the fact that in the steady state

[TABLE]

Working out the commutator in (11), using (10) and the Langreth theorem [105], it is straightforward to derive the following ‘self-consistency’ condition [49]:

[TABLE]

Let us underline again that the expression (12) is exact under the proviso that the couplings to the leads are energy independent, $\Gamma_{\sigma}^{\lambda}(E)\equiv\Gamma^{\lambda}_{\sigma}=\mathrm{const}$. If this condition is violated, as it might be the case in graphene [106, 107, 108], hybrid systems with one (or both) of the electrodes being a superconductor, e.g., d-wave [109], other approaches are needed. For models with energy dependent couplings one still has to rely on approximate relations between the lesser self-energy $\Sigma^{<}_{\sigma}(E)$ and the retarded one, $\Sigma^{r}_{\sigma}(E)$, making use of Ng’s approximation [60], a generalisation of the fluctuation-dissipation theorem [28, 40], or using the equation of motion for the lesser GF [110] and suitable decoupling. For recent attempts including energy-dependent couplings, see the paper [75].

With the above (exact in the steady state and for constant $\Gamma$’s) formulation, the charge current across the two-terminal system can be written in terms of the retarded GF only:

[TABLE]

where $\tilde{\Gamma}_{\sigma}=\Gamma_{\sigma}^{L}\Gamma_{\sigma}^{R}/(\Gamma_{\sigma}^{L}+\Gamma_{\sigma}^{R})$. We stress that the above formula is an exact representation of the current in terms of the retarded GF of the interacting Hamiltonian, which have to be calculated either numerically or analytically. Here we adopt the latter approach, though obviously calculating the GF analytically requires some approximations, see the details given in the Appendix. In brief, the decouplings we are using ensure that the GFs are formally exact up to second order in the couplings. Additionally, we correct the result by phenomenologically introducing the lifetimes [49] of various states, which is important to capture the Kondo effect even in the particle-hole symmetric model.

The expression (13) is analogous to the well-known Meir-Wingreen formula [111]. A direct calculation leads to $I=I_{L}=-I_{R}$, which expresses the current conservation. Similar expressions can be derived for the heat current flowing from the left,

[TABLE]

and right electrode:

[TABLE]

It can be verified that

[TABLE]

in agreement with energy conservation. Here $I=I_{L}=-I_{R}$ is the charge current, and $\dot{Q}=J_{L}+J_{R}$ the total heat current leaving the leads. For later use we define the power, $P=(\mu_{R}-\mu_{L})I/e$, and the voltage bias, $V=(\mu_{R}-\mu_{L})/e$, across the system.

Three-terminal quantum dot devices have been proposed [34, 27, 12, 20] as efficient, easy to control heat engines, and analysed in equilibrium and beyond, both for non-interacting and interacting quantum dots. Our general formulas for the charge and energy currents (8) flowing out of an arbitrary electrode $\lambda$ allow for an easy application to the three-terminal system like the one shown in Fig. 1(b). In the notation of the figure the heat is flowing from the hot to two cold electrodes; assuming that no charge current flows out of or into the hot electrode. The flow of charge is dictated by the gate bias of the two quantum dots, which in the discussed system act as energy filters. The total heat current in the system is

[TABLE]

while the charge current (assuming $I_{H}=0$) reads

[TABLE]

Introducing a ‘load’ between the two cold electrodes, in the figure shown as an external voltage $V$, one defines the power delivered by the system as

[TABLE]

and the system energy harvesting efficiency as

[TABLE]

To calculate currents we only require the retarded Green function, which is still a complicated issue with no exact analytic solution available. There exist a few numerically exact approaches: for example, among the many techniques applied to study the single-impurity Anderson model [52] the quantum Monte Carlo [112] and the numerical renormalization group method [113] should be mentioned. However, our aim is to use the analytic expressions which have been proven to be quantitatively correct [49], and valid far from equilibrium; in addition, in our opinion these are physically more transparent compared to purely numerical results, as discussed below.

One finds the following final expression for the on-dot GF:

[TABLE]

where

[TABLE]

Details and definitions are given in the Appendix. The expression for the Green function (21) agrees with that obtained earlier [49].

It has been proposed by Lavagna [49] that for the correct description of the Kondo effect one has to supplement the above set of equations by two ingredients, see Appendix. First, one should introduce the finite lifetimes of singly and doubly occupied states on the dot. Second, it is important to take care of the many-body renormalization of the dot energy level $\varepsilon_{d}$. The first improvement consists in replacing the $E+i0$ terms in the definitions of various Green and correlation functions by $E+i\tilde{\gamma}_{\alpha}$, where the inverse lifetimes $\tilde{\gamma}_{\alpha}$ of the excited states $\alpha=|\sigma\rangle,|\uparrow,\downarrow\rangle$ are due to higher order processes, and can be calculated up to the desired order via the generalized Fermi rule as

[TABLE]

with $T(E)=\hat{V}+\hat{V}g(E)\hat{V}+\cdots$ being the scattering matrix, and $\hat{V}$ denoting the part of the Hamiltonian describing the coupling between quantum dot and reservoirs. The second improvement amounts to the replacement of $\varepsilon_{d}$ by $\tilde{\varepsilon}_{d}$, to be calculated self-consistently from

[TABLE]

The GF (21) has been shown to fulfil the unitarity limit and describe quantitatively correctly the Kondo effect [49] even of particle-hole symmetric systems in equilibrium and out-of-equilibrium.

3 Kondo effect in equilibrium and far from equilibrium

For illustrative purposes and in order to introduce the framework, we discuss in this section the properties of an interacting quantum dot between two normal electrodes, as illustrated in Fig. 1(a). We start the presentation of the results by showing the density of states (DOS) of the quantum dot in various regimes. From now on, we neglect the spin dependence of the couplings, and slightly change the notation: $\Gamma^{L}_{\uparrow}=\Gamma^{L}_{\downarrow}=\Gamma_{L}$, $\Gamma^{R}_{\uparrow}=\Gamma^{R}_{\downarrow}=\Gamma_{R}$. In addition, we measure all energies in units of $\Gamma_{0}\equiv\Gamma_{L}$. The particle-hole symmetric case is of special importance as it is well known [49] that all previous attempts to use the EOM method failed in this case [40, 65, 70]. In Fig. 2 we present the DOS in the particle-hole symmetric situation with $\varepsilon_{d}=-4\Gamma_{0}$ and $U=8\Gamma_{0}$. Both lower and upper Hubbard bands centered in energy around $\varepsilon_{d}$ and $\varepsilon_{d}+U$ are clearly visible. At the same time, the Abrikosov-Suhl, sometimes also called Kondo resonance develops at the chemical potential.

The external bias shifts the chemical potentials of the leads $\mu_{L/R}=\mu\pm eV/2$ to new positions, and the Kondo peak splits into two with each resonance pinned to the chemical potential of the lead. In particular, Fig. 2 shows the evolution of the two peaks with temperature. At $T=0.001\Gamma_{0}$, they are very sharp, while at $T=0.1\Gamma_{0}$ they are reduced but still clearly visible. The changes of the Kondo resonance with bias and temperature are crucial to understand the behaviour of the thermally induced current and the (non-linear) thermoelectric power as discussed in the next sections. At still higher temperatures (not shown) both Kondo peaks vanish altogether, and only lower and upper Hubbard bands survive. The Kondo resonance is due to spin flip processes on the dot while the Hubbard sub-bands are related to charge fluctuations. This explains the relative robustness of Hubbard sub-bands with increasing temperature, and the fragility of the Kondo correlations. However, voltage and temperature affect the Kondo resonance in a different way. While a voltage $eV>T_{K}$ leads to a splitting of the resonance in two sub-peaks with concomitant decrease of their maximum, the finite temperature only lowers the height of the peak and broadens it. All these features are well reproduced by the used decoupling scheme. Outside the particle-hole symmetric point, the voltage-split Kondo resonances are not symmetric anymore. This is mainly due to the closer proximity of one of the resonances to the lower Hubbard band. The individual resonances are pinned to the Fermi levels of the respective electrodes. This is well visible in Fig. 2(b), and also in Fig. 3(a).

We now add a thermal bias to one of the electrodes, with focus on the question of how the density of states evolves with temperature difference. We start with the split Kondo resonance as shown in Fig. 3(a). The temperature of the right lead is kept constant, at $T_{R}=T$, while we gradually increase the temperature of the left electrode, $T_{L}=T+\Delta T$. The panel (a) of the Fig. 3 shows the voltage-split Kondo peaks for $\Delta T=0$, for easy comparison, while the panel (b) demonstrates the evolution of both peaks with increasing $T_{L}$. The peak pinned to the chemical potential of the left electrode (appearing at $E=\mu_{L}$) is strongly affected by the temperature bias, it broadens and vanishes with increasing temperature. One observes only small changes of the other Kondo peak. Thus under voltage and temperature bias the peaks’ heights and widths depend mainly on the conditions ($V$, $\Delta T$) at the electrode with the chemical potential to which it is pinned. For diminishing external voltage bias both peaks overlap in energy, and one observes a single feature in the density of states corresponding to the average temperature $(T_{L}+T_{R})/2$ of the system.

4 Thermocurrents and their rectification

Recently there has been great interest in the control of heat transport [114] and its use to process information in analogy to electronics. To this end, one typically uses physical systems which are driven by a temperature bias. The field develops fast with many theoretical proposals and related experimental realisations [91, 92, 115, 116, 117, 118, 119, 120]. For a recent overview of the experimental advances in thermal rectification with particular attention to nanoscale devices, the reader is invited to consult the recent review [121], and references therein.

Motivated by the importance of thermal rectification and the fact that it is typically induced by the temperature difference across the nanodevice, we first calculate the thermo-currents, i.e., the thermally induced currents, and later focus on rectification phenomena in systems with broken symmetry. In Fig. 4 we show the current $I(T,\Delta T)=I(T,\Delta T,V=0)$ for the interacting system with $U=10\Gamma_{0}$, $T=0.01\Gamma_{0}$, and for a few values of $\varepsilon_{d}$. As a function of temperature bias $\Delta T$ the thermo-current shows an interesting set of zero values. The number of zeros depends on the interaction strength and the gate voltage, i.e., $\varepsilon_{d}$. To understand this behaviour let us consider a particular value of $\varepsilon_{d}$. For $\varepsilon_{d}=-\Gamma_{0}$ and the assumed $U$, the Kondo temperature is $T_{K}\approx 0.09\Gamma_{0}$. In such a situation the on-dot level is slightly below the chemical potential $\mu=\mu_{L}=\mu_{R}$, and the doubly occupied state $2\varepsilon_{d}+U=8\Gamma_{0}$ is far above it, so the Kondo resonance develops. The schematic energy diagram (without Kondo resonance) is illustrated in Fig. 4(b). For a small average temperature, the temperature difference $\Delta T$ induces the current flowing mainly via the Kondo resonance (located at $E\approx\mu=0$) out of the left electrode, and $I(\Delta T)>0$. For a temperature difference around $\Delta T\approx\Gamma_{0}>T_{K}$, the contribution from the smeared Kondo resonance gets relatively smaller than the contribution from the other states. For the considered, relatively large temperature difference it is mainly holes which flow from the left electrode or electrons mainly flowing from the right to the left due to the difference of the Fermi functions. For still higher $\Delta T$, the thermally excited electrons at the left electrode start to flow from left to right through the upper doubly occupied level, and the current is again positive. In a similar manner one can understand the behaviour of the current at other gate voltages [40], and other $\varepsilon_{d}$. The effect has been thoroughly discussed recently [36], and measured experimentally [95].

According to Eqs. (13) and (14), the charge and heat currents through the quantum dot directly depend on the on-dot Green function. In the linear approximation the matrix of kinetic coefficients is symmetric, according to Onsager’s relations. Outside linearity one expects the violation of these symmetries, and a non-linear dependence of the currents on the voltage and/or thermal bias. An even more interesting behaviour is expected if the device shows an internal asymmetry. The asymmetry of a device induces current rectification effects as discussed earlier in the context of transport via molecular junctions in the pair-tunnelling regime [122]. In the two-terminal structure with the same spectrum of the left and right electrode the only source of asymmetry are the couplings to the leads $\Gamma_{L(R)}$.

To show the thermally induced charge and heat current rectification in the Kondo regime, we assume an asymmetry in the couplings with $\Gamma_{R}/\Gamma_{L}=2,3$. We apply first the temperature bias $\Delta T$ to the left electrode and calculate the current across the system, $I(T,\Delta T)=I_{L}(T_{L})-I_{R}(T_{R})$, and later apply the same temperature difference to the right electrode, keeping the average temperature of the system constant. For the studied temperature range we are deep in the Kondo regime since $T_{\rm{av}}=(T_{L}+T_{R})/2\ll T_{K}$ ($\approx\Gamma_{0}$). In Fig. 5(a) the charge current calculated for $T_{L}=T+\Delta T$, $T_{R}=T$ is compared to the one flowing in the system with the same temperature bias but applied to the right electrode. In both cases the base temperature of the device is the same and equals $T$, so the difference between the curves of different colours directly demonstrates the rectification effect. Thin (thick) curves correspond to an asymmetry of 2 and 3, respectively. The current flow is in the direction consistent with temperature difference, and for sufficiently large values of $\Delta T$ their values are visibly different. The rectification factor amounts to about 30% and depends on $\Delta T$.

Figure 5(b) shows a similar rectification of the heat current for the same model parameters and the same asymmetry, $\Gamma_{R}/\Gamma_{L}=2,3$. However, for clarity the curves calculated for $\Gamma_{R}/\Gamma_{L}=3$ have been shifted vertically by $+2\cdot 10^{-3}$ for better visibility. The heat current rectification effect for this set of parameters is at least an order of magnitude smaller than that for the charge and equals 1% to 2% in the parameter range shown in the figure. This is better visible in the inset which shows the $\Delta T>0.12\Gamma_{0}$ part of the main figure. Calculating the charge and heat currents we have neglected the possible dependence of the chemical potentials of the electrode with temperature.

5 Linear and non-linear thermopower and the role of asymmetry

Non-linear effects are expected to be important in all nano-devices whatever the applied bias [6, 17]. The suitable generalisation of the definition (1) of the Seebeck coefficient to non-linear situations reads

[TABLE]

where $\Delta T$ is the temperature difference applied to the nanostructure. The traditional way to measure the Seebeck coefficient assumes the application of the temperature bias $\Delta T$, and finding such voltage $V$ where the current across the device vanishes, $I(\Delta T,V)=0$. In this formulation the only source of non-linearity is directly given by the value of $\Delta T$, assumed to be large enough to preclude the linear expansion of the charge current $I(\Delta T,V)$ in first powers of $V$ and $\Delta T$.

On the contrary, when the linear expansion is valid ($V\to 0$ and $\Delta T\to 0$), we have

[TABLE]

and the corresponding thermopower is given by the ratio of kinetic coefficients

[TABLE]

Expanding Eq. (13) for small $V$ and $\Delta T$ up to linear order, one finds an explicit expression for $S$.

Sometimes the definition (25) of the non-linear Seebeck coefficient $S_{n}$ is extended to the differential one, $S_{d}$, calculated formally for constant current flowing as a result of the external voltage $V$. This $S_{d}$ measures the response of the system with current flow to the change in temperature difference $\partial\Delta T$. At the applied external voltage $V$ and temperature difference $\Delta T$, $S_{d}$ is defined [33, 76] as the derivative

[TABLE]

As argued earlier, $S_{d}$ should also be accessible experimentally [33]. In that paper the response $\left({\partial I}/{\partial\Delta T}\right)$ has been calculated under the additional approximation that, at a given but arbitrary external voltage, a very small temperature bias $\Delta T$ is applied to the system. The definition (28) can be viewed as analogous to the differential conductance $G(V)={\partial I}/{\partial V}$ often used to characterise the conductance in the non-linear regime. Similarly to $G(V)$ also $S_{d}$ should be accessible experimentally in appropriate ac circuits. We remark that $S_{d}$ has been analysed in Refs. [33, 76], and the authors claim that this non-linear Seebeck coefficient allows for a better understanding of decoherence processes at finite voltage. It may have potential applications in nanoscale temperature sensors.

Here we are interested in the generalisation of $S_{d}$ to arbitrary temperature differences $\Delta T$ and arbitrary bias as well as its comparison to $S_{n}$ and $S$. Of course, for infinitesimally small $\Delta T$ and $V$ all definitions of Seebeck coefficients are equivalent and lead to the same result, $S=S_{n}=S_{d}$. For arbitrary $V$ and $T$ but vanishingly small $\Delta T$, the two non-linear Seebeck coefficients are equal, $S_{n}=S_{d}$. The formula (28) extends the standard definition (25) towards the non-linear regime, caused by both a large $\Delta T$ and additionally a finite (large) externally applied bias voltage $V$. In this regime, $S_{n}$ and $S_{d}$ attain different values.

To gain some feeling about the relative values and the behaviour of the three, in principle different Seebeck coefficients we show in Fig. 6 the linear $S$ (where appropriate), non-linear $S_{n}$, and differential $S_{d}$ coefficients as calculated for the two-terminal device. The three panels in Fig. 6 illustrate their dependence on temperature $T$, temperature difference $\Delta T$, and voltage $V$. We assumed $T_{R}=T$, $T_{L}=T+\Delta T$, $\mu_{L}=\mu+eV/2$, and $\mu_{R}=\mu-eV/2$, and performed the calculations for $U=8\Gamma_{0}$, $\varepsilon_{d}=-5\Gamma_{0}$, with other parameters as given in the figure. We remark that for these parameters the Haldane formula for the Kondo temperature (with $(\Gamma_{L}+\Gamma_{R})/2=\Gamma_{N}$),

[TABLE]

leads to $T_{K}\approx 0.071\Gamma_{0}$.

Figure 6(a) illustrates the $T$-dependence of all Seebeck coefficients for three values of the voltage $eV=-0.1$, [math], and $0.1$, and for a very small temperature difference, $\Delta T\rightarrow 0$. The value $V=0$ in fact denotes a very small voltage, $V\rightarrow 0$. All parameters are given again in energy units of $\Gamma_{0}$. For the actual calculations, we have used $\delta T=10^{-9}$, and to calculate the voltage derivative in (28) we have used $\delta V=10^{-9}$. All coefficients have the same value $S=S_{n}=S_{d}$ if $V=0$. For a relatively large value of $V$ the calculation of $S$ is meaningless; the figure shows that, independent of $T$ and for both values of $V$ one has $S_{n}=S_{d}$, the equality being due to the smallness of the applied temperature difference.

The situation is different if $\Delta T$ is arbitrary. The results for $V=0$ are presented in Fig. 6(b); it is apparent that all coefficients assume different values, except for small $\Delta T$ (the linear case). The differences increase with increasing $\Delta T$, with $S_{d}$ lying below $S_{n}$ (and $S_{n}$ below $S$) for all $T$ and a given set of parameters.

For non-zero voltages the relative order of $S_{n}$ and $S_{d}$ may be different, as is visible from the panel (c) of the figure. The Seebeck coefficient $S_{n}$ provides a generalisation of the standard definition to the non-linear regime. On the other hand, the differential Seebeck coefficient $S_{d}$ characterises the response to the temperature change of the voltage biased system. The data presented in Fig. 6(c) have been obtained for $T=0.01$, and two values $\Delta T=0.005$ and $0.05$. For small temperature bias the curves corresponding to $S_{n}$ and $S_{d}$ are rather close to each other. However, for larger $\Delta T$ the non-linear dependence of the Fermi functions and the on-dot density of states on $V$ and $\Delta T$ lead to various contributions to both $S_{n}$ and $S_{d}$.

The non-zero value of both Seebeck coefficients $S_{n}$ and $S_{d}$ for $V=0$ can be understood by taking into account the lack of particle-hole symmetry, $2\varepsilon_{d}+U\neq 0$, for the set of parameters used. For these parameters the density of states is similar to that shown in the panel (a) of Fig. 2. The differences between the curves $S_{n}(V)$ and $S_{d}(V)$ for the same $\Delta T$ are of the same character as those observed in the middle panel of the figure, while the global similarities between the two sets of curves calculated for $\Delta T=0.005$ and $\Delta T=0.05$ can be traced back to a larger contribution to $(\partial I/\partial\Delta T)$ obtained from Eq. (13):

[TABLE]

Note that the difference

[TABLE]

is an odd function of the voltage, which implies that the resulting curves are nearly anti-symmetric with respect to $V=0$ (and that the ordinates are equal, $S_{n}(0)\approx S_{d}(0)$). The smaller contribution proportional to $\left(-{\partial f_{L}(E)}/{\partial\Delta T}\right)$ depends on $V$ in a non-universal way. This causes the crossing of two coefficients at various voltages and for $\Delta T$ well beyond the validity of the linear approximation.

As a final remark we note that the asymmetry of the couplings $\Gamma_{L}\neq\Gamma_{R}$ also affects the Seebeck coefficients. To see this we assume $V=0$ and a relatively small (but slightly beyond the validity of the linear approximation) temperature difference $\Delta T=0.03$. Figure 7 shows the three Seebeck coefficients $S$, $S_{n}$, and $S_{d}$, calculated for isotropic coupling $\Gamma_{R}/\Gamma_{L}=1$ (thin lines) as well as for an anisotropic system $\Gamma_{R}/\Gamma_{L}=2$ (thick lines). Interestingly, the largest decrease of the Seebeck coefficient is observed at low temperatures, well below the Kondo temperature $T_{K}$, in agreement with recent findings [76] based on the non-crossing approximation. These authors already noted that the effect is largest in the Kondo regime. The asymmetry of the couplings is an important experimental fact [123] which has to be taken into account in any calculation aiming at a comparison with experiment [124]. Experimentally, it has been found [123] that, $e.g.$, the asymmetry shifts the cut-off of the noise emission from quantum dots from lower towards higher frequencies.

6 Three-terminal heat engine: the role of strong Coulomb interactions

As already noted a thermoelectric device can serve as waste heat to electricity converter. The energy harvesting is an important issue at large and small scales. As discussed in the Introduction, the usefulness of the bulk material as heat engine is well characterised by the figure of merit, $ZT$. This parameter ceases to be a good quality indicator in nanostructures [15, 17], and the calculations of both the power and the efficiency are required. The quantity of interest then is the maximum power, and the efficiency at maximum power.

In this section we consider the effects of strong interactions between the on-dot electrons on the characteristics of the three-terminal heat engine consisting of two cold terminals connected via two quantum dots to the hot one, cf. Fig. 1. Such a three-terminal heat engine, for the non-interacting case but well outside equilibrium, has been studied earlier [34, 35]. The work has been later extended [53] to take screening effects into account, treating the interactions within a mean-field approximation. The main conclusion was that albeit the screening modifies the parameters at which the engine is optimal, it does neither change the maximal power nor the efficiency at the maximum power. A similar heat engine has also been optimized with respect to the transmission function [125].

Here we shall pursue the analysis assuming non-linear working conditions and taking strong interactions into account. The calculations include the Kondo regime. To this end, we use the general expressions (7) for the charge and (8) for the heat current flowing out of the $\lambda$ lead, the energy conservation (16), and the general expression (21) for the on-dot Green function. From the charge current flowing from the left to the right electrode, and the heat current flowing from the hot to the cold electrodes, we calculate the performance of the engine as quantified by the maximum power $P$ and the efficiency $\eta$ at maximum power.

Previously we have found [53] that the three-terminal heat engine at optimized conditions attains an efficiency slightly above 0.2 in units of the Carnot efficiency, and that this value is roughly the same as without screening effects. The optimization involved the effective coupling between the dots and the leads, $\Gamma/T_{\rm{av}}$, the temperature difference between the electrodes, $\Delta T/T_{\rm{av}}$, and the load voltage, $V$; $T_{\rm{av}}$ is the average temperature of the system. The calculations have shown that the optimal coupling, i.e., the one leading to the maximal power if other parameters remain fixed, is of the order of the average temperature, $\Gamma/T_{\rm{av}}\approx 1$. (The coupling $\Gamma$ introduced here refers to the symmetric engine with $\Gamma_{L}=\Gamma_{R}=\Gamma_{H}\equiv\Gamma$.) The efficiency has been found to depend rather weakly, and the power strongly on the temperature difference $\Delta T$ between the hot ($H$), $T_{H}=T+\Delta T$, and the two cold ($R,L$) electrodes, $T_{R}=T_{L}=T$.

In order to demonstrate the effect of the Coulomb interaction on the performance of the engine, we show in Fig. 8 the efficiency $\eta$ $vs.$ power $P$. The efficiency is measured in units of the Carnot efficiency, $\eta_{C}={\Delta T}/{T_{R}}=0.4$, and the power is normalized by $(k_{B}T_{\rm{av}})^{2}$. The plot shown in panel (a) of the figure, valid for a three-terminal system, has been obtained by calculating the heat and charge currents as well as the optimal voltage (and the power) for a given difference of the dot’s energy levels, $\Delta E=\varepsilon_{R}-\varepsilon_{L}$, where $\varepsilon_{R}$ ($\varepsilon_{L}$) refers to the energy level of the right (left) dot in the system (see Fig. 1). It has to be noted that for appropriate values of $\varepsilon_{R}$ ($\varepsilon_{L}$) and low temperature the respective dot may show Kondo behaviour. This has an important effect on the $\eta$–$P$ plot.

In panel (a) of Fig. 8 one sees that the $U=0$ curve essentially encompasses all curves obtained for the interacting case. Coulomb interactions generally suppress the performance characteristics – at least so under the assumptions of the present approach, including the wide band approximation. The temperature of all electrodes are such that this system (Fig. 8(a)) is always outside the Kondo regime.

To illustrate the behaviour of the engine working in the Kondo regime, we show in Fig. 8(b) a similar plot, but obtained for a two-terminal system and a much lower value of the average temperature. For the interaction $U=8\Gamma_{0}$ and temperatures $T_{L}=0.15\Gamma_{0}$, $T_{R}=0.25\Gamma_{0}$, the dot enters the Kondo regime for a range of gate voltages (or $\varepsilon_{d}$). The Kondo effect existing for some values of $\varepsilon_{d}$ results in the appearance of the new performance branch on the $\eta$–$P$ plane. As visible in the panel (b) of the figure, this region is characterized by low efficiencies and powers when the system works as a heat engine.

The three-terminal quantum dot working as a refrigerator at very low temperature and in the Kondo regime requires a more careful examination, which is outside the scope of the present paper.

The results shown in Fig. 8 have been obtained under the assumption that the couplings to the leads equals the average temperature, $\Gamma=k_{B}T_{\rm{av}}$, the value which has been found to be optimal (i.e., leading to the maximal power) for the non-interacting system. This facilitates the comparison with the previous results [53]. The red curve in Fig. 8(a), corresponding to the non-interacting case, agrees with previously obtained results for the same set of parameters: $\Gamma_{L}=\Gamma_{R}=\Gamma_{H}=T$, $\Delta T/T_{\rm{av}}=0.5$. Analogous to Ref. [53], changing the dots energy levels we fix the difference $\Delta E=\varepsilon_{R}-\varepsilon_{L}$ and calculate the power and efficiency at each point for the optimized value of the load voltage $V$. The curves in Fig. 8 are thus parametrised by $\Delta E$. As is apparent from the figure, the maximal efficiency is slightly higher than 0.2$\eta_{C}$, which nicely agrees with recent measurements on similar systems [104].

7 Summary and conclusions

We have studied the thermoelectric transport properties of two- and three-terminal systems with quantum dots, paying attention to strong interactions of electrons on the dot(s) and far-from-equilibrium conditions. Of particular importance in such nanostructures is the non-linear regime at large external voltage bias $V$ and large temperature difference $\Delta T$. As discussed in the Introduction, the linear approximation is hardly ever valid in nanostructures. This has been again confirmed here, and is visible as the detrimental effect of strong correlations on the performance of the three-terminal optimized heat engine. Calculating the maximum power $P$ and the efficiency $\eta$ of the optimized heat engine for various interactions $U$, we have observed that except at very low temperatures the curves calculated for $U\neq 0$ are encompassed by the curve obtained for the non-interacting system. Obviously the filtering properties are affected by the interactions which broaden or split the density of states rendering the filter less effective. This agrees with previous attempts at optimizing heat engines by engineering the transmission function [27, 125]. On the other hand, strong interactions are responsible for the second leaf on the efficiency vs. power plane, visible in Fig. 8(b) and appearing at low temperatures when the system enters the Kondo regime, albeit for a relatively small range of gate voltages. This novel feature appearing in the Kondo regime requires further analysis, especially for a three-terminal quantum dot working as a refrigerator at very low temperatures.

The strongly non-linear regime requires proper definitions of the thermopower. The generalisation of the standard linear-response definition of the Seebeck coefficient to the strongly non-linear regime has been considered. This results in two different formulas for two coefficients $S_{n}$ and $S_{d}$. Both are, in principle, valid for arbitrary values of $V$ and $\Delta T$, albeit the first coefficient ($S_{n}$) is easier to handle in systems with zero external voltage but arbitrarily large temperature difference between a given pair of electrodes. The second, called differential Seebeck coefficient, has earlier been applied [33, 76] to systems with finite current flow and a small temperature difference. It has been generalized and studied here for arbitrary $\Delta T$ in the presence of an external voltage bias $V$. Interestingly the asymmetry of the couplings to the external leads has qualitatively similar influence on both Seebeck coefficients. They develop a minimum for temperatures well below the Kondo temperature. The observed quantitative differences between $S_{n}$ and $S_{d}$ are expected to be important for temperature sensing by means of thermopower measurements [126]. The non-symmetric device with different couplings between left and right electrode can, in principle, work as a thermal diode, albeit with a rather small heat rectification factor of the order of 1%. The analogous rectification factor for a charge flow is more than an order of magnitude higher, up to about 30%.

The work reported here has been supported by the M. Curie-Skłodowska University, National Science Center grant DEC-2017/27/B/ST3/01911 (Poland), the Deutsche Forschungsgemeinschaft (DFG, project no. 107745057, TRR 80), and the University of Augsburg.

Appendix A Details of the calculations

There exist many techniques which provide a correct description of the physics leading to the single-impurity Kondo effect [52]. However, most of them are based on numerical approaches very often restricted to equilibrium situations. A quantitatively correct analytical description of the Kondo effect valid both in and out of equilibrium, however, is important as it provides additional insights into the underlying physics [71, 64]. Here we shall show that the appropriately generalized [49] EOM technique supplemented with an analytic calculation of the lifetimes provides an easy, physically transparent, and numerically correct approach to study the Kondo regime.

To calculate the retarded GF in the stationary state, it is convenient to use the EOM for two-time functions [50] in the energy representation. The equation of motion for an arbitrary Green function $\langle\langle A|B\rangle\rangle_{E}$ composed of the fermionic operators $A$ and $B$ reads:

[TABLE]

However, the calculation of lesser Green functions are more conveniently performed on the time contour and continued to real times by means of Langreth’s theorem [105]. This approach has been used to obtain the average number of electrons on the dot, Eq. (12).

To determine the retarded GF, one starts with the equation of motion for $g^{r}_{\sigma}(E)=\langle\langle d_{\sigma}|d^{\dagger}_{\sigma}\rangle\rangle_{E+i0}$, which reads

[TABLE]

where $H$ denotes the Hamiltonian (2) of the system.

Calculating the commutators one finds

[TABLE]

Here the symbol $\bar{\sigma}\equiv-\sigma$ is introduced, and two new GFs show up on the r.h.s. In turn, we calculate both of them, again using the equation of motion (32). The simplest one is given by

[TABLE]

In Eq. (34) we need this function multiplied by $V^{*}_{\lambda k\sigma}$ and summed over $\lambda k$. The result is:

[TABLE]

The factor in front of the GF on the r.h.s. defines the self-energy:

[TABLE]

In the wide band limit, one approximates (37) by its imaginary part:

[TABLE]

and neglects the possible energy dependence of $\Gamma_{\sigma}^{\lambda}(E)=\Gamma_{\sigma}^{\lambda}$. As discussed in Sec. 2, the wide band approximation is essential in the derivation of the exact formula (12).

In principle it is tempting to decouple the GF on the r.h.s. of (34) multiplying $U$ as $\langle\langle n_{\bar{\sigma}}d_{\sigma}|d^{\dagger}_{\sigma}\rangle\rangle_{E}\approx\langle n_{\bar{\sigma}}\rangle\langle\langle d_{\sigma}|d^{\dagger}_{\sigma}\rangle\rangle_{E}$, but it turns out that the equation of motion for this GF introduces important new functions, describing fluctuations of opposite-spin ($\bar{\sigma}$) electrons. In order to obtain a qualitatively correct description of Kondo correlations, this function has to be calculated exactly [56, 55] at this level. In the next step, we obtain:

[TABLE]

Interestingly, the GFs containing one lead operator enter the above equations via the sums

[TABLE]

which are calculated using again the EOM for the new GFs. They read

[TABLE]

[TABLE]

[TABLE]

To close the infinite hierarchy of consecutive equations, we have to employ an approximation at a certain level. The common approximation is to project those new GFs appearing in the above set which contain two lead operators onto the lower order ones, $e.g.$:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

The motivation for the decoupling comes from the knowledge that the many-body Kondo singlet observed between the electron localized on the dot and the conduction electrons is due to anti-ferromagnetic spin flip processes. Thus performing this decoupling, we concentrate on the spin, say up, particle moving in a self-consistently (see below) determined dynamic field of the spin down particles tunneling between the dot and electrodes. These approximations are known as Lacroix decoupling [56]. They are valid up to the second order in the coupling $V_{\lambda k\sigma}$ to the leads. The above approximations result in the appearance of the various self-energies given as follows [49]:

[TABLE]

and

[TABLE]

They can easily be expressed in terms of the retarded ($r$) and advanced ($a$) Green functions $\langle\langle d_{\sigma}|d^{\dagger}_{\sigma}\rangle\rangle^{r,a}_{E}$ of the opposite spin GFs. In these expressions we have introduced two important generalisations as proposed earlier [49]. They are suggested by the work of Van Roermund et al. who have extended the EOM technique and systematically studied the Anderson model up to the fourth order [69], $i.e.$, up to $|V_{\lambda\vec{k}\sigma}|^{4}$. One the main findings of that paper is that the primary effect of higher order processes is to provide a finite lifetime of excited states on the dot, and a renormalization of the on-dot energy level $\varepsilon_{d}$. These two effects are taken into account by replacing the imaginary part in the GFs, $i0$, by $i\tilde{\gamma}_{\alpha}$, and $\varepsilon_{d}$ by $\tilde{\varepsilon}_{d}$. Finally we obtain the following expressions:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

In the above we have introduced $\varepsilon_{1}=\tilde{\varepsilon}_{\sigma}-\tilde{\varepsilon}_{\bar{\sigma}}$, and $\varepsilon_{2}=\tilde{\varepsilon}_{\sigma}+\tilde{\varepsilon}_{\bar{\sigma}}+U$. The subscripts 1 and 2 refer to the excited 1- and 2-electron states of the dot, respectively. The self-energies are calculated iteratively and self-consistently with the GF (21). The remaining two self-energies $\Sigma_{\sigma}^{(1,2)}$ are equal to $\Sigma_{0\sigma}$ for $i\tilde{\gamma}^{\alpha}_{1,2}=i0^{+}$; however, for arbitrary values of $i\tilde{\gamma}^{\alpha}_{1,2}$ they have to be calculated directly from the definition,

[TABLE]