Optimization of a three-terminal non-linear heat nano-engine

Charge and heat transport in nano-structures is dominated by non-equilibrium effects which strongly influence their behaviour. These effects are studied in a setup consisting of three external leads, one of which is considered as a heat reservoir and is tunnel-coupled to two cold electrodes via two independently controlled quantum dots. The energy flow from the hot electrode together with energy filtering provided by quantum dots leads to a voltage bias between the cold electrodes. The heat and charge currents in the device effectively flow in mutually perpendicular directions, allowing for their independent control. The non-equilibrium screening changes the values of the system parameters needed for its optimal performance but leaves the maximal output power and efficiency unchanged. Our results are important from the theoretical point of view as well as for the practical implementation and the control of the proposed heat engine.


Introduction
The efficient harvesting of waste heat is one of the most important challenges of modern technology both at large and small scales. Combining the effective cooling of the integrated circuits with simultaneous conversion of a part of the released heat into electric power could revolutionize electronics. On the road towards effective heat nano-engines a number of important findings have to be noticed. Among them, the observation of the importance of energy filtering [1], and the independent control of heat and charge flow [2,3,4,5,6] are of great interest.
These observations are at the heart of our approach as we use quantum dots as efficient energy filters and the three-terminal setup for independent control of heat and charge flow. The setup we are considering is shown in Fig. 1. It consists of two independently tuned quantum dots and three external terminals. The left and right current junctions contain a quantum dot. The central terminal is the hot one. It can be considered as a cavity [4] connected to an external heat bath. The temperature of the hot reservoir equals T H , while the two other reservoirs are assumed to have lower temperatures T L and T R , respectively. The chemical potentials of the electrodes µ L,R,H may be changed by the external voltage or as a result of the temperature difference between hot and cold electrodes.
Experimentally similar systems [7] have been already produced and implemented in electronic refrigerators [8], proved to be successful in cooling in the mK temperature range, and recently shown to be efficient heat harvesters [9,10,11]. Thermoelectric nano-engines with quantum dots tunnel-coupled to external electrodes in two-and threeterminal geometry have been proposed as effective heat to electricity converters [12,13,14,15,16,17,18,19,20,21,22,23,24,25,26]. Note that besides quantum dots also molecules [27] and nano-wires [28] are useful elements for efficient energy harvesting at the nano-scale. The field of thermoelectric energy harvesting with quantum dots has been been recently reviewed [29,30], while a more general discussion related to energy harvesting can be found in [31].
Consider the system shown in Fig. 1 with quantum dots energy levels differing by ∆E. Assume for a while that the tunnelling via quantum dots is possible only at sharp values of on-dot energies ε L and ε R . In such a case an electron from the left lead with energy ε L can tunnel into the H-lead, and an electron with energy equal to ε R can tunnel from H to the right electrode. Thus each electron transferred between L and R electrode gains an energy ∆E = ε R − ε L from the H electrode. This process is possible if the temperature of the H electrode is the highest (hence hot electrode). The charge is transferred between L and R electrodes at the cost of the heat from the H electrode. The dots act here as efficient energy filters, and charge effectively flows in direction perpendicular to the heat flow. In other words, the electron flow being a result of the temperature difference between hot and cold electrodes gives rise to a voltage bias between the two cold electrodes. One can invert the perspective and argue that in the presence of the (not too large) voltage bias (load) between two cold electrodes the electron flow at the cost of heat from the hot electrode performs work against the bias. The value of the bias at which the charge currents stop to flow is called stopping bias, and is denoted by V stop . The device operates as a heat engine in the voltage range (0, V stop ). The possibility of independent control of heat and charge flow is a main motivation to consider the three terminal geometry.
A similar heat nano-engine has been recently optimised [4] for maximum power. The optimization involved the coupling strength of the dots to external electrodes, the "energy gain" ∆E and the voltage load between left and right electrodes. However, the authors [4] have not considered the non-equilibrium effects related to charge redistribution and screening, being of importance outside the linear transport approximation.
Indeed, recent measurements of the thermoelectric voltage clearly show [32,33] that the observed non-linear effects are related to heating-induced renormalisation of the dot energy levels. The other source of non-linearities, namely, the energy dependence of the transmission function has been found to play a negligible role. The renormalisations of the dot energy states by electric and thermal gradients in a similar system have been recently studied theoretically within the scattering approach [34,35,47].
In other words, beyond the linear approximation the non-equilibrium screening potentials start to play an importnat role. As a result, e.g. the optimal set of parameters of an engine differs from that obtained in the theory which does not take such effects into account. For finite voltage or temperature bias, the charges pile up in the electrodes and quantum dots. Due to the long-range Coulomb interaction they screen other charges and change the injection rates of particles from the electrodes [36]. This observation is especially important for large temperature differences and large load voltages well beyond the validity of linear response.
At the nano-scale the issue of linear response is a tricky one. In principle it is even not well defined. In bulk diffusive systems the small temperature difference between the far ends of the sample allows well defined local temperatures and an average temperature gradient. Similarly a small bias usually leads to a small gradient of the electrochemical potential. In nano-structures even small biases do not imply validity of linear response. To capture non-linearity we shall use the non-equilibrium Green function approach to derive equations for heat and charge currents and consider the effects of non-equilibrium screening of charges and their piling up in the electrodes [37].
Working as an energy harvester the system converts the heat current J into power, P = IV , where I is the charge current flowing between left and right (cold) electrodes. The voltage V is used to power an external device (the load). The efficiency is defined as the ratio between the useful power and the heat current flowing into the system, η = P/J. The efficiency calculated in this way can be contrasted with the ideal Carnot value η C = 1 − T R /T H , and the Curzon-Ahlborn efficiency [38] expected at maximum power, The organization of the paper is as follows. In the next section we present the microscopic model of the system at hand, and calculate the charge and heat currents using the non-equilibrium Green function technique. The non-linear effects in transport are discussed in Sec. III. The results of the optimization of the engine working well outside the linear regime are presented in Sec. IV. We end with summary and conclusions (Sec. V).

