Tunneling conductance in semiconductor-superconductor hybrid structures

TL;DR

This paper investigates how tunneling conductance in semiconductor-superconductor hybrid structures reveals the effects of disorder, magnetic fields, and coupling strength on Majorana zero modes, highlighting challenges in realizing topological quantum states.

Contribution

It introduces a comprehensive model treating the superconductor as an active component, showing how sub-gap states affect conductance features and the stability of Majorana modes.

Findings

Sub-gap states cause particle-hole asymmetry in conductance.

Disorder and magnetic fields lead to breakdown of zero-bias peak quantization.

Coupling strength influences the low-energy conductance behavior.

Abstract

We study the differential conductance for charge tunneling into a semiconductor wire--superconductor hybrid structure, which is actively investigated as a possible scheme for realizing topological superconductivity and Majorana zero modes. The calculations are done based on a tight-binding model of the heterostructure using both a Blonder-Tinkham-Klapwijk approach and a Keldysh non-equilibrium Green's function method. The dependence of various tunneling conductance features on the coupling strength between the semiconductor and the superconductor, the tunnel barrier height, and temperature is systematically investigated. We find that treating the parent superconductor as an active component of the system, rather than a passive source of Cooper pairs, has qualitative consequences regarding the low-energy behavior of the differential conductance. In particular, the presence of sub-gap states in the parent superconductor, due to disorder and finite magnetic fields, leads to characteristic particle-hole asymmetric features and to the breakdown of the quantization of the zero-bias peak associated with the presence of Majorana zero modes localized at the ends of the wire. The implications of these findings for the effort toward the realization of Majorana bound states with true non-Abelian properties are discussed.

Equations

$$H=H_m+H_{sm}+H_{sc}+T_m+T_{sc}+H_{ext},$$

$$H_m=\sum _{i,\delta }{t}_m^{\delta }{\mathrm{a}}_i^{†}{\mathrm{a}}_{i+\delta }+{\mu }_m\sum _i{\mathrm{a}}_i^{†}{\mathrm{a}}_i,$$

$$H_{sm}+H_{ext}$$

$$H_{sc}$$

$$T_m={\tilde{t}}_{m-sm}\sum _{i_y=1}^{N_y}\left({\mathrm{a}}_{(N_m,i_y)}^{†}{c}_{(1,i_y)}+h.c.\right),$$

$$T_{sc}=\tilde{t}\sum _{i_x=1}^{N_{sm}}\sum _{i_y=1}^{N_y}\left(c_i^{†}{a}_i+a_i^{†}{c}_i\right),$$

$${\Sigma }_{sc}(\omega )={\tilde{t}}^2{G}_{sc}(\omega ),$$

$${\Sigma }_{sc}(\omega )=-|\tilde{t}{|}^2{\nu }_{sc}\left(\frac{\omega {\hat{\tau }}_0+{\Delta }_0{\hat{\tau }}_x}{\sqrt{{\Delta }_0^2-{\omega }^2}}+\zeta {\hat{\tau }}_z\right),$$

$$\gamma =|\tilde{t}{|}^2{\nu }_{sc}.$$

$${\Sigma }_{sc}=-\gamma {\tau }_x.$$

$$\bar{\rho }(i,\omega )$$

$$\rho (\omega )$$

$$\sum _{j_x=0}^{N_{tot}+1}\sum _{j_y=1}^{N_y}\sum _{{\sigma }^{\prime }}({\mathcal{H}}_{i\sigma ,j{\sigma }^{\prime }}-\omega {\delta }_{i,j}{\delta }_{\sigma ,{\sigma }^{\prime }}){\Psi }_{j{\sigma }^{\prime }}=0$$

$${\Psi }_{(0,j_y)}^{(\nu ,\sigma )}$$

$${\Psi }_{(1,j_y)}^{(\nu ,\sigma )}$$

$$\displaystyle~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}+\sum_{\nu^{\prime}}\phi_{\nu^{\prime}}(j_{y})\left(\begin{array}[]{c}{[r_{N}]_{\nu\sigma,\nu^{\prime}\uparrow}~{}e^{-ik_{e}^{\nu^{\prime}}a}}\\ {[r_{N}]_{\nu\sigma,\nu^{\prime}\downarrow}~{}e^{-ik_{e}^{\nu^{\prime}}a}}\\ {[r_{A}]_{\nu\sigma,\nu^{\prime}\uparrow}~{}e^{ik_{h}^{\nu^{\prime}}a}}\\ {[r_{A}]_{\nu\sigma,\nu^{\prime}\downarrow}~{}e^{ik_{h}^{\nu^{\prime}}a}}\end{array}\right),$$

$$k_{e(h)}^{\nu }(\omega )={\cos }^{-1}\left(-\frac{{\mu }_m+{ϵ}_{\nu }\pm \omega }{2t_m}\right),$$

$${\Psi }_{(N_{tot}+1,j_y)}^{(\nu ,\sigma )}$$

$${\Psi }_{(N_{tot},j_y)}^{(\nu ,\sigma )}$$

$$q_{e(h)}^{\nu }={\cos }^{-1}\left(-\frac{{\mu }_{sc}+{ϵ}_{\nu }\pm \sqrt{{\omega }^2-{\Delta }_0^2}}{2t_{sc}}\right),$$

$${\tilde{\Psi }}^{(\nu \sigma )}=(\mathcal{H}+Q(\omega )-\omega {)}^{-1}{J}^{(\nu \sigma )}(\omega ),$$

$${\tilde{\Psi }}_{1,j_y}^{(\nu ,\sigma )}$$

$${\tilde{\Psi }}_{N_{tot},j_y}^{(\nu ,\sigma )}$$

$$\times$$

$$u^2(\omega )=1-v^2(\omega )=\frac{1}{2}\left(1+\frac{\sqrt{{\omega }^2-{\Delta }^2}}{\omega }\right).$$

$$Q_{(1,j_y;1,j_y^{\prime })}^{(\nu ,{\nu }^{\prime })}(\omega )=t_m{\delta }_{\nu ,{\nu }^{\prime }}{\delta }_{j_y,j_y^{\prime }}$$

$$\times \left(\begin{array}{ccccccccc}{e}^{i{k}_e^{\nu }a}& & 0& & 0& & 0& & \\ 0& & e^{i{k}_e^{\nu }a}& & 0& & 0& & \\ 0& & 0& & -e^{-i{k}_h^{\nu }a}& & 0& & \\ 0& & 0& & 0& & 0& & -e^{-i{k}_h^{\nu }a}\end{array}\right),$$

$$Q_{(N_{tot};j_y,N_{tot},j_y^{\prime })}^{(\nu ,{\nu }^{\prime })}(\omega )=\left(\begin{array}{ccc}{Q}_{ee}^{(\nu ,{\nu }^{\prime })}(\omega )& & Q_{eh}^{(\nu ,{\nu }^{\prime })}(\omega )\\ -Q_{eh}^{(\nu ,{\nu }^{\prime })}(\omega )& & Q_{hh}^{(\nu ,{\nu }^{\prime })}(\omega )\end{array}\right),$$

$$Q_{ee}^{(\nu ,{\nu }^{\prime })}(\omega )=\frac{t_{sc}}{u^2-v^2}{\delta }_{\nu ,{\nu }^{\prime }}{\delta }_{j_y,j_y^{\prime }}$$

$$\times \left(\begin{array}{ccc}{u}^2{e}^{i{q}_e^{\nu }a}-v^2{e}^{-i{q}_h^{\nu }a}& & 0\\ 0& & u^2{e}^{i{q}_e^{\nu }a}-v^2{e}^{-i{q}_h^{\nu }a}\end{array}\right),$$

$$Q_{hh}^{(\nu ,{\nu }^{\prime })}(\omega )=\frac{t_{sc}}{u^2-v^2}{\delta }_{\nu ,{\nu }^{\prime }}{\delta }_{j_y,j_y^{\prime }}$$

$$\times \left(\begin{array}{ccc}{u}^2{e}^{-i{q}_h^{\nu }a}-v^2{e}^{i{q}_e^{\nu }a}& & 0\\ 0& & u^2{e}^{-i{q}_h^{\nu }a}-v^2{e}^{i{q}_e^{\nu }a}\end{array}\right),$$

