# Single Charge Current in a Normal Mesoscopic Region Attached to Superconductor Leads via a Coupled Poisson Nonequilibrium Green Function Formalism

Authors: David Verrilli; F. P. Marin; Rafael Rangel
Journal: The Scientific World Journal
 (2014)

---

## Abstract

We study theI-V characteristic of mesoscopic systems or quantum dot (QD) attached to a pair of superconducting leads. Interaction effects in the QD are considered through the charging energy of the QD; that is, the treatment of current transport under a voltage bias is performed within a coupled Poisson nonequilibrium Green function (PNEGF) formalism. We derive the expression for the current in full generality but consider only the regime where transport occurs only via a single particle current. We show for this case and for various charging energies values U0 and associated capacitances of the QD the effect on the I-V characteristic. Also the influence of the coupling constants on the I-V characteristic is investigated. Our approach puts forward a novel interpretation of experiments in the strong Coulomb regime.


---

## Body

## 1. Introduction

The overall shape of theI-V characteristic of a variety of systems (metals, semiconductors, and molecular conductors) in the nanometer scale sandwiched between metallic or superconductors leads has been recently a matter of study (see [1, 2] and references therein). In these systems, the energy level discreteness is quite important since level spacing is comparable with other energy scales [3, 4]. Indeed, the coupling with the bath modifies drastically the properties of an otherwise uncoupled nanometer system in a sharp contrast with similar nonequilibrium macroscopic systems [5–13]. They constitute hybrid systems. Theoretical studies [14–20] as well as experimental measurements have been done by many research groups [1–3] on such systems mostly at low enough temperature with negligible thermal and non-equilibrium fluctuations.All the systems mentioned above underlay universal common features with the hybrid superconductor quantum dot devices we want to address in this work [2, 4, 21, 22]: (i) broadened energy levels of the quantum dot due to hybridization with the leads; (ii) spatial potential profile. (iii) a charging energy U0 due to the potential profile. An insight behind these issues has been highlighted recently [23, 24] for molecular dots. The device we study in this work is shown in Figure 1. It constitutes a spin degenerated quantum dot level, which is coupled to a pair of biased superconductors contacts or leads (source and drain). When a source-drain voltage Vd is applied, an electric current flows between the leads and across the quantum dot. The biasing defines a non-equilibrium steady state situation. Such situation is coming from the frustration to establish simultaneously an equilibrium configuration with both leads under a given bias. In addition, a gate voltage Vg sets the quantum dot spectrum. However, the charge energy can modify it whenever the density of states is significant. In response to the applied voltages, an actual potential develops inside the dot; that is, an effective electrostatic profile potential inside the mesoscopic region exists in such a way, that it couples to both the electronic non-equilibrium state population and the non-equilibrium electric current. That approach, as introduced by Datta [4], links the electrostatic profile to the electronic population of the quantum dot [4, 25] via the non-equilibrium Keldysh formalism (NEGF) [26, 27]. The whole system is modeled by coupling capacitances which represents the drain, source, and gate contributions to the self-consistent electrostatic problem. Incoming electrons have to overcome an energy barrier (Coulomb blockade). On the other hand, gate or source-drain voltage can lower or increase this energy barrier. These source, drain, and gate electrodes capacitances (see Figure 2) constitute a simple capacitive model (in experiments [2, 28], these capacitances are measured) from which Uℒ, the Laplacian part of the potential, can be obtained. In addition, the charge in dot can be expressed as the sum of the charges in the coupling capacitances. It yields the Poisson contribution UP to the total potential U, as a function of the dot population. In other words, we solve the self-consistency (SC) of the total electrostatic potential U=Uℒ+UP together with the dot population. After that, the electric current is evaluated.Figure 1
Set-up: single level quantum dot connected with two superconducting leads via coupling constantsΓs.Figure 2
Equivalent capacitive circuit with coupling capacitancesCs, Cg, and Cd, corresponding to the capacitances in the source, gate and drain, respectively.Previous to the self-consistent program, the non-equilibrium current through the dot and electronic occupation in the dot are worked out. We emphasize that the calculation is carried out in a general framework. However, we confine our attention to the single particle current contribution. We adapt the SC to two different approximation regimes. In Section4, the equivalent capacitive circuit (Figure 2) is introduced; the spatial potential profile U is calculated within the capacitive model. The SC scheme is applied to two cases [29, 30]. First, the so-called restricted case, where the gap is the bigger energy scale and the coupling QD-Leads is of the order of the charging energy (Δ≫ΓL,R≃U). In this case, quantitative results are expected to be accurate. We also make calculations for the so-called unrestricted case, where the charging energy is the dominant energy scale Δ≃U≫ΓL,R. In this case the results are quantitatively less accurate. The experiments of Ralph et al. [28] were done in this regime. Their I-V characteristic shows that the spacing of the energy levels is subjected to strong fluctuations. According to our model, the fluctuations are due to complex multilevel charging effects. Our hybrid S/QD/S system has been studied in previous theoretical works [16–19]. However, to our knowledge, the coupled SC scheme which describes charging effects has not been considered so far. This is an important step; then, gauge invariant independence of the results as well independence of the zero reference voltage is fulfilled [31, 32]. Our model uses experimental values of the equivalent capacitances [4]. To this respect, pioneering work is done by Meir et al. [33, 34] for N/QD/N systems, considering the interatomic Coulomb term Un↑n↓ as a measure of the charging energy e2/C. Their purpose was to find the main object of the non-equilibrium formalism, namely, the QD Green-Keldysh function, in which the influence of the leads on the QD is taken into account. Due to the presence of the Coulomb term, its equation of motion generates a two-particle Green-Keldysh function. By ignoring correlations with the leads, the equation of motion for the QD Green function closes after truncation of higher order equations of motion. This solution (their Equation (8)) has two resonances, one at the energy level weighted by the probability that the other spin degenerate level (raised by U) is vacant and another one at the energy level raised by U weighted by the probability that the level is occupied. It is correct for temperatures higher than the Kondo temperature and is exact in the noninteracting limit (U=0) and the isolated limit. Analogously, for S/QD/S hybrid systems Kang [16] has obtained an expression for the current through the QD (his Equation (8)), which is evaluated in the U→∞ limit (his Equation (13)). The QD Green function from the very beginning does not contain off-diagonal terms that involve superconducting pairing, which excludes the possibility of Andreev reflection processes. The presence in the equation for the current (his Equation (14)) of terms proportional to (1-〈n-σ〉) affects the contribution to the current of the considered level. In order to complete the outlined program one has to calculate 〈n-σ〉 self consistently which is not carried out. Instead, Kang calculate the current (his Equation (8)) where the spectral function is calculated in the limit of zero coupling with the leads via a model taken from literature (his Reference [18]) and without taking into account the dependence of the contribution of one level to the current on the occupancy of the other. The point of view which neglects the unavoidable influence of the bath (the leads) on the small system (the QD) is accomplished by factorizing the density matrix (ρ(t)=ρQD⊗ρBaths) and integrating out the leads degrees of freedom which simplify the Liouville-von Neumann equation (Equation (3.140) in [35]). This program is carried out by Kosov et al. for S/QD/S system [36]. In this way, a Markovian master equation is obtained and an expression for the current is calculated. In their Figure 2, they show the I-V characteristic of a nondegenerated QD for a given set of parameters. In this case, the Cooper pair density in the QD is zero [37]. For the sake of comparison, we restrict our calculations to this case. A similar but not identical approach was done by Pfaller et al. [38]. Also, the approach of both Kosov et al. and Pfaller et al. misses the energy levels broadening as discussed in the introduction. This lack of broadening is a general deficit of quantum Markov approach [39]. In particular, Pfaller et al. [38] introduce a phenomenological broadening while our approach derives it from first principles. In fact, within the Keldysh formalism, this broadening appears naturally (see (58) below). Levy Yeyati et al. [17] writes an expression for the current (his Equation (2) and Figure 2). They use that expression to explain the experimental results of Ralph et al. [28]. Their calculation was done in the U→∞ limit. In addition, they include charging effects, although they do not say explicitly in which way these effects are included. In this respect, one has to realize that U has important contributions to the QD mesoscopic charging effect. In t→-∞, the leads and the QD maintain independent thermal equilibrium, that is, are uncoupled systems. When they become coupled, the Keldysh formalism yields the general behavior of the system. After a long enough time, this particular system reaches a steady state.Our point of view is taken from the fact that the charging of the QD is the origin of the Coulomb repulsion between two electron is occupying a two-fold degenerate level. Therefore, we study the behavior of a noninteracting QD att→-∞ where exact expressions are found. In this way, we obtain a formally similar expression (see (61) below) for the current as Equation (12) in the work of Meir and Wingreen [33]. Later on, Coulomb repulsion is introduced via a self-consistent field (SCF) that depends dynamically on the applied bias (HQD+USCF) and, in consequence, on the actual number of electrons in the QD. This approach constitutes the coupled Poisson NEGF formalism that has been discussed in the context of molecular conductors by Datta et al. [4, 24]. We use a capacitive model in Section 4 to calculate USCF and as discussed above, a numerical procedure is used to evaluate the current. Our approach has the known disadvantage of ignoring correlations in the QD (as pointed out in [39]). In that sense, there is a proposal by Datta (Equation (3.4.9) in [4]) that improves the SCF method and permits more accurate quantitative results. In Section 4, we apply this improvement for the case when the Coulomb charging is greater than the value of the coupling constants. We discuss possible improvements of our approach in Section 6.