Model and approach
The Hamiltonian of the system is written aŝ where n λ kσ = c † λ kσ c λ kσ and n iσ = d † iσ d iσ denote particle number operators for the leads and the dots, respectively. The operators c † λkσ (d † iσ ) create electrons in respective states λ kσ (iσ) in the leads (on the dots). Symbols i = 1, 2 refer to the left and right dot, and λ = L, R, H denote the left, right, and hot electrode, respectively. The dot energy levels ε L ≡ ε 1 and ε R ≡ ε 2 can be easily tuned by the gate voltages. They are renormalised by the potentials U i which account for the electron-electron repulsion. This interaction is considered here at the mean-field level.
The bare tunneling amplitudes between dot i and electrode λ are denoted by V iλ k . We introduce the symbol which takes into account that left and right leads are coupled, respectively to the first 1 (= L) and the second 2 (= R) quantum dot, while both dots are coupled to the H electrode.
The current in the electrode λ is calculated as time derivative of the average charge in that electrode N λ = kσ n λ kσ : where the symbol ... denotes the statistical average. Calculation of the heat flux follows that of the charge. From thermodynamics we know the relation between the (internal) energy E, heat Q, and work. Assuming that the only work is related to the flow of mass µdN, one writes The dot over symbols denotes time derivative. Applying this equation to one of the electrodes, say λ, allows to write the heat flux as where H λ = k,σ ε λ kσ n λ kσ is the energy operator for the electrode λ. Evaluating the commutators and defining Keldysh "lesser" functions: one gets The final expressions for the stationary currents can easily be written in the general form [39] denotes the matrix of effective couplings of dots (i, j) to the lead λ. The factor 2 in the formulas for currents stems from the summation over spins. The heat current (11) can be written as a difference between the energy current J E λ and the charge current I λ : To calculate lesser Green function [39] we use the equation of motion method [40]. From now on, we shall work in units with Planck constanth = 1 and Boltzmann constant k B = 1. It is convenient to define the frequency dependent dot matrix Green function G < (ω) with elements For the definition of ...|... , see Ref. [40]. After some algebra one findŝ where the matrix lesser self-energy is given bŷ The equation for the retarded matrix Green function [41] G r ij (ω) can be written in explicit form as with the retarded self-energy The matrix inversion gives all the components of the required retarded function. In the wide band limit approximation one replaces the retarded self-energy by its imaginary part only: and neglects the frequency dependence of Γ λ ij (ω) = Γ λ ij .