$$Q_{eh}^{(\nu ,{\nu }^{\prime })}(\omega )=\frac{t_{sc}}{u^2-v^2}{\delta }_{\nu ,{\nu }^{\prime }}{\delta }_{j_y,j_y^{\prime }}$$

$$\times \left(\begin{array}{ccc}0& & uv(e^{-i{q}_h^{\nu }a}-e^{i{q}_e^{\nu }a})\\ uv(e^{-i{q}_h^{\nu }a}-e^{i{q}_e^{\nu }a})& & 0\end{array}\right).$$

$$J_{(1,j_y)}^{(\nu \sigma )}(\omega )=$$

$$J_{(2,j_y)}^{(\nu \sigma )}(\omega )=$$

$${\tilde{\Psi }}^{\prime (\nu \sigma )}=({\mathcal{H}}_m+{\mathcal{H}}_{sm}+{\mathcal{T}}_m+{\Sigma }_Q+Q_m-\omega {)}^{-1}{J}^{\prime (\nu \sigma )}(\omega ),$$

Full text

I Introduction

The search for quantum states of matter characterized by nontrivial topological properties has gained significant momentum in recent years motivated, in part, by the discovery of topological insulators Hasan and Kane (2010); Qi and Zhang (2011); Bernevig (2013) and by the exciting prospect of using topological objects, such as the Majorana zero modes, Alicea (2012); Leijnse and Flensberg (2012); Beenakker (2013); Stanescu and Tewari (2013) as a platform for fault tolerant quantum computation. Nayak et al. (2008); Stanescu (2017) The proposal for engineering topological superconductors capable of supporting Majorana zero modes that has generated the most promising experimental results involves a semiconductor nanowire with strong spin-orbit coupling proximity-coupled to a standard s-wave superconductor in the presence of an applied magnetic field.Lutchyn et al. (2010); Oreg et al. (2010) Upon increasing the magnetic field, the system undergoes a topological quantum phase transitionRead and Green (2000); Sau et al. (2010) from a topologically-trivial superconductor to a topological superconductor with mid-gap Majorana zero modes localized near the ends of the wire. The presence of the Majorana bound states can be probed by tunneling charge into the end of the wire, which is predicted to generate a zero-bias conductance peak that is quantized ($2e^{2}/h$) at zero temperature.Sengupta et al. (2001); Law et al. (2009); Flensberg (2010); Wimmer et al. (2011); Fidkowski et al. (2012) These theoretical proposals and predictions have led to an impressive experimental effortMourik et al. (2012); Deng et al. (2012); Das et al. (2012); Finck et al. (2013); Chang et al. (2015); Albrecht et al. (2016); Zhang et al. (2016); Deng et al. (2016); Chen et al. (2016) for realizing Majorana zero modes in semiconductor-superconductor hybrid structures. The results of this effort are promising, but a conclusive demonstration of the realization of Majorana bound states in the laboratory remains elusive.

A serious problem is the persistence of significant discrepancies between theory and experiment. While the theoretical possibility of realizing topological superconductivity and Majorana zero modes in semiconductor-superconductor structures that satisfy certain requirements is beyond any reasonable doubt, understanding how well real devices realized in the laboratory satisfy these requirements remains a serious challenge. In addition to careful systematic measurements, overcoming this challenge calls for theoretical studies that go beyond idealized conditions and incorporate all relevant details associated with real devices. Since charge tunneling is, so far, the preferred method of detecting the presence of Majorana zero modes, significant theoretical effort was devoted to numerical and analytical studies of the tunneling differential conductance using both normal metalStanescu et al. (2011); Pientka et al. (2012); Prada et al. (2012); Lin et al. (2012); Rainis et al. (2013); Stanescu et al. (2014); Yan and Wan (2014); Setiawan et al. (2015); Liu et al. (2017) and superconductorPeng et al. (2015); Sharma and Tewari (2016); Chevallier and Klinovaja (2016); Setiawan et al. (2017) leads. Most of these studies, however, do not consider explicitly the parent s-wave superconductor; instead, its role is reduced to that of a passive source of Cooper pairs for the semiconductor nanowire. Recently, it has been emphasized that the parent superconductor should be treated as an active ingredient in the physics of the heterostructure.Reeg and Maslov (2017); Stanescu and Sarma (2017) This treatment of the semiconductor wire and the parent superconductor on equal footing leads to a renormalization of the low-energy modes in the hybrid systemStanescu and Sarma (2017) and the emergence of additional conductance peaks at energies corresponding to the bulk gap of the parent superconductor.Reeg and Maslov (2017)

In this work we investigate theoretically the tunneling conductance of a normal-metal – semiconductor-wire – superconductor hybrid structure and demonstrate that the explicit treatment of the parent superconductor as an active component of the system has major consequences not only at energies on the order of the bulk superconducting gap, but all the way down to zero energy. Our findings are based on a systematic numerical analysis of a tight-binding model of the heterostructure. The conductance calculations are carried out within both the Blonder-Tinkham-Klapwijk (BTK) theoryBlonder et al. (1982) and the Keldysh non-equilibrium Green’s function formalismRammer and Smith (1986). To gain a better physical understanding of the low-energy features revealed by the differential conductance, we compare these features with the local density of states. We establish the qualitative features that characterize the dependence of tunneling conductance on the coupling strength between the semiconductor wire and the parent superconductor, the strength of the tunnel barrier, and temperature. We also consider the possibility of sub-gap states being present in the parent superconductor and identify the consequences of these states hybridizing with low-energy states in the wire. We find that the presence of sub-gap states leads to characteristic particle-hole asymmetric features in the low-energy differential conductance and to a breakdown of the zero-bias peak quantization at zero temperature. Our findings suggest that a detailed modeling of the parent superconductor, including the possible presence of disorder-induced sub-gap states in the presence of a finite magnetic field, is a strong requirement for any attempt to quantitatively (and, at this stage, even qualitatively) reproduce the experimentally measured differential conductance.

The remainder of this paper is organized as follows. In Sec. II we introduce a tight-binding model for the normal-metal – semiconductor-wire – superconductor hybrid structure (Sec. II.1) and briefly describe the BTK-type formalism (Sec. II.2) and the Keldysh approach (Sec. II.3) for calculating the differential conductance. The results of our numerical analysis are presented in Sec. III. First, we compare the two approaches for calculating the tunneling conductance and establish the correspondence between the differential conductance (in the low tunneling regime) and the local density of states (Sec. III.1). The dependence of the low-energy features on the transparency of the semiconductor-superconductor interface and the tunnel barrier height is discussed in Sec. III.2. We show that additional features associated with charge transport via quasiparticles in the parent superconductor occur at energies larger than the bulk superconducting gap. The effects of disorder-induced sub-gap states in the parent superconductor are investigated in Sec.III.3 and the dependence of the Majorana-induced zero-bias conductance peak on the relevant parameters (i.e., tunnel barrier height, temperature, and sub-gap states density) is summarized in Sec. III.4. Our conclusions are presented in Sec. IV.

II Theoretical model

In this section we introduce the tight-binding model for the normal-metal – semiconductor-wire – superconductor hybrid structure and briefly describe the BTK and Keldysh methods for calculating the differential conductance. Both approaches explicitly incorporate the parent superconductor.

II.1 Tight-binding Hamiltonian

We consider charge tunneling from a normal metal lead into a semiconductor wire - superconductor hybrid structure through a controllable potential barrier. The Hamiltonian describing the system has the following generic form:

[TABLE]

where $H_{m}$, $H_{sm}$, and $H_{sc}$ describe the normal metal, semiconductor wire, and superconductor components, respectively, the next two terms describe the coupling between the semiconductor wire and the metallic lead ($T_{m}$) and the coupling between the wire and the parent superconductor ($T_{sc}$), while $H_{ext}$ includes contributions from external fields, such as magnetic fields and electrostatic potentials. The specific form of each term in Eq. (1) depends on on the degree of complexity that we want to incorporate into the model. In this study we consider a simple tight-binding description of the hybrid system involving a lattice consisting of coupled parallel chains. Specifically, we have

[TABLE]