## 2. Single Level QD-Model: Derivation of Nonequilibrium Currents

In macroscopic systems, the task of deriving transport equations or generalized Ginzburg-Landau equations relies on quasiclassical Green functions [7]. In addition, recently non-equilibrium transport in dirty Aluminium quasi-one-dimensional nanowires coupled with normal reservoirs [11] was studied experimentally and theoretically with quasi-classical Green functions [13]. As we want to include the possibility of particle interference effects, we do no resort to such objects. This point of view has been discussed in [14]. Instead, we use the equation of motion method (EOM) technique of Keldysh formalism for generating non-equilibrium states (see [8–10, 26]). We consider a spin degenerated single orbital as a quantum dot connected to superconductors leads. The Hamiltonian which describes this system is a generalized Anderson model [40]. It reads
(1)H=HS+HQD+HT,
where HS, HQD, and HT stand for the superconducting leads, the dot, and the tunneling term, respectively. Also HS=∑η‍Hη=HL+HR, where HL and HR are the left and right lead Hamiltonians, respectively. They are given, within the BCS model [41], by
(2)HS=∑ηk→σ‍Ψηk→σ†Hηk→0Ψηk→σ,
with
(3)Hηk→0=(εηk→Δηk→Δηk→*-εηk→),
where εηk→ is the conduction electron energy and Δηk→ is the superconductor gap of the lead η=L,R. Ψηk→σ† and Ψηk→σ are the Nambu spinors:
(4)Ψηk→σ†=(aηk→σ†aη,-k→,-σ),Ψηk→σ=(aηk→σaη,-k→,-σ†).
Here aηk→σ†(aηk→σ) denotes the creation (annihilation) operator for a conduction electron with wave vector k→ and spin σ in the η=L,R superconductor lead.H
QD is the hamiltonian for the single-level quantum dot of energy E0:
(5)HQD=∑σ‍ϕσ†HQDϕσ,
with
(6)HQD=(E000-E0).The model QD does not contain the Hubbard Coulomb repulsion interaction term. As explained in the introduction, Coulomb repulsion is modeled by means of the inclusion of capacitances, which are taken independent of the charge in the QD. The model also ignores possible superconducting correlations in the QD. For sufficiently small QDs, the discreteness of the single energy levels suppresses these correlations [37]. The position of the energy level will be treated first as fixed by the gate potential with respect to the left lead, while the effect of the applied voltage is taken into account by the coupled Poisson scheme. The tunneling hamiltonian HT is given by
(7)HT=∑ηk→σ‍Ψηk→σ†Hηk→Iϕσ,
with
(8)Hηk→I=(Vηk→00-Vηk→).HT connects the dot to the biased superconducting leads and it allows the electric charge flow. Vηk→ is the hybridization matrix element between a conduction electron in the η=L,R superconductor lead and a localized electron on the dot with energy E0. ϕσ† and ϕσ are the dot spinors:
(9)ϕσ†=(dσ†d-σ),ϕσ=(dσd-σ†);
here, dσ†(dσ) is the creation (annihilation) operator for an electron on the dot.The flow of electric charge from the terminalη is given by
(10)Iη(t)=(-e)[-d〈Nη(t)〉dt]=ieℏ〈[HT(t),Nη(t)]〉,
where -e is the electron charge. 〈⋯〉 is the thermodynamical average over the biased L and R leads at the temperature T, taken at time t0→-∞, as indicated in the Keldysh contour in Appendix A:
(11)〈⋯〉≡Tr⁡(ρ(t0)⋯),ρ(t0)≡e-β(H-μN)Tr⁡(e-β(H-μN)),
and Nη=aηk→σ†aηk→σ is the “number of particles” operator. Book-keeping calculation using (10) leads to
(12)Iη(t)=2eℏVηℜ∑k→σFηk→σ<(t,t).Fηk→σ<(t,t′)=i〈dσ†(t′)aηk→σ(t)〉 is the lesser Keldysh Green function,
(13)Fηk→σ(t,t′)≡-i〈TKaηk→σ(t)dσ†(t′)〉≡-iΘ(t,t′)〈aηk→σ(t)dσ†(t′)〉+iΘ(t′,t)×〈dσ†(t′)aηk→σ(t)〉≡Θ(t,t′)Fηk→σ>(t,t′)+Θ(t′,t)Fηk→σ<(t,t′),
and TK is the time-ordering operator, the action of which is to rearrange product of operators, such that operator with later times, on the Keldysh contour are placed to the left of the product. Hereafter, for simplicity, we replace Vηk→ by an average Vη at the Fermi surfaces (Vηk→≡〈|Vηk→|2〉FS) of the leads L and R. Using the scheme given in Appendix A for the rate of change of (13), we proceed to obtain the equation of motion:
(14)i∂Fηk→σ(t,t′)∂t=δ(t,t′)〈{aηk→σ(t),dσ†(t)}〉-i〈TK[aηk→σ(t),H]dσ†(t′)〉,
which leads to
(15)(i∂∂t-ϵηk→)Fηk→σ(t,t′)=-σΔηℱηk→σ(t,t′)+VηGσ(t,t′),
where
(16)ℱηk→σ(t,t′)=-i〈TKaηk→,-σ†(t)dσ†(t′)〉,(17)Gσ(t,t′)=-i〈TKdσ(t)dσ†(t′)〉.
Note that Gσ(t,t′) is the QD single-particle Green function. Similarly, ℱηk→σ(t,t′) satisfies the equation of motion:
(18)(i∂∂t+ϵηk→)ℱηk→σ(t,t′)=-σΔηFηk→σ(t,t′)-Vη𝒢σ(t,t′),
where
(19)𝒢σ(t,t′)=-i〈TKd-σ†(t)dσ†(t′)〉.
Here 𝒢σ(t,t′) is the QD of two-particle Green’s function.Equations (15) and (18) can be written in a compact form as follows (see Appendix B):
(20)(i∂∂t-ϵηk→σΔησΔηi∂∂t+ϵηk→)(Fηk→σ(t,t′)ℱ~ηk→σ(t,t′)ℱηk→σ(t,t′)F~ηk→σ(t,t′))=Vησz(Gσ(t,t′)𝒢~ηk→σ(t,t′)𝒢σ(t,t′)G~ηk→σ(t,t′)).
We introduce the tilde Keldysh-Green functions:
(21)ℱ~ηk→σ(t,t′)=-i〈TKaηk→σ(t)d-σ(t′)〉,F~ηk→σ(t,t′)=-i〈TKaη-k→-σ†(t)d-σ(t′)〉,𝒢~σ(t,t′)=-i〈TKdσ(t)d-σ(t′)〉,G~σ(t,t′)=-i〈TKd-σ†(t)d-σ(t′)〉.
Consider the following 2×2 matrix whose elements are the unperturbed Green-Keldysh functions, that is, defined for Vη=0:
(22)(gηk→σ(t,t′)fηk→σ(t,t′)fηk→σf~ηk→σ(t,t′)g~ηk→σ(t,t′)),
where
(23)gηk→σ(t,t′)≡-i〈TKaηk→σ(t)aηk→σ†(t′)〉0f~ηk→σ(t,t′)≡-i〈TKaηk→σ(t)aη-k→-σ(t′)〉0fηk→σ(t,t′)≡-i〈TKaη-k→-σ†(t)aηk→σ†(t′)〉0g~ηk→σ(t,t′)≡-i〈TKaη-k→-σ†(t)aη-k→-σ(t′)〉0.According to AppendixA, their equations of motion are given by
(24)(i∂∂t-ϵηk→)gηk→σ(t,t′)+σΔηfηk→σ(t,t′)=δ(t,t′),(i∂∂t-ϵηk→)f~ηk→σ(t,t′)+σΔηg~ηk→σ(t,t′)=0,(i∂∂t+ϵηk→)fηk→σ(t,t′)+σΔηgηk→σ(t,t′)=0,(i∂∂t+ϵηk→)g~ηk→σ(t,t′)+σΔηf~ηk→σ(t,t′)=δ(t,t′).
These equations can be written in matrix form as follows:
(25)(i∂∂t-ϵηk→σΔησΔηi∂∂t+ϵηk→)(gηk→σ(t,t′)f~ηk→σ(t,t′)fηk→σ(t,t′)g~ηk→σ(t,t′))=(δ(t,t′)00δ(t,t′)).Equation (20) can be written as an integral along the Keldysh contour CK (for an explanation see Appendix B)
(26)(Fηk→σ(t,t′)ℱ~ηk→σ(t,t′)ℱηk→σ(t,t′)F~ηk→σ(t,t′))=∫CKdt′′(gηk→σ(t,t′′)f~ηk→σ(t,t′′)fηk→σ(t,t′′)g~ηk→σ(t,t′′))×Vησz(Gσ(t′′,t′)𝒢~σ(t′′,t′)𝒢σ(t′′,t′)G~σ(t′′,t′)).From the last expression one can read forFηk→σ(t,t′) the following equation:
(27)Fηk→σ(t,t′)=Vη∫CK‍dt′′[gηk→σ(t,t′′)Gσ(t′′,t′)-f~ηk→σ(t,t′′)𝒢σ(t′′,t′)].We now apply the procedure explained in AppendixC; in order to obtain the Fηk→σ(t,t′) lesser component, we obtain
(28)Fηk→σ<(t,t′)=Vη{∫-∞∞dt′′[gηk→σ(r)(t,t′′)Gσ<(t′′,t′)-f~ηk→σ(r)(t,t′′)𝒢σ<(t′′,t′)]+∫-∞∞dt′′[gηk→σ<(t,t′′)Gσ(a)(t′′,t′)-f~ηk→σ<(t,t′′)𝒢σ(a)(t′′,t′)]∫-∞∞‍}.
Furthermore, the superscripts (<), (>), (r), and (a) correspond to lesser, greater, retarded, and advanced Green’s functions, respectively.Therefore, from (12), Iη(t) can be written as follows:
(29)Iη(t)=Iη(1)(t)+Iη(2)(t),
with
(30)Iη(1)(t)=2eℏℜ∑σ‍∫-∞∞dt′{[Vη2∑k→gηk→σ(r)(t,t′)]Gσ<(t′,t)+[Vη2∑k→gηk→σ<(t,t′)]Gσ(a)(t′,t)},(31)Iη(2)(t)=-2eℏℜ∑σ‍∫-∞∞dt′{[Vη2∑k→f~ηk→σ(r)(t,t′)]𝒢σ<(t′,t)+[Vη2∑k→f~ηk→σ<(t,t′)]𝒢σ(a)(t′,t)}.When applying the Fourier transformations, (30) and (31) can be expressed as follows:
(32)Iη(1)(t)=2ehℜ∑σ‍∫-∞∞‍dω∫-∞∞‍dω′2πe-i(ω-ω′)t×[Ση(r)(ω)Gσ<(ω,ω′)+Ση<(ω)Gσ(a)(ω,ω′)],Iη(2)(t)=-2ehℜ∑σ{e-2iμηt∫-∞∞dω×∫-∞∞dω′2πe-i(ω-ω′)t×[Ξ~η(r)(ω)σ𝒢σ<(ω,ω′)+Ξ~η<(ω)σ𝒢σ(a)(ω,ω′)Ξ~η(r)]∫-∞∞dω∫-∞∞dω′2πe-i(ω-ω′)t},
with
(2)Vη2∑k→gηk→σ(r)(t,t′)≡∫-∞∞dω2πe-iω(t-t′)Ση(r)(ω),Vη2∑k→gηk→σ<(t,t′)≡∫-∞∞dω2πe-iω(t-t′)Ση<(ω),Vη2∑k→‍f~ηk→σ(r)(t,t′)≡∫-∞∞e-2iμηtσdω2πe-iω(t-t′)Ξ~η(r)(ω),Vη2∑k→‍f~ηk→σ<(t,t′)≡∫-∞∞e-2iμηtσdω2πe-iω(t-t′)Ξ~η<(ω).In AppendicesD to G, we evaluate the unperturbed Green’s functions gηk→σ(r)(t,t′), gηk→σ<(t,t′), f~ηk→σ(r)(t,t′), and f~ηk→σ<(t,t′) in the wide band limit.We summarize these results as follows:(34)Ση(r)(ω)=-Γη[ω-μηΔηζ(Δη,ω-μη)+iζ(ω-μη,Δη)],Ση<(ω)=2iΓηζ(ω-μη,Δη)f(ω-μη),Ξ~η(r)(ω)=Γη[ζ(Δη,ω+μη)+iΔηω+μηζ(ω+μη,Δη)],Ξ~η<(ω)=-2iΓηΔηω+μηζ(ω+μη,Δη)f(ω+μη),ζ(ω,ω′)≡Θ(|ω|-|ω′|)|ω|ω2-ω′2.All these expressions will used below.


