Josephson.html Josephson.mws

Josephson Junction

We investigate a simple Hamiltonian that is used to describe the conduction properties of a Josephson junction. The discussion follows D.V. Averin's article in Scalable Quantum Computers  , S.L Braunstein and H.-K. Lo (eds), Wiley-VCH (Berlin 2001), ISBN 3-527-40321-3.

A Josephson junction is a remarkable device: it is of a mesoscopic scale, i.e., involves many particles, yet it displays quantum behaviour. This is a feature common to Bose condensed states (recently Bose-Einstein condensation has become possible for ensembles of atoms which are super-cooled in magneto-optical traps). This is remarkable, because normally ensembles of quantum systems behave classically on macroscopic scales. Here we can observe (almost) the entire ensemble being in a pure quantum state.

The junction Hamiltonian is in a number representation with respect to Cooper pairs and has tridiagonal form (nearest-neighbour coupling). It is controlled by three parameters:

q = C*V_e/(2*e)   is the charge induced on the junction capacitance C  by the potential V_e ;

E_J  is the energy of Cooper pair tunneling = Josephson coupling energy;

E_C = (2e)^2/(2*C)  is the charging energy of a single Cooper pair.

H = E_C*(n-q)^2 - 0.5*E_J*(|n><n+1| + |n+1><n|)   is the Hamiltonian; in the ket-bra product (an operator in the occupation number representation) summation over the n -index is implied.

The first part of the Hamiltonian is called the 'number' part: it represents the Coulomb charging of the junction capacitance by the Cooper pairs (hence the 2e  for the charge). The second part represents the 'phase' term in the energy associated with the tunneling of a pair. The quantum dynamics of the junction depends on the ratio of the two energy scales: for E_J >> E_C  the junction can behave almost clasically. The capacitance C  of the junction has to be small so that the charging energy of a single Cooper pair is appreciable. This can be achieved in a mesoscopic system.

The origin of the number part of the Hamiltonian: it is an electrostatic energy associated with (a) the charging of the capacitance by the n  Cooper pairs, and (b) also by an external electric field which induces a potential difference between two adjacent superconducting islands, V_e . The polarization charge q=C*V_e/(2*e)  can be varied continuously through V_e . The charge is in units of (2e)  so that it can be compared with the charge of a Cooper pair.

The number of Cooper pairs n  is quantized, and it is this characteristic difference which results in interesting dynamics. The form of the coupling term derives from the fact that the tunneling of a Cooper pair from one superconducting island to the other results in the change of n  by +/- 1 .

One wants to find the ground state energy in units of E_J  as a function of q  for a fixed ratio E_C/E_J . Also of interest is the average number of Cooper pairs under this condition.

We pick the unit:

>    restart; E_J:=1;

E_J := 1

>    E_C:=10*E_J;

E_C := 10

From Fig.1 in the article we read off that we should get for the ground state approximately E0 = 2  when q=1/2, and slightly less than zero when q=0 .

We truncate the number representation at some finite value N  and will carry out a convergence analysis.

>    with(LinearAlgebra):

>    q:='q': N:=15;

N := 15

>    HM:=Matrix(N+1):

We make room for the n=0  Cooper pairs state.

>    for i from 1 to N+1 do: n:=i-1: HM[i,i]:=E_C*(n-q)^2: od:

What do we get when we sandwich the Hamiltonian between number states <i|  and |j>  ? Assume that the states are normalized.

<i|n><n+1|j> + <i|n+1><n|j>  is the off-diagonal coupling matrix element. The number states are orthogonal and normalized, so we get:

delta(i,n)*delta(n+1,j) + delta(i,n+1)*delta(n,j)   where delta  is the Kronecker delta.

>    KD:=(mu::nonnegint,nu::nonnegint)->if mu=nu then 1 else 0 fi;

KD := proc (mu::nonnegint, nu::nonnegint) options operator, arrow; if mu = nu then 1 else 0 end if end proc

>    KD(2,2),KD(3,0),KD(0,0);

1, 0, 1

>    KD(1.2,2.1);

Error, invalid input: KD expects its 1st argument, mu, to be of type nonnegint, but received 1.2

>    KD(-2,0),KD(-2,-2);

Error, invalid input: KD expects its 1st argument, mu, to be of type nonnegint, but received -2

We carry out a brute-force calculation. The pair number parameter n  is being summed over in the ket-bra part of the Hamiltonian.

>    for i from 1 to N+1 do: for j from 1 to N+1 do: HM[i,j]:=HM[i,j]-1/2*E_J*add(KD(i,n)*KD(n+1,j)+KD(i,n+1)*KD(n,j),n=0..N); od: od:

>    evalm(SubMatrix(HM,1..5,1..5));

matrix([[10*q^2, -1/2, 0, 0, 0], [-1/2, 10*(1-q)^2, -1/2, 0, 0], [0, -1/2, 10*(2-q)^2, -1/2, 0], [0, 0, -1/2, 10*(3-q)^2, -1/2], [0, 0, 0, -1/2, 10*(4-q)^2]])

The truncated Hamiltonian matrix has a tridiagonal structure. We compute the eigenvalues.

>    EVs:=Eigenvalues(subs(q=0.5,HM)):

>    print(seq(Re(EVs[i]),i=1..N+1));