where $i=(i_{x},i_{y})$ labels the position, with $1\leq i_{y}\leq N_{y}$ designating the chain and $1\leq i_{x}\leq N_{m}$ the position along the chain, $\delta=(\pm 1,0)$ or $(0,\pm 1)$ is a nearest neighbor vector, $t_{m}^{\delta}=(t_{m}^{x},t_{m}^{y})$ is the hopping matrix element, and $\mu_{m}$ the chemical potential of the metallic lead. Using spinor notation, the electron creation operator on site $i$ is ${\rm a}_{i}^{\dagger}=({\rm a}_{i\uparrow}^{\dagger}~{}{\rm a}_{i\downarrow}^{\dagger})$. The semiconductor wire (including the applied external fields) is described by the tight-binding Hamiltonian

[TABLE]

where $c_{i}^{\dagger}=(c_{i\uparrow}^{\dagger}~{}c_{i\downarrow}^{\dagger})$ is the electron creation operator on the site $i=(i_{x},i_{y})$ of the semiconductor wire, with $1\leq i_{x}\leq N_{sm}$ and $1\leq i_{y}\leq N_{y}$, $t_{sm}^{\delta}=(t_{sm}^{x},t_{sm}^{y})$ is the nearest-neighbor hopping, $\delta=(\pm 1,0)$ or $(0,\pm 1)$ is a nearest neighbor vector, while $\delta_{x}=(1,0)$ and $\delta_{y}=(0,1)$ are (positive) unit vectors in the x- and y-direction, respectively, $\mu_{sm}$ is the chemical potential, and $(\alpha_{x},\alpha_{y})$ are the Rashba spin-orbit coupling coefficients. The position-dependent potential $V_{i}$ describes the tunneling barrier and additional back-gate-induced potentials, while $\Gamma$ represents the Zeeman splitting. The matrices $\hat{\sigma}_{\mu}$, with $\mu=x,y,z$, represent Pauli matrices associated with the spin degree of freedom. The parent superconductor is described within the Bogoliubov-de Gennes (BdG) formalism by

[TABLE]

where $a_{i}^{\dagger}=(a_{i\uparrow}^{\dagger}~{}a_{i\downarrow}^{\dagger})$ is the electron creation operator on the site $i=(i_{x},i_{y})$ of the superconducting chains, with $1\leq i_{x}\leq N_{sc}$ and $1\leq i_{y}\leq N_{y}$, $t_{sc}^{\delta}=(t_{sc}^{x},t_{sc}^{y})$ is the nearest-neighbor hopping matrix element, $\mu_{sc}$ the chemical potential of the superconductor, and $\Delta_{0}$ is the superconducting gap. We assume that the superconductor is equal to or longer than the semiconductor wire, $N_{sc}\geq N_{sm}$. The coupling between the semiconductor wire and the metallic lead is described by the term

[TABLE]

where the nearest-neighbor hopping $\tilde{t}_{m-sm}$ quantifies the coupling strength, while $(N_{m},i_{y})$ and $(1,i_{y})$ label the sites at the right end of the metallic lead and the left end of the semiconductor wire, respectively. Finally, the coupling between the semiconductor and the superconductor is described by the term

[TABLE]

where $i=(i_{x},i_{y})$, $\tilde{t}$ is the hopping across the semiconductor-superconductor (SM-SC) interface, while $c_{i}$ and $a_{i}$ designate annihilation operators on the SM and SC sides, respectively.

It may be convenient to use a Green’s function formalism and calculate the effective Green function for the semiconductor wire by integrating out the degrees of freedom of the parent superconductor.Stanescu and Tewari (2013) The proximity effect induced by the superconductor is captured by a self-energy term,

[TABLE]

where $G_{sc}(\omega)$ is the Green function of the parent superconductor at the SM-SC interface. Assuming that the parent superconductor is a truly bulk system (i.e., wide enough and thick enough), the self-energy becomes local,Stanescu and Tewari (2013) $\Sigma_{sc}(\omega;i,j)=\delta_{i,j}\Sigma_{sc}(\omega)$, with

[TABLE]

where $\nu_{sc}=\sqrt{4t_{sc}\mu_{sc}-\mu_{sc}^{2}}/({2t_{sc}^{2}})$ is the surface density of states of the bulk superconductor at the chemical potential and $\zeta=(2t_{sc}\mu_{sc}-4t_{sc}^{2})/(\mu_{sc}^{2}-4t_{sc}\mu_{sc})$ is a proximity-induced shift of the SM chemical potential. The matrices $\hat{\tau}_{\mu}$, with $\mu=x,y,z$, represent Pauli matrices associated with the particle-hole degree of freedom. We define the effective SM-SC coupling as

[TABLE]

Note that here we do not discuss the physics associated with the finite wire thickness in the direction perpendicular to the SM-SC interface (i.e., the $z$-direction). This is critical when considering electrostatic effects, such as those generated by back-gate potentials. In general, the effective SM-SC coupling $\gamma$ is a band-dependent quantity proportional to the amplitude squared of the wave-functions (corresponding to a given confinement-induced band) at the SM-SC interface and may involve band off-diagonal contributions.Stanescu and Tewari (2013) In the weak-coupling limit, $\gamma/\Delta_{0}\rightarrow 0$, the self-energy (8) can be approximated at low energy by the anomalous contribution

[TABLE]

Note that in this approximation the proximity-induced gap (defined as the minimum quasiparticle gap in the absence of an applied magnetic field) is $\Delta_{ind}=\gamma$ and the minimum critical Zeeman field associated with the topological quantum phase transition is $\Gamma_{c}=\gamma\equiv\Delta_{ind}$. However, in general the critical field and the induced gap have different values.Sticlet et al. (2016); Stanescu and Sarma (2017) To describe the low-energy spectral properties of the hybrid structure it is convenient to use the density of states (DOS) $\rho(\omega)$ and the local density of states (LDOS) $\bar{\rho}(i,\omega)$ at the end of the wire. In terms of the effective wire Green’s function $G(\omega)=[\omega-H_{sm}-\Sigma_{sc}(\omega)]^{-1}$ we have

[TABLE]

II.2 BTK Formalism

The calculation of the differential conductance using the Blonder-Tinkham-Klawijk (BTK) formalismBlonder et al. (1982) involves calculating the reflection and transmission coefficients for incoming (and outgoing) plane waves by solving the Bogoliubov-de Gennes (BdG) equation for the total Hamiltonian given by Eq. (1) with appropriate boundary conditions. The Schrödinger equation that determines the reflection and transmission coefficients is

[TABLE]

where ${\cal H}$ is the first quantized BdG Hamiltonian corresponding to Eq. (1) and $N_{tot}=N_{m}+N_{sm}+N_{sc}$. For convenience, we define the transverse modes $\phi_{\nu}$, with $1\leq\nu\leq N_{y}$, characterized by the wave functions $\phi_{\nu}(i_{y})=\sqrt{2/(N_{y}+1)}\sin[i_{y}\nu\pi/(N_{y}+1)]$. Each transport channel is labeled by the pair of quantum numbers $(\nu~{}\!\sigma)$ corresponding to a given transverse mode and spin orientation and the reflection and transmission coefficients are $2N_{y}\times 2N_{y}$ matrices with matrix elements indexed by these channel labels.

Consider now an incoming plane wave in channel $(\nu,\sigma)$. To ensure that the wave function entering the normal lead is a (propagating) plane wave, we impose (propagating) boundary conditions involving the first two sites of each chain, i.e., $j_{x}=0$ and $j_{x}=1$. Note that similar propagating boundary conditions involving two sites at the rightmost end of the superconductor, i.e., $j_{x}=N_{tot}$ and $j_{x}=N_{tot}+1$, correspond to an outgoing plane wave. Specifically, the (incoming) boundary conditions at the left end of the metallic lead are

[TABLE]

where $r_{N}$ and $r_{A}$ are the normal and anomalous reflection coefficients, respectively, $a$ is the lattice constant, while $k_{e}^{\nu}$ and $k_{h}^{\nu}$ are wave vectors in the normal lead at energies $\omega$ and $-\omega$, respectively,

[TABLE]

with $\epsilon_{\nu}=2t_{m}\cos[\nu\pi/(N_{y}+1)]$. Similarly, the outgoing plane wave dictates the boundary conditions at the rightmost end of the bulk superconductor chains,

[TABLE]

where $t_{N}$ and $t_{A}$ are the normal and anomalous transmission coefficients, respectively, while $q_{e}^{\nu}$ and $q_{h}^{\nu}$ are wave vectors in the bulk superconductor

[TABLE]