Non-linear effects
The long-range nature of the Coulomb interactions is responsible for the back-reaction of the non-equilibrium charge distribution onto the transport properties of the device.
In the Hamiltonian (1) this is represented by the screening potentials U i . Their values depend on the thermoelectric configuration, i.e. voltages and temperatures of all electrodes. This effect has been considered in mesoscopic normal systems first by Altshuler and Khmelnitskii [37], and later by Büttiker and coworkers [36,42], and others [43,44,45,34,46,47]. It has been also explored in metal-superconductor twoterminal [48,49] and three-terminal junctions [50].
Here we follow Ref. [48] and others [49,50], assuming that the long-range interactions modify the on-dot energies ε i , changing them into ε i − eU i . In equilibrium the potentials U i have constant values (independent of V λ and ∆T λ ), which we denote by U i,eq . In the presence of applied voltages V λ and temperature biases ∆T λ , the deviations δU i = U i − U i,eq , in lowest order, are linear functions [48] of V λ and ∆T λ .
Thus we write for the potential on each dot i = 1, 2: where the subscript zero indicates that the partial derivatives have to be evaluated with all V λ , ∆T λ set to zero, and the dots denote higher order terms. The charge densities on the dots n i also depend on the temperature and voltage bias as well as on the potentials U i . Expanding to lowest order in these parameters we get The above equation defines the Lindhard matrix function as The derivatives can be easily calculated by noting that and using equation (15). Another relation between charges and potentials defines the capacitance matrixĈ of the system: The equations (21) and (24) are easily solved, and one finds explicit expressions for the characteristic potentials The knowledge of the characteristic potentials allows to calculate how the temperature difference between hot and cold electrodes and voltages modify the potentials U i of the dots. These changes, in turn, affect the heat and charge currents flowing in the system. For the explicit results presented below, we will assume the small capacitance limit, C ≈ 0. Due to the approximation in Eq. (20) our approach is called "weakly non-linear", as we do not consider higher order corrections.

Results
We assume that the left and right electrodes of the system (cf. Fig. 1) have the same temperature, T L = T R . The temperature of the hot electrode, which is kept fixed from now on, is denoted by T H (i.e. ∆T H = 0). In addition, ∆T R = ∆T L = T H − T R ≡ ∆T , and the average temperature of the system is T = (T H + T R )/2. The current does not flow into or out of the hot electrode which is grounded (I H = 0). This means that charge conservation written in the form I L +I R = 0 serves as a condition for the chemical potential µ H of the hot electrode. The energy conservation J E H + J E L + J E R = 0 may serve as a condition for the actual temperature of that electrode. We shall take another point of view and assume that the hot electrode serves as an energy reservoir characterised by the constant temperature T H . The heat current J H flows out of it towards the L and R electrodes. For ∆E = ε R − ε L = 0 the electrons entering the left electrode at energy ε L leave the right one at energy ε R . As a result the voltage V appears between both electrodes.
To facilitate the calculations we impose additional conditions. The couplings of the quantum dots to external leads fulfil Γ L ij = Γ L δ i1 δ j1 , Γ R ij = Γ R δ i2 δ j2 , and we assume the matrix Γ H ij to be symmetric with elements Γ H 11 , Γ H 12 = Γ H 21 and Γ H 22 . We tune the positions of the dots' energy levels symmetrically with respect to the chemical potential µ H . Also the voltages are assumed to be symmetrical with respect to µ H , i.e. µ L = µ H −eV /2 and µ R = µ H + eV /2. Such choice of parameters assures that I H = 0, and I = I L = −I R . For an arbitrary set of parameters not fulfilling the above symmetries, µ H has to be calculated from the condition I H = 0.