## 3. QD Green Function

We need to evaluate the most important objet for calculations, namely, QD Green’s functions given by (17) and (19), as well as their respective tilde functions:
(35)G~σ(t,t′)=-i〈TKd-σ†(t)d-σ(t′)〉,𝒢~σ(t,t′)=-i〈TKdσ(t)d-σ(t′)〉.Again using the scheme given in AppendixA, their equation of motion is(36)i∂∂t(Gσ(t,t′)𝒢~σ(t,t′)𝒢σ(t,t′)G~σ(t,t′))=(δ(t,t′)-i〈TK[dσ(t),H]dσ†(t′)〉-i〈TK[dσ(t),H]d-σ(t′)〉-i〈TK[d-σ†(t),H]dσ†(t′)〉δ(t,t′)-i〈TK[d-σ†(t),H]d-σ(t′)〉),which develops to
(37)i∂∂tGσ(t,t′)=δ(t,t′)-iE0〈TKdσ(t)dσ†(t′)〉-i∑ηk→‍Vη〈TKaηk→σ(t)dσ†(t′)〉=δ(t,t′)+E0Gσ(t,t′)+∑ηk→‍VηFηk→σ(t,t′),i∂∂t𝒢σ(t,t′)=iE0〈TKd-σ†(t)dσ†(t′)〉+i∑ηk→‍Vη〈TKaη-k→-σ†(t)dσ†(t′)〉=-E0𝒢σ(t,t′)-∑ηk→‍Vηℱηk→σ(t,t′),i∂∂t𝒢~σ(t,t′)=-iE0〈TKdσ(t)d-σ(t′)〉-i∑ηk→Vη〈TKaηk→σ(t)d-σ(t′)〉=E0𝒢~σ(t,t′)+∑ηk→Vηℱ~ηk→σ(t,t′),i∂∂tG~σ(t,t′)=δ(t,t′)+iE0〈TKd-σ†(t)d-σ(t′)〉+i∑ηk→Vη〈TKaηk→σ†(t)d-σ(t′)〉=δ(t,t′)-E0G~σ(t,t′)-∑ηk→Vηℱ~ηk→σ(t,t′).
This can be written as follows:
(38)(i∂∂t-E000i∂∂t+E0)(Gσ(t,t′)𝒢~σ(t,t′)𝒢σ(t,t′)G~σ(t,t′))=(δ(t,t′)00δ(t,t′))+∑ηk→Vησz(Fηk→σ(t,t′)ℱ~ηk→σ(t,t′)ℱηk→σ(t,t′)F~ηk→σ(t,t′)).
When Vη=0, one has
(39)(i∂∂t-E000i∂∂t+E0)(G0(t,t′)00G~0(t,t′))=(δ(t,t′)00δ(t,t′)),
with
(40)G0(t,t′)=-i〈TKdσ(t)dσ†(t′)〉0,G~0(t,t′)=-i〈TKd-σ†(t)d-σ(t′)〉0,G0(t,t′)≡Gσ(t,t′)|Vη=0,G~0(t,t′)≡G~σ(t,t′)|Vη=0.The last two equations can be written as follows:(41)(i∂∂t-E000i∂∂t+E0)×(Gσ(t,t′)-G0(t,t′)𝒢~σ(t,t′)𝒢σ(t,t′)G~σ(t,t′)-G~0(t,t′))=∑ηk→Vησz(Fηk→σ(t,t′)ℱ~ηk→σ(t,t′)ℱηk→σ(t,t′)F~ηk→σ(t,t′)).We write the last equation in its equivalent convolution integral along the Keldysh contour (see AppendixB):
(42)(Gσ(t,t′)-G0(t,t′)𝒢~σ(t,t′)𝒢σ(t,t′)G~σ(t,t′)-G~0(t,t′))=∫CKdt′′(G0(t,t′′)00G~0(t,t′′))×∑ηk→Vησz(Fηk→σ(t′′,t′)ℱ~ηk→σ(t′′,t′)ℱηk→σ(t′′,t′)F~ηk→σ(t′′,t′)).An equivalent way to write the last equation (using (25)) as a convolution of Σσ(t,t′) and Gσ(t,t′) is
(43)Gσ(t,t′)=G0(t,t′)+∫CKdt′′G0(t,t′)Σσ(t′,t′′)Gσ(t′′,t′),
with
(44)Gσ(t,t′)≡(Gσ(t,t′)𝒢~σ(t,t′)𝒢σ(t,t′)G~σ(t,t′)),G0(t,t′)≡(G0(t,t′)00G~0(t,t′)),Σσ(t,t′)≡∫CKdt′′(Vη2∑ηk→gηk→σ(t′′,t′)-Vη2∑ηk→f~ηk→σ(t′′,t′)-Vη2∑ηk→fηk→σ(t′′,t′)Vη2∑ηk→g~ηk→σ(t′′,t′)).We are interested in two regimes: a first regime in whichU0~Γ<Δ and the Coulomb blockade effects are neglected because in this case the couplings to the leads are not extremely small and the dot capacitance is large enough, a second regime for U0~Δ>Γ where Coulomb blockade effects must be taken into account. For both regimes and from now on, we are interested in the case eV>Δ, where multiple Andreev reflection [42] processes are strongly suppressed. Therefore only the single particle current (SP) has to be considered ISP. From the above considerations, we have that the Keldysh Green function 𝒢σ(ω), which carries information of the quantum dot two-particle Green's function, can be neglected and all relevant information is contained in Gσ(ω).The Keldysh Green function becomes spin independent;Gσ(ω)≡G(ω). Element 11 of (43) is given by
(45)G(t,t′)=G0(t,t′)+∫CKdt′′∫CKdt′′′G0(t,t′′)×Σ(t′′,t′′′)G(t′′′,t′).Again, using the recipe given in AppendixC, we obtain for G<(t,t′) and G(a)(t,t′) the following:
(46)G<(t,t′)=G0<(t,t′)+[∫-∞∞dt′′∫-∞∞‍dt′′′G0(r)(t,t′′)×Σ(r)(t′′,t′′′)G<(t′′′,t′)+G0(r)(t,t′′)Σ<(t′′,t′′′)×G(a)(t′′′,t′)+G0<(t,t′′)×Σ(a)(t′′,t′′′)G(a)(t′′′,t′)∫-∞∞dt′′∫-∞∞‍dt′′′G0(r)(t,t′′)Σ(r)],G(a)(t,t′)=G0(a)(t,t′)+∫-∞∞dt′′∫-∞∞‍dt′′′G0(a)(t,t′′)Σ(a)(t′′,t′′′)G(a)(t′′′,t′).
Taking the Fourier transform of (46) results in a set of algebraic equations:
(47)G<(ω,ω′)=2πδ(ω-ω′)G0<(ω)+G0(r)(ω)Σ(r)(ω)G<(ω,ω′)+G0(r)(ω)Σ<(ω)G(a)(ω,ω′)+G0<(ω)Σ(a)(ω)×G(a)(ω,ω′),G(a)(ω,ω′)=2πδ(ω-ω′)G0(a)(ω)+G0(a)(ω)Σ(a)(ω)G(a)(ω,ω′).
Dot Keldysh Green’s functions Gσ<(ω,ω′) and Gσ(a)(ω,ω′) are below straightforward evaluated. In this regime, quantities such as currents are independent of time. Therefore, we have
(48)G<(ω,ω′)=2πδ(ω-ω′)G<(ω),G(a)(ω,ω′)=2πδ(ω-ω′)G(a)(ω).
Therefore (47) result in
(49)G<(ω)=G0<(ω)+G0(r)(ω)Σ(r)(ω)G<(ω)+G0(r)(ω)Σ<(ω)G(a)(ω)+G0<(ω)Σ(a)(ω)G(a)(ω),(50)G(a)(ω)=G0(a)(ω)+G0(a)(ω)Σ(a)(ω)G(a)(ω).
Solving (50),
(51)G(a)(ω)=1G0(a)(ω)-1-Σ(a)(ω)=1ω-E0-Σ(a)(ω)=G(r)(ω)*.
Moreover, we know that
(52)G0<(ω)∝δ(ω-E0),G0(a)(ω)=(ω-E0-i0+)-1,
Resulting in
(53)G0<(ω)Σ(a)(ω)G(a)(ω)=-G0<(ω).
Equation (49) is reduced to
(54)G<(ω)=G0(r)(ω)Σ(r)(ω)G<(ω)+G0(r)(ω)Σ<(ω)G(a)(ω),(55)G<(ω)=Σ<(ω)G(a)(ω)G0(r)(ω)-1-Σ(r)(ω)=Σ<(ω)|G(r)(ω)|2=πΣ<(ω)-ℑG(r)(ω)/πℑ(G(r)(ω))-1=πΣ<(ω)-ℑΣ(r)(ω)ρ(ω).
Here ρ(ω) is the so-called quantum dot spectral function which is given in terms of the imaginary part (ℑ) of the retarded Keldysh Green function G(r)(ω):
(56)ρ(ω)=-1πℑG(r)(ω)=-1πℑΣ(r)(ω)ω-E0-ℜΣ(r)(ω)2+ℑΣ(r)(ω)2.
From (30), the single particle current (ISP) results in
(57)IηSP(V,E0)=4ehℜ∫-∞∞dω[Ση(r)(ω)G<(ω)+Ση<(ω)G(a)(ω)Ση(r)].
Substituting (51) and (55) in (57),
(58)IηSP(V,E0)=4eh∫-∞∞dω[πℑΣη(r)(ω)ℑΣ<(ω)ℑΣ(r)(ω)ρ(ω)+ℑΣη<(ω)ℑG(r)(ω)ℑΣ<(ω)ℑΣ(r)(ω)]=4πeh∫-∞∞dωρ(ω)[ℑΣη(r)(ω)ℑΣ<(ω)ℑΣ(r)(ω)ρ(ω)-ℑΣη<(ω)ℑΣ<(ω)ℑΣ(r)(ω)]=4πeh∫-∞∞dωρ(ω)Γ(ω)[Γη(ω)ℑΣ<(ω)-Γ(ω)ℑΣη<(ω)],
with Γη(ω)=-ℑΣη(r)(ω)=Γηζ(ω,Δη) and Γ(ω)=∑η‍Γη(ω). In our regime, eV>Δ; therefore, ℜΣ(r)(ω) in the above equations is zero. We use the expression for Σ(r)(ω) from Appendix D and obtain the single particle current ISP≡(IR,SP-IL,SP)/2:
(59)ISP(V,E0)=8πeh∫-∞∞‍dωΓL(ω)ΓR(ω+eV)ΓL(ω)+ΓR(ω+eV)×ρ(ω)[f(ω)-f(ω+eV)].-eV=μL-μR corresponds to the applied voltage between the superconductors electrodes with chemical potential μη. In the following, we fix the chemical potential μL=0 and use eV as a measure of μR. In addition, the QD energy E0 is measured with respect to μL. On the other hand, the limits of integration are given by the functions ΓL(ω) and ΓR(ω+eV). The extra 2π factor arises from the dot Keldysh Green functions. ρ(ω) and Γ(ω) are given by
(60)ρ(ω)=Γ(ω)/π(ω-E0)2+Γ2(ω),(61)Γ(ω)=ΓL(ω)+ΓR(ω+eV).At steady state there is no net flow into or out of the mesoscopic channel or quantum dot which yields a stationary particle number in it. The population numberN, at the dot, is given by
(62)N=2[-iG<(t,t)]=2∫-∞∞dω2πiG<(ω),
which becomes a weighted average over the L and R contacts:
(63)N=2∫-∞∞dωρ(ω)[ΓL(ω)Γ(ω)f(ω)+ΓR(ω+eV)Γ(ω)f(ω+eV)].
For the N/QD/N case, ΓR,L are just constants. This case was studied in the context of the generalized quantum master approach (section IV in [39]). That approach permits the inclusion of broadening in a natural way. They obtained Equations similar to (59)–(63).


