Tunneling transport in NSN junctions made of Majorana nanowires across the topological quantum phase transition

We theoretically consider transport properties of a normal metal (N)- superconducting semiconductor nanowire (S)-normal metal (N) structure (NSN) in the context of the possible existence of Majorana bound states in disordered semiconductor-superconductor hybrid systems in the presence of spin-orbit coupling and Zeeman splitting induced by an external magnetic field. We study in details the transport signatures of the topological quantum phase transition as well as the existence of the Majorana bound states in the electrical transport properties of the NSN structure. Our theory includes the realistic nonperturbative effects of disorder, which is detrimental to the topological phase (eventually suppressing the superconducting gap completely), and the effects of the tunneling barriers (or the transparency at the tunneling NS contacts), which affect (and suppress) the zero bias conductance peak associated with the zero energy Majorana bound states. We show that in the presence of generic disorder and barrier transparency the interpretation of the zero bias peak as being associated with the Majorana bound state is problematic since the nonlocal correlations between the two NS contacts at two ends may not manifest themselves in the tunneling conductance through the whole NSN structure. We establish that a simple modification of the standard transport measurements using conductance differences (rather than the conductance itself as in a single NS junction) as the measured quantity can allow direct observation of the nonlocal correlations inherent in the Majorana bound states and enables the mapping out of the topological phase diagram (even in the presence of considerable disorder) by precisely detecting the topological quantum phase transition point.