Linear transport parameters, power factor, and efficiency
From charge and heat currrents we calculate linear transport characteristics of the device including charge (G) and thermal (κ) conductances, and Seebeck coefficient S. For convenience we change the notation and denote the charge current as I 1 and the heat current as I 2 . Expanding the currents to linear order in bias and temperature forces X 1 = eV /T and X 2 = ∆T /T 2 , we write the fluxes in standard notation: The transport coefficients are given by the parameters L ij . In accordance with standard definitions one finds the conductance G = (I 1 /V ) ∆T =0 = L 11 /T , the Seebeck coefficient which we ocasionally refer to as thermopower S = −(V /∆T ) I 1 =0 = L 12 /(T L 11 ), and the thermal conductance κ = (I 2 /∆T ) I 1 =0 = (L 22 − L 12 L 21 /L 11 )/T 2 . The Peltier coefficent defined as Π = (I 2 /I 1 ) ∆T =0 is given by L 21 /L 11 . The combination of these parameters defines the thermoelectric figure of merit ZT = GS 2 T /κ, which enters the expression for the efficiency η lin of the thermoelectric heat engine [30]: The expected efficiency based on linear coefficients will serve as a reference below. In Fig. 2 we show the dependence of the linear conductance G, thermal conductance κ, Seebeck S and Peltier Π coefficients as well as power factor GS 2 vs. ∆E/k B T for the The linear transport coefficients and the thermoelectric figure of merit ZT depend on ∆E and γ. While both, the conductance G and thermal conductance κ monotonically increase, the thermopower decreases with increasing coupling γ. The increase of G and κ with γ is related to the fact that the currents (8), (9) are proportional to γ. To understand the decrease of the thermopower with γ it is useful to recall that in nano-structures S is directly related to the slope of the density of states at the Fermi energy [51], which is higher for smaller γ. In the linear approach the thermoelectric figure of merit yields the efficiency of the engine. Taking into account the behaviour of ZT on γ shown in the inset in the right panel of Fig. 2 and the formula (29), we conclude that maximum efficiency, approaching the Carnot value η C , is obtained for vanishingly small γ. However, the power of the engine approaches zero value, rendering the engine practically useless.
We are interested in achieving the maximum power, which can be realised by optimising ∆E and γ. In the linear theory the appropriate parameter characterising the obtained power is the power factor defined as GS 2 . Its dependence on ∆E for a few values of γ is shown in Fig. 3. As we shall see in the next section the power factor shows a dependence on ∆E and γ qualitatively similar to that found in the exact non-linear approach, however, with a few important differences to be discussed later.