## 4. Coupled Poisson Nonequilibrium Green Function Scheme: The Capacitive Model

So far, we are not including the side effects of a potential profile inside the mesoscopic channel. On the one hand, its inclusion takes in order zero or Hartree approximation the electron-electron interaction in the QD. Its inclusion also guarantees current independence from the choice of zero potential [32]. Such potential is induced by the action of source, drain and gate applied voltages. In principle, we have to couple the number of population equations. Equation (62), with electric field U. However, since the number of quantum levels in the channel is small and the particle number variation is negligible, the potential profile variation inside the channel is negligible. Then it is appropriate to visualize the channel as an equivalent circuit framework (Figure 2). In this framework, we associate capacitances Cd, Cs, and Cg with the drain, source and gate, respectively. Whenever drain, source, and gate bias potentials Vd, Vs, and Vg, respectively, are present, there is an electrostatic potential VQD inside the QD, which induces an energy shift of the QD energy level U=-e(VQD-V0), V0 are channel electrostatic potential before we apply the source and drain biases, respectively.The electronic populations before and after we apply the biases mentioned above are given by(64)-eN0=CdV0+CsV0+CgV0,-eN=Cd(VQD-Vd)+Cs(VQD-Vs)+Cg(VQD-Vg),
respectively. It leads us to
(65)-eΔN≡-e(N-N0)=CE(VQD-V0)-CdVd-CsVs-CgVg,
where CE=Cd+Cs+Cg. Therefore, the energy shift U is given by
(66)U=Uℒ+e2CEΔN,
where
(67)Uℒ≡CdCE(-eVd)+CsCE(-eVs)+CgCE(-eVg).
In the expression for U, Uℒ represents a uniform shift for all levels, whereas the second term (the Poisson contribution denoted by UP in the introduction) represents a level of repulsion which is proportional to the averaged occupation of the QD level denoted by N0 and proportional to the charging energy U0=e2/CE.On the other hand, one hasΔN from (62) and (65) given by
(68)ΔN=2∫-∞∞‍dω2πi[G<(ω,U)-G<(ω,-eV0)].
In the expression for G<(ω,U) (see (55)), the energy level shifts only (E0⇒(E0+U)) in the expression for the QD spectral function ρ(ω-U). Equations (66) and (68) are coupled nonlinear equations with unknowns U and ΔN. We solve the coupled equations via an iteration procedure. First we guest a value for ΔN plug this value in U, and then we calculate ΔN with the following equation:
(69)ΔN=2∫-∞∞dωρ(ω-U)×ΓL(ω)f(ω)+ΓR(ω+eV)f(ω+eV)ΓL(ω)+ΓR(ω+eV),
and so on until convergence is achieved. With the final value of U obtained for a given bias voltage V, ISP is calculated via the equation:
(70)ISP(V,U)=8πeh∫-∞∞dωΓL(ω)ΓR(ω+eV)ΓL(ω)+ΓR(ω+eV)×ρ(ω-U)[f(ω)-f(ω+eV)].In summary, the procedure for computingI consists of the following steps. (i) Determine the spectral density. (ii) Specify Vg, Vd, Vs, and coupling constants. (iii) Iteratively solve (69) and (66). (iv) Evaluate the current from (70) for the Vg, Vd, and Vs. Once a converged U has been found, the current is finally evaluated.The way we consider electron-electron interactions imposes restrictions on the possible values of the charging energyU0. For the self-consistent scheme to be valid, we have to assume that Δ≫ΓL,R≃U0. However, less precisely quantitative results, although qualitative correct results can be obtained if Δ≃U0≫ΓL,R, when the so-called Coulomb Blockade energy dominates over the coupling constants. For this case, we use the improvement of the SCF method discussed in the introduction [4, 29, 30]. The self-consistent generalizes to
(71)U↑=Uℒ+e2CE(N↓-N0),(72)U↓=Uℒ-e2CE(N↑-N0),
where the up-spin level feels a potential due to the down-spin electrons and viceversa. Notice the different signs which reflects the Coulomb repulsion between otherwise degenerate levels:
(73)N↑=∫-∞∞dωρ(ω-U↑)×ΓL(ω)f(ω)+ΓR(ω+eV)f(ω+eV)ΓL(ω)+ΓR(ω+eV),(74)N↓=∫-∞∞dωρ(ω-U↓)×ΓL(ω)f(ω)+ΓR(ω+eV)f(ω+eV)ΓL(ω)+ΓR(ω+eV),(75)N=N↑+N↓.
Here, N↑ and N↓ are the population of the spin-up and spin-down levels. Once the values of U↑ and U↓ are calculated, ISP is calculated from
(76)ISP(V,U↑,U↓)=4πeh∫-∞∞dωΓL(ω)ΓR(ω+eV)ΓL(ω)+ΓR(ω+eV)×{ρ(ω-U↑)+ρ(ω-U↓)}×[f(ω)-f(ω+eV)].
As Datta has pointed out [4], the approach described above (called unrestricted SCF) can lead to a better quantitative agreement in comparison with a conceptually correct multi-level Master equation calculation.