with $\epsilon_{\nu}=2t_{sc}\cos[\nu\pi/(N_{y}+1)]$ and $|\omega|>\Delta_{0}$. To solve the BdG equation (II.2) with the boundary conditions (27) and (38) it is convenient to express the $16N_{y}$ wave function components $\Psi_{(j_{x},j_{y})}^{(\nu,\sigma)}$ corresponding to $j_{x}=0,1,N_{tot},N_{tot}+1$ in terms of the $4N_{y}$ reflection and $4N_{y}$ transmission coefficients corresponding to the incoming channel $(\nu,\sigma)$. Thus, Eq. (II.2) reduces to a system of $4N_{tot}N_{y}$ linear equations with $4N_{tot}N_{y}$ unknown coefficients that include the $4(N_{tot}-2)N_{y}$ wave function components $\Psi_{(j_{x},j_{y})}^{(\nu,\sigma)}$, with $2\leq j_{x}\leq N_{tot}-1$, and the $8N_{y}$ reflection and transmission coefficients. Writing theses unknown coefficients as a vector $\widetilde{\Psi}$ with $4N_{tot}N_{y}$ components, we have

[TABLE]

Explicitly, we define the components of the new vector $\widetilde{\Psi}^{(\nu,\sigma)}$ to be the same as the wave function components, $\widetilde{\Psi}^{(\nu,\sigma)}_{(j_{x},j_{y})}\equiv\Psi^{(\nu,\sigma)}_{(j_{x},j_{y})}$, for $2\leq j_{x}\leq N_{tot}-1$, while

[TABLE]

where $u$ and $v$ are BCS coherence factors:

[TABLE]

The matrix $Q$ introduces frequency-dependent correction at the boundaries,

[TABLE]

with

[TABLE]

All other elements of $Q$ are zero, $Q^{(\nu,\nu^{\prime})}_{(j_{x},j_{y};~{}\!j_{x}^{\prime},j_{y}^{\prime})}(\omega)=0$ if $j_{x}\neq j_{x}^{\prime}$ or $j_{x}=j_{x}^{\prime}\neq 1,N_{tot}$. Finally, the vector $J$ contains the “free terms” from the boundary conditions (i.e., those that do not depend on the reflection and transmission coefficients), which are associated with the incoming plane wave. Explicitly, we have

[TABLE]

and $J^{(\nu\sigma)}_{(j_{x},j_{y})}(\omega)=0$ for $j_{x}\geq 3$. Note that all the information about the boundary conditions is contained in the matrix $Q(\omega)$ and the vector $J^{(\nu\sigma)}(\omega)$. The coefficients $[r_{N(A)}]_{\nu\sigma,\nu^{\prime}\sigma^{\prime}}$ and $[t_{N(A)}]_{\nu\sigma,\nu^{\prime}\sigma^{\prime}}$ are obtained by solving Eq. (45) $2N_{y}$ times, once for each incoming channel.

An interesting observation is that the matrix in Eq. (45) has the structure of a Green’s function, $G_{Q}(\omega)=({\cal H}+Q(\omega)-\omega)^{-1}$. Generalizing the approach used in Sec. II.1 to calculate the effective Green’s function for the semiconductor wire, we integrate out the degrees of freedom of the bulk superconductor and write an equation for the reflection coefficients involving only states in the semiconductor and the metallic lead:

[TABLE]

where ${\cal H}_{m}$, ${\cal H}_{sm}$, and ${\cal T}_{m}$ are the matrices corresponding to the first quantized Hamiltonians for the metallic lead, the semiconductor wire, and the lead-wire coupling, respectively, $Q_{m}$ is the matrix corresponding to the boundary condition at the left end of the normal lead and $\Sigma_{Q}$ is a “self-energy” contribution that incorporates the effect of the parent superconductor, including the the outgoing boundary conditions. Note that below the bulk gap there is no normal current (i.e., current carried by quasiparticles) in the superconductor; hence the outgoing boundary conditions become trivial and we have $\Sigma_{Q}(\omega)=\Sigma_{sc}(\omega)$. In general, since the boundary conditions are diagonal in real space while $H_{sc}$ is diagonal in mode space, the self energy $\Sigma_{Q}$ has to be calculated numerically. Nonetheless, using Eq. (51) instead of Eq. (45) provides significant advantages, as one can calculate the self-energy $\Sigma_{Q}$ once, then solve the reduced equation repeatedly for different values of relevant parameters, such as potential barrier height and Zeeman field, instead of solving the full equation (45) for each set of parameters.

The differential conductance for tunneling into the hybrid structure can be written in terms of reflection coefficients asBlonder et al. (1982)

[TABLE]

where $V$ is the bias voltage. To include the effect of finite temperature, the conductance is broadened by convolving it with the Fermi function. Explicitly, we have

[TABLE]

where $G_{0}(\epsilon)=dI/dV$ is the zero temperature conductance given by Eq. (52) at voltage bias $\epsilon$.

II.3 Keldysh Formalism

Consider a normal metal-semiconductor junction. The current through the junction can be expressed in terms of the number operator for the component to the right of the junction as

[TABLE]

where $\ket{\Psi_{0}}$ is the equilibrium quantum state of the composite system, $H$ is the total Hamiltonian, and $N_{R}=\sum_{r}c^{\dagger}(r)c(r)$ is the number operator for the right component. The label $r$ indicates both position (on the right side of the junction) and spin. The number operator $N_{R}$ can only change if some electrons cross the junction. Hence, one can rewrite the current asBerthod and Giamarchi (2011)

[TABLE]

where $T_{m}$ is the coupling between the metal and the semiconductor given by Eq. (6). The differential conductance can be extracted by taking the derivative of the current with respect to the bias voltage in the lead, $\frac{dI}{dV}$. Notice that the propagator $iG^{<}(l,r;0^{+})=\bra{\Psi_{0}}\rm{a}^{\dagger}(l)c(r)\ket{\Psi_{0}}$ is the so-called lesser Green’s function at time $t=0^{+}$.Rammer and Smith (1986) Using the relation between $G^{<}(l,r;t)$ and the Keldysh Green’s function, Eq. (55) becomes

[TABLE]

where we have used the property $\int d\omega G^{K}(\omega)=G^{K}(l,r;0^{+})$. Expanding $G^{K}$ in terms of Green’s functions for the left and right sub-systems, we have

[TABLE]

The Keldysh Green’s function becomes

[TABLE]

where $G^{+}$ and $G^{-}$ are the retarded and advanced Green’s functions, respectively. Next, we take the derivative with respect to the bias voltage. We assume the voltage drops at the junction, so we only take the derivative in the Green’s function of the lead, $G_{L}$ (i.e., we assume that $G_{R}$ is independent of $V$). The main contribution comes from the derivative of the Fermi function with respect to the voltage, hence from $G^{K}_{L}(\omega)=[1-2f_{L}(\omega)][G^{+}_{L}(\omega)-G^{-}_{L}(\omega)]$, where $f_{L}(\omega)=f(\omega-eV)$. Explicitly, the differential conductance is given by

[TABLE]

Finite temperature can be included in the same way as in the BTK method.

III Results

In this section we discuss the results of our numerical analysis focusing on features in the calculated differential conductance that are only present if the parent superconductor is explicitly treated as an active component of the system and on the dependence of the Majorana zero bias conductance peak on the relevant system parameters.

The numerical results described below are obtained by setting the chemical potentials for the metal and the superconductor ($\mu_{m}$ and $\mu_{sc}$, respectively) to the middle of their respective bands. The other model parameters are chosen as follows: The lattice constant is $10~{}$nm along the wire (i.e., in the $x$ direction) and $33~{}$nm across the wire (i.e., in the $y$ direction); if the wire is $1~{}\mu$m long, there are 100 sites in each SM chain. The hopping parameters are $t_{sm}^{x}=9.5~{}$meV, $t_{sm}^{y}=0.95~{}$meV, $t_{m}^{x}==3.8~{}$meV, $t_{m}^{y}=0.38~{}$meV, $t_{sc}^{x}=1.5~{}$meV, $t_{sc}^{y}=0.15~{}$meV, and $\tilde{t}_{m-sm}=3.0~{}$meV. The spin orbit coupling parameters are $\alpha_{x}=2~{}$meV and $\alpha_{y}=0.6~{}$meV. The bias potential is applied to the metallic lead, which is $40$ sites long (the shortest length that is consistent with plane wave boundary conditions), while the superconductor is grounded.