Non-linear transport: maximum power and efficiency of the heat nano-engine
We start the presentation of the results obtained with non-equilibrium screening potentials taken into account by showing the parameters u i,λ and z i,λ and their dependence on the couplings Γ ij and the energy difference ∆E = ε R − ε L . These together with the voltage load V are the main optimization parameters.
In Fig. 4 we show the dependence of the characteristic potentials u 1,L , u 2,L and z 1,L , z 1,R on ∆E. The behaviour of the other parameters u 1,R , u 2,R , z 2,L and z 2,R can be inferred from their symmetries. For the symmetric system we are dealing with, the parameters u i,λ are even functions of ∆E and fulfil the relations: u 1,L (∆E) = u 2,R (∆E), u 1,R (∆E) = u 2,L (∆E). The z parameters are antisymmetric functions of ∆E, e.g. z 1,L (∆E) = −z 1,L (−∆E) and are related as: z 1,L (∆E) = −z 2,R (∆E), and z 1,R (∆E) = −z 2,L (∆E). The "off-diagonal" characteristic potentials are much smaller than the "diagonal" ones. This is especially true for the thermal potentials, z 1,R and If the couplings Γ H 12 , Γ H 21 are assumed to vanish, then also the off-diagonal characteristic potentials u 1,R , u 2,L , and z 1,R , z 2,L vanish. This fact, however, has only a small effect on the performance of the engine.
The most important parameters of the engine are the maximum output power P max and the efficiency η max at maximum power. The dependence of the maximum power, scaled by (k B T ) 2 , on ∆E/k B T with all couplings equal to γ is shown in the left panel of Fig. 5. Different curves correspond to different γ in units of k B T , and for each of them the power has been optimized with respect to the applied voltage load. The calculations were performed for ∆T /T = 0.3 with the non-linear effects taken into account. Interestingly, for ∆E/k B T up to about 12 the largest power is obtained for γ = k B T , but for ∆E/k B T > 12, (typically) higher values of γ lead to larger power. This non-monotonic dependence of the maximal power on the effective width of the resonance can be traced back to the strong ∆E dependence of the characteristic potentials, which in turn renormalise the dots energy levels ε L and ε R (and thus ∆E).
The right panel of Fig. 5 shows the efficiency of the engine corresponding to maximum power η max /η C , measured in units of the Carnot efficiency. The maximum value of η max /η C strongly depends on the coupling γ. For the optimal value of γ = k B T and for ∆T /T = 0.3, it exceeds 20% of the Carnot value. The efficiency as well as the maximum power are increasing functions of temperature difference. The linear approximation for the on-dot potentials presumably precludes reliable results for ∆T /T > 0.3. For a given value of ∆T /T the maximal efficiency increases with decreasing γ, tending to the Carnot limit when γ → 0. At the same time the power diminishes towards zero. This agrees with our analysis in the linear approximation carried out in Section 4.1, and the recent analytical treatment of the two-terminal system [26]. In Fig. 6 the comparison of the maximum power calculated with full renormalization of the on-dot energy levels (lines) with that obtained in the approach neglecting screening (symbols) is presented. Two important features are apparent. First, the non-linearities strongly affect the value of the dots energy difference, ∆E, for which the system performs with maximum power. While without screening effects the optimal ∆E corresponds to ≈ 6k B T , the screening shifts the optimal value to ≈ 9k B T . In order to understand this behaviour, one has to note that for a given bias V = V R −V L and temperature difference ∆T , one may define the "effective" energy difference (∆E) eff ≈ ∆E −eu 1,L V −2ez 1,L ∆T ; here we used the symmetries of the u i 's and z i 's, as well as the smallness of the offdiagonal potentials, cf. Sec. 4.2. As a result, the upward shift of the optimal ∆E is such that the effective energy difference roughly agrees with the optimal value found without screening effects. In fact, the shift due to screening evaluated at ∆E/k B T = 6 is about 2k B T , to be compared with the difference of ∆E-values mentioned above, namely ≈ 3k B T , for the considered parameters (∆T = T , γ = k B T ). The remaining difference, ≈ k B T , can be understood by observing that linear terms appear in the denominators of the Green functions which are integrated over energies with the Fermi functions in the nominator. Second, the value of the power itself and the efficiency remain unchanged. Other differences are less important for the question of the maximum power but should be noted. For example, one observes much stronger asymmetries of the P max (∆E) curves with the non-linear screening effects taken into account.
Comparison of the optimal power, Fig. 5, with the power factor, Fig. 3, shows that there exist important differences between these quantities. In particular, the energy difference ∆E at which the power is maximal markedly differs from that leading to the maximum power factor: this is related to the fact that in the full theory the voltage serves as an additional optimisation parameter.

Summary and conclusions
The three-terminal heat nano-engine has been analysed in the linear and the non-linear approximation. In the latter case it has been optimised with non-equilibrium screening effects taken into account. In the linear limit the Onsager symmetry relations are fulfilled due to the floating character of the hot electrode. In this limit we have calculated charge and thermal conductances, the Seebeck coefficient and the thermoelectric figure of merit ZT . The expected efficiency calculated from the figure of merit is compared with that calculated self-consistently in the weakly non-linear limit in the right panel of Fig. 5. For the optimal value of the coupling (γ = k B T ) the linear efficiency surprisingly well describes the performance of the system [4] calculated exactly.
Our calculations show that the optimal value of the coupling is essentially unchanged by the non-equilibrium screening effects, and equals γ = k B T . These effects do not change the maximum value of the output power and the efficiency obtained for the optimal value of the coupling constant and for a given value of ∆T /T . The optimal distance between the dots energy levels, ∆E, changes as a result of the screening effects. The change is directly related to modifications of the on-dot energy levels by the potentials U i . The system efficiency at maximum power in units of the ideal Carnot efficiency exceeds 20% for ∆T /T = 0.3.
In agreement with other studies of heat nano-engines [52,53,35,54,55], we have found that the large value of the thermoelectric figure of merit ZT does not necessarily imply the usefulness of the device as an efficient energy harvester. Interestingly, for this particular device, we have found an overall agreement between the efficiency obtained within linear approximation and in the full theory (cf. two sets of data for γ/k B T = 1 in the right panel of Fig. 5). Hence the three-terminal system under study markedly differs from a typical two-terminal engine.