## 5. Numerical Results and Remarks on Experiments

### 5.1. First Case:Δ≫ΓL,R≃U0

In this regime, there are the multiple Andreev reflections [42] for voltages such that eV<Δ (MAR). Also, there is the possibility for quasiparticle co-tunneling current for energy levels far from μL. These cases will be considered in a future work that involves the whole expressions we have derived for the currents (see (30) and (31)) and eventually more accurate Green-Keldysh functions and the use of master equations [17]. This case was studied experimentally in [43]. For given values of capacitances and source voltage, we iterate (66) and (68) in order to find the potential U. Then the single particle current ISP(V,U) is evaluated (see (70)). We put the charge before biasing N0=0, such that Coulomb repulsion with the QD-energy level is absent. Anyhow, in this regime, the effect of the second term in (74) is negligible. Consequently, the Laplace term Uℒ essentially positions the QD degenerate energy level (with respect to μL=0). In Figure 3, we show I-V characteristic for gate voltage values Vg=0 and KBT=Δ, whereas Figure 4 shows the occupation number ΔN. These curves are symmetric, due to the assumed equality of the coupling capacitances (Cd/CE=0.5). Otherwise, the I-V shifts to right or to the left for Cd/CE=0.5>0.5 or <0.5 respectively. For this case, we show in Figure 7 the spectral density ρ(ω) for eVd=-6Δ. Notice that the position of the energy level is essentially -4.5Δ, that is, just the sum of E0+Uℒ. Qualitatively, these results are similar to Levy Yeyati et al.’s [17]. Characteristic is the broadening of the BSC singularity. The effect of bigger values of ΓR,L is a more pronounced round off the BCS-type singularity. We discuss this issue below. For large enough bias, the current approaches the normal saturation value ISat.Figure 3
Zero temperatureI-V characteristics for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=0.005Δ, Cd/CE=0.5, and ΓL=ΓR=0.005Δ.Figure 4
Zero temperature number of electronseVd/Δ graph for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=0.005Δ, Cd/CE=0.5, and ΓL=ΓR=0.005Δ.