1.99386586289294598, 2.99363153646673118, 62.5020833901770132, 22.5062525355293630, 202.500625001306332, 122.501041672996450, 422.500297619183300, 302.500416667044931, 562.500223214342668, 722.50017361...
1.99386586289294598, 2.99363153646673118, 62.5020833901770132, 22.5062525355293630, 202.500625001306332, 122.501041672996450, 422.500297619183300, 302.500416667044931, 562.500223214342668, 722.50017361...
1.99386586289294598, 2.99363153646673118, 62.5020833901770132, 22.5062525355293630, 202.500625001306332, 122.501041672996450, 422.500297619183300, 302.500416667044931, 562.500223214342668, 722.50017361...
1.99386586289294598, 2.99363153646673118, 62.5020833901770132, 22.5062525355293630, 202.500625001306332, 122.501041672996450, 422.500297619183300, 302.500416667044931, 562.500223214342668, 722.50017361...

>   

The eigenvalues come out sorted in value. The fround state is close to what can be read off from Fig. 1a in the reference.

Exercise 1:

Observe the convergence behaviour of the two lowest eigenenergies with truncation size N .

Two issues emerge:

1) A question: can this matrix be diagonalized exactly?

2) The calculation of a tridiagonal matrix should be done more efficiently by an operator method which utilizes the sparseness of the matrix.

So let us construct a graph of the behaviour of the fround and first excited states with the charge value q . We place the eigenvalue calculation into a loop for a sequence of q -values, and sort the list of eigenvalues in order to extract the first two as the lowest states.

>    ET0:=[]: ET1:=[]: for i_q from 1 to 80 do: q0:=evalf(0+(i_q-1/2)/40); EVs:=convert(map(Re,Eigenvalues(subs(q=q0,HM))),list): EVs:=sort(EVs);

>    E0:=evalf(EVs[1]):  ET0:=[op(ET0),[q0,E0]]: E1:=evalf(EVs[2]):  ET1:=[op(ET1),[q0,E1]]: od:

>    PL0:=plot(ET0,color=red): PL1:=plot(ET1,color=blue): plots[display](PL0,PL1,title="energy gap vs charge q",labels=[q,E]);

[Maple Plot]

We observe periodic behaviour in the energy values as a function of q .

How do we get the average Cooper pair number in the ground state? Averin refers to the junction free energy, which can be studied as a function of temperature. We do the simple quantum mechanics calculation at zero temperature.

We use information from the eigenstates given that we should know the composition of the ground state in terms of the number states.

The average Cooper pair number can be calculated by using the eigenvectors from the matrix diagonalization. Their interpretation is that they provide us with the expansion coefficients for the eigenstates in terms of the number state basis. The magnitude squared of the coefficient gives the occupation probability, and therefore we can add these probabilities weighted with the particle number to calculate the average Cooper pair number in the energy eigenstate as a function of q .

On a Maple technical note: the sorting is more involved in this case, as we need to re-sort not just a list of eigenvalues, but the entire construct associated with the eigenvectors. For a first look at the ground-state behaviour it appears to be sufficient to just track the first eigenvalue/eigenvector result; this is misleading, however, for certain values of q  of does not get away without sorting!

>    AN:=[]: for i_q from 1 to 90 do: q0:=evalf((i_q-1/2)/80); EVs:=Eigenvectors(subs(q=q0,HM)):

>    j0:=1: E0:=Re(EVs[1][1]): for j from 2 to N do: if Re(EVs[1][j])<E0 then j0:=j: E0:=Re(EVs[1][j]): fi: od: AN:=[op(AN),[q0,add(Re(EVs[2][j,j0])^2*(j-1),j=1..N)]]; od:

>    plot(AN,labels=[q,"<n>"]);

[Maple Plot]

The step-like behaviour is a manifestation of the quantization of Cooper pairs in the mesoscopic junction. For smaller E_C/E_J values the step behaviour is somewhat smoothed out. This reflects the physical situation that for large values of the tunneling energy of a pair E_J  the junction electrodes are so strongly coupled that the junction behaves classically and the average pair number <n>  just tracks the induced polarization charge, i.e., is proportional to q . Another way how this can be phrased is: large quantum fluctuations of n  induced by the tunneling term in the Hamiltonian drown out the quantization effects observed in the step like behaviour in this case.

Exercise 2:

Repeat the calculation for different values of E_C/E_J and observe the approach towards the continuous result <n>=q  as the ratio becomes small.

Now we repeat the findings for the eigenenergy vs q , and for the average pair number vs q  for the first excited state and compare the finding with the ground-state result. We also extend the q -range.

>    AN0:=[]: AN1:=[]: for i_q from 1 to 240 do: q0:=evalf((i_q-1/2)/80); EVs:=Eigenvectors(subs(q=q0,HM)): j0:=1: E0:=Re(EVs[1][1]): for j from 2 to N do: if Re(EVs[1][j])<E0 then j0:=j: E0:=Re(EVs[1][j]): fi: od: j1:=N: E1:=Re(EVs[1][N]): for j from 1 to N-1 do: if j<>j0 and Re(EVs[1][j])<E1 then j1:=j: E1:=Re(EVs[1][j]): fi: od:

>    AN0:=[op(AN0),[q0,add(Re(EVs[2][j,j0])^2*(j-1),j=1..N)]]; AN1:=[op(AN1),[q0,add(Re(EVs[2][j,j1])^2*(j-1),j=1..N)]]; od:

>    PL3:=plot(AN0,color=red): PL4:=plot(AN1,color=blue): plots[display]({PL3,PL4},labels=[q,"<n>"],title="average Cooper pair number vs polarization charge q");

[Maple Plot]

This result is cute. It shows that as one moves across the point q=0.5*(2e)  for the induced polarization charge the ground and first excited states trade status as far as Cooper pair number is concerned.

Due to the mesoscopic quantum behaviour of JJs they are considered a candidate for the implementation of quantum computing.

>