Majorana Zero Modes Protected by Hopf Invariant in Topologically Trivial Superconductors

TL;DR

This paper demonstrates that Majorana zero modes can exist in topologically trivial superconductors through a model based on the Hopf invariant, challenging the conventional association with topological phases.

Contribution

It introduces a minimal single-band model showing Majorana zero modes in trivial superconductors using the Hopf map and invariant.

Findings

Majorana zero modes can exist without topological order

The Hopf invariant protects zero modes in trivial phases

New avenues for realizing Majorana modes in superconductors

Abstract

Majorana zero modes are usually attributed to topological superconductors. We study a class of two-dimensional topologically trivial superconductors without chiral edge modes, which nevertheless host robust Majorana zero modes in topological defects. The construction of this minimal single-band model is facilitated by the Hopf map and the Hopf invariant. This work will stimulate investigations of Majorana zero modes in superconductors in the topologically trivial regime.

Equations

$$\displaystyle H({\bf k},\lambda)=\left(\begin{array}[]{cc}\xi_{\bf k}(\lambda)&\Delta_{\bf k}(\lambda)\\ \Delta_{\bf k}^{*}(\lambda)&-\xi_{-{\bf k}}(\lambda)\\ \end{array}\right)$$

$$H(\mathbf{k},\lambda )=\sum _{i=x,y,z}{d}_i(\mathbf{k},\lambda ){\tau }_i,$$

$$d_i=z^{†}{\tau }_{i}z,$$

$$z_1$$

$$z_2$$

$$H_{\eta }(\mathbf{k},\lambda )=\eta (d_x{\tau }_x+d_y{\tau }_y)+d_z{\tau }_z,$$

$$\lambda =n\theta ,$$

$$N_h=-\int d^{2}kd\theta \mathbf{j}\cdot \mathbf{A}.$$

$$H_{\mathrm{eff}}(k_y,\lambda )=-v{k}_y{\sigma }_x+M\lambda {\sigma }_y,$$

$$[iv{\sigma }_x{\partial }_y+(M/R)y{\sigma }_y]\psi =E\psi ,$$

$$E_{\mathrm{outer}}^2(s_z,n)=2(n+1/2)vM/R-s_{z}vM/R,$$

$$E_{\mathrm{inner}}^2(s_z,n)=2(n+1/2)vM/R^{\prime }+s_{z}vM/R^{\prime },$$

$${\psi }_{\mathrm{inner}}\sim \exp (-M{y}^2/2v{R}^{\prime })|s_z=-1⟩,$$