We model the parent superconductor as a long, thick slab having $1000$ sites in the $x$ direction (i.e., $10~{}\mu$m) and $300$ sites in the $z$ direction. To handle the problem concerning the large number of degrees of freedom associated with the parent superconductor, we calculate the (surface) Green’s function of a large slab once and store the information. Subsequently, we “attach” this large superconductor to semiconductor wires (of various lengths) and calculate the corresponding interface self-energy, which completely contains the effect of the superconductor. More specifically, we first calculate numerically the superconductor Green’s function, which determines both the self energy $\Sigma_{sc}$ given by (7) and the boundary-condition-dependent “self-energy” $\Sigma_{Q}$ from Eq. (51), in the transverse mode space, where the problem is diagonal and can be solved separately for each of the modes. Then, we calculate the surface Green’s function for frequencies up to $3~{}$meV taking $100$ equally spaced values in the range $2$ to $3~{}$meV. If necessary, we extrapolate for frequencies in-between these these discrete values. Note that below the bulk gap (i.e., below $2~{}$meV), we can use the analytical solution given by Eq. (8). These calculated values of the surface Green’s function are stored for later use in the calculation of the self-energies corresponding to wires of different lengths.

III.1 Method comparison and differential conductance vs. LDOS

The BTKBlonder et al. (1982) and KeldyshRammer and Smith (1986) formalisms are two possible methods for calculating the differential conductance. Since the two approaches are formally rather different, the first question that we want to address is whether or not the predictions based on these methods are the same. In addition, when interpreting charge transport measurements on semiconductor -superconductor (SM-SC) Majorana hybrid structures it is essential to understand the connection between the measured differential conductance and the low-energy spectral properties of the hybrid system. To address these questions, we calculate the density of states (DOS) and the local density of states (LDOS) at the end of the wire for a nanowire of length $L=1~{}\mu$m proximity-coupled to a $\Delta_{0}=2~{}$meV superconductor so that the induced gap is $\Delta_{ind}=0.25~{}$meV. We compare these quantities with the differential conductance for tunneling into the end of the same wire, which is calculated using both methods described in the previous section. The results are shown in Fig. 1.

Panels (a) and (b) show the DOS and LDOS, respectively, as functions of the chemical potential and energy. The bottom panels in Fig. 1, i.e., (c) and (d), show the differential conductance calculated using the BTK and Keldysh methods. The “stripy” features are due to the quantization of states in a finite wire and reflect the variation of the eigenvalues of the BdG Hamiltonian with the chemical potential. In our model the energy difference between successive states becomes larger as the chemical potential increases. Note that $\mu_{sm}=0$ corresponds to the bottom of the band, which is clearly revealed by the large DOS in panel (a). By contrast, there is no clear signature associated with the bottom of the band in the LDOS [panel (b)] or in the differential conductance [panels (c) and (d)]. This is due to the fact that the LDOS is calculated at the end of the proximitized wire adjacent to the potential barrier, where states from the bottom of the band have very low amplitude. These states are weakly coupled to the lead and, consequently, do not contribute significantly to the differential conductance, as shown in panels (c) and (d). Another way of understanding this behavior is in terms of velocity. Penetration of the barrier is determined by the group velocity. An electron occupying a state with higher velocity has a greater chance to get through the barrier and couple to the lead. Since the velocity vanishes at the bottom of the band, the corresponding states will have a small contribution to the conductance. As the chemical potential increases, the wavelength of the states is shortened and their velocity becomes larger, which results in higher values of the conductance.

Based on this analysis we conclude that (i) the BTK and Keldysh formalisms predict similar values for the differential conductance and (ii) the dependence of the differential conductance on the chemical potential and bias voltage is qualitatively similar to that of the LDOS at the end of the wire adjacent to the tunnel barrier. Of course, the second property holds in the weak tunneling regime,Berthod and Giamarchi (2011) but breaks down when the transparency of the barrier is large. In Fig. 1 the only notable difference between the differential conductance and the LDOS is the decrease of $dI/dV$ with increasing bias voltage $|V_{Bias}|$, which is not reflected in the dependence of the LDOS on energy. To understand this behavior, we note that for voltage values below the gap of the parent superconductor charge can only travel through the superconductor as a supercurrent. Tunneling electrons into the wire involves anomalous (Andreev) reflection processesBlonder et al. (1982); He et al. (2014) (an incoming electron returns as a hole) and the probability for such a process depends on the properties of the quasiparticle states of the wire-superconductor hybrid system. Near the induced gap ($\Delta_{ind}\approx 0.25$ meV in Fig. 1), where every state is an almost equal blend of particles and holes, the probability of anomalous reflection is large, but it decreases at higher (in absolute value) energies, where the quasiparticles have a dominant particle or hole character.

We have verified our conclusions for a broad range of parameters, including bias voltage values larger than the parent superconductor gap, and we found them to be rather generic properties. Below we calculate the conductance using the most convenient method for a specific setup or physical aspect that we want to investigate. Since we are mostly interested in the weak tunneling regime, we will interpret the results as (approximately) representing the local density of states near the end of the wire-superconductor system adjacent to the tunnel barrier.

III.2 Dependence on the semiconductor-superconductor coupling strength and tunnel barrier height

Proximity coupling a semiconductor wire to a conventional superconductor results in the opening of an induced gap in the quasiparticle spectrum at the chemical potential. At zero momentum, $k=0$, and in the absence of a magnetic field the induced gap $\Delta_{ind}$ is determined by the gap $\Delta_{0}$ of the parent superconductor and by the effective semiconductor-superconductor (SM-SC) coupling constant $\gamma$. The simplest way to account for the emergence of an induced gap is to introduce a pairing potential in the effective BdG Hamiltonian for the wire. However, this approach fails to describe accurately the low-energy physics of the hybrid structure (i.e., the induced ”bulk” quasiparticle gaps and the in-gap bound states), with the exception of the weak coupling limit $\gamma/\Delta_{0}\rightarrow 0$. However, experimentally-relevant devices are not in this limit and the explicit treatment of the parent superconductor, e.g., using Eq. (8), is necessary.Stanescu and Sarma (2017) What features of the differential conductance are directly related to the parent superconductor being an active component of the hybrid structure?

First, we consider the dependence of the low-energy features of a finite wire in the topological regime on the effective SM-SC coupling strength. The differential conductance as function the bias voltage for three different ratios $\gamma/\Delta_{0}$ is shown in Fig. 2. The calculations are done for a wire of length $L=1~{}\mu$m in the presence of a Zeeman field $\Gamma=0.5~{}$meV. Note that the induced gap (i.e., the quasiparticle gap at $\Gamma=0$) is the same in all three cases, $\Delta_{ind}=0.25~{}$meV. However, the critical Zeeman splitting corresponding to the topological quantum phase transition (TQPT) is different, as it is determined by the effective coupling strength $\gamma$, rather than the induced gap.Stanescu and Sarma (2017); Sticlet et al. (2016) For the Zeeman field value used in the calculations, $\Gamma=0.5~{}$meV, all three systems are in the topological regime, as signaled by the presence of a Majorana-induced zero-bias conductance peak. The differential conductance in Fig. 2(a), which corresponds to $\gamma\ll\Delta_{0}$, is the same as that predicted by the effective pairing approximation defined by Eq. 10. However, increasing the effective coupling [panels (b) and (c)] results in several effects that cannot be captured by this approximation. First, the weight of the Majorana-induced zero-bias peak (ZBP) decreases with increasing coupling. The main reason for this behavior is the fact that the spectral weight of the Majorana mode within the wire decreases with increasing coupling. In other words, as $\gamma$ increases, the amplitude of the Majorana wave function within the wire (i.e., in the vicinity of the tunnel barrier) is reduced, while its amplitude within the parent superconductor (i.e., farther away from the tunnel barrier) is enhanced.Stanescu and Sarma (2017) This results in a lower “visibility” of the Majorana mode, as revealed by the ZBP. The second effect is the renormalization of the quasiparticle energies and the reduction of the topological gap.Stanescu and Sarma (2017) In essence, these effects are due to the fact that all low-energy states reside both in the wire and in the parent superconductor, with more and more spectral weight being transfered to the superconductor as the coupling strength increases.

