Helical Majorana surface states of strongly disordered topological superconductors with time-reversal symmetry

Noncentrosymmetric superconductors with strong spin-orbit coupling and the B phase of ${}^3\mathrm{He}$ are possible realizations of topological superconductors with time-reversal symmetry. The nontrivial topology of these time-reversal invariant superconductors manifests itself at the material surface in the form of helical Majorana modes. In this paper, using extensive numerical simulations, we investigate the stability and properties of these Majorana states under strong surface disorder, which influences both bulk and surface states. To characterize the effects of strong disorder, we compute the level spacing statistics and the local density of states of both two- and three-dimensional topological superconductors. The Majorana surface states, which are located in the outermost layers of the superconductor, are protected against weak disorder due to their topological characteristic. Sufficiently strong disorder, on the other hand, partially localizes the surface layers, with a more pronounced effect on states with energies close to the gap than on those with energies close to zero. In particular, we observe that for all disorder strengths and configurations there always exist two extended states at zero energy that can carry thermal current. At the crossover from weak to strong disorder the surface state wave functions and the local density of states show signs of critical delocalization. We find that at this crossover the edge density of states of two-dimensional topological superconductors exhibits a zero-energy divergence, reminiscent of the Dyson singularity of quasi-one-dimensional dirty superconductors.

Figure 1

Spectral function $A(\omega ,{\mathbf{k}}_{\parallel })$, Eq. (8b), on a log scale as a function of energy $\omega$ and surface momentum $k_x$ with $k_y=0$ for the first three outermost layers at the (001) surface of a three-dimensional topological NCS. The strength $\gamma$ of the Gaussian distributed surface disorder increases from $\gamma =0t$ in (a) to $\gamma =30t$ in (e).

Figure 2

(a)–(h) Layer and energy-resolved probability distribution $P\left[{\tilde{\rho }}_l\left(\omega \right)\right]$ of the normalized local density of states ${\tilde{\rho }}_l\left(\omega \right)$ for an ensemble of one hundred disordered three-dimensional topological NCSs, with disorder strengths $\gamma =5t$ [panels (a)–(d)] and $\gamma =30t$ [panels (e)–(h)]. Each column shows the probability distribution $P\left[{\tilde{\rho }}_l\left(\omega \right)\right]$ for a different layer, with $l=1$ the surface layer [panels (a) and (e)] and $l=4$ the fourth inward layer [panels (d) and (h)]. The energy dependence is indicated by the color scale, with blue representing $\omega =0$ and red $\omega =\Delta /3$. Panels (i) and (j) show the wave function probability density $|{\mathit{\phi }}_m{(l,\mathbf{r})|}^2$ in the first $\left(l=1\right)$ and second inward layer $\left(l=2\right)$, respectively, for the three lowest positive energy wave functions ${\mathit{\phi }}_0,{\mathit{\phi }}_1$, and ${\mathit{\phi }}_2$ with successive eigenvalues $E_0=0<E_1<E_2$ (red, green, blue), at the disorder strength $\gamma =5t$. We find that at this disorder strength the wave functions show signs of critical delocalization. The amplitude in panel (j) has been multiplied by a factor of 5 for clarity. ${\mathit{\phi }}_0$ (red) is the helical Majorana state at exactly zero energy.

Figure 3

(a) Disorder-averaged total density of states ${\rho }_{\mathrm{tot}}\left(\omega \right)$ for a three-dimensional topological NCS with disorder strength $\gamma =5t$ (blue solid line) and $\gamma =30t$ (red solid line). The average is taken over 1000 disorder configurations. For comparison, the green dashed line displays the density of states of a clean Majorana cone. (b) Level spacing distribution function $P\left(S\right)$ for in-gap states with energies within the interval $\left|\omega \right|<\Delta /3$ with two different disorder strengths $\gamma =5t$ (blue circles) and $\gamma =30t$ (red triangles). The black dashed line is the generalized Wigner surmise for class DIII (i.e., $\alpha =1$ and $\beta =4$). The black dotted line represents the Poisson distribution. (c) Difference between the numerical data and the class DIII level statistics.

Figure 4

Spectral function $A(\omega ,k_x)$ on a log scale as a function of edge momentum $k_x$ for the first three outermost layers at the (01) edge of a two-dimensional topological NCS on a square lattice of size $80\times 40$. The strength of the edge disorder increases from $\gamma =0t$ in (a) to $\gamma =30t$ in (e).

Figure 5

(a)–(h) Layer and energy-resolved probability distribution $P\left[{\tilde{\rho }}_l\left(\omega \right)\right]$ of the normalized local density of states ${\tilde{\rho }}_l\left(\omega \right)$ for an ensemble of one hundred disordered two-dimensional topological NCSs with disorder strengths $\gamma =5t$ and $\gamma =30t$. The energy dependence is indicated by the color scale, with blue representing $\omega =0$ and red $\omega =\Delta /3$. Each column shows the probability distribution $P\left[{\tilde{\rho }}_l\left(\omega \right)\right]$ for a different layer.

Figure 6

Disorder-averaged total density of states ${\rho }_{\mathrm{tot}}\left(\omega \right)$ for a two-dimensional topological NCS with disorder strengths $\gamma =0.25t$ (green), $\gamma =5t$ (blue), and $\gamma =30t$ (red). The average is taken over 100 disorder configurations.

Figure 7

(a)–(c) Spectral function $A(\omega ,k_x)$ on a log scale as a function of surface momentum $k_x$ with $k_y=0$ for the first three outermost layers at the (001) face of a three-dimensional topological NCS in the presence of magnetic surface disorder with $\gamma =1t$. In (a), (b), and (c) the impurity spins are polarized along the $x,y$, and $z$ axes, respectively. (d)–(f) show the same as (a)–(c) but for the (01) edge of a two-dimensional topological NCS on a square lattice of size $80\times 40$.