$$\displaystyle\lambda(\theta)=\left\{\begin{array}[]{cc}\lambda_{1},&\theta\in[0,\pi/2]\\ \lambda_{2},&\theta\in[\pi/2,\pi]\\ \lambda_{3},&\theta\in[\pi,2\pi].\end{array}\right.$$

$$z_1$$

$$z_2$$

$$d_x$$

$$d_y$$

$$d_z$$

$$d_x$$

$$d_y$$

$$d_z$$

$${\xi }_{\mathbf{k}}$$

$${\Delta }_{\mathbf{k}}$$

$${\hat{H}}_0$$

$${\hat{H}}_{\Delta }$$

$${\Delta }_{(x,y);(x+1,y)}=\eta (\frac{3}{2}-\cos \lambda -i\sin \lambda ),$$

$${\Delta }_{(x,y);(x,y+1)}=\eta (\sin \lambda -i\cos \lambda +\frac{3}{2}i),$$

$${\Delta }_{(x,y);(x,y+1)}=i{\Delta }_{(x,y);(x+1,y)},$$

$${\theta }_{x,y}=\arctan \frac{y-y_0}{x-x_0},$$

$${\lambda }_{x,y}=n{\theta }_{x,y}.$$

$${\hat{H}}_0$$

$${\hat{H}}_{\Delta }$$

$$\hat{H}=\sum _{x,y;x^{\prime },y^{\prime }}{\Psi }^{†}(x,y)H_{x,y;x^{\prime },y^{\prime }}\Psi (x^{\prime },y^{\prime }),$$

Full text

Majorana zero modes (MZMs) or Majorana bound states are exotic excitations predicted to exist in the vortex coresRead and Green (2000); Volovik (1999) of two-dimensional (2D) topological superconductorsHasan and Kane (2010); Qi and Zhang (2011); Bansil et al. (2016); Bernevig and Hughes (2013); Shen (2013) and at the ends of 1D topological superconductorsKitaev (2001). Spatially separated MZMs give rise to degenerate ground states, which encode qubits immune to local dechoerenceKitaev (2001); Nayak et al. (2008). Furthermore, unitary transformations among the ground states can be implemented by braidingMoore and Read (1991); Wen (1991); Ivanov (2001); Nayak and Wilczek (1996); Das Sarma et al. (2005) or measurementsVijay and Fu (2016); Bonderson et al. (2008) of these modes, indicating that such qubits may become building blocks in topological quantum computation and informationAasen et al. (2016); Karzig et al. (2016, 2016); Van Heck et al. (2012); Landau et al. (2016); Deng and Duan (2013). Therefore, MZMs have been vigorously pursued in condensed matter physicsAlicea (2012); Beenakker (2013); Stanescu and Tewari (2013); Leijnse and Flensberg (2012); Elliott and Franz (2015); Sarma et al. (2015); Sato and Fujimoto (2016).

There have been a great variety of proposals for topological superconductors, including 2D semiconductor heterostructuresSau et al. (2010); Alicea (2010), topological insulator-superconductor proximityFu and Kane (2008); Qi et al. (2010); Chung et al. (2011); Law et al. (2009); Akhmerov et al. (2009), 1D spin-orbit-coupled quantum wiresOreg et al. (2010); Lutchyn et al. (2010); Alicea et al. (2011); Lutchyn et al. (2011); Stanescu et al. (2011); Potter and Lee (2010); Rainis et al. (2013); Prada et al. (2012); Das Sarma et al. (2012), spiral magnetic chains on superconductorsChoy et al. (2011); Martin and Morpurgo (2012); Nadj-Perge et al. (2013); Klinovaja et al. (2013); Vazifeh and Franz (2013), Shockley mechanismWimmer et al. (2010); Deng et al. (2015), and cold atom systems in 2DSato et al. (2009); Zhang et al. (2008); Tewari et al. (2007a); Liu et al. (2014) and 1DJiang et al. (2011); Diehl et al. (2011), etc. Experimentally, suggestive signatures of MZMs in both 1DMourik et al. (2012); Nadj-Perge et al. (2014); Rokhinson et al. (2012); Deng et al. (2012); Das et al. (2012); Finck et al. (2013); Churchill et al. (2013); Albrecht et al. (2016); Deng et al. (2016); Franz (2013); Pawlak et al. (2016) and 2DXu et al. (2015); Sun et al. (2016); Lv et al. (2016); Wang et al. (2012, 2016); He et al. (2016) topological superconductors have been found.

It is often implicitly assumed that topological superconductivity is a prerequisite for MZMs, accordingly, the chiral edge states go hand in hand with the vortex zero modes in 2D superconductors. In this Letter we show that certain topological defectsTeo and Kane (2010); Chiu et al. (2016); Shiozaki and Sato (2014); Lee et al. (2007); Qi et al. (2008a); Ran et al. (2009) in 2D topologically trivial superconductors can support robust MZMs. Somewhat surprisingly, single-band superconductors suffice this purpose. The model Hamiltonian is related to the Hopf maps, which originally refer to nontrivial mappings from a 3D sphere $S^{3}$ to a 2D sphere $S^{2}$, characterized by the integer Hopf invariantWilczek and Zee (1983); Nakahara (2003). Mappings from a 3D torus $T^{3}$ to $S^{2}$ inherit the nontrivial topology from the mappings $S^{3}\rightarrow S^{2}$. The Hopf invariant has found interesting applications in nonlinear $\sigma$ models and spin systemsWilczek and Zee (1983); Fradkin (2013), Hopf insulatorsMoore et al. (2008); Deng et al. (2013, 2014, 2016); Kennedy (2016); Liu et al. (2016), liquid crystalsAckerman and Smalyukh (2017), and quench dynamics of Chern insulatorsWang et al. (2016); Fläschner et al. (2016).

Our model describes topologically trivial superconductors with zero Chern number and no chiral edge state. Nevertheless, a topological point defect is characterized by a Hopf invariant defined in the $(k_{x},k_{y},\theta)$ space, where $k_{x},k_{y}$ are crystal momenta and $\theta$ is the polar angle111Converting a momentum variable to an angle variable was adopted in Refs.Qi et al. (2008b)(Fig.1a). The parity (even/odd) of Hopf invariant determines the presence (absence) of robust MZMs, though the superconductor for every fixed $\theta$ is topologically trivial. Stimulated by this mechanism, which significantly differs from the magnetic-vortex origin of zero mode in topological $p$-wave superconductorRead and Green (2000); Volovik (1999), we design trivial-superconductor-based (and vortex-free) T-junctions harboring MZMs.

*Zero modes.–*Before studying topological defects, we consider spatially uniform 2D single-band Bogoliubov-de Gennes (BdG) Hamiltonians parameterized by $\lambda$:

[TABLE]

where ${\bf k}=(k_{x},k_{y})$, $\xi_{\bf k}=E_{\bf k}-\mu$, $E_{\bf k}$ and $\mu$ is the energy and chemical potential, respectively, and $\Delta_{\bf k}$ is the Cooper pairing. It describes single-band spinless (or spin-fully-polarized) superconductors. This Hamiltonian can be written in terms of the Pauli matrices $\tau_{i}$ as

[TABLE]

with $d_{x}={\rm Re}\Delta_{\bf k}$, $d_{y}=-{\rm Im}\Delta_{\bf k}$, $d_{z}=\xi_{{\bf k}}$ (we have $\xi_{\bf k}=\xi_{-{\bf k}}$ in our model). For reason to become clear shortly, we take

[TABLE]

where $z=(z_{1},z_{2})^{T}$ and

[TABLE]

with $m_{0}=\frac{3}{2}$. We can check that $\Delta_{-{\bf k}}=-\Delta_{\bf k}$ and $\Delta_{(-k_{y},k_{x})}=i\Delta_{(k_{x},k_{y})}$, thus the pairing is $p$-wave. Given Eq.(6), the pairing $\Delta_{\bf k}$ is of the same order as the hopping $\xi_{\bf k}$. To describe weakly-pairing superconductors, one may consider

[TABLE]

with a small but nonzero $\eta$. Nevertheless, tuning the value of $\eta$ does not close the energy gap, hence it does not qualitatively change the results. Thus we will simply take $\eta=1$ below. The mathematical form of Eq.(6) has been introduced for 3D Hopf insulatorsMoore et al. (2008); Deng et al. (2013, 2014, 2016); Kennedy (2016); Liu et al. (2016), with $\lambda$ replaced by the third momentum $k_{z}$. The physical system we will study is nevertheless not directly related to Hopf insulator.

The familiar Chern numberThouless et al. (1982); Read and Green (2000); Qi et al. (2010) $C$ that characterizes 2D topological superconductors can be obtained by a straightforward numerical calculation, which yields $C(\lambda)=0$ for every $\lambda$.

A topological defect can be generated if the parameter $\lambda$ depends on spatial coordinates, in a manner that the configuration cannot be smoothly deformed to a spatially uniform one. Let us focus on the defects with $\lambda$ depending on the polar angle $\theta$ (Fig.1a) as

[TABLE]

where $n$ is an integer. These configurations are topologically nontrivial due to a nonzero Hopf invariant, as we now explain. The unit vector ${\bf\hat{d}}({\bf k},\theta)\equiv\frac{1}{\sqrt{\eta^{2}(d_{x}^{2}+d_{y}^{2})+d_{z}^{2}}}(\eta d_{x},\eta d_{y},d_{z})$ maps the 3D torus $T^{3}$ ($k_{x},k_{y},\theta$ are defined modulo $2\pi$) to the 2D unit sphere $S^{2}$. For nonzero $n$, the inverse-image circles of two points on $S^{2}$ are linkedsup . To quantify such linking, the Hopf invariant can be definedWilczek and Zee (1983); Moore et al. (2008): $N_{h}=-\frac{1}{4\pi^{2}}\int d\theta d^{2}k\epsilon^{\mu\nu\rho}a_{\mu}\partial_{\nu}a_{\rho}$, where the integrating range is the Brillouin zone for ${\bf k}$ and $[0,2\pi]$ for $\theta$, $a_{\mu}=-i\langle\psi({\bf k},\theta)|\partial_{\mu}|\psi({\bf k},\theta)\rangle$, with $\mu,\nu,\rho=k_{x},k_{y},\theta$, and $|\psi\rangle$ is the negative-energy eigenfunction of $H_{\eta}({\bf k},\theta)$. Alternatively, we can define $A_{\mu}=a_{\mu}/2\pi$, $j^{\mu}=\epsilon^{\mu\nu\rho}\partial_{\nu}A_{\rho}=(1/8\pi)\epsilon^{\mu\nu\rho}\hat{{\bf d}}\cdot(\partial_{\nu}\hat{{\bf d}}\times\partial_{\rho}\hat{{\bf d}})$, thenWilczek and Zee (1983); Moore et al. (2008)

[TABLE]

It can be calculated numerically by discretizing the Brillouin zonesup . The numerical result for $n=1$ is shown in Fig.1b. More generally, we have $N_{h}=n$. We will call the topological defects defined by Eq.(8) as Hopf defects.

To obtain energy spectra, we Fourier-transform the BdG Hamiltonian to real-space lattice, then numerically solve the Hamiltoniansup . For $n=1$, two zero modes are found (Fig.2a), one of which is sharply localized around the defect, the other is localized at the sample boundary. The profile of particle component ($\tau_{z}=1$ component, denoted by $u$) and hole component ($\tau_{z}=-1$ component, denoted by $v$) is shown in the main figure and the inset, respectively. It is apparent that the zero modes are equal-weight superpositions of particle and hole components, which is a feature of MZMs (inspection of the wavefunction confirms that $u=v^{*}$).

To check the robustness of MZMs, we add an impurity potential $Uc^{\dagger}_{A}c_{A}$ ($c$ is the fermion operator) at a single site $A$ (specified in Fig.1a), which amounts to adding a $U\tau_{z}$ term at site $A$ in the real-space BdG Hamiltonian. The numerical result for $U=1.0$ is shown in Fig.2b. The energies of MZMs remain pinned to $E=0$, though the mode profile is changed compared to Fig.2a.

For $n=2$, we find two localized modes in the defect and two at the boundary (near $\theta=0$ and $\theta=\pi$). Unlike the $n=1$ case, the energies of defect modes are not pinned to zero. Higher $n$’s are also calculated, and the results support the conclusion that there is a single robust MZM for odd-integer $n$, and no robust MZM for even-integer $n$, therefore, it is the Z2 Hopf invariant (even/odd) that determines the existence of MZM in the defect. One may notice that the Hopf invariant takes the form of a Chern-Simons invariantTeo and Kane (2010); Qi et al. (2008b), which is not accidental, because the latter is indeed a general topological invariant, nevertheless, it has been appliedTeo and Kane (2010); Chiu et al. (2016) only to topologically nontrivial superconductors, for which it is just the product of the Chern number and the vorticity of pairing phase. Our model shows that nonzero Chern number is not a necessary condition for MZM.

*Edge theory.–*It is desirable to have an intuitive understanding of the MZM from the perspective of edge theory, which, as we will show, differs significantly from that of the chiral topological superconductorRead and Green (2000); Fendley et al. (2007); Stone and Roy (2004). First, we numerically solved the edge states of open-boundary systems for various values of $\lambda$, and found that gapless edge modes exist only for $\lambda=0$. In Fig.3a, we show the energy bands for a ribbon along $y$ direction. The gapless edge modes for $\lambda=0$ are shown as the solid blue lines. They are non-chiral, and are immediately gapped out when $\lambda$ is tuned away from 0 (edge modes of $\lambda=\pi/20$ are shown in dashed curves), in other words, the edge modes are not topologically robust. This is consistent with the vanishing of Chern number.

This numerical observation is confirmed by analytic solutions. We consider a semi-infinite geometry with the sample occupying $x<0$ region, $k_{y}$ being a good quantum number. For $\lambda=0$, we obtain two degenerate edge modes at $k_{y}=0$, both of which are eigenfunctions of $\tau_{x}$ with eigenvalue $-1$sup , thus they are equal-weight superpositions of particle and hole components. We introduce Pauli matrices $\sigma_{x,y,z}$ (unrelated to the $\tau_{x,y,z}$ matrices) in this two-dimensional space, so that the two eigenfunctions have $\sigma_{z}=\pm 1$, respectively. Including small $k_{y}$ and $\lambda$ as perturbations, we derive an effective theorysup :

[TABLE]

where the effective parameters $v,M$ are found to be both $3/4$ in our specific modelsup . Thus the edge-state spectra are $E_{\pm}(k_{y})=\pm\sqrt{v^{2}k_{y}^{2}+M^{2}\lambda^{2}}$. It is immediately clear that the edge states become gapped when $\lambda$ moves away from 0, which is consistent with the numerical finding in Fig.3a. As a comparison, we note that the edge spectrum of a chiral topological superconductorRead and Green (2000), $E(k_{y})=vk_{y}$, cannot be gapped out.

Based on this effective edge theory, we proceed to study a hollow disk with polar-angle-dependent parameter, $\lambda=\theta$ (Fig.3b). We are only concerned with low-energy modes, therefore, we focus on the neighborhood of $\theta=0$. Suppose that both the outer and inner radiuses $R,R^{\prime}$ are large. On the outer boundary, we have $\lambda=\theta=y/R$, thus the edge state spectra are given by solving $H_{\rm eff}(k_{y}\rightarrow-i\partial_{y},\lambda\rightarrow y/R)\psi=E\psi$. More explicitly, it reads

[TABLE]

which squares to $(-v^{2}\partial_{y}^{2}+M^{2}y^{2}/R^{2}-vM/R\sigma_{z})\psi=E^{2}\psi$. This equation resembles the Schrödinger equation of harmonic oscillators, though $E$ is replaced by $E^{2}$, and there is a crucial additional $-vM/R\sigma_{z}$ term. The eigenfunctions are $\sigma_{z}$-eigenvectors (eigenvalues are denoted by $s_{z}=\pm 1$), with energies given by

[TABLE]

where $n=0,1,\cdots$. There is a MZM in the $s_{z}=+1$ sector, with $n=0$, which is illustrated as the blue bump in Fig.3b. Since this mode is the eigenfunction of $\tau_{x}$, it is an equal-weight superposition of particle and hole components.

For a semi-infinite geometry with sample occupying the $x>0$ region, the effective edge theory is almost the same as Eq.(10), except that the sign of the first term reversedsup . On the inner boundary of the hollow disk (again near $\theta=0$), we have $\lambda=\theta=y/R^{\prime}$, thus the edge-mode spectrum can be obtained from $[-iv\sigma_{x}\partial_{y}+(M/R^{\prime})y\sigma_{y}]\psi=E\psi$, analogous to Eq.(11). The energies are given by

[TABLE]

which features a MZM in the $s_{z}=-1$ sector. The zero-mode wavefunction is

[TABLE]

which is exponentially localized near $y=0$, namely $\theta=0$ (illustrated by the green bump in Fig.3b). All nonzero energies grow as $1/R^{\prime}$ as $R^{\prime}$ is decreased, while the MZM remains at zero energy, evolving to the defect mode shown in Fig.2. For a hollow disk with $\lambda=2\theta$, there are two MZMs on the inner boundary for large $R^{\prime}$, near $\theta=0$ and $\theta=\pi$, respectively. Shrinking $R^{\prime}$ causes overlapping between them, which splits the two zero energies to nonzero values. This is consistent with the absence of MZM in the $n=2$ defect.

It is useful to compare our systems with the chiral topological superconductor, for which a magnetic vortex with $\pi$-flux hosts a MZMRead and Green (2000); Volovik (1999); Fu and Kane (2008); Tewari et al. (2007b). In a hollow-disk geometry, this MZM comes from the chiral edge states on the boundary circleRead and Green (2000); Fendley et al. (2007); Stone and Roy (2004). The MZM wavefunction is evenly distributed on the circle, which implies its sensitiveness to the magnetic flux. In contrast to this picture, the MZM in our model is not derived from chiral edge state, which is simply absent here, moreover, the MZM is exponentially localized near $\theta=0$ (Fig.3b), thus it is insensitive if a magnetic flux is inserted.

T-*junctions.–*So far, we have only studied configurations with $\lambda$ continuously varied. It is conceivable that the smooth Hopf defect defined by Eq.(8) can be imitated by a discontinuous one, for instance, we may consider a T-junction:

[TABLE]

$\lambda_{1,2,3}$ being three unequal constants. Such T-junctions will presumably be easier to realize than configurations with $\lambda$ smoothly varying in space.

We will study the simpler superconductor-superconductor-vacuum T-junction by replacing the $\lambda_{3}$-region by the vacuum. Since the value of $\lambda$ in the vacuum is not well defined, this replacement is not fully justified in advance. Nevertheless, the numerical results thus obtained indicate that MZMs do exist in such T-junctions, as shown in Fig.4a for $\lambda_{1}=0.1\pi$, $\lambda_{2}=\pi$. There is certain arbitrariness in choosing the hopping at the boundary between the $\lambda_{1,2}$ regions, for which we keep only the nearest-neighbor hopping (discarding the next-nearest-neighbor hopping and the pairing)222Keeping only the nearest-neighbor hopping means that (in real space) $H_{{\bf r}{\bf s}}=(1.5-\cos\lambda)\tau_{z}$ for the boundary sites ${\bf r}=(L/2,y)$ and ${\bf s}=(L/2+1,y)$, where we take $\lambda$ as $\lambda_{2}=\pi$ (Another choice can be $\lambda=\lambda_{1}$. The choice is not unique).. The energy eigenvalues near zero are shown in Fig.4c (we show 12 of them), from which it is clear that the zero-mode levels are separate from all other energy levels by a finite gap in the $L\rightarrow\infty$ limit.

We have also studied superconductor-insulator-vacuum T-junctions. To this end, we consider the $H_{\eta}$ in Eq.(7), in which taking $\eta=0$ amounts to removing the Cooper pairing. We notice that $H_{\eta=0}({\bf k},\lambda=\pi)$ describes an insulator without any Fermi surface [In contrast, $H_{\eta=0}({\bf k},\lambda=0.1\pi)$ describes a metal]. Now we can design superconductor-insulator-vacuum T-junctions by taking $\eta=1$ in the $\lambda_{1}=0.1\pi$ region, and $\eta=0$ in the $\lambda_{2}=\pi$ region. We find one MZM for each T-junction (shown in Fig.4b), and the energy gap between the zero-mode levels and other energy levels is apparent in Fig.4d. We emphasize that each region by itself is topologically trivial, in particular, the superconductor (with $\lambda=0.1\pi$) is a topologically trivial one without gapless edge state.

*Conclusions.–*We have investigated the intriguing possibility of creating MZMs in 2D topologically trivial superconductors. The Hopf defect is constructed as a minimal model for this purpose. Furthermore, we studied the more accessible T-junctions constructed from topologically trivial superconductors. Hopefully, the trivial-superconductor-based approach will broaden the scope of searching MZMs in various superconductors. In particular, absence of chiral Majorana edge state in a 2D superconducting sample does not necessarily imply absence of robust MZM in its point defects.

We conclude with several remarks. First, we have focused on a single-band model, while many materials have multi-bands. We emphasize that the MZMs found here are nevertheless robust to small mixing with other bands, because a single localized MZM cannot move away from zero energy, as required by the intrinsic particle-hole symmetry of the BdG HamiltonianHasan and Kane (2010). Second, we have taken a simple model BdG Hamiltonians as our starting point (like Ref.Read and Green (2000)). More realistic Hamiltonians should be adopted when dealing with real materials, for instance, the semiconductor-superconductor heterostructuresSau et al. (2010); Alicea (2010), for which our theory implies that robust MZMs can exist in certain defects (e.g. judiciously constructed T-junctions), without requiring the uniform system being tuned to the topologically nontrivial regime. This will be left for future works.

*Acknowledgements.–*We would like to thank Suk Bum Chung for helpful suggestions on the manuscript. This work is supported by NSFC (No. 11674189). Z.Y. is supported in part by China Postdoctoral Science Foundation (No. 2016M590082).

I I. Explicit expressions of $d_{x,y,z}$

In the main article, we have taken $d_{i}=z^{{\dagger}}\tau_{i}z$ with $z=(z_{1},z_{2})^{T}$ and

[TABLE]

For general $m_{0}$, the explicit form of these $d_{i}$’s read

[TABLE]

For the special case $m_{0}=\frac{3}{2}$, they become

[TABLE]

II II. Explicit expressions of BdG Hamiltonian in the real space: uniform system

In the main article, the Hamiltonian is written in a compact form in the ${\bf k}$ space. More explicit (less compact) expressions for the kinetic energy $\xi_{\bf k}$ and the pairing gap $\Delta_{\bf k}$ are

[TABLE]

For many purposes, it is useful to do a Fourier transformation to the real space. For a uniform system (i.e. spatially independent $\lambda$) , the BdG Hamiltonian is given by $\hat{H}=\hat{H}_{0}+\hat{H}_{\Delta}$ with the kinetic energy term

[TABLE]

and the pairing term

[TABLE]

where $x,y$ are real-space coordinates taking integer values. As explained in the main article, $\lambda$ is simply a parameter of the Hamiltonian. We can see that the pairing between the nearest-neighbor sites $(x,y)$ and $(x+1,y)$ is

[TABLE]

and the pairing between sites $(x,y)$ and $(x,y+1)$ is

[TABLE]

thus

[TABLE]

showing a $p$-wave character in real space. Since the Chern number is zero (for any value of parameter $\lambda$), the superconductor is topologically trivial, as we have explained in the main article.

III III. Explicit expressions of BdG Hamiltonian in the real space: with a defect

In the previous section, we have focused on uniform systems with spatially independent $\lambda$. As explained in the main article, a topological defect can be created if the Hamiltonian parameter $\lambda$ depends on the polar angle $\theta$, which is measured from the defect center $(x_{0},y_{0})$, i.e.

[TABLE]

such that the configuration cannot be smoothly deformed to a uniform one. The spatial dependence of the parameter $\lambda$ in our topological defects is given by

[TABLE]

Before proceeding, it is useful to point out that, since the hopping terms and the pairing terms are defined on the links instead of sites, it is natural to take the middle point of the link to define $\theta$, for instance, $\cos\theta\,c^{\dagger}_{x,y}c_{x+1,y}$ is taken as $\cos\theta_{x+1/2,y}\,c^{\dagger}_{x,y}c_{x+1,y}$ with $\theta_{x+1/2,y}=\arctan[(y-y_{0})/(x+1/2-x_{0})]$. It is also viable to take a different convention, say $\cos\theta_{x,y}\,c^{\dagger}_{x,y}c_{x+1,y}$, which does not affect the existence of Majorana zero mode, as verified in our numerical calculations. The reason is that $|\theta_{x,y}-\theta_{x+1/2,y}|\rightarrow 0$ as $(x-x_{0})^{2}+(y-y_{0})^{2}\rightarrow\infty$, thus different choices merely differ in the near-defect-core region. The topological character of the defect is fixed by the Hamiltonian far from the defect center, which is not affected by modifying the near-defect-core region.

With the above technical aspect explained, we are ready to write down the real-space BdG Hamiltonian in the presence of a defect:

[TABLE]

which simply replace the parameter $\lambda$ in Eq.(23) and Eq.(24) by $n\theta$.

It is a conventional step to define $\Psi(x,y)=(c_{x,y},c_{x,y}^{{\dagger}})^{T}$, and the real-space BdG Hamiltonian becomes

[TABLE]

with nonzero elements of $H_{x,y;x^{\prime},y^{\prime}}$ given by

[TABLE]

The rank of the matrix $H_{x,y;x^{\prime},y^{\prime}}$ is $2L^{2}$, $L$ being the linear size of a square sample.

IV IV. Numerical calculation of Hopf invariant

We start from the ${\bf k}$-space BdG Hamiltonian:

[TABLE]

For $\eta\ll 1$, $H_{\eta}$ describes a superconductor with weak pairing (in the main article, we focused on the $\eta=1$ case).

As we explained in the main article, the unit vector

[TABLE]

maps the 3D torus $T^{3}$ ($k_{x},k_{y},\theta$ are defined modulo $2\pi$) to the two-dimensional unit sphere $S^{2}$. For nonzero $n$, the inverse-image circles of two points on $S^{2}$ are linked, which is illustrated in Fig.5 for the $n=1$ case.

The Hopf invariant, which characterizes the topological property of the Hamiltonian, is given by

[TABLE]

where ${\bf j}({\bf k},\theta)=(j^{x},j^{y},j^{\theta})$ (Note that the superscripts $x,y$ stand for $k_{x},k_{y}$) takes the form of

[TABLE]

and ${\bf A}=(A^{x},A^{y},A^{\theta})$ (again, the superscripts $x,y$ stand for $k_{x},k_{y}$) is the gauge potential satisfying $\nabla\times{\bf A}={\bf j}$, where $\nabla=(\partial/\partial k_{x},\partial/\partial k_{y},\partial/\partial\theta)$. It is difficult to do the integration in Eq.(53) analytically. To do it numerically, we rewrite Eq.(53) into a discrete form:

[TABLE]

where $N$ is the number of lattice sites (the lattice constant is set to unit) and $k_{x,y},\theta$ take discrete values in $\{-\pi,-\pi+\frac{2\pi}{N},...,\pi-\frac{2\pi}{N}\}$.

To obtain the expression of ${\bf A}$ from ${\bf j}$, we do the following Fourier transformation:

[TABLE]

where ${\bf q}=(q_{x},q_{y},q_{\theta})$, whose components $q_{x,y,\theta}$ take discrete values in $\{-\frac{N}{2},-\frac{N}{2}+1,...,\frac{N}{2}-1\}$. Under the gauge ${\bf q}\cdot{\bf A}=0$, it is readily found that

[TABLE]

thus the Hopf invariant becomes

[TABLE]

The numerical integration converges quite rapidly as we increase $N$. For $\lambda=n\theta$, we find $N_{h}=n$. The $n=1$ case is shown in Fig.1b of the main article.

V V. More supporting data on the effects of impurity potential

In the main article, we have shown that the Majorana zero mode for $n=1$ is robust to an impurity potential, and remarked that robust Majorana zero mode is absent for even-integer Hopf invariant. Here, more data on this is shown in Fig.6. As explained in the main article, an impurity potential $Uc^{\dagger}_{A}c_{A}$ at a single site $A$ (indicated in Fig.1a in the main article) is added to test the robustness of Majorana zero modes. For odd-integer Hopf invariant, Fig.6a shows that the impurity potential cannot move the zero mode energies away from $E=0$ (in consistent with the results given in the main article). For even-integer Hopf invariant, as shown in Fig.6b and Fig.6c, zero modes are generally absent in the defect in the presence of impurity potential (the blue lines in Fig.6b and Fig.6c are the energies for modes near the boundary, not near the defect at the system center). Note that several nonzero energy eigenvalues are almost independent of $U$ in Fig.6a, Fig.6b, and Fig.6c, because their eigenfunctions are quite far from the impurity site $A$.

VI VI. Edge theory for various sample geometries

VI.1 Sample occupying $y>0$

First, we study a semi-infinite geometry with the sample occupying the $y>0$ region (illustrated in Fig.7a). The momentum $k_{x}$ remains a good quantum number, while $k_{y}$ is not. In the $y$ direction, we use the real-space coordinate, which is an integer in our lattice model. The wave functions and energy eigenvalues can be obtained by solving the following equation

[TABLE]

where the subscript $j$ of $\psi_{j}$ is the $y$ coordinate, each $\psi_{j}$ is two-component, and

[TABLE]

Eq.(74) can be written compactly as $H(k_{x},\lambda)|\Psi\rangle=E|\Psi\rangle$, $H$ denoting the large matrix at the left-hand-side of Eq.(74). Let us first find the zero-energy ($E=0$) solutions, if any, at $k_{x}=0$ and $\lambda=0$. When $k_{x}=0$ and $\lambda=0$, $H_{i}$’s simplify to

[TABLE]

We can define $\tau_{\pm}=(\tau_{z}\pm i\tau_{x})/2$, which is the raising/lowering operator of $\tau_{y}$, thus $H_{0}=-(\tau_{+}+\tau_{-})/4$, $H_{1}=H_{2}=-\tau_{+}$, and $H_{1}^{\dagger}=H_{2}^{\dagger}=-\tau_{-}$. Operating on the $\tau_{y}$ eigenvectors $|\chi_{\pm}\rangle=(1,\pm i)^{T}/\sqrt{2}$, they produce

[TABLE]

We can check that the zero-energy wavefunctions take the form of

[TABLE]

where $|j\rangle$ is localized on the $j$-site:

[TABLE]

with the coefficients $a_{j}$ satisfying the following iteration relation:

[TABLE]

We can see that $\psi_{j}=a_{j}|\chi_{-}\rangle$. It is straightforward to find two normalizable solutions for $a_{j}$’s, which we will denote as $a_{j}=\alpha_{j}$ ($j=1,2,3,\cdots$) and $a_{j}=\beta_{j}$ ($j=1,2,3,\cdots$). The explicit expressions of $\alpha_{j}$ and $\beta_{j}$ are

[TABLE]

and

[TABLE]

where $C_{\alpha}$ and $C_{\beta}$ are two normalization constants such that $\sum_{j}|\alpha_{j}|^{2}=1$ and $\sum_{j}|\beta_{j}|^{2}=1$. It is straightforward to find that $C_{\alpha}=\sqrt{3}/2$ and $C_{\beta}=\sqrt{27/80}$.

We can check that the two wave functions $|\Psi_{1}\rangle=\sum_{j}\alpha_{j}|j\rangle\otimes|\chi_{-}\rangle$ and $|\Psi_{2}\rangle=\sum_{j}\beta_{j}|j\rangle\otimes|\chi_{-}\rangle$ are not orthogonal. The orthogonalization can be achieved by the Gram-Schmidt orthogonalization. We have

[TABLE]

where $\mathcal{N}=\sqrt{5}/2$ is a normalization constant and the superscript “o” stands for “orthogonalization”.

So far, we have focused on $k_{x}=0$, $\lambda=0$. In the neighbourhood of $k_{x}=0$, $\lambda=0$, we can expand the Hamiltonian $H$ to the first order of $k_{x}$ and $\lambda$, namely, $H_{i}(k_{x},\lambda)=H_{i}(0,0)+\Delta H_{i}(k_{x},\lambda)$ with

[TABLE]

It is readily seen that

[TABLE]

Straightforward calculations lead to

[TABLE]

and

[TABLE]

Noting the mathematical identities

[TABLE]

we can see that

[TABLE]

Similarly, we find that

[TABLE]

To summarize the above calculations, we have the following low-energy effective Hamiltonian for the edge states near $k_{x}=0$, $\lambda=0$:

[TABLE]

or more compactly,

[TABLE]

VI.2 Sample occupying $y<0$

Now we study the sample occupying the $y<0$ region (Fig.7b), then the edge modes can be obtained by solving the following eigenvalue problem,

[TABLE]

Following the same procedures in the previous section, we first find the solutions for $k_{x}=0$ and $\lambda=0$. We find two wave functions with $E=0$, one of which is $|\Psi_{1}^{\prime}\rangle=\sum_{j}\alpha_{j}|j\rangle\otimes|\chi_{+}\rangle$, and the other is $|\Psi_{2}^{\prime}\rangle=\sum_{j}\beta_{j}|j\rangle\otimes|\chi_{+}\rangle$. The orthogonalzation of the two wavefunctions are achieved by the Gram-Schmidt orthogonalization,

[TABLE]

It is readily checked that

[TABLE]

Straightforward calculations yield

[TABLE]

thus, the low-energy effective Hamiltonian for the edge state reads

[TABLE]

VI.3 Sample occupying $x>0$

For sample occupying the $x>0$ region (Fig.7c), the edge modes can be obtained by solving the following eigenvalue problem,

[TABLE]

where

[TABLE]

When $k_{y}=0$ and $\lambda=0$,

[TABLE]

It is not difficult to see that

[TABLE]

where $|\xi_{\pm}\rangle=(1,\pm 1)^{T}/\sqrt{2}$ are the two eigenvectors of $\tau_{x}$. Following similar procedures as previous sections, we find two solutions at $E=0$ for $k_{y}=\lambda=0$. One is $|\tilde{\Psi}_{1}\rangle=\sum_{j}\alpha_{j}|j\rangle\otimes|\xi_{+}\rangle$, and the other is $|\tilde{\Psi}_{2}\rangle=\sum_{j}\beta_{j}|j\rangle\otimes|\xi_{+}\rangle$. Again we adopt the Gram-Schmidt orthogonalization to define

[TABLE]

In the neighbourhood of $k_{y}=0$ and $\lambda=0$, $\tilde{H}_{i=0,1}$ can be expanded to the first order of $k_{y}$ and $\lambda$:

[TABLE]

with

[TABLE]

It is straightforward to check that

[TABLE]

which lead to

[TABLE]

Therefore, the low-energy effective Hamiltonian for the edge state takes the form of

[TABLE]

VI.4 Sample occupying $x<0$

Finally, we study the sample occupying the $x<0$ region (Fig.7d). We need to solve the eigenvalue equation:

[TABLE]

Following the same steps as in previous sections, we find two $E=0$ modes at $k_{y}=0$ and $\lambda=0$, one of which is $|\tilde{\Psi}_{1}^{\prime}\rangle=\sum_{j}\alpha_{j}|j\rangle\otimes|\xi_{-}\rangle$, and the other is $|\tilde{\Psi}_{2}^{\prime}\rangle=\sum_{j}\beta_{j}|j\rangle\otimes|\xi_{-}\rangle$, where we continue to use $|\xi_{\pm}\rangle=(1,\pm 1)^{T}/\sqrt{2}$ to denote the two eigenvectors of $\tau_{x}$. After orthogonalzation, the wave functions take the form of

[TABLE]

In the neighborhood of $k_{y}=\lambda=0$, we expand the Hamiltonian to the first order of $k_{y}$ and $\lambda$: $\tilde{H}_{i}(k_{y},\lambda)=\tilde{H}_{i}(0,0)+\Delta\tilde{H}_{i}(k_{y},\lambda)$, which satisfy

[TABLE]

therefore, we have

[TABLE]

thus the low-energy effective Hamiltonian for the edge state is given by

[TABLE]