### 5.2. Second Case:Δ≃U0≫ΓL,R

In this regime, the charging energy acts effectively in lifting the degeneracy of the otherwise single degenerate QD energy level. For this regime, we use the couple system defined by (71)–(75) and calculate the current according to (76). This is the unrestricted SCF method mentioned in the introduction. The transport begins through one level as long as there is in average less than one electron in it. For the given parameters, the onset of current is similar to the first case (no interaction with residual charge in the QD is considered). However, when the average occupation exceeds one, the other degenerate levels float according to the resulting values U↑ and U↓. These values push down the position of this second level and push up the already occupied energy level. In Figures 5-6, we show I-V and the number of electrons. In Figure 8, it is shown that the spectral density for eVd=-8Δ. In this case, E0+Uℒ=-5.5Δ. The values obtained from the SCF calculation position the energy levels to -4.846Δ and to -6.396Δ (see (72)).Figure 5
Zero temperatureI-V characteristic showing the Coulomb blockade for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓL=ΓR=0.01Δ.Figure 6
Zero temperature number of electrons-eVd/Δ for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓL=ΓR=0.01Δ.Figure 7
Spectral density of the quantum dot-ω for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, eVd=6.0Δ, U0=0.01Δ, Cd/CE=0.5, and ΓL=ΓR=0.05Δ.Figure 8
Spectral density of the quantum dot-ω for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, eVd=8.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓL=ΓR=0.20Δ.The experimental work of Ralph et al. [28] corresponds to this second case. In their Figure 3, they show the current for single state level peaks I~5 pA at a bias V~2.4 mV. According to Figures 5-6, for this sample at 30 mK, there should be another peak at V~4.4 mV with a current value ISP≃10 pA. The charging energy for this sample is U0≃2.0 mV. However, in their samples with radios ~2.5 nm or greater, the level spacing of the energy levels is such, that ΔE<Ec; therefore another current signal may occur before and a quantitative proper description would be a multilevel QD-model. One notice, however, the strong fluctuations in the spacing in Figure 2 [28], indicating complex charging for many levels of phenomena. Notice that the theoretical explanation of Levy Yeyati et al. [17, 18] does not contain our prediction.