I. INTRODUCTION
The subject of topological superconductors (SCs) hosting non-Abelian quasiparticles has become one of the most intensively investigated topics in condensed matter physics. 1,2 In particular, one-dimensional topological superconductors have been predicted to support zeroenergy particle-hole symmetric non-Abelian Majorana bound-states (MBS) localized at the ends. 2 Beyond their intrinsic fundamental interest, MBS have attracted attention for their potential use in fault-tolerant topological quantum computation schemes. 3 Far from being a subject of purely theoretical interest, concrete experimental proposals to realize these exotic states of matter have been put forward recently, 4-8 some of which have been implemented experimentally. [9][10][11][12][13] In particular, Refs. [6][7][8] showed that a one-dimensional semiconductor nanowire in proximity to a bulk s-wave superconductor, and subjected to strong Rashba spin-orbit coupling can be driven into a topologically non-trivial phase with MBS localized at the ends, upon the application of an external Zeeman magnetic field. In this topologically non-trivial phase, the nanowire becomes effectively a helical spinless p-wave superconductor, realizing an idea originally proposed by Kitaev for the localization of isolated MBS in a physical system. 2 Other experimental setups involving arrays of magnetic atoms on s-wave SCs, 14 or coldatomic systems 15 have also been proposed and are currently under experimental consideration. It is important to mention here that the real significance of course is the creation of isolated zero-energy MBS at the ends of the nanowire which are well-separated from each other so that they can be considered topologically protected.
On the experimental side, one of the most relevant questions is how to establish the presence of "true" MBS in a real experiment. In principle, the tunneling conductance at the end of the topological SC nanowire should reveal an MBS as a quantized zero-bias peak (ZBP) of magnitude 2e 2 /h in the conductance at zero temperature, which is a direct manifestation of the perfect Andreev reflection associated with the MBS. 6,[16][17][18] Recent experiments implementing the proposal in Refs. [6][7][8] have shown an intriguing ZBP, in apparent agreement with theoretical predictions for the existence of MBS, which appears upon application of a Zeeman field, providing compelling preliminary evidence of the Majorana scenario. [9][10][11] However, the interpretation of these experiments seems to be considerably more complex than the ideal models originally proposed and show several deviations from the predicted behavior, among which we mention the important ones: a) the smallness of the ZBP in comparison to the ideally theoretical value of 2e 2 /h (i.e., 0.1 − 0.2 e 2 /h in the low temperature limit), b) the presence of a continuum of fermionic excitations in the subgap region (i.e., the so-called "soft-gap" feature) instead of a well-defined SC gap, and c) the lack of evidence for the closing and then reopening of this SC softgap upon increasing the Zeeman field across the putative critical field V c . We stress that this last feature is a crucial prerrequisite for the existence of MBS, which would be indicative of a topological quantum phase transition (TQPT) taking place in the sample where the gap must vanish.
Contrasting with the interpretation that the recent nanowire experimental observations are indeed evidence for the isolated existence of MBS in a topological SC system, it has been pointed out that other ZBPs (or near-ZBPs) sharing similar features with the MBS are generically allowed in spin-orbit-coupled nanowires subject to a magnetic field in the presence of disorder or smooth confining potentials, both in the topologically trivial and non-trivial phases, a fact that would hinder the observation of bona fide Majorana-type excitations. [19][20][21][22][23] In particular, disorder is known to have strong detrimental effects in p-wave SCs. [24][25][26][27][28][29][30][31][32][33][34][35][36][37] Motrunich et al. showed more than a decade ago that Andreev subgap states induced by disorder tend to proliferate in one-dimensional systems described by Bogoliubov-de Gennes Hamiltonians with broken time and spin-rotational symmetry (symmetry class D, like the nanowires in Refs. [6][7][8], and render the system gapless. 24 These authors predicted that for weak disorder an infinite system realizes a topologically non-trivial phase with two degenerate zero-energy MBS localized at the ends of the wire. In a finite-length system of size L w , this degeneracy is lifted by an exponential splitting ∆ε ∼ e −Lw/ξ , where ξ is the superconducting coherence length. Increasing the amount of disorder generates low-energy Andreev bound states, and the (averaged) scaling of the splitting energy changes to ∆ε ∼ e −Lw/ξ+Lw/(2 e ) , where e is the elastic mean-free path of the system. 27 Beyond a critical disorder amount, defined by the condition e = ξ/2, the system experiences a TQPT induced by disorder and enters a nontopological insulating phase with no end-MBS. 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 intimate connection between localization and topology in disordered topological superconductors has been stressed in a series of the-oretical works. 24,30,31,35 The interplay among disorder, superconductivity, and possible Majorana zero modes is still very much an important open problem in the subject, and whether the experimentally observed ZBP is indeed the manifestation of the theoretically predicted MBS can only be sorted out definitively by accurately understanding the precise role of disorder in the experimental systems. In particular, a key question is the effect of disorder on the TQPT itself, which is a central topic of the current work.
Concerning the rather ubiquitous presence of in-gap states (or "soft gap") in the experiment, it is important to note the lack of evidence of a well-defined superconducting gap even in absence of an applied magnetic field, when the time-reversal symmetry is not broken (i.e., symmetry class DIII or BDI). The soft gap associated with the proximity effect in the nanowire induced by the swave superconducting substrate is one of the most puzzling problems in this subject which transcends the issue of TQPT or MBS since it exists even in zero magnetic field. The soft gap may or may not be directly tied to disorder (i.e. disorder may not be necessary) although strong disorder will undoubtedly produce some soft gap type features feature in the proximity-induced superconductivity. From the theoretical point of view, it is extremely important to understand the physical origin of this soft gap for the correct modeling and description of the experiment (as well as to help produce hard gap systems for future Majorana experiments). Since the topological protection of the MBS is directly provided by the existence of the SC gap, it is of obvious importance to understand the soft gap so that future experiments can realize hard gap nanowire proximity effect. Stanescu et al. have suggested recently that intrinsic quasiparticle broadening effects due to the hybridization of the SC with the normal metallic lead, could explain this feature. 38 Indeed, it is well-known that a highly transparent NS barrier can produce large subgap conductance 39 [i.e. Blonder-Tinkham-Klapwijk (or BTK) barrier parameter Z → 0], and therefore could also induce a large broadening of Bogoliubov quasiparticles, and hence a soft gap through this "inverse proximity effect" of the normal metallic lead on the SC nanowire. However, this does not seem to be the complete explanation of the experiment. Recent experiments where the transparency of the NS contact was systematically reduced have shown that the soft gap persists even in the lowtransparency limit (i.e., "pinching off" the quantum point contact when the inverse proximity effect should be exponentially suppressed). 13 An alternative explanation for the soft gap, valid in the limit of low transparency (i.e., large BTK barrier parameter Z → ∞), was proposed in Ref. 40. Among the many different pair-breaking mechanisms that might be operative in Majorana nanowires as considered in Ref. 40 (e.g., finite temperature T , presence of magnetic impurities, quasiparticle broadening, etc.) realistic parameter considerations point to the predominance of a special kind of inhomogeneity, which was not considered before in the present context: the spatial fluctuations in the proximity-induced pair potential ∆ (x). 40 Physically, spatial fluctuations in ∆ (x) are likely to be introduced by disorder or inhomogeneities at the SC/semiconductor contact. ? By improving the quality of the semiconductor/SC interface using molecular beam-epitaxy methods, recent experiments have reported much harder gaps, 41 following the suggestion of the theoretical explanation given in Ref. 40. This constitutes a qualitative improvement in the fabrication of Majorana nanowires, and hopefully a new generation of experiments where disorder effects are dramatically reduced will be soon available with hard proximity gaps (i.e. no subgap fermionic excitations) and well-defined MBS. We incorporate this aspect of the soft gap physics in the current work through a simple model approximation which mimics the spatial variation in the proximity-induced superconducting pair potential arising from the inhomogeneities at the superconductor-nanowire interface [see Eq. (4) below and the associated discussion].
The above discussion describes the rather complex situation faced in the experiments in order to detect "true" MBS in the topological phase. In this article we focus on a specific configuration, the normal-topological superconductor-normal (NSN) configuration, which is currently under experimental study. The SC part of this NSN (i.e. the S-part) junction is the semiconductor nanowire which has proximity-induced superconductivity from an underlying ordinary s-wave SC system. Many of the recent experiments have focused specifically on just the simple NS junction, but NSN junctions are essentially "equally easy" to study, and they have been studied also. We believe that NSN junctions have some intrinsic advantages over the minimal NS junction transport for studying MBS physics and the associated TQPT. We provide a comprehensive theoretical analysis of its transport properties taking into account the effects of disorder, inhomogeneities and temperature. As noted in previous works, the NSN configuration allows to extract the same information as in the simpler NS contacts, but contains additional interesting new physics arising from non-local correlations. [42][43][44][45] In this work, in addition to describing the transport properties of NSN disordered Majorana wires, we expand the idea first suggested in Ref. 35, on how to extract information about the TQPT in a realistic disordered sample. While this is not a "smoking-gun" experiment, it might be an extremely useful experimental tool providing information about the topological phase diagram of the system, complementary to non-local shot noise correlations. [42][43][44][45] Observation of non-local correlations as well as studying the TQPT itself using our suggested transport techniques in NSN junctions taken together may in fact serve as the smoking gun evidence for the existence of Majorana modes in nanowire systems. An associated significant advantage of the NSN junctions over the much-studied NS junctions in the context of Majorana physics in nanowires, which should be obvious from the above discussion and is emphasized throughout this work, is that transport in NSN junctions potentially studies both the end MBS and the bulk topological SC phase whereas NS junction tunneling properties may very well be dominated by the end MBS so as to suppress the manifestation of the bulk TQPT and the nonlocal correlations between the two end MBS which must go through the bulk nanowire. This is the key reason for our promoting NSN junction transport studies as an important tool for the Majorana investigation. In order to illustrate the main motivation of this article, let us first consider a "dirty" proximity-induced SC Majorana nanowire, electrically connected to ground and attached to normal contacts in a NSN configuration, as shown schematically in Fig. 1(a). Here we consider a generic situation where inhomogeneities are present both in the form of spatial fluctuations of the (proximity-induced) pairing potential, and in the form of quenched disorder in the on-site chemical potential fluctuations [Figs. 1(b) and (c)]. We also assume an external Zeeman field applied in the direction parallel to the nanowire, which allows to drive the system across the TQPT. A relevant experimental quantity is the differential conductance matrix, defined as where I i and V j are, respectively, the current and voltage applied in the {i, j} = {L, R} normal contact. In ideal conditions, the local conductances G LL and G RR should reveal the presence of end-MBS as a quantized ZBP peak of magnitude 2e 2 /h at T = 0. 6,[16][17][18] In practice, however, disorder, finite temperatures, quasiparticle poissoning, etc., might hinder or even destroy the purported topological phases and, therefore, the MBS. Since we are motivated by the current experiments, we start by showing a typical example of our numerical simulations of tunneling transport in Figs. 2(a) and 2(b), and leave the explanation of the theoretical details for Secs. II,III and IV. In these plots have computed the local conductances G LL and G RR for a disordered wire at a finite temperature as a function of the local bias voltages V L and V R , respectively, and for different values of the applied Zeeman field. In contrast to the ideal case 46 (i.e., clean system and T = 0), where a vanishing singleparticle excitation gap signals the TQPT across the critical Zeeman field, with the ZBP emerging on the topological side at higher magnetic field, here the presence of the abovementioned non-idealities renders the situation much less clear to determine the TQPT and the nature of the ZBPs. In other words, the information about the ZBP has been "washed out" by a combination of thermal effects, disorder and quasiparticle broadening, although the conductance results in Fig. 2 are explicitly obtained theoretically in a system where the MBS definitively exists in the ideal situation. (As an aside, we mention that the theoretical conductance results depicted numerically in Fig. 1  The article is divided as follows. In Sec. II we present the theoretical framework, the model and the main approximations. In Sec. III we describe the method used to determine theoretically the topological phase diagram of a disordered Majorana wire. In Sec. IV we present the theoretical technique to describe the differential conductance of a generic disordered Majorana wire in the NSN configuration and analyze the physical content in the analytical expressions. In Sec. V we describe in detail a proposal to extract information about the TQPT and to assess the topological stability of MBS. Sec. VI is intented to provide a simple intuitive theoretical understanding the physics underlying our proposal, in Sec. VII we present a summary and our conclusions, and finally in Appendix A we give a detailed derivation of the Eqs. (11)- (14) for the conductance matrix in the SNS configuration.

II. THEORETICAL MODEL
In accordance with previous works on Majorana wires, 7,8,46 we consider the following Hamiltonian describing a disordered semiconductor nanowire of length L w , subjected to Rashba spin-orbit coupling and a Zeeman field, Here, ψ † σ (x) creates a fermion with spin projection σ, andσ i (with i = x, y, z) are the Pauli matrices acting on spin space. The parameter α R is the Rashba spin-orbit coupling strength and V x is the Zeeman field along the wire, and summation over repeated indices σ is implied. The term H ∆ represents the effect of a proximate bulk s-wave SC on the nanowire [not shown in Fig. 1(a)], which induces a mean-field SC pairing potential ∆ (x) through the proximity effect. For simplicity, we have assumed single-channel occupancy in the nanowire with no loss of generality. As we will explain later, our results are generic and this single-channel (or single-subband) assumption does not affect the main conclusions in the case of many occupied subbands (as long as an odd number of subbands are occupied which is a necessary condition for the existence of the MBS for many occupied subbands). We recall that H NW is only an effective onedimensional model describing the system at low temperatures. A more realistic model should involve an explicit coupling t ⊥ to the proximate bulk SC, which is the source of superconducting correlations, and a self-consistent determination of ∆ (x). However, this task is beyond the scope of this work and does not change our results qualitatively since all we need in our model is the existence of a pairing potential in the nanowire. For more details, we refer the reader to Refs. 46,47 where a deeper discussion on this issue is provided, which is not particularly germane for our consideration in the current work where we are interested in the realistic manifestation of the MBS themselves rather the issue of proximity effect.
Disorder and inhomogeneities enter in the above model through two physically different mechanisms: a) Local fluctuations of the chemical potential µ (x) = µ 0 +δµ (x), with µ 0 a uniform value which in principle can be controlled by external gates, and the fluctuations δµ (x) are physically related to the presence of impurities, vacancies, etc. in the environment (both the nanowire itself and the surrounding). We assume δµ (x) to be a Gaussian random variable fully characterized by δµ (x) = 0 and δµ (x) δµ (y) = υ 2 µ δ (x − y), with the standard deviation υ µ representing the "strength" of disorder [see the horizontal axis in Fig. 2(c)]. For one single realization of disorder, once the NW is deposited and electrically contacted, we assume this parameter to be fixed throughout the experiment. b) Local variations in the (induced) pair potential ∆ (x), which for concreteness (and numerical convenience) here we model as for 0 < x < L w , i.e., a smooth profile that vanishes at the ends of the nanowire. Here ∆ 0 is the value in the bulk (i.e., right next or beneath the bulk SC), and d ∆ is an adjustable parameter that controls the slope of the profile. As mentioned above, a more rigorous treatment of this mean-field Hamiltonian should involve a self-consistent determination of this profile, but for our present purposes this simplification is well justified. In contrast to Ref. 40, here we only consider the deterministic profile Eq. (4) and we neglect other random inhomogeneities in ∆ (x) introduced by disorder. More details on disorderinduced SC pairing potential fluctuations can be found in Ref. 40.
In the absence of disorder and in the uniform case (i.e., limit υ µ = d ∆ = 0), the Hamiltonian H in the limit L w → ∞, can be easily diagonalized in momentum space k. In that case, the dispersion relation for the Bogoliubov For given values µ 0 , ∆ 0 and α R , this model has a TQPT as a function of magnetic field V x (i.e., the Zeeman spin splitting) from a topologically trivial phase to a nontrivial phase with the appearance of MBS localized at the ends of the nanowire at the critical Zeeman field value V x,c = ∆ 2 0 + µ 2 0 , as originally shown by Sau et al. 5 In the presence of disorder and other spatial fluctuations of the parameters in the model, the critical field V x,c typically shifts to larger values and its value depends on the precise details of the disorder realization. 28,31,32,34,35 The determination of the critical field defining the TQPT is then non-trivial and has to be done numerically for a given disorder realization. This is the subject of the next section.

