Vortex lattices in the superconducting phases of doped topological insulators and heterostructures

Majorana fermions are predicted to play a crucial role in condensed matter realizations of topological quantum computation. These heretofore undiscovered quasiparticles have been predicted to exist at the cores of vortex excitations in topological superconductors and in heterostructures of superconductors and materials with strong spin-orbit coupling. In this work, we examine topological insulators with bulk $s$-wave superconductivity in the presence of a vortex lattice generated by a perpendicular magnetic field. Using self-consistent Bogoliubov–de Gennes calculations, we confirm that beyond the semiclassical, weak-pairing limit the Majorana vortex states appear as the chemical potential is tuned from either side of the band edge so long as the density of states is sufficient for superconductivity to form. Further, we demonstrate that the previously predicted vortex phase transition survives beyond the semiclassical limit. At chemical potential values smaller than the critical chemical potential, the vortex lattice modes hybridize within the top and bottom surfaces, giving rise to a dispersive low-energy mid-gap band. As the chemical potential is increased, the Majorana states become more localized within a single surface but spread into the bulk toward the opposite surface. Eventually, when the chemical potential is sufficiently high in the bulk bands, the Majorana modes can tunnel between surfaces and eventually a critical point is reached at which modes on opposite surfaces can freely tunnel and annihilate leading to the topological phase transition previously studied in the work of Hosur et al. [Phys. Rev. Lett. 107, 097001 (2011)].

Figure 1

Intraorbital pairing order parameters ${\Delta }_A$ (solid symbols) and ${\Delta }_B$ (hollow symbols) vs $\left|\mu \right|$ at different $\left|U\right|$. $\left|\mu \right|$ and $\left|U\right|$ are in units of eV. Via our definition for pairing in Eq. (10), $\Delta$ remains unitless. The mass term $\mathbb{M}=-1.5$ eV. The system contains periodic boundary conditions in $x$, $y$, and $z$ directions. The simulations are performed on a lattice grid of size $80a\times 80a\times 10a$.

Figure 2

A $4\times 4$ vortex lattice. In this example, the number of magnetic unit cells is $N_x\times N_y=4\times 2$. Each black solid circle denotes a vortex location and each magnetic unit cell contains two vortices.

Figure 3

The energy spectra of the low-energy states vs $\left|\mu \right|$ for periodic boundary conditions along the $z$ axis (vortex rings) for (a) $\left|U\right|=2$ eV, (b) $\left|U\right|=2.5$ eV, and (c) $\left|U\right|=2.8$ eV. The systems have a $N_x\times N_y=10\times 5$ (100 vortices) vortex lattice and the size of the magnetic unit cell is $l_x\times l_y\times l_z=12a\times 24a\times 10a.$ We note that in this work, we assume that the Zeeman splitting due to the applied magnetic field is negligible.

Figure 4

The spatial pairing order parameter profiles $|{\Delta }_{A,\vec{r}}|$ within a unit cell at $\left|U\right|=2.8$ eV and at (a) $\left|\mu \right|=0.96$ eV, (b) $\left|\mu \right|=0.98$ eV, and (c) $\left|\mu \right|=1.3$ eV. The vertical axis represents the pairing order parameter magnitudes and the horizontal plane is the $xy$ plane. Note the variation of the pairing order parameter magnitudes at the unit-cell center ${\Delta }_{\text{center}}={\Delta }_{A,\vec{r}\in \text{unit}\text{cell}\text{center}}$ for different $\mu$: (a) ${\Delta }_{\text{center}}\simeq 0$, (b) ${\Delta }_{\text{center}}\simeq 0.12$, and (c) ${\Delta }_{\text{center}}\simeq 0.28$. ${\Delta }_{\text{center}}\simeq 0$ indicates that the vortex core stays onsite, whereas ${\Delta }_{\text{center}}\ne 0$ means that the vortex core moves off the lattice vertices and into the plaquette. We note that the profiles for $|{\Delta }_{B,\vec{r}}|$ look similar.

Figure 5

The quasiparticle band structure for open vortex lines in the bulk superconducting TI. The interaction strength is chosen at $\left|U\right|=2.8$ eV and the chemical potentials are (a) $\left|\mu \right|=0.6$ eV, (b) $\left|\mu \right|=0.9$ eV, (c) $\left|\mu \right|=1.0$ eV, and (d) $\left|\mu \right|=1.3$ eV. The inset in (c) denotes the magnetic Brillouin zone for the square vortex lattice. The magnetic unit-cell sizes are $l_x\times l_y\times l_z=8a\times 16a\times 6a$.

Figure 6

The energy splitting $\delta E$ vs thickness $l_z$. The magnetic unit cells are $8\times 16\times l_z$ at $\left|\mu \right|=1$ and $1.1$ eV. The Hubbard interaction is $\left|U\right|=2.8$ eV.

Figure 7

Spatial slices (side views in the $yz$ plane at $x=\pm 2a$) of the probability density for the Majorana modes and order parameter density in the bulk superconducting TI at different $\mu$. Upper panels: (a)–(d) show the evolution of the Majorana mode distributions. Brighter regions represent the higher probability density. Lower panels: (e)–(h) show the distribution of pairing order parameters ($|{\Delta }_A+{\Delta }_B|$). The chemical potentials are $\left|\mu \right|=0.6$ eV in (a) and (e), $\left|\mu \right|=0.9$ eV in (b) and (f), $\left|\mu \right|=1$ eV in (c) and (g), and $\left|\mu \right|=1.3$ eV in (d) and (h), respectively.

Figure 8

The spatial side view of the pairing order parameter distribution of a six-layer $s$-wave/TI heterostructure. The interaction strength is chosen at $\left|U\right|=2.8$ eV, and the surface and bulk chemical potentials are $|{\mu }_S|=1$ eV and $|{\mu }_B|=0.55$ eV, respectively. The dark blue tubes indicate the region that pairing is suppressed and form vortex lines. The magnetic unit-cell sizes are $8a\times 16a\times 6a$.

Figure 9

(a) The quasiparticle band spectrum for the six-layer $s$-wave/TI heterostructure. (b) The spatial side view of the Majorana modes. The interaction strength is chosen at $\left|U\right|=2.8$ eV and $|{\mu }_S|=1$ eV and $|{\mu }_B|=0.55$ eV.

Figure 10

The comparison between the mid-gap sizes ${\delta }_m\left(T\right)$ and $k_{B}T$. The intersection occurs at $T=T_0\sim 0.025$ K, indicating that as $T<T_0$, the Majorana modes can stably exist on the surface of the doped topological insulators.