Electrical Detection of Topological Quantum Phase Transitions in Disordered Majorana Nanowires
Electrical Detection of Topological Quantum Phase Transitions in Disordered Majorana Nanowires
Electrical Detection of Topological Quantum Phase Transitions in Disordered Majorana Nanowires
nanowires
Benjamin M. Fregoso,1, 2 Alejandro M. Lobos,2 and S. Das Sarma2
an exponential splitting eL/ , where is the superconducting coherence length. Increasing the strength of
disorder induces a proliferation of low-energy Andreev
bound states (i.e., quantum Griffiths effect), and the
splitting scales as eL/+L/(2`e ) , where `e is the elastic mean-free path of the system.26 Beyond a critical
disorder strength (defined by the condition `e = /2),
the system enters a nontopological insulating phase with
no end-MBSs. At both sides of the TQPT, the system
is localized at zero energy, and exactly at the critical
point separating these phases, the wave functions become delocalized and the smallest Lyapunov exponent
(i.e., the inverse of the localization length of the system)
vanishes. This key observation links the physics of localization and topological properties of a disordered D-class
SC wire2832 . For example, it has been shown that the
Q2M
quantity Q = sign (Det r) = sign( n=1 tanh n ), where
r is the reflection matrix and n is the Lyapunov exponent related to the n-th transmission channel, is a suitable topological invariant for a disordered class D superconductor, which changes sign whenever a Lyapunov exponent crosses zero28,33 . It was suggested that the delocalized nature of the wave functions, with their non-local
correlations that appear at the topological critical point,
could be observed in the quantization of the thermal con2
ductance Gth /G0 = 1 (with G0 = 2 kB
T /6h at temperature T ) or the onset of quantized non-local current-noise
correlations, constituting evidence for the TQPT. Unfortunately, the highly challenging requirements of these
experiments have hindered further progress.
In this work we propose a different yet simple, electrical transport experiment to detect the TQPT, measuring
directly the non-local correlations in the bulk appearing exactly at the critical point. We study a disordered
topological SC coupled to left and right normal leads in
2
a normal-superconducting-normal (NSN) device as depicted in Fig. 1. Instead of computing the left-right conductance GLR , our proposal consists in calculating the
local conductance at one end of the NW, while tuning
the coupling to the opposite lead. As we show below, this
procedure allows to extract information about the nonlocal correlations, which in turn could be used to identify the TQPT in the bulk of the wire. We stress that
this is different from measuring GLR , which vanishes,
or from measuring the non-local correlations in the shot
noise34 . Our method can be immediately implemented in
on-going experiments looking for zero bias conductance
peak in the Majorana nanowires.
Theoretical model. We firstmotivate our results by
studying a 1D solvable model of a disordered D-class SC
consisting of spinless Dirac fermions with a random pwave gap (x). In the Majorana basis the Hamiltonian,
from x = 0 to x = L, is26,28 Hw = i~vF z x + y (x),
where ( x , y , z ) is the vector of Pauli matrices acting on the space of right- and left-moving Majorana
fields. At zero energy, this Hamiltonian has localized
Majorana modes
x at the ends of the wire, e.g., (x) =
exp[(1/~vF ) 0 dx0 (x0 ) x ] (0) is a localized mode
at the left end (i.e., x = 0). The reflection r1 =
1
tanh(L/~v
(L/~v
F ) and transmission t1 = cosh
F)
amplitudes are obtained by imposing (0) = (1, r1 )T ,
cosh2 (L/~v
F ), of width equal to the Thouless energy
of the system, i.e., ~vF /L. This result is a consequence
of the particular reflection-less boundary condition imposed at the right end of the wire, x = L. However,
one can assume a more general situation introducing a
barrier at the end of the wire, described by a generic
scatterer with reflection and transmission amplitudes r2
and t2 , respectively (subject to the unitarity constraint
2
2
|r2 | + |t2 | = 1). This would correspond to an imperfect
coupling to the right lead or to any backscatterer which is
external to the wire itself. The total reflection amplitude
becomes
r2
r = r1 + t1
t1 .
(1)
1 + r2 r1
Intuitively, the last term in Eq. (1) represents processes
in which the right-moving Majorana mode is transmitted to the right end of the wire with amplitude t1 , and
is reflected back with amplitude r2 as a left-mover. This
result means that the quantity r in general contains nonlocal contributions from the scattering occurring at the
right end35 . Note that for a non-vanishing r2 the point
r = 0 is shifted with respect to the topological transition
= 0. This shift is of the order
in the bulk of the wire
of the intrinsic width ~vF /L of the TQPT, and hence
IL
(a)
tR
GL L= dIL /dVL
VL
IS
IL
(b)
N
tR + t R
N
IS
VL
FIG. 1. (Color online) Schematic representation of the NS-N system under consideration. The basic idea consists of
measuring the conductance at one end of the NW (e.g. left
end) for different values of R , the coupling to the opposite
(i.e., right) lead. The difference of conductances GLL [see
Eq. (6)] contains information about the non-local correlations
between the end-points of the NW, which can be used to
detect the TQPT.
does not affect its experimental detection. Let us now assume that r2 is another external tunable parameter in the
system, in which case a small variation r2 around r2 = 0
allows to extract the non-local contributions in Eq. (1),
X
hlmi,s
X
cl,s cm,s
x
cl,s (l VZ ss
0 ) cl,s0
l,s
y
i cl,s ss
, (2)
0 cl+1,s0 + 0 cl cl + H.c.
l,s
with effective hopping parameter t and lattice parameter a. Here cl,s creates an electron with spin projection
s = {, } at site l in the tight-binding chain, is the
Rashba SOC parameter, VZ is the Zeeman energy due to
an external magnetic field along x, and 0 is the induced
s-wave gap which must be calculated self-consistently.
Since the precise numerical value of 0 does not modify
our conclusions, here we make the simplifying assump-
3
tion that 0 already satisfies the self-consistent SC gap
equation.
Short-ranged nonmagnetic static disorder in the semiconductor NW is included through a fluctuating chemical potential l = 0 + l about the average 0 . For
simplicity we assume l to be a delta-correlated random variable with Gaussian distribution, hl m i =
v02 lm . Hamiltonian in Eq. 2 is another particular example of a disordered Dclass SC23 . As a function of
the external Zeeman field VZ , and in absence of disorder this model has a TQPT from a topologically trivial
phase topnon-trivial phase with end MBS at the value
VZ,c = 20 + 20 as shown originally by Sau et al.37 .
In the presence of disorder, the critical field VZ,c shifts
to higher values, and its precise value depends on the
particular details of the disorder26,27,32 .
We describe the coupling to the external leads
P
tR fRk,s
cN,s )
e2
r
r
L R g1s,N
L R f1s,N
s0 () and [teh ()]s,s0
s0 (),
r
r
where gls,ms0 () and fls,ms0 () are the normal and
anomalous retarded Greens functions in the chain (see
Appendix A).
The topological phase occurring for VZ > VZ,c is
characterized by a quantized zero-bias peak at GLL =
2e2 /h, which is a direct consequence of an MBS localized at the left end of the NW4043 . However, the
proliferation of disorder-induced subgap Andreev bound
states near zero energy results in a power-law singular
ity h1, ()idis 1/ || in the disorder-averaged density of states, and complicates the interpretation of this
zero-bias peak15,18,2427 . On the other hand, the experimental detection of the predicted delocalization TQPT
GLL =
is hindered by the fact that the non-local electrical conductance GLR (0) vanishes, since Tr [tee tee ] = Tr [teh teh ]
at the transition28 . In addition, the predicted quantized
thermal conductance Gth /G0 = 2Tr [tee tee + teh teh ] is
experimentally very difficult to observe. This calls for
alternative methods to detect the TQPT.
In analogy to Eq. (1), Eq. (3), despite being a local quantity computed at the left lead, contains information about the non-local correlations in the NW. To
see this, we make use of the Greens function identity44
G r () G a () = G r () [r () a ()] G a (), where
G r(a) () = [ HBdG
r(a) ]1 is the Greens function
w
matrix, defined in terms of the 2 2 Nambu blocks
r
g () f r ()
,
G r(a) () =
fr () gr ()
HBdG
is the BdG Hamiltonian corresponding to Eq. (2),
w
and r(a) () = (i/2) (L l,1 + R l,N ) s,s0 is the retarded (advanced) self-energy due to the coupling Hmix .
This allows us to express Eq. (3) in a more suggestive
form
GLL =
e2
h {2Tr
[reh reh
] + Tr [tee tee + teh teh ]},
(5)
(6)
4
sic property of Hw ), the maximum of GLL is shifted
FIG. 2. (Color online). (a) Transmission probability T1N between the end-points (sites 1 and N ) of the isolated NW and
topological phase diagram as a function of disorder strength
v0 and Zeeman field VZ . At the topological
critical point,
Q
the topological invariant Q = sign( 4n=1 tanh n ), where n
are the Lyapunov exponents, changes sign and the transmission probability T1N becomes 1 [see also Fig. 3(a)]. The
bold dashed line follows the TQPT. (b) Difference of left
conductances GLL [see Eq. (6)]. We have chosen parameters L = 0.2t, R = 105 t, R = 103 t, N = 300,
= 0.07t, = 0.37t and 0 = 2t. In the clean case and at
Vz = 0.45t = 6.4 we obtain the maximal ratio of L/ 10
where L is the length of the NW.
is given by M =
Ml =
QN
l=1
tl +(Vz 0 )
t2 +2
l +t(Vz 0 )
t2 +2
l +t(Vz +0 )
t2 +2
tl (Vz +0 )
t2 +2
t
t2 +2
t2 +2
t2 +2
t
t2 +2
0
0
0
0
FIG. 3. (Color online). (a) Topological invariant for an isolated NW. The localization -delocalization transition is evident form the appearance of a quantized peak in the transmission probability, T1N . (b) difference in conductance in the
left lead for different couplings to the right lead in the regime
R R . The main point is that these two quantities follow each other which can be used to detect the topological
phase transition. Here we have chosen fixed disorder amplitude v0 = 3.14, see dashed vertical lines in Fig. 2.
5
Appendix A: Calculation of the conductance matrix in a SNS contact
Here we provide the details of the calculations of the conductance through a generic NSN system. Our derivation is
standard and makes use of the so-called Hamiltonian formalism45 , which is equivalent to the the more frequently-used
scattering or BTK formalism38,39 , provided the Green functions are calculated to all orders in the coupling across the
SN interface. Our model Hamiltonian in the main text is
H = Hw + Hmix + Hlead, L + Hlead, R ,
X
X
X
y
x
i cl,s ss
Hw = t
cl,s cm,s
cl,s (l VZ ss
,
0 ) cl,s0 +
0 cl+1,s0 + 0 cl cl + H.c.
hlmi,s
Hmix =
l,s
(A1)
(A2)
l,s
(tL fLk,s
c1,s + tR fRk,s
cN,s ) + H.c.,
(A3)
Hlead,j =
k fj,k,s
fj,k,s .
(A4)
k,s
We assume that each lead is in equilibrium at a chemical potential j = eVj controlled by external voltages, where
j = {L, R}, and that the SC NW is grounded, i.e., S = 0 (see Fig. 1). The expression for the electric current
calculated through the contacts is Ij = ehdNj /dti = ieh[H, Nj ]i/~ = ieh[Hmix , Nj ]i/~, which can be written in terms
of the Green function at the contacts as45,46
E D
Ei
ie X hD
tL cL, c1, c1, cL, ,
~
E D
Ei
ie X hD
IR =
tL cR, cN, cN, cR, .
~
(A5)
IL =
(A6)
With these definitions, note that the currents are positive if particles move into the leads (i.e., exit the SC), and
negative otherwise. On the other hand, charge conservation demands that IL + IR + IS = 0, where IS is the excess
current that flows to earth through the SC. Within the Baym-Kadanoff-Keldysh formalism44 we define the lesser
Green function
D
E
<
ci, cj, (t) ,
gi,j
0 (t) ie
(A7)
e X d <
<
tL
gL,1 () g1,L
() ,
~
2
X
e
d <
<
IR = tR
gR,N () gN
,R () ,
~
2
(A8)
IL =
(A9)
Using equations of motion, we can express Eqs. A8 and A9 in terms of local Greens functions as45,46
h
i
e X
0,<
0,>
>
<
IL = t2L
d gL,L
() g1,1
() gL,L
() g1,1
() ,
h
h
i
X
e
0,<
0,>
>
<
IR = t2R
d gR,R
() gN
,N () gR,R () gN ,N () .
h
(A10)
(A11)
0,
Our first step to obtain the expression of the currents is to specify the unperturbed Greens functions gj,j () in
the leads, with j = {L, R}:
0,<
gj,j
() = 2i0j, () nj () ,
0,>
gj,j
() = 2i0j, () [nj () 1] ,
6
where nj () = nF ( + j ) are the Fermi distribution functions at the leads. Substituting these expressions gives
X
r
ie
a
<
d 0L, () nL () g1,1
() g1,1
() + g1,1
() ,
2t2L
h
X
r
ie
a
<
IR = 2t2R
d 0R, () nR () gN
,N () gN ,N () + gN ,N () ,
h
IL =
(A12)
(A13)
where we have used the identity44,47 g > () g < () = g r () g a (). Obtaining an explicit expression for the
currents IL and IR is quite cumbersome. Since we will be interested only in the conductance, we note that there is
an enormous simplification if we compute directly the conductance matrix by deriving the currents with respect to
the voltages VL , VR . Then
GLL
GLR
GRL
GRR
X
dIL
ie2
= 2t2L
dVL
h
X
dIL
ie2
= 2t2L
dVR
h
X
dIR
ie2
= 2t2R
dVL
h
X
dIR
ie2
= 2t2R
dVR
h
0L,
()
d 0L, ()
<
dg1,1
()
dnL () r
a
g1,1 () g1,1
() +
d (eVL )
d (eVL )
<
dgN
,N ()
d (eVL )
0R,
(A14)
<
dg1,1
()
,
d (eVR )
d 0R, ()
()
(A15)
,
(A16)
<
dgN
dnR () r
,N ()
a
gN ,N () gN
()
+
,N
d (eVR )
d (eVR )
)
,
(A17)
Therefore, we see that the problem is reduced to finding the Greens functions in the superconducting system. In
a non-interacting system, the full Greens function verifies the Dysons equation in Nambu space45
G () = [1 + G r () (T L + T R )] G 0, () [1 + (T L + T R ) G a ()] ,
G
(r,a)
() = G
0,(r,a)
() + G
0,(r,a)
() (T L + T R ) G
(r,a)
() ,
(A18)
(A19)
gi,j
0 (z) fi,j 0 (z)
,
(A20)
tj 0
0 tj
.
(A21)
0
G 0,<
i,j 0 () = 2ii,j 0 () nF () ,
(A22)
G 0,>
i,j 0
(A23)
2i0i,j0
() [nF () 1] ,
i 0
0
1
0,r
0
i,j 0 () i,j 0 ()
i,j0 () = Im G i,j0 () =
,
0
i,j
0i,j0 ()
0 ()
() =
We only need the derivative with respect to the voltages, which are only in the leads. This gives,
(A24)
X dnL
dg1,1
d
nL 0 r
2
0 r
a
a
= 2itL
g
g
+
f
f
,
d (eVL )
d (eVL ) L 1,1s 1s,1 d (eVL ) L 1,1s 1s,1
s
X dnR
dg1,1
d
nR 0 r
r
a
a
= 2it2R
0R g1,N
g
+
f
f
,
s N s,1
d (eVR )
d (eVR )
d (eVR ) R 1,N s N s,1
s
X dnL
dgN ,N
d
nL 0 r
r
a
a
= 2it2L
0L gN
g
+
f
f
,
,1s 1s,N
d (eVL )
d (eVL )
d (eVL ) L N ,1s 1s,N
s
X dnR
dgN ,N
d
nR 0 r
2
0 r
a
a
= 2itR
g
g
+
f
f
.
d (eVR )
d (eVR ) R N ,N s N s,N d (eVR ) R N ,N s N s,N
s
r
a
Substituting into Eqs. A14-A17, and using the result gj,j
() gj,j
() = 2ij (), where we have defined the
local density of states j () j,j (), yields
GLL
GLR
GRL
GRR
"
#
X dnL
X d
e2 X
dnL
nL
r
a
r
a
d L ()
21
L g1,1s g1s,1
L f1,1s f1s,1 ,(A25)
=
h
d (eVL )
d (eVL )
d (eVL )
s
s
d
n
dnR r
e2 X
R
r
a
d L R
g
ga
+
L R f1,N
,
(A26)
=
s fN s,1
h ,s
d (eVR ) 1,N s N s,1 d (eVR )
dnL ()
d
nL
e2 X
r
a
r
a
d
R L gN ,1s g1s,N +
R L fN ,1s f1s,N ,
=
(A27)
h ,s
d (eVL )
d (eVL )
X dnR
dnR
e2 X
r
a
2N
R gN
d R ()
=
,N s gN s,N
h
d (eVR )
d
(eV
)
R
s
X d
nR
r
a
R fN ,N s fN s,N ,
(A28)
d (eVR )
2t2j 0j
(A29)
() .
(A30)
In particular at T = 0 and zero-bias, and assuming electron-hole symmetry in the leads (i.e., j = j ), we obtain
GLL
GLR
GRL
GRR
"
#
X
2 X 2 r 2
e2 X
2 r
L f1,1s
,
=
2L 1
L g1,1s +
h
s
s
=0
h
2 r
i
e2 X
r
f1,N s 2
=
L R g1,N
,
s
h ,s
=0
h
2 r
i
e2 X
r
fN ,1s 2
=
R L gN
,
,1s
h ,s
=0
#
"
X
2
2 X 2 r
e2 X
2 r
2R N
R gN ,N s +
R fN ,N s
=
h
s
s
(A31)
(A32)
(A33)
(A34)
=0
38,39
To make contact
, we can express these results in a more standard form by recalling
that ML =
P with BTK theory
P
2Tr [L 1 ] = 2L 1 () is the number ofmodes in the lead
L,
and
M
=
2Tr
[
]
=
2
R
R
R N
N
(),
L(R)
0
L(N ), ()
0
where we have defined the matrices L(R) =
, and 1(N ) = 2
(see
0
L(R)
0
L(N ), ()
8
Ref. 48). On the other hand, defining the matrices
r
r
L g1,1
L g1,1
rLL
=
r
ee
gr
L g1,1
L 1,1
r
r
R gN ,N R gN
,N
rRR
=
ee
gr
gr
R N ,N r R N,N r
L R g1,N
L R g1,N
LR
tee =
r
r
L R gN
L R g1,N
,N
r
r
L f1,1
L f1,1
r
fr
L f1,1
L 1,1
r
r
R fN ,N R fN
,N
=
fr
fr
R N ,N r R N,N r
L R f1,N
L R f1,N
=
r
r
L R fN
L R f1,N
,N
rLL
eh =
rRR
eh
tLR
eh
GLL =
(A35)
GLR
(A36)
GRL
GRR
(A37)
(A38)
In order to make explicit the non-local terms in these expressions we make use of the identity48
G r () G a () = G r () [r () a ()] G a () ,
(A39)
(A43)
and hence, substituting into Eqs. A25-A28, we obtain
GLL =
e2 X 2 r
a
r
a
r
a
2L f1,1s f1s,1
+ L R g1,N
,
s gN s,1 + f1,N s fN s,1
h ,s
(A44)
GLR =
r
e2 X
a
r
a
L R g1,N
s gN s,1 f1,N s fN s,1 ,
h ,s
(A45)
GRL =
r
e2 X
a
r
a
R L gN
,1s g1s,N fN ,1s f1s,N ,
h ,s
(A46)
GRR =
e2 X 2 r
a
r
a
r
a
,
2R fN ,N s fN
s,N + R L gN ,1s g1s,N + fN ,1s f1s,N
h ,s
(A47)
Appendix B: Transmission probability and topological phase diagram for a closed system obtained via the
Transfer Matrix method
The equations of motion for the fermionic operators cn, , cn, in the isolated N -site NW (see Eq. 2) are
d
cn, = t(cn+1, + cn1, ) n cn, + Vz,n cn, + (cn+1, cn1, ) + n cn,
dt
d
i cn, = t(cn+1, + cn1, ) n cn, + Vz,n cn, (cn+1, cn1, ) n cn,
dt
i
(B1)
(B2)
where we have included possible inhomogeneity in the chemical potential, paring potential and magnetic field (random
hopping could also be easily incorporated). In the Majorana basis cn, = (an + ibn )/2, cn, = (an ibn )/2, cn, =
(
an + ibn )/2, c = (
an ibn )/2 the equations of motion are
n,
d
an
dt
d
bn
dt
d
a
n
dt
d
bn
dt
(B3)
(B4)
(B5)
= t(
an+1 + a
n1 ) n a
n + Vz,n an (an+1 an1 ) n an
(B6)
We are
of the NW which we assume are linear combinations of Majorana operators,
P interested in the normal modes
Q = n (n an + n a
n + in bn + i
nbn ). For clarity we suppress a label indexing the modes. The coefficients s and
s are determined by requiring that the operator be an eigenmode of energy E, i.e., idQ/dt = EQ. Using the fact the
Majorana operators are complete and matching like terms we obtain the discrete form of the Schrodinger equation
n
t
=
n
t
E n =
n
n+1
t
+
t
n+1
n+1
t
+
t
n+1
n1
n
Vz,n + n
n
+
,
t
n1
Vz,n n
n
n
n1
n
Vz,n n
n
+
,
t
n1
Vz,n + n
n
n
(B7)
(B8)
(B9)
(B10)
P
P
At zero energy we can define two independent Majorana operators Q1 n (n an + n a
n ) and Q2 n (n bn + nbn )
which contain the same information about the localization properties of the system. Focusing on the transfer matrix
for the an , a
n modes, we obtain at zero energy,
1
n un 1
n
Mn =
(B11)
n
0
n
n =
(B12)
n
t
n =
(B13)
t
n
Vz,n + n
un =
(B14)
Vz,n n
n
and hence the Mn matrix is a 4 4 matrix,
10
tn +(Vz,n n )
Mn =
t2 +2
n +t(Vz,n n )
t2 +2
n +t(Vz,n +n )
t2 +2
tn (Vz,n +n )
t2 +2
t
t2 +2
t2 +2
0
0
t2 +2
t
t2 +2 .
0
0
(B15)
In the presence of disorder there is no translational invariance and the transfer matrices Mn will site dependent. The
topological invariant can be constructed from the eigenvalues of the full transfer matrix
M=
N
Y
Mn .
(B16)
n=1
1
f 0 (z)
nf =
dz
,
(B17)
2i |z|=1
f (z)
should be odd. The above considerations give a concrete way to find the phase boundary between topological and
non-topological regions in a closed system. In the clean case, where the transfer matrices Mn are all equal, the physics
of localization is determined by anyp
of the Mn matrices, and the
p well-known Pfaffian criterion for a topological phase
transition in an isolated NW2 , i.e., (2t + )2 + 2 < Vz < (2t )2 + 2 , is recovered.
From the full transfer matrix we obtain the transmission matrix tt using the identity49
[2 + M M + (M M )1 ]1 =
1
4
tt 0
.
0 t0 t0
(B18)
10
11
12
13
14
15
16
17
18
19
20
21
22
23
11
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49