III. THERMAL TRANSPORT AND TOPOLOGICAL PHASE DIAGRAM OF A DIRTY MAJORANA NANOWIRE
Let us now focus on the topological phase diagram of the disordered Majorana nanowire. In order to make progress, we have discretized the Hamiltonian in Eqs. (2) and (3), and obtained a N −site tight-binding model with the lattice parameter a (see Ref. 46) where c † l,s , µ l and ∆ l are the discrete versions of ψ † σ (x), µ (x) and ∆ (x), respectively, and t = 2 /2m e a 2 is the effective hopping parameter. Here α = mα 2 R /2 is the corresponding Rashba coupling parameter in the tightbinding model. The first site at the left end corresponds to l = 1 and the final site at the right is l = N .
We consider a single distribution of µ l (disorder realization), and systematically vary its dispersion υ µ around the mean value µ 0 . As mentioned above, υ µ is not an experimentally tunable parameter, but it is useful and instructive to visualize the topogical phase diagram as a function of varying disorder. Presumably, a fixed disorder realization is closer to the experiment, where the semiconductor NW is in the mesoscopic regime, and it is not clear that disorder necessarily self-averages at the very low experimental temperatures. We mention that whether the experimental temperatures are low enough so that the system is not self-averaging (so that mesoscopic fluctuations are important as one goes from one sample) is currently not known for the Majorana experiments, and the issue of whether to ensemble average over disorder realizations or not for quantitative comparison with experiments remains open at this stage.
We compute the topological phase diagram of the isolated nanowire (i.e., in absence of the normal contacts) using the transfer-matrix approach 31,35 for the model Hamiltonian Eq. (5), as a function of the disorder strength υ µ and the external Zeeman field V x . Physically, the transfer matrix relates states in the left end to states in the right end of the wire. This statement can be made more precise introducing the Majorana basis c l,σ = (a σ,l + ib σ,l ) / √ 2, where the Majorana operators obey the anticommutation relations {a σ,l , a s,m } = {b σ,l , b s,m } = δ l,m δ σ,s and zero otherwise. In terms of these operators, a generic eigenmode Ψ of H NW satisfying the eigenvalue equation H NW Ψ = EΨ can be written as with real coefficients γ σ,l and η σ,l . At E = 0, defining the matrix κ = t −α α t and the vector of coefficients ψ l = (γ ↑,l , γ ↓,l ) T , the above eigenvalue equation can be written as 0 = κ † ψ l−1 + κ ψ l+1 + µ l ψ l , and from here we obtain the transfer equation where is the l−th transfer matrix relating the vectors ψ l+1 and ψ l−1 . Then, the full transfer matrix of the nanowire, from site l = 1 to site l = N , is simply given by M = N l=1 M l . The eigenvalues of M can be written as e ±N λn , where λ n are the (dimensionless) "Lyapunov exponents" of the system, 48 which represent the inverse of the localization length. The connection to localization properties are better understood recalling that the transmission probability the transmission eigenvalue corresponding to the n−th channel . 48 The connection between the localization and the topological properties of a "dirty" class D nanowire was made explicit by Akhmerov et al , 30 who obtained the topological invariant Q = sign( 2M n=1 tanh λ n ). These authors have shown that in the clean case this topological invariant actually reduces to the one derived by Kitaev which is given in terms of the Pfaffian of the Hamiltonian in momentum space. 2 Here we see explicitly that Q changes sign when one of the Lyapunov exponents vanishes and changes sign. This signals the TQPT. As discussed in Refs. 24,30, the TQPT of a class D SC corresponds to a delocalization point for zero-energy particles, i.e., one of the Lyapunov exponents λ n vanishes and changes sign at the TQPT inducing a "perfect" transmission probability T n = 1. Everywhere else in the parameter space the system is localized at zero energy, i.e., all λ n are finite. This crucial result will be addressed in detail in Section VI. For the moment, we can check that this idea also works in the clean case: for a clean NW, sufficiently close to the TQPT on the topological side, the MBS wavefunctions are localized within the SC correlation length is the effective SC quasiparticle gap controlled by the Zeeman field. The TQPT is reached at the critical field V x,c = ∆ 2 0 + µ 2 0 , where the quasiparticle gap ∆ (V x,c ) → 0 and the localization length ξ clean → ∞. When ξ clean L w , the MBS localized at opposite ends can "see" each other and overlap forming a Majorana "channel" that connects the left ing that the coupling to the leads can be controlled in situ experimentally [e.g., using pinch off gates (not shown here)], the system can be effectively disconnected from the right lead and turned effectively into an NS junction. and the right end. Therefore, for a clean system near the TQPT the smallest Lyapunov exponent is λ clean ∝ ξ −1 clean . Since the Majorana channel has equal contributions of electrons and holes at E = 0, the current sustained by electron-like states exactly cancels the current of hole-like states, and the total electric current vanish. Therefore, the perfectly quantized transmission coefficient T n = 1 occurring at the TQPT is physically related to the thermal conductance (and not to the electrical conductance). We will return to this point in Sec. IV.
In Fig. 2(c) we show a 2D color map of the thermal transmission coefficient T 1N for a dirty wire as a function of disorder "strength" υ µ and applied Zeeman field V x , fixing all other parameters (chemical potential, pairpotential profile, etc.) according to Table I. These parameters correspond exactly to those used in Figs. 2(a) and 2(b) for the same configuration of disorder potential, and each dot corresponds to each one of the curves in those figures. The blue regions indicate the points for which the transmission coefficient is close to the maximal value T 1N = 1, and therefore indicate the approximate location of the TQPT. Therefore, Fig. 2(c) allows to determine the phase boundary separating the topological from the non-topological region. Note that this boundary has an intrinsic width which scales as ∝ 1/N ∝ 1/L w . More precisely, as we will see in Sec. VI, the width corresponds to the Thouless energy v F /L w . For the parameters in Table I and in the absence of disorder, we estimate an upper bound L w /ξ 15, where we have used the estima-tion for the minimal value of the SC correlation length ξ clean = v F /∆ 0 20 a and L w = 300 a.
Contrasting Figs. 2(a), 2(b) and 2(c), we note that the last four curves correspond to dots in (c) which are closer to the topological phase boundary, where the topological protection is expected to be more fragile. This seems to be in agreement with the fact that Fig. 2(a) shows a splitting in the ZBP. The ZBP appearing in the corresponding curves in Figs. 2(b) and the apparent inconsistency with Fig. 2(a) (i.e., peaks not correlated) will be addressed and discussed in the next section. In constrast, the four central curves in both figure (a) and (b) show a more robust ZBP, which is consistent with the corresponding points in Fig. 2(c) located further from the boundary. As we see, this analysis showing all curves "side-to-side" is potentially helpful to interpret the experimental transport results. A natural question arises: is it possible to access the information in (c) experimentally? This will be the subject of the next sections.