The bias voltage range shown in Fig. 2 corresponds to energies below the gap of the parent superconductor. In this regime there is no normal current flowing through the superconductor. In the BTK approach, this means that it is sufficient to impose boundary conditions at the end of the normal lead. However, if we want to explore an energy window larger than the bulk gap, we have to supplement this with boundary conditions at the end of the parent superconductor. These boundary conditions account for the normal current carried by the quasiparticles in the parent superconductor. Note that these boundary conditions are diagonal in real space, while the superconductor Hamiltonian is diagonal in the mode space. Consequently, finding an analytical expression is rather difficult; instead, we calculate the self-energy $\Sigma_{sc}(\omega)$ (which incorporates the boundary conditions) numerically.

Including the normal transmission is required at bias voltages above the bulk gap, as demonstrated by the data shown in Fig. 3. Here we calculate the differential conductance for a short nanowire of length $L=0.3~{}\mu$m coupled to a superconductor with $\Delta_{0}=2~{}$meV. The effective semiconductor-superconductor coupling is $\gamma=1.5~{}$meV. In the top panel, the differential conductance is calculated without including the contribution from the normal transmission through the parent SC. All the features in this panel represent anomalous current contributions, which become rather small at energies larger than the bulk gap. Note that there is no feature associated with the bulk gap edge. By contrast, in the bottom panel, which shows the full calculation that includes the normal current contribution, one can clearly observe an additional feature associated with the parent superconductor gap (orange/light gray area). The emergence of an additional conductance peak at the energy corresponding to the gap $\Delta_{0}$ of the parent superconductor, which was recently discussed in Ref. Reeg and Maslov, 2017, confirms the the predictions of Ref. Cole et al., 2015, which are based on a local density of states analysis. The presence of features associated with the parent superconductor may generate certain difficulties when interpreting $dI/dV$ measurements for hybrid systems in the intermediate and strong coupling regimes. In these cases, the induced gap $\Delta_{ind}$ and the bulk gap of the parent superconductor $\Delta_{0}$ have comparable values and, therefore, it is important to be able to distinguish the corresponding features. To this end, we analyze their dependence on the effective SM-SC coupling strength and on the tunnel barrier height.

As discussed in the relation to Fig. 2, increasing the coupling between the semiconductor wire and the superconductor results in a renormalization of the low-energy spectrum and a reduction of the spectral weight within the wire. In addition, varying $\gamma$ has a strong impact on the induced gap $\Delta_{ind}$. In a multi-band system the specific dependence of the induced gap on the SM-SC coupling strength is rather complicated (except in the weak coupling regime) and involves additional parameters, such as the inter-band spacing and chemical potential.Stanescu and Sarma (2017) On the other hand, the feature in the differential conductance associated with the bulk superconductor has a relatively weak dependence on the SM-SC coupling, as shown in Fig. 4. In this figure, we calculate again at the differential conductance for a nanowire coupled to a superconductor with $\Delta_{0}=2~{}$meV, but the length of the nanowire is now $L=1~{}\mu$m; the effective SM-SC coupling takes three different values, $\gamma=0.5,~{}2.0,$ and $4.0~{}$meV. Note, that the peak associated with the bulk gap edge becomes stronger with increasing $\gamma$. The features corresponding to low-energy semiconductor states (blue/darker gray) reveal the dependence of the induced gap on the SM-SC coupling and the proximity-induced energy renormalization discussed above. In particular, the energy spacing between successive low-energy states decreases with $\gamma$, so that a finite system characterized by discrete spectral features in the low-coupling regime [panel (a)] appears as having a continuous spectrum for a high-enough coupling strength [panel (c)].

The dependence of the differential conductance on the tunnel barrier height is illustrated by the data shown in Fig. 5. We calculate $dI/dV$ as a function of the bias voltage for different barrier heights and two different wire lengths $L=0.3~{}\mu$m and $L=0.6~{}\mu$m. As before, the parent SC gap is $\Delta_{0}=2~{}$meV, while the SM-SC coupling is $\gamma=1.5~{}$meV. First, we note that the sub-gap features associated with the presence of the proximitized semiconductor wire (blue/darker gray) are strongly dependent on the length of the structure. By contrast, the parent superconductor contribution (orange/light gray) is practically independent on $L$. In very short wires [panels (a), (c), and (e)] finite size quantization results in a discrete spectrum that generates a set of well-separated peaks in the differential conductance. Increasing the length of the wire [panels (b), (d), and (f)] reduces the intervals between successive peaks and, eventually, the spectrum becomes continuous [see, for example, Fig. 4 (c) and Fig. 6]. Note that the energy separation required for resolving discrete peaks depends on temperature, barrier transparency, and strength of the disorder present in the system. Above the bulk gap of the parent superconductor the semiconductor states hybridize with the superconductor states, which are overwhelmingly more numerous and, practically, determine the conductance response. Consequently, above $\Delta_{0}$ the conductance is nearly independent of the wire length. Note, however, that the semiconductor plays an important role here. One can view all the hybrid states above the gap as superconductor states that have acquired a small “tail” that extends into the wire, including the region adjacent to the tunnel barrier. While the spectral weight carried by the tails is negligible, they are essential for ensuring the coupling between the normal lead and the “superconductor states.” In other words, the orange/light gray component of the differential conductance corresponds to hybrid states that have most of their spectral weight inside the parent superconductor, but have small tails extending into the barrier region.

Varying the strength of the barrier potential affects the differential conductance both below and above the bulk gap (see Fig. 5). However, the differential conductance below the bulk gap (blue/darker gray) decreases faster than the conductance above the gap (orange/light gray) with increasing barrier potential. In fact, it is possible that for a large enough barrier potential the sub-gap conductance is “completely” suppressed (for a certain finite resolution), while the parent superconductor contribution is still measurable. A rigorous proof of this statement would imply performing a (self-consistent) calculation of the barrier potential profile, i.e., finding a (self-consistent) solution of the Schrödinger-Poisson equation. Intuitively, one can understand the effect of the barrier potential as resulting from the huge difference between the carrier densities in the semiconductor and the superconductor. Increasing the potential of the tunneling gate results in the potential barrier penetrating into the proximitized segment of the nanowire (i.e., under the superconductor) and “pushing away” the low-energy states that reside inside this component of the hybrid structure. On the other hand, the effect on the superconductor is negligible due to the screening provided by the large carrier density. To illustrate this effect, we calculate the differential conductance for a nanowire with $L=0.6~{}\mu$m, similar to that corresponding to the right panels in Fig. 5, and a barrier potential profile that penetrates more into the wire. The results are shown in Fig. 6. Note the large difference between the conductance above the (bulk SC) gap and the sub-gap conductance, which appears as a small tail of the coherence peak. Remarkably, this shape looks quite similar to recent experimental results on epitaxial Al-InAs superconductor-semiconductor nanowiresChang et al. (2015) (more specifically, Fig. 5(c) for the half-shell device). To experimentally disentangle the induced and bulk contributions in the strong SM-SC coupling limit (where they have similar energy scales) we propose an analysis of the dependence of various features on the barrier potential. Starting, for example, with an extremely high barrier will only reveal bulk-type contributions.The induced features will progressively emerge as the potential barrier is lowered.

III.3 Parent superconductor with sub-gap states

So far we have assumed an “ideal” parent superconductor characterized by a perfectly clean gap, even in the presence of finite magnetic fields. However, real superconductors may have sub-gap states induced by disorder and finite external fields. Nonetheless, if a “dirty” superconductor has a few localized sub-gap states, we expect the system to behave qualitatively similar to a clean system. Indeed, the additional sub-gap states may hybridize with low-energy states from the semiconductor wire, which will result in the hybrid structure having additional low-energy modes and smaller quasiparticle gaps. However, the structure of the low-energy spectrum will be qualitatively similar to that of a clean system. Moreover, if the sub-gap states are localized far from the ends of the wire the Majorana modes will not be affected. The situation is completely different if the parent superconductor has a small but finite density of sub-gap states. In this situation, the entire gap is populated with hybrid states having most of their spectral weight inside the parent superconductor and some of it inside the semiconductor nanowire.Stanescu and Sarma (2017) The question that we want to address is the following: what are the characteristic signature in the tunneling conductance of a hybrid structure with a parent superconductor having a finite density of sub-gap states?

