Projectwork
Andreevstates in Josephson junctions
supervisor:
Dr. Jonas Fransson
Institutet fysik och astronomi
Susanne Wagner, BSc
11258888
Master Technische Physik
E 066 461
Contents
1 Introduction 1
2 Technical Details on Green function Formalism 1
3 Calculating Green functions for different systems 4
3.1 First Application. BCS Hamiltonian . . . 4
3.2 Andreev States in Josephson Junction . . . 7
3.3 Andreev States in Josephson Junction with Quantum Dot . . . 10
4 Analyzing electronic density of states starting from the Green function 12
4.1 BCS superconductor . . . 13
4.2 Josephson Junction . . . 14
4.3 Quantum Dot in a Josephson Junction . . . 19
5 Conclusion and Outlook 21
6 Appendix 23
A Calculation of Commutators with BCS Hamiltonian . . . 23
B Calculation of Commutators with JJ-Hamiltonian . . . 23
1
Introduction
The aim of this project work is to derive and understand the electronic density of states of a BCS superconductor when being isolated and when interacting with another superconductor or another superconductor and a quantum dot. What is observed is that through this interaction the density of states is different from the isolated system. The method that was used to derive the density of states is the Green function formalism.
In the first problem I want to derive the electronic density of states in a superconductor when interacting with another superconductor through a thin non-superconducting layer. This device is called a Josephson junction. It will be shown that new discrete states in the energy spectrum of the electrons will emerge. These new states are called Andreev bound states.
Josephson junctions are electronic devices of two superconductors which are parted by a thin layer of either an insulating material or a normal metal [?]. Cooper pairs which are the carriers of su-percurrent are able to tunnel through this thin layer and therefore the susu-percurrent is not stopped by the junction. The effect that enables the tunneling of the Cooper pairs is called Proximity effect. Due to the tunneling also new energystates in both superconductors will emerge. These so called Andreev states are of high importance because they may carry a high spectral weight and are mainly responsible for the flow of the supercurrent through the junction [2, S.2].
For the second problem there will be a quantum dot placed in between the superconduct-ing leads of the Josephson junction. Quantum dots are most commonly semiconductors formsuperconduct-ing a nano crystal or particle with diameters in the micro or nano meter regime. Due to the confinement of the electrons in these small structures the electronic states are discrete. The spacing between the states depends on the size and the material of the quantum dot. But not only semiconductors can be used as quantum dots. Also molecules, mostly carbon nanotubes are used as quantum dots in experiments.
Now the density of states for the quantum dot is of interest. When coupled to electrodes, either metallic or superconducting ones, the electronic structure of the quantum dot will change. For quantum dots coupled to superconducting electrodes again the Proximity effect will cause leaking of Cooper pairs into the quantum dot. This will lead to an alteration of the quantum dot’s density of states which also includes new discrete states, the Andreev states [2, S.1].
2
Technical Details on Green function Formalism
The Green function formalism is a very practical tool when describing a system in which the many-body-wavefunction or equivalently the partition function cannot be found analytically anymore. This might be the case because of interaction terms in the Hamiltonian. Green functions then help to calculate approximate expectation values of physically measurable quantities without even having to know the partition function. One can choose the degree to which one wants to approach the exact solution.
The central object is the Green function. The causal one particle Green function is defined as follows:
Gkk0σσ0(t, t0) := (−i)hT ckσ(t)c†
k0σ0(t0)i (2.1)
Here ckσand c†kσare the annihilation- and creation operators for any particle in Heisenberg picture.
ci(t), cj(t0)† += δ(t − t 0)δ ij and c†i(t), c † j(t 0) +=ci(t), cj(t 0) += 0 (2.2)
The average is taken over the grand canonical ensemble [3, S.133].
T is Wick’s time ordering operator and defined according to (2.3). The operators a and b are assumed to be fermionic operators, Wick’s time ordering operator includes a minus sign for each permutation of fermionic operators. For bosonic operators the sign would stay the same [3, 132ff].
T a(t)b(t0) = (
a(t)b(t0) if t > t0
−b(t0)a(t) if t < t0 (2.3)
The Green function as defined in (2.1) describes the propagation properties of a single electron and is called normal Green function. The corresponding Green function for the propagation properties of holes is simply the adjoint of the electronic Green function.
G†kk0σσ0(t, t0) := (−i)hT c
†
kσ(t)ck0σ0(t0)i (2.4)
The Green functions is a probability amplitude for a single electron (hole) to change it’s properties
from momentum k and spin σ at time t to momentum k0 and σ0 at time t0.
The aim of this formalism is to be able to calculate physically measurable quantities without having to know the partition function of the system. The central formula which allows this is the Spectral Theorem (2.5). From it one can connect correlation functions to the spectral density, which again can be derived from the Green function. The Spectral Theorem for two operators a and b states [3, S.142]: hb(t0)a(t)i = Z ∞ −∞ dωSab(ω) eβω+ 1e (−iω(t−t0)) (2.5)
In this formula Sab(ω) is the so called spectral density belonging to the Green function Gaab(t, t0) =
(−i)hT a(t)b(t0)i. The spectral density can be obtained from the causal Green function G(ω) in
energy space by analytically continuing the Green function and extracting the imaginary part. This is done by this formula:
Sab(ω) = i 2π Gab(ω + i0 + ) − Gab(ω − i0 + ) (2.6)
For a Green function which is already defined on the imaginary axis as well, the way to get the spectral density is [3, formula (3.188), S.149] given by:
S(ω, p) = −
1 π
eβω+ 1
eβω− 1Im[G(ω, p)] (2.7)
In both cases the Green function G(ω) is just the Fourier transform (2.13) of the Green function G(t, t0). For determining the analytical continuation Gab(ω ± i0
+
1
x − x0± i0+ = P
1
x − x0 ∓ iπδ(x − x0) (2.8)
Now one can see that from the Spectral Theorem any expectation value of the form hb(t)a(t0)i can
be obtained just by defining suitable Green functions Gab. In the case of the one electron Green
function (2.1), where a = ckσ(t) and b = c
†
k0σ0(t0) one obtains the occupation number operator
nkσ= c†kσckσ by putting all primed indices equal the not primed ones.
Now that Green functions were introduced and also the way in which they are useful, what is still left to understand is how one can actually calculate them. Therefore the equation of motion for the Green function will be derived by taking the total time derivative of the Green function as defined in (2.1). i∂tGkk0σσ0(t, t0) = ∂t Θ(t − t0)hckσ(t)c†k0σ0(t0)i − Θ(t0− t)hc † k0σ0(t0)ckσ(t)i (2.9)
With integration by parts, follows:
i∂tGkk0σσ0(t, t0) = δ(t − t0) hckσ(t)c†k0σ0(t0)i + hc † k0σ0(t0)ckσ(t)i +Θ(t − t0)h∂tckσ(t)c†k0σ0(t0)i − Θ(t0− t)hc † k0σ0(t0)∂tckσ(t)i (2.10)
In the next step Wicks time ordering operator is reintroduced.
i∂tGkk0σσ0(t, t0) = δ(t − t0)hckσ(t), c†k0σ0(t0)
+i + (−i)hT i∂tckσ(t)c
†
k0σ0(t0)i (2.11)
The term i∂tckσ(t) is replaced by the commutator of ckσ(t) and the Hamilton operator, by using the
Heisenberg equation of motion for an operator in the Heisenberg picture. Also the commutation relations for fermionic particles (2.2) are used in the next step.
i∂tGkk0σσ0(t, t0) = δ(t − t0)δkk0δσσ0+ (−i)hTckσ(t), H(t) −c
†
k0σ0(t0)i (2.12)
This is the final form for the equation of motion which is used in the following calculations. The commutation relations in the last term of (2.12) is usually proportional to one ore more electronic creation or annihilation operators depending on the form of the Hamiltonian. Therefore, formally the last term in the equation of motion (2.12) is again a Green function. If the Hamiltonian con-tains two electron interactions the resulting Green function is then no longer a one electron but a two electron Green function and is said to be of higher order than the first one. In the equation of motion for this Green function then the four electron Green function will appear next which describes the interaction of two electrons intermediated by two others. This way, a infinite chain of Green functions comes about where each Green function depends on the next higher order Green function. At some point one has to neglect the higher orders of interaction in order to be able to close the system of equations and calculate the Green functions.
Finding Gkk0σσ0 from equation (2.12) can be done easier after a Fourier transform. This way
Gkk0σσ0(ω) = Z +∞ −∞ d(t − t0)Gkk0σσ0(t − t0)eiω(t−t 0) (2.13) Gkk0σσ0(t − t0) = 1 2π Z +∞ −∞ dωGkk0σσ0(ω)e−iω(t−t 0) (2.14)
Applying the Fourier transform to the equation of motion (2.12) the formula which is found is this. ωGkk0σσ0(ω) = δkk0δσσ0+ F T (−i)hTckσ(t), H(t) −c † k0σ0(t0)i (2.15)
In energy space the equation of motion yields a coupled system of simple algebraic equations.
3
Calculating Green functions for different systems
3.1
First Application. BCS Hamiltonian
As a first application I want to show how to calculate the Green function for a single supercon-ductor. The theory from which we take the Hamiltonian to describe superconductivity is BCS theory, named after Bardeen, Cooper and Schrieffer who proposed this theory in 1957 [1]. It was the first microscopic theory on superconductivity.
With greek indices giving the defining quantum numbers of a Bloch state φ and with cα and c†α
creating or annihilating an electron in this Bloch state, in second quantization the most general form for the operator of the potential energy that considers two particle interaction is given by [3, S.27]: ˆ V =1 2 X α,β,γ,δ hφαφβ| ˆV |φδφγic†αc † βcγcδ (3.1)
In BCS theory additional assumptions and procedures will still change the appearance of this term in the Hamiltonian. I will not go deep into the derivation of the model Hamiltonian, with which we deal in the end, still I want to give some underlying assumptions and aspects of BCS theory One basic assumption of BCS theory is that Fermi liquid theory, which was developed to describe metals, also applies for superconductors only with one additional interaction, which is the
attrac-tive potential between electrons mediated by phonons. The interaction term ˆV therefore does not
only contain Coulomb-repulsion but also an attractive term. This enables electrons (and holes) to build so called Cooper pairs, the quasi particles which are responsible for the flux of supercurrent [4, S.55]. In Fermi liquid theory the important step was made from describing strongly interacting electrons to quasi particles which are independent from each other [4, S.57].
The greek indices contain information about momentum and spin direction of the electron. Ac-cording to BCS theory it is most likely that electrons or holes with different spin direction and opposite momentum form Cooper pairs. Therefor:
The interaction term also enables the Cooper pairs to scatter from k → k0 [4, S.58]. Anyway, here
we neglect the scattering between Cooper pairs and therefore put k ≡ k0 [4, S.60].
This interaction term can be further simplified by applying an effective field theory. The reason why one can apply an effective field theory here is that we are dealing with a large number of electrons (actually quasi particles). Therefore the amplitude (or probability) for exciting a Cooper-pair from a given state does not depend very much on the state itself anymore, as long as it is either the ground state or a low lying excitation state [4, S.63]. Therefore one can formally
substitute the operator c−k↓ck↑ by its thermodynamic average and only making a small error of
c−k↓ck↑− hc−k↓ck↑i. Using this substitution for both operators, c−k↓ck↑ and c†−k↓c†k↑ one finds
this when neglecting second order terms:
c†−k↓c†k↑c−k0↓ck0↑= hc† −k↓c † k↑ic−k0↓ck0↑+ c† −k↓c † k↑hc−k0↓ck0↑i + hc† −k↓c † k↑ihc−k0↓ck0↑i (3.3)
After applying this, the resulting full Hamiltonian takes the form:
HBCS = X p,σ pcpσ† cpσ+ X p ∆pcp↑† c†−p↓+ ∆∗pc−p↓cp↑ (3.4)
The new parameter ∆pis dependent on the interaction potential as well as on the thermodynamic
average hc−p↓cp↑i. Due to this procedure the BCS Hamiltonian formally does not contain two or
more electron interactions. So, the chain of higher order electron interactions does not appear in the equations of motion for the Green functions anymore.
Due to Fermi liquid theory the interacting electrons in the superconductor are replaced by nonin-teracting quasi particles which contain the interaction potential in their eigenenergies. Therefore the particles cannot scatter and stick with their initial momentum k and their initial spin direction. This means, that Gkk0σσ0 ≡ δkk0δσσ0Gkk0σσ0. Since the momentum indices are equal for simplicity
they will be omitted from now on. For the spin indices there are still two possibilities for the
electronic Green function, which are G↑↑ and G↓↓. They are completely equal since the system is
spin degenerate. Therefore we concentrate on the description of only one Green function G↑↑.
When calculating the normal electronic Green function, in the first step one inserts the
Hamilto-nian (3.4) into the equation of motion (2.12) for the normal Green function Gσσ and expresses
the commutator between the creation operator and Hamiltonian. The full calculation of the com-mutator is done in Appendix A.
ckσ(t), HBCS(t)
− = −kckσ(t) − ∆kδ↑σc
†
−k↓(t) + ∆−kδ↓σc†−k↑(t) (3.5)
The commutator (3.5) plugged into the equation of motion (2.15) yields for G↑↑:
(i∂t+ k)G↑↑(t − t0) = δ(t − t0) − ∆k(−i)hT c†−k↓(t)c†k↑(t
0)i (3.6)
Similarly the equation of motion for G†↓↓ yields:
(i∂t− k)G†↓↓(t − t0) = δ(t − t0) − ∆∗k(−i)hT c−k↑(t)ck↓(t0)i (3.7)
The commutatorc†kσ(t), HBCS(t)
− used for the Green function G
†
by taking its adjoint (see Appendix A).
Here one sees that two new Green functions (−i)hT c†−k↓(t)c†k↑(t0)i and (−i)hT c−k↑(t)ck↓(t0)i
ap-peared. Taking the second one, it describes the propagation of a single electron which previous to time t0 built a Cooper-pair together with an electron with momentum k and spin ↓. After this
elec-tron has been removed at time t0 the propagation up to the time t of the remaining electron with
momentum −k and spin ↑ is described by the Greenfuction called F↑↓ := (−i)hT c−k↑(t)ck↓(t0)i.
This Green function together with F↓↑† := (−i)hT c†k↓(t)c†k↑(t0)i, which describes the corresponding
mechanism for holes are called anomalous Green functions, whereas G↑↑and G†↓↓are called normal
Green functions [2, S.18].
The equations of motion for the anomalous Green functions depend again on the normal Green functions. They are:
(i∂t− k)F↑↓(t − t0) = −∆k(−i)hT c†k↓(t)ck↓(t0)i (3.8)
(i∂t+ k)F↓↑† (t − t0) = ∆∗k(−i)hT ck↑(t)c†k↑(t0)i (3.9)
We have reached a point where we deal with four different types of Green functions (Gσσ0, G†σσ0,
Fσσ0 and F†
σσ0) and further each of them is a 2×2 matrix in spin space.In total this gives 16
different Green functions. But as one writes down all the equations it can be easily seen that half of them are zero and that due to spin degeneracy the rest of them splits into two equal systems of equations. This is why Nambu space supplies a practical tool here. Instead of treating all 16 Green functions, in Nambu space one only looks at four Green functions, which are exactly the Green functions that couple to each other when deriving the equations of motion [?].
In Nambu space the Green function looks like this:
Gkkσσ0 =
Gkk↑↑ Fkk↑↓
Fkk↓↑† G†kk↓↓
(3.10)
For simplicity we drop the momentum- and spin indices, so Gkk↑↑≡ G.
G = GF† GF† := (−i) hTck ↑c † k ↑i hTc−k ↑c † k ↓i hTc−k ↓ck ↑† i hTck ↓c†k ↓i ! (3.11)
After a Fourier transform according to (2.13) of all equations of motion found so far (3.6-3.9) what is left is a linear system of algebraic equations. It can be written in matrixform like this:
ω − k −∆k
∆∗
k −(ω + k)
G = σz (3.12)
The solution for G can then be easily obtained by inverting a 2×2 matrix. The solution is called
the bare Green function gk(ω) in distinction from the dressed Green function which will come up
When calculating the bare Green function it is found that: gk(ω) = 1 ω2− (2 k+ |∆|2) ω + k ∆k ∆∗k ω − k (3.14)
The bare Green function can also be split into two terms using u2
k = (1 + k/Ek)/2 and vk2 =
(1 − k/Ek)/2. Also the electron energy E is introduced by Ek2= 2k+ |∆|2.
gk(ω) = 1 ω − Ek u2k ∆k/2Ek ∆∗k/2Ek v2k + 1 ω + Ek vk2 −∆k/2Ek −∆∗ k/2Ek u2k (3.15)
3.2
Andreev States in Josephson Junction
In the previous chapter the electronic states in a superconductor were calculated starting from the BCS Hamiltonian. Now we are interested in a system of two superconductors divided by a insulating layer, which is called a Josephson junction.
We want to find how the electronic states of the first superconductor are changed by the presence of a second superconductor when the insulating layer is thin enough so that tunneling between the two superconductors is still possible. The insulating layer carries a certain tunneling probability for electrons, which is assumed to be independent of the electronic momentum.
The Hamiltonian (3.16) for this system is built from two BCS Hamiltonians, one for the right
(HR) and one for the left (HL) superconducting lead as well as an additional tunneling term
(HT). To distinguish between left and right superconductor electrons one denotes the momentum
of electrons in the left superconductor with p and momentum of electrons in the right supercon-ductor with q. The probability of an electron with momentum p tunneling from left to right
through the insulating layer ending up with momentum q is τpq. The probability for the reverse
process is τpq∗. These two probabilities are of course only different from each other if there is a voltage applied to the system. Here we do not deal with that case and therefor assume for any two momenta p and q that τpq= τqp.
HJ J = HL+ HR+ HT (3.16) HL= X p,σ pc†pσcpσ+ X p (∆Lc†p↑c†−p↓+ ∆∗Lc−p↓cp↑) (3.17) HR= X q,σ qc†qσcqσ+ X q (∆Rc†q↑c † −q↓+ ∆∗Rc−q↓cq↑) (3.18) HT = X p,q,σ (τpqc†pσcqσ+ τpq∗c†qσcpσ) (3.19)
Gpp= Gpp Fpp Fpp† G†pp := (−i) hTcp↑c † p↑i hTc−p↑c † p↓i hTc−p↓cp↑† i hTcp↓cp↓† i ! (3.20)
The Greenfuction for the right lead Gqq is defined exactly the same way with p and q interchanged.
For the equations of motion for Gpp one finds:
i∂tGpp(t, t0) = Iδ(t − t0) + (−i)
hT cp↑, HJ J −c † p↑i hTc−p↑, HJ J −cp↓i hTc† −p↓, HJ J −c † p↑i hTc † p↓, HJ J −cp↓i ! (3.21)
The relevant commutators were calculated in B and are given by:
ckσ(t), HJ J(t) −= kckσ+ ∆k c†−k↓δ↑σ− c†−k↑δ↓σ +X q0 τkq0cq0σ +X p0 τp∗0kcp0σ (3.22) c†kσ(t), HJ J(t) −= −kc † kσ− ∆ ∗ k c−k↓δ↑σ− c−k↑δ↓σ −X q0 τkq∗0c†q0σ −X p0 τp0kc† p0σ (3.23)
The index k can still be either a q- or a p-momentum. Note that ∆k becomes ∆L when k=p and
∆R when k=q. As well, τkp becomes zero for k being a p-momentum and τkq becomes zero for k
being a q-momentum.
With the commutators (3.22) and (3.23) inserted into (3.21) there will appear more Green
func-tions of the form Gpq and Gqp which describe the tunneling through the insulating layer. The
definition of these Green functions is again completely analogous to (3.20). As a trick one
mul-tiplies the whole equation with the third Pauli matrix σz. A Fourier transform then leads to the
equations of motion in energy space, which are given in matrix notation by:
ω − p −∆L ∆∗L −ω − p Gpp Fpp F† pp G†pp = σz+ X q0 τpq0 0 0 τpq0∗ Gq0p Fq0p Fq†0p G † q0p (3.24)
We assume a constant τ for all possible momentum transfers and introduce a new matrix: ˆτ =
diag(τ, τ∗).
Bringing the matrix on the left side to the right side, the first term is again the bare Green function as it is defined in (3.13). It describes the properties of the left superconductor not experiencing the presence of the right superconductor.
In the same manner as for Gpp one finds the equations of motion for the remaining Greenfunction
Gqpwhich then closes the equation for Gpp.
ω − q −∆R ∆∗R −ω − q Gqp= X p0 ˆ τ∗Gp0p≈ ˆτ∗Gpp (3.25)
As we are dealing with interaction free electrons (quasi particles), the Green functions appearing on the right hand side can be approximated according to Gpp0 ≡ δpp0Gpp0.
In terms of the bare Green function gq for the right superconductor equation (3.25) can also be
Gqp= gqσzτˆ∗Gpp (3.26)
Inserting (3.26) into (3.25) and introducing the self-energy ΣR(ω) =Pqˆτ gq(ω)ˆτ∗σz one finds for
Gpp the following: hω − p −∆L ∆∗L −ω − p − ΣR(ω) i Gpp = σz (3.27)
Now I want to attend to the approximate solution for the self-energy ΣR(ω). Therefore we use the
bare Green function defined as in (3.15). From looking at the matrix elements in the bare Green function, one sees that the upper left and the lower right element are the same. The lower left
and upper right elements are also the same just with ∆R and ∆∗Rinterchanged.
X q gq(ω) =X q 1 ω − Eq u2 q ∆R/2Eq ∆∗R/2Eq v2q + 1 ω + Eq v2 q −∆R/2Eq −∆∗ R/2Eq u2q (3.28)
In the first step we make the well known approximation where one replaces the summation over the q-momenta by an integral over the energies. The density of states which comes forward is assumed
to be constant, NR() ≈ NR. This is concluded from the assumption that the superconducting
leads have flat conduction bands. Because these conduction bands are also symmetric around the Fermi energy, the integration interval for can be set to (−∞, ∞) [2, S.8].
X q u2q ω − Eq + v 2 q ω + Eq ≈ NR Z ∞ −∞ u2q ω − Eq + v 2 q ω + Eq d (3.29)
Here the integrand can be simplified by inserting the expressions for u2q and vq2 and considering that the integration over an odd function gives zero. This way one obtains for the diagonal matrix elements the following:
NR Z ∞ −∞ u2q ω − Eq + v 2 q ω + Eq d =NR 2 Z ∞ −∞ 1 ω − Eq + 1 ω + Eq d (3.30)
The integration variable is now changed to E via d = (E2− ∆2
R)
−1/2
EdE. The fact that E is not defined in the interval [−∆R, ∆R] is taken care of in the additional step function.
NR 2 Z 1 ω − Eq + 1 ω + Eq d = NR 2 Z 1 ω − Eq + 1 ω + Eq E pE2− ∆2 R Θ(|E| − |∆R|)dE (3.31)
In the next step, we make the analytical continuation and write the imaginary variable ω as ω + iδ and use the Dirac identity (3.32).
Z dx F (x) x ± iδ = P Z dxF (x) x ∓ iπ Z dxF (x)δ(x) (3.32)
(3.31) ≈ −iπNR
|ω|
pω2− |∆
R|2
Θ(|ω| − |∆R|) := −i|ω|ΓR (3.33)
For the off-diagonal elements in gq the steps are similar and one obtains:
X q ∆R 2Eq(ω − Eq) − ∆R 2Eq(ω + Eq) ≈ −iπNR ∆R pω2− |∆ R|2 Θ(|ω| − |∆R|) := −i∆RΓR (3.34)
Altogether, the expression one finds for the self-energy:
ΣR(ω) = −iΓR |τ |2|ω| τ2∆ R (τ∗)2∆∗ R |τ |2|ω| σz (3.35)
The solution for Gpp is now called the dressed Green function:
Gpp= hω − p −∆L ∆∗ L −ω − p − ΣR(ω) i−1 σz (3.36)
In contrast to the bare Green function the dressed Green function now contains the influence of the right superconductor through the self-energy ΣR. The full expression for Gpp gives:
Gpp= 1 (ω + iΓR|τ |2|ω|)2− 2 p− |∆L|2− iΓR((τ∗)2∆∗R∆L+ τ2∆R∆∗L) + (ΓR|τ |2|∆R|)2 × ω + p+ iΓR|τ |2|ω| ∆L+ iΓRτ2∆R ∆∗L+ iΓR(τ∗)2∆∗R ω − p+ iΓR|τ |2|ω| (3.37)
When analyzing the poles of the Green function, one obtains information about the positions of the energy levels which can be occupied by electrons. The analysis and physical interpretation is done in the last chapter 4.
3.3
Andreev States in Josephson Junction with Quantum Dot
We introduce new creation and annihilation operators specifically for electrons populating the
quantum dot energy level, dσand d†σ. These operators do not carry any momentum index because
the momentum only depends on the energy of the quantum dot, which is fixed. The full Hamiltonian looks like:
H = HR+ HL+ HQD+ HT (3.38)
with HL and HR as before. Also, we stick to all assumptions made in the previous chapter
concerning the superconducting leads. The quantum dot and the tunneling Hamiltonian are given by: HQD = X σ 0d†σdσ (3.39) HT = X p,σ τLc†pσdσ+ X q,σ τRc†qσdσ+ X p,σ τL∗d†σc−pσ+ X q,σ τR∗d†σc−qσ (3.40)
Defining the Green function of the quantum dot,
GQ:= (−i)
hTd↑d↑†i hTd↑d↓i
hTd↓†d↑†i hTd↓†d↓i
!
, (3.41)
the equation of motion for GQ becomes:
i∂tGQ(t, t0) = Iδ(t − t0) + (−i)
hT d↑, H−d†↑i hTd↑, H−d↓i hTd† ↓, H −d † ↑i hTd † ↓, H −d↓i ! (3.42)
The commutators which are relevant now have been derived in Appendix C and are given by:
dσ, H − = 0dσ+ X p τL∗c−pσ+ X q τR∗c−qσ (3.43) d† σ, H − = −0d † σ− X p τLc†pσ− X q τRc†qσ (3.44)
When multiplying σz to both sides of (3.42), inserting the commutators and performing a Fourier
transform one finds the following equation:
ω − 0 0 0 −ω − 0 GQ= σz+ X p ˆ τL∗GpQ+ X q ˆ τR∗GqQ (3.45)
Here, new Green functions GpQand GpQappear which are the Green functions that describe the
tunneling from the quantum dot either to the left or the right superconductor. Since there is a
summation over all momenta p and q, the indices in the operators like c−pσ can be changed to
The commutators which are needed to express the equations of motion for the new Green functions
GpQhave been derived in Appendix C and are given by:
cpσ, H − = pcpσ− ∆L c † −p↑δσ↓− c†−p↓δσ↑ + τLdσ (3.46) c† pσ, H − = −pc † pσ+ ∆∗L c−p↑δσ↓− c−p↓δσ↑ − τL∗d†σ (3.47)
Using the commutators (3.46) and (3.47) leads to the equation for GpQ.
GpQ= gpσzτˆLGQ (3.48)
The bare Green function gpis already known from the problem with only one superconductor
(Chapter 3.1) and defined as in equation (3.13).
The equation for GqQ is exactly the same with p and q as well as L and R interchanged. With
this, one can close equation (3.45) and finds this expression for GQ, which is:
hω − 0 0 0 −ω − 0 − ΣL− ΣR i GQ= σz (3.49)
Here, ΣL and ΣR are the self-energies corresponding to the left and right superconductor, which
have introduced in the previous section. They are defined like:
ΣL:= X p ˆ τL∗gpτˆLσz (3.50) ΣR:= X q ˆ τR∗gqτˆRσz (3.51)
The full expression for the self energy was derived in the previous chapter and is given by:
ΣL(ω) = −iΓL |τ |2|ω| τ2∆ L (τ∗)2∆∗ L |τ |2|ω| σz (3.52)
For simplicity we only look at the case where ∆L = ∆R, τL = τR and NL= NR. With this the
Green function for the quantum dot becomes:
GQ= 1 ω2− ( 0− i2Γ|τ |2|ω|)2− 4Γ2|τ |4|∆|2 ω + 0− i2Γ|τ |2|ω| −i2Γτ2∆ −i2Γ(τ∗)2∆∗ −ω + 0− i2Γ|τ |2|ω| (3.53)
4
Analyzing electronic density of states starting from the
Green function
function the appearing expectation value is the one for the occupation number hni() = hc†ci(). In this chapter the spectral densities for the introduced systems shall be derived and also the occupation number in the BCS superconductors in the non-interacting case and when interacting with another superconductor. Finally, the spectral density for the quantum dot in between the superconducting leads will be derived.
4.1
BCS superconductor
To get the spectral density one has to analytically continue the Green function (3.37) and then using formula (2.6) and (3.32). For the electronic spectral density in an isolated BCS supercon-ductor one gets:
S(ω) = u2kδ(ω − Ek) + vk2δ(ω + Ek) (4.1)
From the spectral density one sees, that the electron energy is now restricted to the energies Ek = p2k+ |∆|2. Consequently, there are no electron states defined in the interval [−|∆|, |∆|]
which is the reason why |∆| is called the gap parameter. The spectral density is depicted in figure 1.
Figure 1: Spectral density of BCS system as a contour plot (blue=0, white=infinity) Using the spectral theorem (2.5) one can then calculate the occupation number for electrons in the superconductor. hnkσi = Z ∞ −∞ dω S(ω) eβω+ 1 = u2 k eβEk+ 1 + v2 k e−βEk+ 1 (4.2)
the electronic quasi particles is shown in figure 2. If one does not want look at the occupation of electronic quasi particles but rather on the occupation of electronic and hole quasi particles,
one has to replace the energy p by |p| in the spectral density. This way the electronic Green
function G↑↑ then gives for negative energies p the hole Green function G↓↓. An excitation of an
electron with negative energy is now viewed as an excitation of a hole with positive energy. The holes are regarded as the anti-particles of electrons. The occupation number with the described substitution can be seen in figure 3. All plots and calculations were done with Mathematica and can be found in ’spectral density BCS.nb’.
Figure 2: electronic occupation of the energy states in a BCS superconductor
Figure 3: occupation of the energy states by electrons and holes in a BCS superconductor Alternatively to calculating the spectral density, one can also just look at the poles of the Green function (3.14). Putting the denominator equal to zero gives the secular equation, from which the excitation spectrum ω can be found. The secular equation directly gives:
ω = q
2
k+ |∆|2 (4.3)
4.2
Josephson Junction
h ω2− 2 p− |∆L|2− (ΓR)2|τ |4(|ω|2− |∆R|2) i + ih2ωΓR|τ |2|ω| − ΓR(τ∗)2∆∗ R∆L+ τ2∆R∆∗L i = 0 (4.4) Recalling the definition of ΓR (3.34), there are two cases which have to be analyzed. First if
|ω| < |∆R|, ΓR will be zero and therefore the secular equation gives exactly the same results as if
the left superconductor was isolated. This means that as long as there is the right superconduc-tor’s gap opposite to the left superconductor, no electrons are available to tunnel from right to left. Therefor the left side doesn’t notice the presence of the right one.
In the second case when |ω| > |∆R|, ΓR= πNR/pω2− |∆R|2. For the secular equation one then
obtains: h ω2− Ep2− π 2N2 R|τ | 4i+ ih2πNR(|τ |2ω|ω| − ((τ∗)2∆∗R∆L+ τ2∆R∆∗L)) p|ω|2− |∆ R|2 i = 0 (4.5)
The zeros of the real part in this equation give the positions of possible excitations in the energy spectrum, whereas the imaginary part gives the width of the peaks. The zeros for the real part of equation (4.5) are:
ω1/2 = ± q 2 p+ |∆R|2+ π2NR2|τ |4= ± q E2 p+ π2NR2|τ |4 (4.6)
By taking the imaginary part of the causal electronic Green function (3.37) and plugging it into (2.7) one gets the spectral density for the left superconductor. Because of the
defini-tion of ΓR, the spectral density is defined by different analytical functions in the two regions
n
−∞, −|∆R| ∧ |∆R|, ∞
o
and −|∆R|, |∆R|. In the gap region of the right superconductor,
the left superonductor’s spectral density looks exactly like the one for the isolated BCS supercon-ductor. For the region where |ω| > |∆| one can find the analytical form of the spectral density by taking the imaginary part of the electronic Green function and using formula (2.7). For the spectral density one finally gets:
S(ω, p) = −1 π eβω+ 1 eβω− 1Im[G↑↑] if |ω| > |∆R| 0.5(u2 pδ(ω − Ep) + vp2δ(ω + Ep)) if |ω| < |∆R| (4.7)
Here, the imaginary part of the Green function is given by:
Im[G↑↑] = 2Γ Rτ2(ω + p)(ω|ω| − ∆R∆L) + ΓR|τ |2|ω|(ω2− Ep2− NR2π 2|τ |4) (ω2− E2 p− NR2π2|τ |4)2+ 4N 2 Rπ2|τ |4(|ω|ω − ∆R∆L)2/(|ω|2− |∆R|2) (4.8)
Like for the spectral density of the BCS superconductor, we will make the substitution p−→
To get the occupation number one has to compute the integral: n(p) = Z ∞ −∞ dωS(ω, p) eβω+ 1 (4.9)
Due to the definition of spectral density (4.7) this integral splits into three parts:
n(Ep) = − 1 π Z −|∆R| −∞ dωIm[G(ω, p)] eβω− 1 + Z |∆R| −|∆R| dωu 2 pδ(ω − Ep) + vp2δ(ω + Ep) eβω+ 1 − 1 π Z ∞ |∆R| dωIm[G(ω, p)] eβω− 1 (4.10)
The first and third integral cannot be solved analytically. In order to get the numerical solution Mathematica was used. In the files ’spectral density JJ smaller.nb’, ’spectral density JJ bigger.nb’ and ’spectral density JJ equal.nb’ the analytical form of the spectral density was used to numer-ically integrate the integral in 4.10. Although the numerical integration sometimes failed to con-verge in the regions where the spectral density gets very large or has high oscillations, the last approach to the integral was taken.
The second integral in this expression is non zero only if the delta-functions are located inside the intervalh−|∆R|, |∆R|
i
. This is only fulfilled for certain values of p. More precisely it is the case
for all 2p < |∆R|2− |∆L|2. This can of course be only the case if |∆R| > |∆L|. This case as well as
the other two cases where |∆L| = |∆R| and the case where |∆L| < |∆R| shall be treated in more
detail now.
• |∆L| = |∆R|
The spectral density of the BCS superconductor (figure 1) was a perfect Delta-peak at the
energies ω = p2
k+ ∆2. If the electrons and holes can tunnel into the right
superconduc-tor the Delta-peaks experience a broadening due to the finite lifetime of the electron and hole energy states. Also there is a shift of the maximum which can be derived by looking at
the secular equation (4.5). The maxima now lie at the energies ω = ±p2
k+ ∆ 2
Figure 4: Spectral density as a function of p and ω
The electronic density then follows by numerical integration of (4.10). Also the transforma-tion from free electron energies p to the BCS electron energies Ep=
q 2
p+ ∆2L was done.
Figure 5: Occupation number as a function of Ep
• |∆L| < |∆R|
In the interval [−|∆R|, |∆R|] the left superconductor does not get influenced by the right
superconductor. Here the spectral density looks exactly like the one for the isolated
super-conductor (figure 1). Only when ω becomes greater than ∆R the electrons and holes can
Figure 6: Spectral density as a function of p and ω
This can also be seen in the plot for the occupation of energy states by electrons (Ep > 0)
and holes (Ep < 0) in (figure 6). Here it can be seen that as soon as the excitation energy
of an electron or a hole becomes greater than the right superconductor’s gap parameter, the electrons or holes are now able to tunnel into the opposite superconductor. Therefor the
step in the occupation occurs at the energy Ep= ∆R= 2.
Figure 7: Occupation number as a function of Ep
The effect of the left superconductor on the occupation number of the right one can be seen from the next case.
• |∆L| > |∆R|
From the spectral density (figure 7) one can see that the gap approaches the right super-conductor’s gap parameter for high energies. It is exactly the other way round as for the
previous case (figure 6), where the Andreev states go from the gap edge ∆L at low energies
Figure 8: Spectral density as a function of p and ω
Compared to the occupation number plot of an isolated BCS system (figure 3) the electron and hole density in the Josephson junction is lower. This is because the electrons and holes can leave their states through tunneling.
Figure 9: Occupation number as a function of Ep
4.3
Quantum Dot in a Josephson Junction
By taking the same way as before, we calculate the spectral density from the Green function (3.53) of the quantum dot located in a Josephson junction. For simplicity we will now assume that left
and right superconductor are equal. Therefore, NR= NL= N , ∆R= ∆L = ∆ and ΓR= ΓL= Γ.
For the different energy regimes one finds:
Plotting the spectral density in Mathematica using different values for the parameters N, τ and 0 gives the following figures.
Figure 10
Figure 11
Figure 13
What can be seen is that compared to the isolated case, where the quantum dot’s spectral density is a delta peak at the quantum dot energy, this single energy experiences a broadening due to tunneling into the superconductors which results in a finite lifetime of the state. The lifetime or width of the state is dependent on the tunneling probability and on the density of states at the Fermi level N. The more electrons there are in the superconductors to be able to tunnel into the quantum dot and the lower the tunneling barrier is, the more likely an actual tunneling event is and the lifetime of the states in the quantum dot decreases.
Because in both superconductors the gap region cannot be occupied by electrons also tunneling from these energies into the corresponding energies in the quantum dot is not possible. This is why the superconducting gap causes the cut in the quantum dot’s spectral density. At the edge of the gap, where the spectral density for the superconductors themselves diverges, there are also many electrons available to tunnel. Therefor additional peaks occur at the edges of the gap as it can be seen in Figure 12.
The secular equation gives:
ω2− (0− 2iΓ|τ |2|ω|)2− 4Γ2|τ |4|∆|2= 0 (4.12)
The zeros of the real part give the positions of the peaks in the spectral density. They are:
ω1/2= ± q
2
0− 4π2N2|τ |4 (4.13)
5
Conclusion and Outlook
6
Appendix
A
Calculation of Commutators with BCS Hamiltonian
Useful identities for this calculation are the commutator-identity and the commutation relations for fermionic creation and annihilation operators.
A, BC −=B, A+C − BC, A+ (6.1) ckσ, ck0σ0 +=c † kσ, c † k0σ0 += 0 (6.2) ckσ, c†k0σ0 + = δkk 0δσσ0 (6.3)
For this calculation k and σ in HBCSare now denoted as k00and σ00in order not to get confused
with the indices in the definition of the Green function.
ckσ(t), HBCS(t) − = X k00σ00 k00ckσ, c†k00σ00ck00σ00 −+ X k00 ∆k00ckσ, c†k00↑c † −k00↓ −+ ∆ ∗ k00ckσ, c−k00↓ck00↑ − = X k00σ00 k00(c†k00σ00, ckσ +ck 00σ00− c†k00σ00ck00σ00, ckσ +)+ X k00 ∆k00(c† k00↑, ckσ +c † −k00↓− c † k00↑c † −k00↓, ckσ +) + ∆∗k00(c−k00↓, ckσ +ck 00↑− c−k00↓ck00↑, ckσ +) = X k00σ00 k00(δk00kδσσ00ck00σ00) + X k00 ∆k00(δk00kδ↑σc† −k00↓− c † k00↑δ−k00kδ↓σ) = kckσ+ ∆kc†−k↓δ↑σ− ∆−kc†−k↑δ↓σ (6.4)
Because HBCS is hermitian the commutatorc†kσ(t), HBCS(t)
− is easily derived from the first
commutator by taking the adjoint.
c† kσ(t), HBCS(t) −= −ckσ(t), HBCS(t) † − = −kc†kσ− ∆ ∗ kc−k↓δ↑σ+ ∆∗−kc−k↑δ↓σ (6.5)
B
Calculation of Commutators with JJ-Hamiltonian
Again, the indices in HJ J are all replaced double-primed indices. The commutation relations of
the creation and annihilation operator for electrons in the left superconductor with HL and for
electrons in the right superconductor with HR are of course the same as for HBCS. However,
superconductor. ckσ(t), HT(t) −= X p00q00σ00 τp00q00ckσ, cp†00σ00cq00σ00 −+ τ ∗ p00q00ckσ, c†q00σ00cp00σ00 − = X p00q00σ00 τp00q00ckσ, c† p00σ00 +cq 00σ00+ τp∗00q00ckσ, c†q00σ00 +cp 00σ00 = X p00q00σ00 τp00q00δkp00δσσ00cq00σ00+ τp∗00q00δkq00δσσ00cp00σ00 =X q00 τkq00cq00σ +X p00 τp∗00kcp00σ (6.6)
Of course any τpp0 or τqq0 is zero.
The commutatorc†kσ(t), HT(t)
− is again be obtained from the first commutator.
c† kσ(t), HT(t) − = −ckσ(t), HT(t) † − = −X q00 τkq∗00c†q00σ −X p00 τp00kc†p00σ (6.7)
C
Calculation of Commutators with QD-Hamiltonian
This calculations concern the commutators between dσ or ckσand the Hamiltonian.
dσ(t), H(t) − = X σ0 0dσ, d†σ0dσ0 −+ X pσ0 τLdσ, c†pσ0dσ0 −+ X qσ0 τRdσ, c†qσ0dσ0 − +X pσ0 τL∗dσ, d†σ0c−pσ0 −+ X qσ0 τR∗dσ, d†σ0c−qσ0 − =X σ0 0dσ, d†σ0 +dσ 0+ X pσ0 τL∗dσ, d†σ0 +c−pσ 0+ X qσ0 τR∗dσ, d†σ0 +c−qσ 0 =X σ0 0δσσ0dσ0 + X pσ0 τL∗δσσ0c−pσ0+ X qσ0 τR∗δσσ0c−qσ0 = 0dσ+ X p τL∗c−pσ+ X q τR∗c−qσ (6.8) d† σ(t), H(t) −= −0d † σ− X p τLc†−pσ− X q τRc†−qσ (6.9)
For expressing the commutator ckσ, H
− , k can again be either q or p, depending in which
get the solution for k=q, one just has to replace p by q, and L by R. cpσ, H −=cpσ, HL −+cpσ, HT − (6.10) with: cpσ, HL − = pcpσ− ∆L c † −p↑δσ↓− c†−p↓δσ↑ and: cpσ, HT − = X p0σ0 τLcpσ, c†p0σ0dσ0 −+ X p0σ0 τL∗dσ, d†σ0c−pσ0 − =X p0σ0 τLcpσ, c†p0σ0 +dσ 0 = τLdσ resulting in: cpσ, H − = pcpσ− ∆L c † −p↑δσ↓− c†−p↓δσ↑ + τLdσ (6.11) c† pσ, H − = −pc † pσ+ ∆∗L c−p↑δσ↓− c−p↓δσ↑ − τL∗d†σ (6.12)
References
[1] Cooper Bardeen and Schrieffer. Microscopic theory of superconductivity. Physical Review 106, 1957.
[2] Machiel Geert Flokstra. Proximity effects in superconducting spin-valve structures. PhD thesis, Universiteit Leiden, 2010.