### 5.3. Influence of Coupling Constants

The influence of the coupling constant is to broaden the otherwise sharp energy QD level. However, the broadening is not equally strong and depends on the relative values ofΓR/ΓL. Ralph et al. [28] determined for the sample in his Figure 3  ΓR/ΓL≫1. Figure 10 shows the I-V characteristic for the restricted case when both coupling constants are equal. For larger values of the coupling constants, the broadening is stronger. If they are dissimilar in value, the broadening is stronger when ΓR>ΓL. This effect is shown in Figures 9, 10, 11, 12, 13, and 14. This effect is due to the stronger involving of the BSC-DOS singularity of the left lead in the integral expression for the current (see (59)) and the particular choice of the zero bias voltage (see Figure 2).Figure 9
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=0.05Δ, Cd/CE=0.5, and ΓL=ΓR.Figure 10
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=0.05Δ, Cd/CE=0.5, and ΓL=4ΓR.Figure 11
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=0.05Δ, Cd/CE=0.5, and ΓR=4ΓL.Figure 12
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓL=ΓR.Figure 13
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓL=4ΓR.Figure 14
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓR=4ΓL.


## 5.1. First Case:Δ≫ΓL,R≃U0

In this regime, there are the multiple Andreev reflections [42] for voltages such that eV<Δ (MAR). Also, there is the possibility for quasiparticle co-tunneling current for energy levels far from μL. These cases will be considered in a future work that involves the whole expressions we have derived for the currents (see (30) and (31)) and eventually more accurate Green-Keldysh functions and the use of master equations [17]. This case was studied experimentally in [43]. For given values of capacitances and source voltage, we iterate (66) and (68) in order to find the potential U. Then the single particle current ISP(V,U) is evaluated (see (70)). We put the charge before biasing N0=0, such that Coulomb repulsion with the QD-energy level is absent. Anyhow, in this regime, the effect of the second term in (74) is negligible. Consequently, the Laplace term Uℒ essentially positions the QD degenerate energy level (with respect to μL=0). In Figure 3, we show I-V characteristic for gate voltage values Vg=0 and KBT=Δ, whereas Figure 4 shows the occupation number ΔN. These curves are symmetric, due to the assumed equality of the coupling capacitances (Cd/CE=0.5). Otherwise, the I-V shifts to right or to the left for Cd/CE=0.5>0.5 or <0.5 respectively. For this case, we show in Figure 7 the spectral density ρ(ω) for eVd=-6Δ. Notice that the position of the energy level is essentially -4.5Δ, that is, just the sum of E0+Uℒ. Qualitatively, these results are similar to Levy Yeyati et al.’s [17]. Characteristic is the broadening of the BSC singularity. The effect of bigger values of ΓR,L is a more pronounced round off the BCS-type singularity. We discuss this issue below. For large enough bias, the current approaches the normal saturation value ISat.Figure 3
Zero temperatureI-V characteristics for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=0.005Δ, Cd/CE=0.5, and ΓL=ΓR=0.005Δ.Figure 4
Zero temperature number of electronseVd/Δ graph for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=0.005Δ, Cd/CE=0.5, and ΓL=ΓR=0.005Δ.