IV. ELECTRONIC TRANSPORT PROPERTIES IN THE SNS CONFIGURATION
We now turn to quantities with more relevance to current experimental measurements. To that end, we introduce a term in the Hamiltonian describing the coupling to external normal leads [see Fig. 3(a)] where the term where t L(R) is the coupling to the left (right) lead and d † L(R)k,s is the corresponding creation operator for fermions with quantum number k and spin σ. The external leads are modeled as large Fermi liquids with Hamiltonian H lead,j = k,σ k d † j,k,σ d j,k,σ , where j = {L, R}. We assume that each lead is in equilibrium at a chemical potential µ j = eV j controlled by external voltages, and that the SC NW is grounded. The expression for the electric current flowing through the contacts is These formulas are standard 35,39,51,52 and we do not derive them here. The reader will find more details in the aforementioned references and in Appendix A. We have defined the normal reflection and transmission matrices (i.e., with subindex "ee" or "electron-electron") and the Andreev reflection and transmission matrices (i.e., with subindex "eh" or "electron-hole") where g r ls,ms (ω) and f r ls,ms (ω) are the normal and anomalous retarded Green's functions 53 in the NW respectively (see Appendix A). We have also defined the effective couplings to the leads where ρ 0 j (ω) is the density of states in the j−lead. Assuming a large bandwidth in the normal contacts, in the following we set ρ 0 j (ω) = ρ 0 j (0), the value at the Fermi level.
Let us analyze the physical content of Eqs. (11)- (14). We first focus on the "local" conductances Eqs. (11) and (14). In these expressions, the first term corresponds to the local contribution 2Tr r LL eh r LL eh † ω and 2Tr r RR eh r RR eh † ω , i.e., the Andreev reflection probability at the left and right lead, respectively. These terms are the only terms appearing in the case of NS or SN contacts, 17 and they already contain the information about the presence of a MBS localized at the corresponding end. From this perspective, the quantized value of the conductance 2e 2 /h corresponds to a "perfect" Andreev reflection Tr r LL eh r LL eh † ω=0 = 1 at T = 0, due to the presence of the MBS. However, note that in Eqs. (11) and (14) we also encounter a non-local contribution , which physically corresponds to particles that travel from one to the other end of the wire, and return to the original lead with information about the opposite lead. These nonlocal terms are present only in the NSN junctions, and not in the simple NS configuration so far studied extensively in the literature. Our primary motivation for considering NSN junctions (with the 'S' part being the NW carrying MBS under suitable conditions) is to study the effect of these nonlocal terms in the transport experiments, since nonlocality is the key concept underlying MBS in the topological phase. Therefore, this contribution must be proportional to the electron-electron and electron-hole transmission coefficients, Tr t LR ee t LR ee † ω and Tr t LR eh t LR eh † ω respectively, and vanishes if either γ L or γ R vanishes by, e.g., "pinching off" one of the quantum point contacts using underlying gates (i.e., pinch off gates) [see Fig. 3(b)]. Note that the presence of such a non-local contribution is expected in multi-terminal phase-coherent mesoscopic systems. 54 Interestingly, the thermal conductance G th across the wire 54-56 where G th,0 = π 2 k 2 B T /6h is the thermal quantum of conductance, is closely connected to the non-local contribution in Eqs. (11) and (14), and allows to make a link with our previous discussion in Sec. III. The connection with the thermal transmission probability T 1N at zero . In principle, the information about the TQPT is contained in this expression. However, such thermal measurements are in general very challenging experimentally, and we need to come up with a different approach which is experimentally feasible. In particular, it is desirable to use electrical measurements (i.e., electrical conductance) for observing the nonlocal MBS correlations at the TQPT and beyond.
We now briefly discuss Eqs. (12) and (13) (i.e., the transconductances G LR and G RL , which obey G RL = G LR ), where a more explicit difference with respect to the NS geometry appears. Physically, the minus sign in these expressions appears because while electrons contribute with a plus sign to the transport, a hole will contribute a minus sign. As discussed previously in Sec. III, note that if the system is in the topological phase with end-MBS, particle-hole symmetry dictates that the contributions Tr t LR ee t LR ee † ω=0 and Tr t LR eh t LR eh † ω=0 must be identical, and therefore the transconductance must vanish. 30 This might seem to rule out the possibility to see the TQPT via electrical measurements of G LR . In Figs. 4(a) and 4(b) we reproduce the same Figs. 2(a) and 2(b), computed for the parameters in Table I, at a temperature T /∆ 0 = 0.02 (which approximately corresponds to the experimental temperature T exp ≈ 60 mK in Ref. 9), and where we have assumed effective couplings γ L = 0.85 t and γ R = 0.95 t, corresponding to an open wire condition (i.e., "good" electrical contact with the leads). In Fig. 4(c) we present a plot for G LR vs V R for the same parameters. Note that this quantity is rather featureless, and is vanishingly small near zero bias as expected. For comparison, in Figs. 4(d), 4(e) and 4(f) we show G LL , G RR and G LR , respectively, for the same parameters but for a much lower temperature T /∆ 0 = 2×10 −4 and smaller couplings γ L = γ R = 0.1 t . In these conditions the thermal and quasiparticle broad-enings decrease dramatically and we realize that the preliminary information about the ZBPs in Figs. 4(a) and 4(b) is misleading: the system does not have zero-bias excitations and the peaks are actually split (rather than being a single zero energy peak) in Figs. 4(d) and 4(e) at low temperatures and at low transparency of the contacts. This picture is actually consistent with Fig. 2(c), where the dots corresponding to the largest magnetic fields are very close to the topological phase boundary, and therefore we expect the MBS to recombine into Dirac fermions and consequently the peaks to shift away from zero bias voltage. This allows to interpret the uncorrelated ZBPs for G LL and G RR . The results in Fig.  4 shows that already for the simple model of Eq. (5), detection of a "true" Majorana ZBP based only on the information about the local conductances might be very tricky. 20 Therefore, the presence of ZBPs in the local conductances G LL and G RR cannot by itself be considered as a "smoking-gun" evidence of the Majorana scenario without some critical considerations of the correlations in the existence of these ZBPs arising from the two end conductances.