To address this question, we use a simple description of the sub-gap states in terms of an imaginary contribution to the parent superconductor Green’s function that is obtained through the substitution $\omega\longrightarrow\omega+i\delta$ in $G_{sc}(\omega)$, where $\delta$ is a magnetic field-dependent phenomenological parameter that describes the finite density of sub-gap states. We note that developing models for the parent superconductor based on a detailed understanding of the low-energy physics at the microscopic level represents an outstanding problem in this field.Cole et al. (2016); Hui et al. (2015) We also note that a finite density of sub-gap states is equivalent to having an additional (i.e., non-superconducting) equilibrium bath,Martin and Mozyrsky (2014) which generates dissipation. In turn, the presence of dissipation introduces particle-hole asymmetry into the finite energy conductance spectrum.Martin and Mozyrsky (2014); Bauriedl et al. (1981); Yazdani et al. (1997) The importance of this mechanism in understanding certain features of the measured conductance data in semiconductor-superconductor hybrid structures was emphasized recently.Das Sarma et al. (2016); Liu et al. (2017)

Consider now the dependence of the differential conductance of a short wire $L=0.4~{}\mu m$ on the chemical potential and bias voltage as in Fig. 7. Here, $\Delta_{0}=2~{}meV$ and $\gamma=0.3~{}meV$ As shown in Fig. 7 (a), in the case of a “clean” parent superconductor the features are particle-hole symmetric and strongly suppressed at large bias voltages, as discussed in the context of Fig. 1. By contrast, a finite density of sub-gap states in the parent superconductor introduces particle-hole asymmetry, as shown in panel (b) of Fig. 7. Also note that the particle-hole asymmetric “stripy” features have similar slopes at positive and negative values of $V_{Bias}$ and are not suppressed at large bias voltages. These types of features are present in the experimentally measured differential conductance of SM-SC hybrid structuresChen et al. (2016), suggesting that parent superconductors such as NbTiN may have a non-vanishing density of sub-gap states.

We note that these results are consistent with the existence of some generic dissipationLiu et al. (2017). Here, we propose the presence of a small but finite density of sub-gap states in the parent superconductor as the natural mechanism responsible for dissipation effects. Thus, dissipation is just one among several important low-energy effects that are directly linked to the parent superconductor, which can be captured theoretically by explicitly incorporating the bulk SC into the modeling. Two relevant examples – the dependence of the quasiparticle spectral weight within the SM wire on the SM-SC coupling and the renormalization of the low-energy spectrum – are discussed below. Explicitly incorporating the parent SC into the theory provides a unified framework for describing these effects, which otherwise require different ad hoc ingredients. This suggests that a thorough investigation of the parent superconductor (in the presence of disorder and finite magnetic fields) is absolutely necessary for understanding in detail the low-energy physics of the hybrid structure.

Figure 8 shows differential conductance curves in the topological regime for three different values of the ratio between the SM-SC coupling and the parent SC gap. The induced gap is kept constant. Notice that the height of the ZBP decreases with increasing $\gamma/\Delta_{0}$. This dependence is the result of the dissipation effect generated by the non-vanishing density of subgap states in combination with the reduction of the spectral weight within the SM wire associated with the Majorana mode. More specifically, as the (relative) coupling strength $\gamma/\Delta_{0}$ increases, the spectral weight of the low-energy states (including the Majorana mode) is transfered from the SM wire to the parent SC.Stanescu and Sarma (2017) Consequently, for a given value of the tunnel barrier transparency, the weight of the corresponding differential conductance peak decreases as the coupling strength increases, as illustrated in Fig. 8. The shape of the ZBP will be studied in more detail in the next section.

The dependence of the ZBP height on the coupling strength cannot be captured using an effective pairing potential model with an ad hoc dissipation term. To emphasize the importance of explicitly including the parent SC, we compare the dependence of the low-energy spectrum on the applied Zeeman field for i) a 1D wire with an effective pairing potential (model I) and ii) a 1D system that explicitly includes the parent SC (model II). The length of the SM wire is $L=1.5~{}\mu$m and the SM-SC coupling strength is $\gamma=0.447~{}$meV. To better simulate the experimental conditions, we assume that the bulk SC gap $\Delta_{0}$ is suppressed by the applied magnetic field. Explicitly, we have $\Delta_{0}(\Gamma)=\widetilde{\Delta}_{0}\sqrt{1-\Gamma/\Gamma^{*}}$, with $\widetilde{\Delta}_{0}=0.3~{}$meV and $\Gamma^{*}=1.35~{}$meV. The corresponding induced gap (defined as the minimum quasiparticle gap at zero field) isStanescu and Sarma (2017) $\Delta_{\rm ind}=0.2~{}$meV. Since the induced gap is one of the most robust and experimentally accessible features, we consider an effective pairing model characterized by the same value of the induced gap and a field-dependent pairing of the form $\Delta(\Gamma)=\Delta_{\rm ind}\sqrt{1-\Gamma/\Gamma^{*}}$. Finally, we assume that the parent SC has sub-gap states corresponding to $\delta=12~{}\mu$eV. A corresponding dissipation term is included in the effective pairing potential model. We expect $\delta$ to be strongly dependent on the applied magnetic field, but we do not incorporate this effect.

The dependence of the density of states within the wire on the Zeeman field is shown in Fig. 9, while the signature of the Majorana mode in the differential conductance is illustrated in Fig. 10. Several observations are warranted. First, the two models predict different values of the critical field corresponding to the emergence of the Majorana mode, despite the fact that they are characterized by the same value of the induced gap at $\Gamma=0$. Indeed, assuming that the chemical potential is $\mu=0$, the effective pairing approximation (model I) predicts a critical field given by the induced gap, $\Gamma_{c}=\Delta_{\rm ind}$, while incorporating the parent SC (model II) gives a critical field controlled by the SM-SC coupling strength,Stanescu and Sarma (2017) $\Gamma_{c}=\gamma$. Of course, one can construct effective pairing models that give the same value of the critical field as model II, but the corresponding quasiparticle gaps will be significantly overestimated in the intermediate and strong coupling regimes, i.e. for $\gamma>\Delta_{0}$. Second, we note the strong energy renormalization that leads to the suppression of the Majorana splitting oscillations in model II. One can still have Majorana splitting in the intermediate/strong coupling regime, but it requires very short wires. Note that in model I the collapse of the SC gap at $\Gamma^{*}$ has no visible impact on the amplitude of the oscillations, which is increasing with the Zeeman field (see top panels in Figs. 9 and 10). On the other hand, in model II the collapse of the SC gap results in a reduction of the ZBP height. This is another manifestation of the dependence of the quasiparticle spectral weight on $\gamma/\Delta_{0}$. Since $\Delta_{0}$ is field-dependent, the relative SM-SC coupling strength increases with $\Gamma$, which results in a decrease of the spectral weight of the Majorana mode reflected in the conductance. We emphasize that in these calculations the chemical potential is constant, hence the Fermi k-vector and the corresponding Fermi velocity of the Majorana band increase with the applied field. Consequently, the transparency of the barrier increases with $\Gamma$. In model I, this results in a ZBP that becomes stronger with $\Gamma$; i.e., its weight increases with the field and so does its height (if it is below the quantized value at finite $\delta$ and/or $T$) whenever the wire is long-enough (to avoid energy splitting). In model II, on the other hand, the effect of a higher barrier transparency is offset by the loss of spectral weight, the net result being the eventual decrease of the ZBP strength upon approaching $\Gamma^{*}$. Finally, we note that the dissipation-induced particle-hole asymmetry responsible for the “stripy” features in Fig. 7 is present even at very low energies, as illustrated in the top panel of Fig. 10. This effect could explain recent experimental observationsChen et al. (2016) showing similar low-energy features.

III.4 The shape of the zero bias conductance peak