## 5.2. Second Case:Δ≃U0≫ΓL,R

In this regime, the charging energy acts effectively in lifting the degeneracy of the otherwise single degenerate QD energy level. For this regime, we use the couple system defined by (71)–(75) and calculate the current according to (76). This is the unrestricted SCF method mentioned in the introduction. The transport begins through one level as long as there is in average less than one electron in it. For the given parameters, the onset of current is similar to the first case (no interaction with residual charge in the QD is considered). However, when the average occupation exceeds one, the other degenerate levels float according to the resulting values U↑ and U↓. These values push down the position of this second level and push up the already occupied energy level. In Figures 5-6, we show I-V and the number of electrons. In Figure 8, it is shown that the spectral density for eVd=-8Δ. In this case, E0+Uℒ=-5.5Δ. The values obtained from the SCF calculation position the energy levels to -4.846Δ and to -6.396Δ (see (72)).Figure 5
Zero temperatureI-V characteristic showing the Coulomb blockade for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓL=ΓR=0.01Δ.Figure 6
Zero temperature number of electrons-eVd/Δ for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓL=ΓR=0.01Δ.Figure 7
Spectral density of the quantum dot-ω for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, eVd=6.0Δ, U0=0.01Δ, Cd/CE=0.5, and ΓL=ΓR=0.05Δ.Figure 8
Spectral density of the quantum dot-ω for superconductor-quantum-dot-superconductor system, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, eVd=8.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓL=ΓR=0.20Δ.The experimental work of Ralph et al. [28] corresponds to this second case. In their Figure 3, they show the current for single state level peaks I~5 pA at a bias V~2.4 mV. According to Figures 5-6, for this sample at 30 mK, there should be another peak at V~4.4 mV with a current value ISP≃10 pA. The charging energy for this sample is U0≃2.0 mV. However, in their samples with radios ~2.5 nm or greater, the level spacing of the energy levels is such, that ΔE<Ec; therefore another current signal may occur before and a quantitative proper description would be a multilevel QD-model. One notice, however, the strong fluctuations in the spacing in Figure 2 [28], indicating complex charging for many levels of phenomena. Notice that the theoretical explanation of Levy Yeyati et al. [17, 18] does not contain our prediction.


## 5.3. Influence of Coupling Constants

The influence of the coupling constant is to broaden the otherwise sharp energy QD level. However, the broadening is not equally strong and depends on the relative values ofΓR/ΓL. Ralph et al. [28] determined for the sample in his Figure 3  ΓR/ΓL≫1. Figure 10 shows the I-V characteristic for the restricted case when both coupling constants are equal. For larger values of the coupling constants, the broadening is stronger. If they are dissimilar in value, the broadening is stronger when ΓR>ΓL. This effect is shown in Figures 9, 10, 11, 12, 13, and 14. This effect is due to the stronger involving of the BSC-DOS singularity of the left lead in the integral expression for the current (see (59)) and the particular choice of the zero bias voltage (see Figure 2).Figure 9
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=0.05Δ, Cd/CE=0.5, and ΓL=ΓR.Figure 10
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=0.05Δ, Cd/CE=0.5, and ΓL=4ΓR.Figure 11
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=0.05Δ, Cd/CE=0.5, and ΓR=4ΓL.Figure 12
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓL=ΓR.Figure 13
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓL=4ΓR.Figure 14
Zero temperatureI-V characteristic for superconductor-quantum-dot-superconductor system for various values of the coupling Γ, calculated using the self-consistent field (SCF) method, with E0=1.5Δ, eVg=0.0Δ, U0=1.0Δ, Cd/CE=0.5, and ΓR=4ΓL.


## 6. Summary and Perspectives

We have studied the single particle current through a quantum dot coupled with two superconductor leads via a coupled Poisson Nonequilibrium Green function (PNEGF) formalism. In a systematic and self-contained way, we derived the expressions for the current in full generality. In this work we focused only on the weak coupling regime where single particle current is the dominant one. The QD is a single degenerate energy level system modeled via a capacitive circuit. The influence of the potential on the QD and on theI-V characteristic is calculated for relevance values of the coupling and capacitances and the implication of experiments is discussed. This was done in the weak coupling regime and for Δ≫ΓL,R≃U0. A second case when Δ≃U0≫ΓL,R also in the weak coupling regime was analyzed. Admittedly, our model of a hybrid system S/QD/S possesses potentially physical extensions. One important missed point is dephasing. This physical effect due to scattering of transport electrons can be incorporated in the self energy phenomenologically [44, 45], or in a stochastic fashion [46]. Another point is to consider a QD with many energy levels and within the self-consistent scheme, to consider the strong and intermediate regimes and many body correlations due to different kinds of electron-electron interaction. Here we have to notice that it is not just to scale the level spacing by the charging energy [47]. It is a genuine many body problem. But the most important missed point was correlations. As pointed out by Datta [4, chapter III], there has been much effort in order to find a suitable SCF that considers correlations. For example, to modify (66) to consider occupancies probabilities. As discussed in the introduction, Kang [16] and Meir and Wingreen [33] find a solution for the QD Green function that contains this type of correlation. In other words, one could go to scheme where a more accurate Green function for the QD is used together with a multielectron picture and corresponding master equation. These would mean to use the Anderson model with a Coulomb interaction U that is obtained from a SCF. We want to check if this point of view is correct. Work in this direction is in progress.



---

*Source: 721671-2014-04-09.xml*