V. DETECTING TOPOLOGICAL PHASE TRANSISTIONS IN THE SNS CONFIGURATION
As mentioned before, in the case of clean wires, the TQPT should be observed in the closing and reopening of the gap of electronic excitations in the nanowire. This re-organization of the fermionic spectrum is necessary in order to accomodate a new MBS at zero energy. However, the experiments so far have been unable to report any definitive closing of the gap. It has been suggested that this negative results might originate because while the tunneling occurs at the end of the nanowire, the information about the gap-closing is contained in wavefunctions with the majority of the weight in the bulk of the nanowire. Therefore, measurements of the LDOS in the middle of the wire, 57 or capacitive measurements of the total DOS 58 should reveal this gap-closing ocurring at the TQPT, but no experimental evidence of these predictions have been reported so far.
In addition to characterizing the transport properties of disordered NSN Majorana wires, another goal of the present work is to explore experimental proposals to determine the topological phase diagram. We believe that the NSN geometry offers an interesting possibility to achieve this goal, and to provide information about the topological stability of the MBS. In Sec. III we stressed that the TQPT in Majorana wires corresponds to a delocalization point at zero energy, a fact that can be detected in the thermal transmission probability across the system. On the other hand, in Sec. IV we showed that in the NSN geometry, the local conductance of a phasecoherent Majorana NW depends on the non-local transmission probability Tr t LR  exactly corresponds to the (dimensionless) thermal transport at energy ω [see Eq. (16)]. This enables mapping out the topological phase diagram by purely electrical measurements. 35 In this section we provide more details in the way the non-local information could be extracted in the SNS configuration in order to obtain the topological phase diagram.
Following Ref. 35, we define the following quantity i.e., the difference of local zero-bias conductances computed for different values of couplings to the opposite lead, while keeping all other parameters fixed. The zerobias conductance at one end G jj (0) is computed for a given value of γj (with compact notationL = R and R = L) and G jj (0) is computed for a different value γ j . From the experimental point of view, this means using γ R and γ L as tuning parameters, something that could be achieved varying the pinch-off gates underneath the ends of the NW. [9][10][11][12][13]59 This constitutes a new experimental knob which has not been explored so far in the Majorana experiment. Note that this quantity [as defined by Eq. (17)], being a difference, is not quantized and can take either positive or negative values. For this reason in what follows we will take the absolute value.
A priori, it might seem counter-intuitive that the transport through a disordered medium could be influenced by the change of a boundary condition at the far-end. However, this intuition is typically built upon the more usual case of trivial Anderson-localized 1D systems, where any amount of disorder localizes the wavefunctions and therefore any object placed at distances larger than the localization length ξ loc has essentially no effect. The crucial difference with class D conductors is that ξ loc ∝ λ −1 n → ∞ at the TQPT (since the gap closes here), i.e., the de- localization point. Therefore, assuming that L w < L φ , where L φ is the phase-relaxation length, the aforementioned intuition is usually correct except at the TQPT. The physical idea behind using ∆G jj (0) as an indicator of the TQPT can be seen quite simply in the extreme case when γ j = 0. In this case, in Eqs. (11) and (14) for G LL and G RR respectively, one completely suppresses the coupling to the opposite lead and the transmission coefficients vanish. The remaining part (i.e., Andreev reflection 2Tr r jj eh r jj eh † ω ) is a purely local contribution.
Therefore, Eq. (17) must correspond to a non-local contribution which contains information about the TQPT. This statement is not entirely correct because modifying the coupling γj to γ j also modifies the local Andreev reflection coefficient through the local anomalous Green's functions f r 1s,1s (ω), which contains information about the entire system. Only in the perturbative limit where δγj ≡ γ j − γj γj, one can show 35 that at the lowest order in δγj, Eq. (17) becomes i.e., proportional to the thermal transmission.
In Fig. 5 we show |∆G LL (0)| as a function of V x for the same parameters as before (see Table I), for υ µ /∆ 0 = 2 and for γ R = 0.94 t and γ R = 0.01 t, and γ L = 0.85 t. The black solid line corresponds to the experimental temperature T exp /∆ 0 = 0.02, while the gray solid line corre-sponds to a temperature T /∆ 0 = 2 × 10 −3 (one order of magnitude smaller). We also show a vertical cut of T 1N [corresponding to υ µ /∆ 0 = 2 in Fig. 2(c)] as a red dashed line. That curve indicates the location of the TQPT (i.e., whenT 1N ≈ 1) in a theoretically isolated wire, which occurs at V x,c /∆ 0 ≈ 5.2 (indicated by a blue dot in the horizontal axis in Fig. 5). For V x > V x,c , strictly speaking the system remains in the topologically non-trivial phase, but the strong fluctuations of T 1N indicate a very fragile topological protection of the Majorana modes.
Comparing the black and gray lines we can see that thermal effects dramatically decrease the magnitude of ∆G jj (0), as can be seen in the overall reduction of the signal when the temperature increases from T /∆ 0 = 2 × 10 −3 to T /∆ 0 = 2 × 10 −2 . This is particularly true near the critical field V x,c , where the signal drops to |∆G jj (0)| 3 × 10 −3 e 2 /h . Although this value is still experimentally detectable, it would be desirable to have a stronger signal. In practice, this proposal is expected to work best for short wires, where the maximal ratio L w /ξ is not too large (in Fig. 5, using the paremeters in Table I, we estimate an upper bound L w /ξ ≈ 15). This is because the visibility of the electrical signal crucially depends on the width v F /L w of the peak in |∆G LL (0)|. Therefore, a very narrow peak might be hard to detect, or could be washed away by thermal effects or other dissipative mechanisms not considered here. Using a shorter nanowire should yield a stronger signal, albeit at the cost of less resolution. We stress the fact that system should be shorter than the phase-relaxation length L w < L φ for the two end-MBS to hybridize coherently. Despite these limitations, these conditions are very well met in the experiment. The fluctuations at higher fields V x > V x,c are still visible at T /∆ 0 = 2 × 10 −2 , and closely follow the fluctuations in T 1N . Physically, detecting such fluctuations in |∆G LL (0)| in an experiment would be an indication of a fragile topological protection.
Overall, we conclude that |∆G jj (0)| is a bona fide indicator of the TQPT and the topological stability of the MBS at low enough temperatures. Although we have suggested using the pinch off gates as a physical way to effectively tune the coupling to the normal contacts, this is not necessarily the only way to change the parameter γ jj . For instance, schemes using quantum dots (QDs) between the normal contact and the Majorana NW 60-62 (i.e., N-QD-S-QD-N setups) will also serve the same purpose. In this case, it would be relatively easy to modify the transparency of the coupling to the lead by changing the gate voltages in each QD. However, the QDs should be large enough to avoid strong Coulomb effects, which might introduce unwanted effects (e.g., Kondo effect 63 ) complicating the experimental interpretation. We finally mention that in Ref. 30, an alternative method to detect the TQPT based on the measurement of the current shot noise was proposed, which would be a complementary to the idea discussed here.