We conclude our analysis with a summary of the dependence of the Majorana-induced zero-bias conductance peak on the relevant system parameters. In this section we will study the ZBP of a hybrid system consisting of a very short nanowire of length $L=0.2~{}\mu$m coupled to a superconductor with $\Delta_{0}=2~{}$meV, the effective coupling strength being $\gamma=0.25~{}$meV. The Zeeman field is set to $\Gamma=0.3~{}$meV. We note that having a short wire ensures that the Majorana mode is separated by a large gap from finite energy excitations, so that the area associated with the zero-bias conduction peak can be accurately estimated. Figure 11 shows the Majorana-induced zero-bias peak (ZBP) at $T=0$ for two different values of the potential barrier in a system with either a clean parent superconductor (solid lines) or a superconductor with a finite density of sub-gap states (dashed lines). First, we note that in the clean system the ZBPs are quantized at $2e^{2}/h$ and that their widths (and, implicitly, the areas under the corresponding $dI/dV$ curves) are reduced by increasing the barrier potential (i.e., reducing the transparency of the barrier). Second, in the structure with a “dirty” parent superconductor (dashed lines in Fig. 11) the height of the ZBP is significantly lower than the quantized value. In fact, the effect of sub-gap states being present in the superconductor is to introduce a broadening of the Majorana mode. The area under the $dI/dV$ curve is the same as in the clean system (for a given potential barrier), but the width of the peak is intrinsically finite. In other words, increasing the barrier height generates a thin quantized peak in the clean system, but results in a small wide peak in the presence of sub-gap states (see solid and dashed orange/light gray lines in Fig. 11). Intuitively, we can understand the broadening of the ZBP as a result of the Majorana mode hybridizing with the sub-gap states from the parent superconductor. The resulting hybrid states will have energies distributed within a certain $\delta$-dependent energy window about zero. As a result, the LDOS at the end of the wire will be characterized by a finite width peak centered at $E=0$.Stanescu and Sarma (2017) By contrast, in a clean system the corresponding signature of the Majorana mode is a delta function-type contribution to the LDOS. Note that the broadening effect of the sub-gap states is qualitatively similar to the effect of finite temperature.

The dependence of the zero-bias conduction peak on the transparency of the barrier, temperature, and density of sub-gap states is shown in Fig. 12. First, we note that the transparency of the barrier depends on both the height and width of the barrier potential. For simplicity, we model the potential barrier as a Gaussian, $V_{b}(x)=V_{0}\exp[-(x-x_{0})^{2}/2\sigma^{2}]$. Hence, the height of the barrier is parametrized by $V_{0}$ and its width by the standard deviation $\sigma$. Increasing $V_{0}$ or $\sigma$ reduces the transparency of the barrier, which results in a smaller area under the $dI/dV$ curve [see panel (a)]. Note that in a clean system at $T=0$ the height of the ZBP is quantized for all values of $V_{0}$ and $\sigma$. For a given tunnel barrier potential, introducing finite temperature or finite density of sub-gap states (in the parent SC) does not modify the area of the ZBP, but broadens the peak. Consequently, its height gets reduced. The dependence of the ZBP height on $T$ and $\delta$ for different potential barriers is shown in panels (b) and (c), respectively. Note the similarities between the dependence of the ZBP height on temperatureLin et al. (2012); Liu et al. (2017) (panel b) and on the density of subgap states (panel c).

In the high temperature limit, which is barrier-dependent, see Fig. 12 (b), the dependence of the height of the ZBP on temperature is given by the analytical expression $A/(8T)$, where $A$ is the area of the ZBP. Similarly, in the limit of high density of sub-gap states, the dependence of the height of the ZBP on $\delta$ is given by $A\pi/(2\delta)$. As the zero-bias peak approaches its quantized value with lowering $T$ (or $\delta$), the dependence of its height on temperature (or on the density of sub-gap states) becomes more complicated. We note that all the curves in panel (b) can be collapsed into a single curve corresponding to a ZBP of area $A_{0}$ through the scaling transformation $T\rightarrow\frac{A_{0}}{A}T$, where $A$ is the area under the ZBP corresponding to an arbitrary curve. Similarly, any two curves in panel (c) are related to one another by the transformation $(\delta\rightarrow\frac{A_{2}}{A_{1}}\delta)$ and can be collapsed into a single curve. However, the transformation $(T\rightarrow 4\pi\delta)$ mapping a curve with area $A$ in panel (b) to a curve in panel (c) with the same area under the ZBP only works in the high-temperature/high density of sub-gap states regimes.

We note that all transport experiments on semiconductor-superconductor hybrid structures reported so far have observed zero-bias peaks having heights significantly smaller than the (theoretically predicted) quantized value, $2e^{2}/h$. These small values are inconsistent with expectations based on the nominal temperature of the system. Of course, the effective temperature could be much higher. However, the presence of sub-gap states is a natural alternative explanation for the low values of the observed zero-bias conduction peaks.

IV Conclusions

We have calculated the tunneling conductance of a normal metal - semiconductor wire - superconductor structure within the Blonder-Tinkham-Klapwijk and Keldysh formalisms by treating the parent superconductor as an active component of the hybrid system. This treatment ensures that several key effects that control the low-energy physics of the hybrid structure are accurately captured by the theory. These effects include the proximity-induced renormalization of the low-energy states, including the renormalization of the topological gap that protects the Majorana zero modes,Stanescu and Sarma (2017) the emergence of additional conductance peaks at energies corresponding to the gap edge of the parent superconductorCole et al. (2015); Reeg and Maslov (2017); Stanescu and Sarma (2017), and the broadening of the low-energy features due to hybridization of the states within the semiconductor wire with sub-gap states in the parent superconductorStanescu and Sarma (2017).

We find that the differential conductance is a good indicator for the local density of states at the end on the hybrid structure adjacent to the tunnel barrier. This includes both the semiconductor wire and the nearby parent superconductor. A clearly defined feature associated with the parent superconductor emerges at energies corresponding to the bulk gap. In the intermediate/strong coupling regime the features associated with sub-gap states (e.g., the induced gap) and those corresponding to contributions from the parent superconductor have similar characteristic energies. We predict that one can disentangle them by varying the strength of the potential barrier, which has a stronger effect on the sub-gap features. For example, in the high barrier limit one can almost completely suppress the sub-gap features while having a measurable parent superconductor contribution. In the topological regime, varying the transparency of the tunnel barrier changes the area of the zero-bias conductance peak, while its height is quantized ($2e^{2}/h$) at zero temperature. Finite temperature reduces the height of the ZBP but preserves its area.

The most dramatic implication of treating the parent superconductor as an active component of the system occurs if we assume that the superconductor itself has sub-gap states induced by disorder and finite magnetic fields. These states hybridize with the low-energy states in the semiconductor wire significantly altering their properties. We find that in the presence of sub-gap states in the parent superconductor the differential conductance develops characteristic “stripy” features (as function of the chemical potential and bias voltage) that break particle-hole symmetry. In addition, the presence of sub-gap states has an effect similar to finite temperature and destroys the quantization of the Majorana-induced zero-bias conductance peak. More generally, the presence of these states results in a broadening of the low-energy features that characterize the hybrid structure. From the perspective of quantum computation, this is a deadly effect, as it destroys the topological gap that protects the low-energy Majorana subspace. All the signatures identified here as being associated with the presence of sub-gap states in the parent superconductor can be recognized in the transport measurements on semiconductor-superconductor hybrid structures reported so far. Therefore, we find it imperative to understand in detail the low-energy properties of currently-used parent superconductors, such as Al and NbTiN, and to consider a systematic effort to optimize their properties or to engineer better suited parent superconductors.

In summary, we have shown that explicitly treating the parent superconductor as an active component of the semiconductor-superconductor structure generates a number of features in the calculated differential tunneling conductance that would otherwise be absent. Some of these features have not been previously discussed in any detail. One of the new features is the reduction in weight of the differential conductance peaks as the coupling to the parent superconductor increases, as shown in Fig. 2. Since the reduction (and eventual collapse) of the parent SC gap at high magnetic fields implies an increase of the effective coupling (relative to the SC gap), a consequence of this feature is the decrease of the ZBP strength at high Zeeman fields, as illustrated in Fig. 10. Another feature is the presence of a superconducting resonance at the energy scale corresponding to the bulk gap (see Fig. 3) with an amplitude that has a dependence on the tunnel barrier height different from that of sub-gap (induced) peaks, as shown in Fig. 5. In addition, the possible presence of sub-gap states in the parent superconductor results in particle-hole asymmetric “stripy” features in the differential conductance, like those illustrated in Fig. 7, that persist even at very low energy, as shown in Fig. 10. Finally, explicitly considering the parent superconductor can strongly suppress the Majorana energy splitting oscillations and generate a finite density of states inside the topological gap, as shown in Fig. 10. Considering that all these features, some of them new, which were calculated within a unitary framework that explicitly incorporates the superconductor, have been observed experimentally in tunneling measurements on semiconductor-superconductor structures, we conclude that our work provides strong evidence regarding the critical role of the parent superconductor in Majorana devices.

This work is supported by NSF DMR-1414683.