VI. INTUITIVE THEORETICAL PICTURE
In this section we provide a simple theoretical framework to interpret our numerical simulations in previous sections. To that end we focus on a simplified version of the Hamiltonian H NW in Eqs. (2) and (3), which will allow us to obtain an exact solution, therefore providing a valuable physical insight, while retaining at the same time the relevant physics. These simplifications will not modify our main conclusions because they do not depend on the details of H NW itself, but on its symmetry class (i.e., class D in this case) which is a robust feature. Therefore, for the present purposes we assume a uniform chemical potential µ 0 = δµ (x) = 0. In this simplified model, disorder enters only through the inhomogeneous pair potential ∆ (x), which we now assume to be generic and not necessarily of the form (4).
It is simpler to start the analysis from the uniform case with periodic boundary conditions, where the band theory helps to visualize the relevant physics related to the TQPT occurring near the point k = 0, at the intersection of the spin-orbit coupled bands with different spin projection [see Fig. 6]. The modes at finite momentum ±k F are assumed to be gapped by the SC paring interaction (not shown in the picture), and decouple from the relevant sector at k = 0. Projecting the original fermionic operator around this point and linearizing the bands results in a helical liquid model described by the Hamiltonian 21 where ψ R ψ ↑ (x) and ψ L ψ ↓ (x) result from spinmomentum locking around k = 0 due to the spin-orbit interaction. We now introduce the Majorana basis in terms of which (18) splits into two independent modes where we have introduced the Pauli matricesτ i acting on LR space. The emergence of Majorana zero-modes can be easily seen by solving the eigenvalue equation for whose solution is We now define the zero-energy eigenmodes in terms of which the expression for a generic MBS localized at the origin (i.e., the left end of the wire) is Although the form of the eigenmodes (24) is more convenient for our purposes, we mention here that one can easily bring this expression into the more familiar form of the Jackiw-Rebbi soliton solution 64 applying a rotation along theŷ−axisR = e i π 4τ y , which transforms to the usual eigevectors of the operatorτ z . In order to ensure the existence of MBS we need to find normalizable solutions that decrease sufficiently fast as x → ∞ and that satisfy generic boundary conditions. For a wire of length L w , we can define the quantity which in the limit L w → ∞ corresponds to the Lyapunov exponent of the system at zero energy for the channel n. In terms of these quantities, note that there are two possible situations: 2 1. Both λ 1 and λ 2 have the same sign, in which case we need to choose either η (x) = a 1 χ + 1 (x) + a 2 χ + 2 (x) or η (x) = a 1 χ − 1 (x) + a 2 χ − 2 (x) in Eq. (25), the sign depending on which of the modes decays for x > 0. Since there are two decaying contributions allowed, we can satisfy generic boundary conditions at the origin. For instance, if the system is a trivial insulator for x < 0, then the boundary condition η (0) = (0, 0) T must be imposed. This is verified with a 1 + a 2 = 0. Other boundary conditions for open wires will be analyzed later.
2. The Lyapunov exponents λ 1 and λ 2 have different signs. Then Eq. (25) is a linear combination of spinors χ + and χ − . This makes it impossible to satisfy generic boundary conditions, except for accidental situations which are not protected against local perturbations. For instance, in our previous example of a vanishing boundary condition at the origin, the condition η (0) = (0, 0) T implies that a 1 + a 2 = a 1 − a 2 = 0, which can only be satisfied for a 1 = a 2 = 0. Therefore the MBS does not exist.
From this analysis, we conclude that the TQPT occurs when one of the Lyapunov exponents λ n passes through zero and changes sign, making explicit the connection between the localization properties of a D-class nanowire and its topological properties. This is a robust feature which is independent of the details of the microscopic Hamiltonian as it depends only on the symmetry classification. Assuming that the magnetic field V x is such that λ 2 > 0, the condition for the topological phase reduces to where we have defined the average gap∆ ≡

A. Open wires
We now assume that our wire is connected to conducting leads at both ends, and focus on the transport across the NSN configuration at zero energy, as depicted in Fig.  7(a). We define the scattering matrix 30 where ψ ν (x) ≡ (η 1,ν (x) , η 2,ν (x)) T are ν−moving (with ν = {L, R}) scattering Majorana states in the left and right leads (x = 0 and x = L w , respectively). In the Majorana basis, S 0 is a real orthonormal matrix (18), the modes n = {1, 2} are decoupled and independent, the transmission and reflection matrices, t 0 , t 0 and r 0 , r 0 respectively, acquire a diagonal form, an can be diagonalized independently with diagonal elements t 0,n , t 0,n and r 0,n , r 0,n . Without loss of generality, in what follows we assume that the only incident modes are right-moving modes arriving from the left lead. This imposes the boundary conditions [see Fig. 7 which in combination with Eqs. (23) and (28), allow to obtain closed analytical expressions for the reflection and transmission coefficients and we recover Eq. (9) for the transmission probability. Exactly at the TQPT, the determinant of the reflection matrix vanishes and a "Majorana channel" with perfect transmission opens at zero energy. As mentioned in Sec. III, Akhmerov et al. 30 derived a suitable topological invariant for a dirty class D nanowire directly in terms of the reflection matrix as Q = sign Det r 0 = sign Det r 0 = n tanh λ 0,n , and suggested that the TQPT could be observed as a quantized peak in the thermal conductance through the nanowire G th /G 0 = Tr t 0 t † 0 = n cosh −2 (L w λ n ), where we recover the result in Sec. III. This is consistent with the results in Ref. 24, where the authors predicted that the transition from topologically trivial to topologically non-trivial phases should be a delocalization transition, and at both sides of this point the system should be generically localized at zero energy. However, unfortunately a Majorana channel is necessarily neutral (i.e., particles and holes have equal weight in the MBS wavefunction) and therefore cannot support an electrical current. On the other hand, direct thermal transport measurements could provide evidence of the transition, 30 but this remains an experimental challenge.
To understand better our experimental proposal in Sec. V we first note that the form of Eqs. (30) and (31) are a consequence of the particular "reflectionless" boundary conditions (29) at the right end. In other words, in Fig. 7(a) at the right end of the wire, the barrier at x = L w is "transparent", and all right-moving Majorana states ψ R (L w ) that are transmitted to the rightend of the nanowire dissapears in the right lead. As shown in Ref. 35, this is not the most general situation. The generic presence of a barrier V b (x) at the end of the nanowire induces some probability of reflection, and imposes a non-vanishing amplitude ψ L (L w ) [see Fig. 7(b)]. More physically, any potential profile at the end nanowire, or the presence of pinch off gates could play the role of a barrier inducing a non-ideal coupling to the right-lead. For simplicity, let us consider a point-like scatterer sitting at some point x b > L w , as depicted in Fig. 7(b). The crucial point is that, in the presence of this new barrier, the reflectionless boundary conditions (29) are no longer possible. Assuming that the potential barrier V b (x) induces reflection and transmission amplitudes, r b,n and t b,n respectively (subject to the unitarity condition |r b,n | 2 +|t b,n | 2 = 1), the scattering matrix obeys and the new boundary conditions for right-moving Majorana states arriving from the left lead are In combination with Eqs. (28) and (32), we obtain the transmission and reflection amplitudes at the left-end of the complete system (nanowire and barrier): r n = r 0,n + t 0,n r b,n 1 − r b,n r 0,n t 0,n .
In particular, the last term in Eq. (35) physically represents processes in which the right-moving Majorana mode is transmitted to the right-end of the nanowire with amplitude t 0,n and is reflected back as a left-mover with amplitude r b,n . The denominator in r b,n (1 − r b,n r 0,n ) −1 = r b,n + r b,n r 0,n r b,n + . . . , represents an infinite sum of all backward and forward internal reflection processes occurring in the wire. Importantly, even though r n is a local quantity involving the reflection at the left-end of the nanowire, Eq. (35) explicitly contains non-local contributions involving scattering at the right-end, and its form is closely related to Eqs. (11) and (14) for the local differential conductances. This is a milestone result in phase-coherent mesoscopic transport which has been well-known for almost thirty years. 54 The above considerations summarize the main theoretical ideas in this work. Assuming that r b,n is a parameter that can be modified in situ in the experiment (as is the case of the pinch off gates in Ref. 9), Eq. (35) shows that a small variation δr b,n gives rise to a modification δr n ∝ cosh −2 (L w λ n ) δr b,n , which would be nonvanishing precisely at the TQPT and which could be detected in electrical measurements. This is the main idea of our proposal, and the main reason for us to propose experiments in the NSN geometry in order to establish the existence of the TQPT and the MBS in Majorana nanowires.

VII. SUMMARY AND CONCLUSIONS
We have explored the transport properties of disordered Majorana nanowires in the NSN configuration with the nanowire being the superconducting S part and the two N parts are ordinary metallic tunneling contacts with suitable gates controlling their tunnel barriers. This type of geometry is being explored at present by experimental groups studying Majorana bound states, and consequently, our study might be of relevance for the interpretation of these results. The NSN configuration allows to access qualitatively new information about the topological properties of the system through a direct study of nonlocal correlations inherent in the MBS which cannot be done in the NS geometry mostly used in the experimental Majorana measurements so far. Physically, this is possible because in the NSN configuration one can test the bulk properties, in addition to the boundary properties which are the only properties accessible in NS contacts. In our work we have adopted a comprehensive point of view, which links the deep theoretical aspects (i.e., topological invariants, topological classification, topological quantum phase transition and topological phase diagram) with the experimental observables (i.e., tunneling transport). We have also proposed a useful tool, i.e., the difference of local conductances Eq. (17), to detect the TQPT occurring as a function of the applied Zeeman field and to assess the topological protection of a given system experimentally. The experimental signal Eq. (17) is expected to be stronger and more robust to thermal broadening effects for "short" wires with ratio L w /ξ not too large (L w /ξ ≈ 15 in this work) and L w < L φ , i.e., smaller than the phase-relaxation length. We stress that this proposal to detect the TQPT is qualitatively different from the study of non-local correlations in the shot noise measurements. [42][43][44][45] Despite the simplifications assumed in this work, we note that our main ideas do not rely on the details of our model, but on generic symmetry properties of class D Bogoliubovde Gennes Hamiltonians. In particular, the fact that the TQPT correspond to a delocalized point at zero energy is a robust feature in these non-interacting Hamiltonians. In the presence of interactions the theoretical description of transport becomes much harder and remains an open issue. However, we speculate that the main idea behind Eq. (17) should remain valid in that case too. Interestingly, using the framework of Abelian bosonization, in Ref. 29 it was shown that the lowtemperature properties of a disordered class D wire with repulsive short-range electron-electron interactions (i.e., dimensionless Luttinger parameter ? K < 1) are adiabatically connected to those of a non-interacting wire (i.e., with K = 1), provided the system remains in the topological phase as the interaction is adiabatically "turned on". In particular, the delocalized nature of the TQPT in the interacting case can be inferred using an instanton calculation in the presence of disorder, where the equivalent of the localization length (i.e., the exponent of the instanton action) diverges at the critical point. 29 dω ρ 0 L,σ (ω) n L (ω) g r 1σ,1σ (ω) − g a 1σ,1σ (ω) + g < 1σ,1σ (ω) , (A11) dω ρ 0 R,σ (ω) n R (ω) g r N σ,N σ (ω) − g a N σ,N σ (ω) + g < N σ,N σ (ω) , where we have used the identity g > (ω)−g < (ω) = g r (ω)−g a (ω). 65,66 Obtaining an explicit expression for the currents I L and I R in the general case is quite cumbersome. However, since we will be interested only in the conductance, we note that there is a great simplification if we compute directly the conductance matrix by deriving the currents with respect to the voltages V L , V R . Then Replacing these expressions into Eqs. (A13)-(A16), and using the result g r jσ,jσ (ω) − g a jσ,jσ (ω) = −2πiρ jσ (ω), where we have defined the local density of states ρ jσ (ω) ≡ ρ jσ,jσ (ω), yields where for convenience we have ommitted the argument ω inside the brackets. In order to make explicit the non-local terms in these expressions we make use of the identity 68 From here, the following results are obtained g r 1σ,1σ (ω) − g a 1σ,1σ (ω) = −2πiρ 1,σ (ω) = −2πi s (A35) = t 2 L ρ 0 L g r 1σ,1s g a 1s,1σ + t 2 Lρ 0 L f r 1σ,1sf a 1s,1σ + t 2 R ρ 0 R g r 1σ,N s g a N s,1σ + t 2 which correspond to Eqs. (11)- (14).