Introduction

A recent breakthrough in the theory of topological superconductors (Tsc) is that certain crystalline symmetries can protect a rich variety of higher-order Tsc phases1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16, featuring different patterns of Majorana corner or hinge modes (see Fig. 1a for an example). These higher-order Tsc phases can intrinsically exist without utilizing proximity effects and are distinct from the well-known first-order Tsc with Majorana edge or surface modes17,18,19. For instance, two-dimensional (2D) superconductors with rotational9,15 or inversion7,8,10,13 symmetries can host Majorana corner modes, as shown from both mathematical analyses using real-space classification7,15 and numerical observations in a candidate material, monolayer WTe27.

Fig. 1: Majorana boundary patterns for all possible 3D time-reversal superconducting phases in space group No. 2.
figure 1

Schematics that show the correspondence between our derived symmetry indicators (SIs) \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\), \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},i}\), \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},ij}\) and the Majorana boundary patterns (colored by pink) on a cubic geometry for 3D time-reversal superconductors with inversion and translation symmetries. The Majorana patterns and the explicit forms of SIs are obtained in Steps 1 and 2, respectively, and their correspondence is derived from the basis-matching procedure in Step 2. Note that the Majorana surface modes associated with the winding number N3 are not shown here. Nonetheless, since the strong SI \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\) detects the evenness and oddness of N3, the strong phases with odd \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\) still show one copy of the surface modes.

Full size image

However, despite the rapid theoretical development and extensive experimental efforts made over the past decades, unambiguously confirmed materials hosting either first- or higher-order Tsc phases remain extremely rare, especially beyond 1D. This is partly because, unlike topological insulators, the topology in Tsc phases is determined not only by the electronic band structures but also by the pairing symmetry. Progress in Tsc material discovery is therefore hindered by the lack of systematic guiding principles for material search and design that are derived from fundamental understanding.

A natural strategy to design the topology of the superconducting state is to identify the necessary normal state properties and pairing symmetry. In this regard, there is a sharp contrast between first-order and higher-order Tsc: To achieve a first-order Tsc, having the correct pairing symmetry alone can be sufficient regardless of whether the normal state is topological, e.g., an order parameter of Δ = px + ipy is well known to lead to a first-order Tsc with Majorana edge modes17. In contrast, to achieve a higher-order Tsc, past works have shown that the normal state topology can play an essential role on top of the pairing symmetry5,7,9,12,13,20,21.

Here, we propose and demonstrate a four-step strategy to systematically obtain guiding principles for identifying and designing candidate materials of higher-order Tsc. Our strategy is described as follows for superconducting materials from a given space group G. In Step 1, we find all possible first- and higher-order Tsc phases and their Majorana boundary patterns by a real-space classification analysis22,23,24,25,26,27,28 (Fig. 1a). In Step 2, we derive the set of superconducting topological invariants κsc that can diagnose these Majorana boundary patterns from the superconducting band structures, following a protocol developed by some of us13. In Step 3, we identify the set of topological invariants κn that characterize the normal-state topology and fermiology. In Step 4, we obtain the function fΔ that relates the two sets of invariants by

$${\kappa }_{{{{\rm{sc}}}}}={f}_{\Delta }({\kappa }_{{{{\rm{n}}}}})$$
(1)

for all possible pairing symmetries Δ in space group G. From the master equation, Eq. (1), guiding principles for the normal state properties described by κn and pairing symmetry Δ can be systematically obtained for different higher-order Tsc phases labeled by κsc.

To demonstrate our strategy, an ideal material platform is group-VI transition metal dichalcogenides (TMD), many of which exhibit intrinsic superconductivity with highly tunable normal states exhibiting spin-orbit couplings or band topology29,30,31,32,33,34,35,36,37,38,39. For example, pressure-induced superconductivity was recently discovered in 3D MoTe2 in the centrosymmetric \(1{T}^{{\prime} }\) structure40, where the normal state was proposed to possess higher-order band topology41,42. Motivated by the observed superconductivity in 1\({T}^{{\prime} }\)-MoTe2, we will demonstrate our strategy on 3D Tsc protected by the inversion symmetry. Following the recipe derived from Eq. (1), we predict that 1\({T}^{{\prime} }\)-MoTe2 is a plausible Tsc candidate with both first- and higher-order topology that hosts surface and corner modes. We support this prediction by a density functional theory (DFT) calculation and microscopic mean-field analysis using a realistic 44-band tight-binding model.

Results

Step 1: Majorana boundary patterns

Motivated by the superconducting centrosymmetric MoTe240, we consider 3D time-reversal superconductors in the simplest inversion-symmetric space group (space group No. 2), which contains the inversion symmetry \({{{\mathcal{I}}}}\) and the three translation symmetries \({T}_{\hat{r}}\), r = x, y, z. Such superconductors are described by a Bogoliubov de Gennes (BdG) Hamiltonian \({H}_{{{{\rm{BdG}}}}}={\sum }_{{{{\bf{k}}}}}{\Psi }_{{{{\bf{k}}}}}^{{\dagger} }H({{{\bf{k}}}}){\Psi }_{{{{\bf{k}}}}}\), where

$$H({\mathbf{k}}) = \left(\begin{array}{cc} h({\mathbf{k}}) & {\Delta}({\mathbf{k}}) \\ {\Delta}{\left({\mathbf{k}}\right)}^{\dagger} & -h{\left(-{\mathbf{k}}\right)}^{\ast} \end{array}\right),$$
(2)

the Nambu basis \({\Psi}_{\mathbf{k}}=({\hat{c}}_{{\mathbf{k}},\uparrow}, {\hat{c}}_{{\mathbf{k}},\downarrow}, {\hat{c}}^\dagger_{-{\mathbf{k}},\uparrow}, {\hat{c}}^\dagger_{-{\mathbf{k}},\downarrow})^\top\), and the indices for degrees of freedom other than spin s = ↑, ↓ are suppressed. Here, the normal state h(k) is invariant under the inversion operation \({{{\mathcal{I}}}}({{{\bf{k}}}})h({{{\bf{k}}}}){{{\mathcal{I}}}}{({{{\bf{k}}}})}^{-1}=h(-{{{\bf{k}}}})\) with \({{{\mathcal{I}}}}(-{{{\bf{k}}}}){{{\mathcal{I}}}}({{{\bf{k}}}})=1\). For the superconducting gap Δ(k), we focus on the odd-parity cases where \({{{\mathcal{I}}}}({{{\bf{k}}}})\Delta ({{{\bf{k}}}}){{{\mathcal{I}}}}{(-{{{\bf{k}}}})}^{-1}=\eta \Delta (-{{{\bf{k}}}})\) with gap parity η = − 1. Together with the particle-hole and time-reversal symmetries \({{{\mathcal{P}}}}\) and \({{{\mathcal{T}}}}\), the symmetry group of H(k) obeys the following group relations:

$$\begin{array}{ll}\quad\,{{{{\mathcal{T}}}}}^{2}=-1,\,{{{{\mathcal{P}}}}}^{2}=1,\,{{{{\mathcal{I}}}}}_{{{{\rm{BdG}}}}}^{2}=1\\ \left[{{{\mathcal{T}}}},{{{\mathcal{P}}}}\right]=0,\,\left[{{{\mathcal{T}}}},{{{{\mathcal{I}}}}}_{{{{\rm{BdG}}}}}\right]=0,\,\{{{{\mathcal{P}}}},{{{{\mathcal{I}}}}}_{{{{\rm{BdG}}}}}\}=0,\end{array}$$
(3)

where \({{{{\mathcal{I}}}}}_{{{{\rm{BdG}}}}}=\,{{\mbox{diag}}}\,({{{\mathcal{I}}}},\eta {{{\mathcal{I}}}})\) is the BdG inversion operator that acts on the Nambu basis, and translations simply commute with all other symmetries. Importantly, for odd-parity superconductors, \({{{{\mathcal{I}}}}}_{{{{\rm{BdG}}}}}\) and \({{{\mathcal{P}}}}\) anticommute so that the particle hole partners have opposite parities.

To obtain all possible Majorana boundary patterns that a 3D time-reversal centrosymmetric superconductor can support, we first compute the classification group \({{{{\mathcal{C}}}}}_{{{{\rm{r}}}}}\) for crystalline Tsc phases described by HBdG, then examine the boundary signature of each phase. This can be achieved by using a well-developed real-space classification method called the Topological Crystal Approach13,22,23,24,25,26,27,28,43,44,45,46. The key idea is that although it is hard to compute \({{{{\mathcal{C}}}}}_{{{{\rm{r}}}}}\) directly in the presence of nonlocal crystalline symmetries, one can dissect the full 3D superconductor into lower-dimensional building blocks that respect only the local internal symmetries but not the nonlocal crystalline symmetries. Specifically, these building blocks are db-dimensional topological states with 0 ≤ db ≤ 3, where their classification groups and boundary modes have been well studied in the prior literature18 (see “Methods” section). By stacking these building blocks into different configurations that respect all the symmetries and checking various consistency conditions13, we can determine the Majorana boundary signature of each configuration. These topologically distinct configurations with different Majorana signatures are dubbed topological crystal states, where each of them provides a minimal model for each of the Tsc phases. Importantly, any 3D superconducting material that respects a given set of crystalline symmetries can be adiabatically connected to a certain topological crystal state13,23,25. We therefore expect that the Majorana boundary pattern we obtain for a topological crystal state can also be found in a realistic lattice model for a superconductor in the same Tsc phase.

The classification group \({{{{\mathcal{C}}}}}_{{{{\rm{r}}}}}\) for our current case of 3D time-reversal superconductors with inversion and translation symmetries is obtained using this approach as follows. First, we identify that the nontrivial building blocks are the time-reversal 1D, 2D, and 3D Tsc states hosting Majorana end, edge, and surface modes, respectively. These building blocks can be stacked on Wyckoff positions in different symmetry-allowed configurations to form inversion-symmetric 3D superconducting states with different Majorana boundary patterns. By identifying all inequivalent and robust configurations and excluding those that lead to atomic superconductors without Majorana modes, we find that the real-space classification group is given by \({{{{\mathcal{C}}}}}_{{{{\rm{r}}}}}=({\mathbb{Z}}\times {{\mathbb{Z}}}_{4})\times {({{\mathbb{Z}}}_{4})}^{3}\times {({{\mathbb{Z}}}_{2})}^{3}\) (see “Method” section).

Next, for each phase captured in \({{{{\mathcal{C}}}}}_{{{{\rm{r}}}}}\), we now discuss the Majorana signature and the protecting symmetries obtained from its building block configuration and will leave the explicit forms of topological invariants to Step 2. Specifically, we find that the \({\mathbb{Z}}\) factor in \({{{{\mathcal{C}}}}}_{{{{\rm{r}}}}}\) corresponds to first-order strong phases with Majorana surface modes, which can be trivialized by breaking the time-reversal symmetry and is described by a nonzero integer topological invariant N3. The first \({{\mathbb{Z}}}_{4}\) factor corresponds to inversion-protected higher-order strong phases with Majorana hinge and corner modes. While the \({({{\mathbb{Z}}}_{2})}^{3}\) corresponds to weak phases protected by two translation symmetries along the xy, yz, or xz-directions, the \({({{\mathbb{Z}}}_{4})}^{3}\) corresponds to the mixed phases protected by the inversion and/or translation symmetries. We label them by topological invariants \({\kappa }_{{{{\rm{mixed}}}}}^{i}\), i = x, y, z, where the \({\kappa }_{{{{\rm{mixed}}}}}^{i}=1\) phases are purely protected by the translation symmetry along i-direction, the \({\kappa }_{{{{\rm{mixed}}}}}^{i}=2\) phases are protected simultaneously by the inversion and translation symmetries, and the \({\kappa }_{{{{\rm{mixed}}}}}^{i}=3\) phases are the stackings of the former two. Note that if we quotient out the phases with an even topological invariant N3, the classification group of the strong phases \({\mathbb{Z}}\times {{\mathbb{Z}}}_{4}\) becomes \({{\mathbb{Z}}}_{8}\), consistent with the findings in previous works that did not consider N36,8,11. With this adjustment, the real-space classification group becomes \({{{{\mathcal{C}}}}}_{{{{\rm{r}}}}}={{\mathbb{Z}}}_{8}\times {({{\mathbb{Z}}}_{4})}^{3}\times {({{\mathbb{Z}}}_{2})}^{3}\).

The resulting Majorana signatures for these strong, weak, and mixed Tsc phases are summarized in Fig. 1a, which we obtain by systematically checking the robustness and consistency relations when stacking the building blocks13 (see the details for this standard procedure of Topological Crystal Approach in “Methods” section). Among these Tsc phases, there are first-order phases with Majorana surface states, higher-order Tsc phases with Majorana hinge or corner modes, as well as a rich variety of Tsc phases with different dimensional Majorana modes, which is dubbed hybrid-order Tsc phases.

Step 2: superconducting state topological invariants κ sc

We now turn to the momentum space to derive explicit forms of a set of topological invariants \({\kappa }_{{{{\rm{sc}}}}}=\{{N}_{3},{\kappa }_{{{{\rm{sc}}}}}^{s},{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},i},{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},ij}\}\), i, j = x, y, z, which can diagnose the complete Majorana boundary signatures for a given centrosymmetric superconductor (see Fig. 1a). In the following, we show that N3 is the well-known 3D winding number18,47,48,49 for 3D time-reversal Tsc, while \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\), \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},i}\), and \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},ij}\) for the strong, mixed, and weak phases are functions of band symmetry data at the high-symmetry points (TRIMs) only. Such invariants are termed symmetry indicators (SIs)6,8,10,11,14,50,51,52,53,54,55,56,57,58. In particular, we find that the SI for the strong phases \({\kappa }_{{{{\rm{sc}}}}}^{s}\) detects both the higher-order phases and the evenness and oddness of N3.

First, we calculate the momentum space classification group \({{{{\mathcal{C}}}}}_{{{{\rm{k}}}}}\) to determine the nature of each topological invariant. Our calculation is performed using a classification method for topological crystalline phases called Twisted Equivariant K Theory59. We find that the full classification group is given by \({{{\mathcal{K}}}}={{\mathbb{Z}}}^{9}\), which consists of a subgroup \({\mathbb{Z}}\) defined on the entire 3D Brillouin zone (BZ) and a subgroup \({{{\mathcal{{K}}}^{{\prime} }}}={{\mathbb{Z}}}^{8}\) restricted only to the TRIMs (see “Methods” section). Since this subgroup \({\mathbb{Z}}\) originates from 3D class-DIII Tsc18 (see Table 3 in Method), the phases captured by the \({\mathbb{Z}}\) subgroup can be labeled by the well-known integer 3D winding number18,47,48,49 for 3D time-reversal Tsc with Majorana surface modes (see Eq. 36 in ref. 18). On the contrary, the phases in the \({{{\mathcal{{K}}}^{{\prime} }}}\) subgroup are labeled by the SIs \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\), \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},i}\), and \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},ij}\). These phases include strong, weak, and mixed phases with various Majorana signatures, as well as all atomic superconductors, which do not host Majorana boundary modes. After removing the contribution from atomic superconductors {AS}, we find that the remaining Tsc phases with non-trivial band topology at TRIMs are classified by \({{{{\mathcal{C}}}}}_{{{{\rm{k}}}}}={{{\mathcal{{K}}}^{{\prime} }}}/\{\,{{\mbox{AS}}}\,\}={{\mathbb{Z}}}_{8}\times {({{\mathbb{Z}}}_{4})}^{3}\times {({{\mathbb{Z}}}_{2})}^{3}\) (see “Methods” section). Therefore, the momentum-space classification \({{{{\mathcal{C}}}}}_{{{{\rm{k}}}}}\) is consistent with the real-space classification \({{{{\mathcal{C}}}}}_{{{{\rm{r}}}}}\) we find in Step 1.

Having shown that the topological invariants that correspond to \({{{{\mathcal{C}}}}}_{{{{\rm{k}}}}}\) are SIs, we now derive the explicit SI expressions that can diagnose the Majorana signatures of the Tsc phases in \({{{{\mathcal{C}}}}}_{{{{\rm{k}}}}}\) (see Fig. 1a). Specifically, these SIs were proposed to be different linear combinations of a \({\mathbb{Z}}\)-invariant6 defined at each of the eight TRIMs k = Γ, X, Y, Z, U, T, R, S:

$${\kappa }_{{{{\rm{sc}}}}}^{\eta }=\mathop{\sum}\limits_{k}{\alpha }_{k}^{\eta }{n}_{k},\quad {n}_{k}=\frac{1}{2}\left({N}^{+}[H(k)]-{N}^{+}[{H}_{{{{\rm{ref}}}}}(k)]\right),$$
(4)

where N+[HBdG(k)] is the number of even-parity occupied states of Hamiltonian H(k), and Href(k) is a trivial BdG Hamiltonian that serves as a reference point13 (see our choice of Href(k) in “Methods” section).

The coefficients \(\{{\alpha }_{k}^{\eta }\}\) for these SIs are further obtained by performing a basis matching procedure13, where we establish a transparent correspondence between the resulting SIs \({\kappa }_{{{{\rm{sc}}}}}^{\eta }\) and the Majorana patterns shown in Fig. 1a. Specifically, to ensure that the \({{\mathbb{Z}}}_{8}\), \({{\mathbb{Z}}}_{4}\), and \({{\mathbb{Z}}}_{2}\) SIs correspond exclusively to the strong, mixed, and weak phases, respectively, we explicitly check the SI values \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\), \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},i}\), and \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},ij}\) for the real-space minimal models we obtained in Step 1 for each of the strong, mixed, and weak Tsc phases in Fig. 1a (see “Methods” section). These minimal models are models for different topological crystal states built by different building block configurations, where Majorana signatures are evident (see “Methods” section). This procedure is necessary because without it, the \({{\mathbb{Z}}}_{8}\), \({{\mathbb{Z}}}_{4}\), and \({{\mathbb{Z}}}_{2}\) SIs in general would each correspond to some profound mixture of strong, mixed, and weak phases due to a basis ambiguity13. Finally, we arrive at the following SI expressions for time-reversal invariant Tsc phases with inversion and translation symmetries:

$$\begin{array}{l}{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\,=\,\mathop{\sum}\limits_{k\in \,{{\mbox{TRIMs}}}\,}{n}_{k}\,\,{{\mbox{mod}}}\,\,8,\\ {\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},z}\,=\,{n}_{Z}+{n}_{U}+{n}_{T}+{n}_{R}\,\,{{\mbox{mod}}}\,\,4,\\ {\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},x}\,=\,{n}_{X}+{n}_{S}+{n}_{U}+{n}_{R}\,\,{{\mbox{mod}}}\,\,4,\\ {\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},y}\,=\,{n}_{Y}+{n}_{S}+{n}_{T}+{n}_{R}\,\,{{\mbox{mod}}}\,\,4,\\ {\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},xy}\,=\,{n}_{S}+{n}_{R}\,\,{{\mbox{mod}}}\,\,2,\\ {\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},yz}\,=\,{n}_{U}+{n}_{R}\,\,{{\mbox{mod}}}\,\,2,\\ {\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},xz}\,=\,{n}_{T}+{n}_{R}\,\,{{\mbox{mod}}}\,\,2,\end{array}$$
(5)

where the superscripts s, m, w stand for strong, mixed, and weak phases. The SIs for this symmetry class have been reported in previos works ref. 6,11 without performing the basis matching procedure. Therefore, the explicit correspondence between the SIs and the Majorana boundary signatures was not explicitly established in previous works. This set of SIs satisfies the bulk-boundary correspondence so that they can fully distinguish not only all the distinct band topology in the bulk, but also all the Majorana boundary patterns shown in Fig. 1a. We expect that these SIs are applicable to realistic material-based models since these minimal models are adiabatically connected to any lattice model in the same Tsc phases.

Before moving on to Step 3, we point out that the complete Majorana signatures are characterized by not just the SIs we find in Eq. (5), but also the 3D winding number N3. In fact, the winding number N3 and SIs \({\kappa }_{{{{\rm{sc}}}}}^{\eta }\) are not mutually independent. Specifically, the parity of \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\) for strong phases is equal to the winding number N3 modulo 247,48,49. In “Methods” section, we explicitly show that the pair \(({N}_{3},{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}})\) is isomorphic to the group \(({\mathbb{Z}}\times {{\mathbb{Z}}}_{4})\), which agrees with the real space classification. In terms of the winding number N3 and the SI for strong phases \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\), first-order Tsc is indicated by \(({N}_{3} \,>\, 0,{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}=0,1)\), higher-order Tsc is indicated by \(({N}_{3}=0,{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}=2,4)\), and all other cases correspond to various hybrid-order Tsc phases.

Step 3: normal-state invariant κ n

To characterize the normal state, we adapt the topological invariant for time-reversal topological crystalline insulators (TCI) in the same space group, which contains the inversion \({{{\mathcal{I}}}}\) and translation symmetries \({T}_{\hat{r}}\). The invariant for such strong TCI phases was proposed to be a \({{\mathbb{Z}}}_{4}\) integer that depends on the electron band parity data at the TRIMs k only54:

$${\kappa }_{I}^{{{{\rm{s}}}}}=\frac{1}{4}\mathop{\sum}\limits_{k\in \,{{\mbox{TRIMs}}}\,}\left({N}^{+}[{h}_{I}(k)]-{N}^{-}[{h}_{I}(k)]\right)\,\,{{\mbox{mod}}}\,\,4,$$
(6)

where the superscript s stands for strong phases, N±[hI(k)] is the number of even- and odd-parity occupied bands in the Hamiltonian hI(k) for the insulator. The TCI phases with \({\kappa }_{{{{\rm{strong}}}}}^{0}=1,2,3\) exhibit electronic surface modes, hinge modes, and a combination of both, respectively.

To adapt this \({{\mathbb{Z}}}_{4}\) invariant \({\kappa }_{I}^{{{{\rm{s}}}}}\) for characterizing the metallic normal states h(k) in Eq. (2), we now allow it to take both integers and half integers values:

$${\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}={\kappa }_{I}^{{{{\rm{s}}}}}{| }_{{h}_{I}(k)\to h(k)}=0,\frac{1}{2},1,\cdots \,,\frac{7}{2}.$$
(7)

Depending on the normal-state fermiology, there are two cases: When all the Fermi surfaces are away from TRIMs, \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}\) remains a \({{\mathbb{Z}}}_{4}\) integer and the normal state can be viewed as a doped TCI that carries the same band topology as the underlying TCI state. When the Fermi surfaces circle at least one TRIM, \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}\) may be a half integer or integer, depending on the number of Fermi pockets circling TRIMs and the topology of the fully occupied bands. For instance, doping a higher-order TCI with hinge modes will lead to a normal state of \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}=2\) if Fermi pockets are away from TRIMs. In contrast, a doped trivial insulator with an even-parity electron pocket at k = Γ is characterized by \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}=\frac{1}{2}\). Note that each band is two-fold degenerate due to the time-reversal and inversion symmetries. Note that instead of rigorously describing the topology of the normal state, \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}\) should be viewed as a computational device for obtaining the recipes in Step 4. Moreover, \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}\) is not restricted to simple doped topological insulators, but applicable to general cases where the Fermi level cuts through multiple connected bands since \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}\) is not sensitive to the order of band parities.

Step 4: recipes for higher-order Tsc states

Equipped with the superconducting and normal state strong invariants \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\) and \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}\), which are \({{\mathbb{Z}}}_{8}\) and \({{\mathbb{Z}}}_{4}\) numbers respectively, we are now ready to obtain the master equation in Eq. (1) that relate the two. Although \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\) and \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}\) in Eqs. (5) and (7) are written in terms of the BdG and normal bands, respectively, the relation fΔ between them can be found in the weak-pairing limit where \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\) can be expressed in terms of the normal band parities for a given pairing symmetry Δ. The weak-pairing limit is the limit where the gap is smaller than the normal band spacings, which is mostly true for superconductors with low Tc. This is done by expressing the \({\mathbb{Z}}\)-invariant nk at each TRIM k in Eq. (4) in terms of the normal band parities as6,13

$${n}_{k}=\frac{1}{2}\left({N}^{+}[h(k)]-{N}^{-}[h(k)]\right).$$
(8)

Given Eqs. (8), (5), and (7), we find that the relating function fΔ has a simple form

$${\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}=2{\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}},$$
(9)

when the superconducting gap Δ is parity-odd. In contrast, when the pairing gap Δ is parity-even, the classification is trivial such that all superconducting invariants \({\kappa }_{{{{\rm{sc}}}}}^{i}\) vanish57. This indicates that an even-parity time-reversal nodeless gap always leads to a topologically trivial superconductor without Majorana modes even when the normal state is topological.

From the relation in Eq. (9), we can deduce recipes for higher-order Tsc phases that consist of conditions on the normal state \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}\) in the presence of odd-parity pairing gap Δ. Here, we discuss two example recipes. First, if the normal state is a doped strong TI labeled by \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}=1\) whose Fermi surfaces are away from TRIMs, introducing an odd-parity gap will drive the system into a second-order Tsc with Majorana hinge modes since \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}=2\). Physically speaking, these Majorana hinge modes are leftover normal-state surface states that cannot be gapped out by the superconducting gap due to the odd-parity nature. Second, if the normal state is a doped higher-order TI featuring inversion-protected hinge modes (\({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}=2\)), we expect an exotic third-order Tsc with Majorana corner modes (\({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}=4\)) when the system develops an odd-parity pairing gap. For the metallic normal state to have \({\kappa }_{{{{\rm{n}}}}}^{{{{\rm{s}}}}}=2\), the doping-induced Fermi pockets can either be away from TRIMs, or there can be pairs of Fermi pockets that have opposite band parities at TRIMs. The latter case is relevant to the MoTe2 case, as shown below. For both recipes, the superconducting gap has to be not only odd-parity but also time-reversal symmetric. When the Fermi pockets are away from TRIMs, an example gap is a spin-triplet px-wave gap whose nodal line does not intersect with the Fermi pockets. When the Fermi pockets circle TRIMs, an example gap is a 3He-B-phase-like Balian-Werthammer (BW) gap with winding number N3 = 1: \({\Delta }_{{{{\rm{BW}}}}}({{{\bf{k}}}})\propto {k}_{x}\left\vert \uparrow \uparrow -\downarrow \downarrow \right\rangle +i{k}_{y}\left\vert \uparrow \uparrow +\downarrow \downarrow \right\rangle +{k}_{z}\left\vert \uparrow \downarrow +\downarrow \uparrow \right\rangle\)18,19. Nonetheless, on top of the Majorana corner modes indicated by \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}=4\), we also expect Majorana surface modes indicated by the non-zero 3D winding number N3. In a realistic superconductor, depending on the actual hopping parameters and on-site potentials, coexisting Majorana corner and surface modes could experience various levels of hybridization effects.

Symmetry indicators in MoTe2

The second recipe suggests that superconducting MoTe2 in the centrosymmetric lattice structure is a plausible candidate for such a third-order Tsc with Majorana corner modes. This is because previous DFT calculations on \(1{T}^{{\prime} }\)-MoTe2 have reported higher-order band topology along with Fermi pockets located at TRIMs41,42. In the following, we will numerically obtain the Majorana boundary signatures in superconducting MoTe2 to examine our prediction made from the last recipe. To this end, we need to numerically compute the full set of SIs \(\{{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}},{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},i},{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},ij}\}\), i, j = x, y, z since the existence of mixed and weak phases is also important for determining the full Majorana boundary signatures. Specifically, we perform a DFT calculation on centrosymmetric MoTe2 in an experimentally relevant geometry60 without the spin-orbit coupling, using Vienna Ab-initio Simulation Package (VASP)61,62. A monoclinic primitive unit cell contains 4 Mo atoms and 8 Te atoms (Fig. 2a). Calculation details are described in “Methods” section. We add an on-site Coulomb repulsion (Hubbard U) term of 3.0 eV for the Mo 4d orbitals, within the DFT+U method63, since it was shown that with this addition, the calculated band structure agree well with the experimental angle-resolved photoemission spectrum (ARPES)40.

Fig. 2: Crystal structure, Fermi surface, and band structure of 1T’-MoTe2.
figure 2

a Crystal structure of MoTe2 (Mo: blue, Te: gray) with the lattice vectors. b Fermi surfaces at chemical potential μ = − 46 meV below the Fermi level in the BZ. c DFT Band structure without SOC where the μ value used in b is marked with an arrow and a maroon solid line.

Full size image

When the superconducting gap is small, fully gapped, and parity-odd, which can be energetically favored in the presence of nearest-neighbor attractions, we can obtain the set of SIs \(\{{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}},{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},i},{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},ij}\}\) following Eq. (5). This is equivalent to obtaining the BdG band parities from the BdG Hamiltonian consisting of a normal state constructed by the DFT bands and a small pairing gap that is fully gapped and parity-odd. We find that the full set of SIs is given by

$${\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}=4,\,\,{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},i}={\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},ij}=0\,\,\,{{\mbox{for}}}\,\,i,j=x,y,z$$
(10)

at a chemical potential μ = − 46 meV below the Fermi level. Table 1 shows the numbers of bands with positive parity and negative parity at the 8 time-reversal invariant momentum (TRIM) points at μ = − 46 meV marked by an arrow in Fig. 2c. We have checked that the SIs do not change without the U value or with spin-orbit coupling. Besides the computed SIs, there is also a non-zero winding number N3 given the considered BW pairing gap ΔBW. We expect that N3 is even since the computed strong SI \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}=4\) is even and there are two Fermi pockets at this chemical potential μ (see Fig. 2b) that each develops an N3 = 1 BW gap.

Table 1 The numbers of occupied bands Ne, positive-parity bands n+, and negative-parity bands n− at the 8 TRIM points at μ = − 46 meV for MoTe2
Full size table

Corner modes in MoTe2

According to the SIs in Eq. (10) and the winding number N3, we expect that centrosymmetric MoTe2 with odd-parity pairing is a hybrid-order Tsc with Majorana corner modes (see the \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}=4\) phase in Fig. 1a) on top of Majorana surface modes. To numerically verify this expectation, we construct a 44-band tight-binding model for MoTe2 based on 44 Wannier functions obtained using WANNIER9064,65,66. The 44 Wannier functions consist of dxy, dyz, dxz, \({d}_{{x}^{2}-{y}^{2}}\), and \({d}_{{z}^{2}}\) orbitals of all four Mo sites and px, py, and pz orbitals of all eight Te sites in the unit cell. The constructed tight-binding model reproduces the DFT band structure in the energy range of [-6,3] eV relative to the Fermi level (Fig. 5a). Details can be found in “Methods” section.

Using this tight-binding model as the normal state, we construct a BdG Hamiltonian with the superconducting gap being the BW gap ΔBW given in the second recipe in Step 4. We consider the BW gap since the normal state has Fermi surfaces at TRIMs (Fig. 2b) such that this is the simplest gap structure under which the superconductor is parity-odd but also fully gapped. By diagonalizing the constructed BdG Hamiltonian on a 3D open geometry with a system size of L = 9 (We are limited to a small system size by the computing power and the not-so-sparse realistic model), we find near-zero-energy eigenstates that are localized at a pair of inversion-related corners (Fig. 3), where the choice of which corners is likely determined by the microscopic positions of Wannier orbitals.

Fig. 3: Near-zero-energy corner modes in superconducting 1T′-MoTe2 with an odd-parity gap.
figure 3

a The BdG spectrum for MoTe2 at a chemical potential labeled in Fig. 2c and with a BW gap ΔBW on a finite lattice of 9 × 9 × 9 unit cells computed by Krylov method. The gap between the blue states is due to both the finite-size effect and a likely small hybridization between the Majorana corner and surface modes. b The spatial probability distribution ∣ψm∣2 of the near-zero-energy BdG eigenstates ψm labeled in blue in a, demonstrating the existence of corner modes when the hybridization is small. The geometry preserves the inversion symmetry, and ax, ay, az are the lattice constants in x, y, z directions, respectively.

Full size image

Importantly, instead of exact zero-energy modes well-separated from other finite-energy quasiparticle states, these low-energy corner modes are burried in a gapless spectrum (see Fig. 3a). This is expected from the surface modes indicated by the non-zero 3D winding number N3, but not captured by the SIs. Consequently, the hybridization between the surface and corner modes can open a small gap and induce some degree of delocalization of the corner modes into the surfaces. This effect is visually not evident in Fig. 3b since we choose a chemical potential μ at which the hybridization strength is likely small to demonstrate the existence of corner modes, but is evident in the gapless spectrum in Fig. 3a. Therefore, our numerically observed near-zero energy corner modes in Fig. 3 supports our prediction that 3D centrosymmetric MoTe2 is a hybrid-order Tsc with \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}=4\) and a non-zero winding number N3 at the chemical potential μ with an odd-parity superconducting gap ΔBW.

Discussions

The effect of hybridization between different dimensional Majorana boundary modes has not been systematically investigated before, but prior studies on topological insulators have shown analytically and numerically that boundary modes can survive the hybridization with higher-dimensional states67,68. In our hybrid-order Tsc case, our numerical results in Fig. 3 also support the robustness of corner modes in \(1{T}^{{\prime} }\)-MoTe2 with any odd-parity superconducting gap. In the actual material, whether the corner modes appear at zero- or near-zero energies depends on the coexistence of surface modes and the hybridization strength. Although the pairing symmetry is yet to be determined experimentally, given the observed superconductivity and our prediction of possible hybrid-order Tsc, we urge further experimental efforts in investigating the superconducting properties of \(1{T}^{{\prime} }\)-MoTe2. We expect that the predicted boundary signatures can be detected by scanning tunneling microscope or through transport measurements.

Methods

Topological crystal approach

In this section, we explain how we perform the Topological Crystal Approach in Step 1. First, we need to perform a real-space cell-decomposition to break the full unit cell down to 0D, 1D, 2D, and 3D building blocks that do not respect any non-local crystalline symmetries (see Fig. 4). We denote the dimension of the building block as db. The db = 3 building block is the 3d time-reversal invariant topological superconductor in AZ class DIII (3d TSC). Decorating the 3-cells with this building block simply gives the usual 3d TSC with a non-trivial strong indicator κstrong = 1. The db = 2 building block is the 2d time-reversal invariant topological superconductor in AZ class DIII (2d TSC). 2d TSCs will be decorated on the 2-cells. The db = 1 building block is the 1d time-reversal invariant Kitaev chain (1d TSC), which will be decorated on the 1-cells. Finally, we also have db = 0 building block, which is described by a 0d BdG Hamiltonian. The resulting topological crystals and the superconductors that are adiabatically connected to these states are regarded as atomic superconductors (ASC), which are superconducting analog of atomic insulators. We view such ASC as topologically trivial because they do not host topologically protected boundary zero modes on open geometries. By quotienting out the ASC, we find that the classification of topological superconductors with non-trivial boundary modes are given by \(({\mathbb{Z}}\times {{\mathbb{Z}}}_{4})\times {{\mathbb{Z}}}_{4}^{3}\times {{\mathbb{Z}}}_{2}^{3}\). Table 2 shows the decoration patterns of the topological crystals and their corresponding symmetry indicators. The \(({\mathbb{Z}}\times {{\mathbb{Z}}}_{4})\) factor contains the strong first, second and third order TSCs. We denote the strong first, second, and third order TSCs by (a, b, c), where a is an integer corresponding to the strong first-order TCS only protected by the internal symmetry in class DIII, b and c are \({{\mathbb{Z}}}_{2}\) number corresponding to the strong second and third order TCSs (the first two entries in Table 2). Physically, a is characterized by the 3D winding number. They satisfy the following non-trivial stacking relations57:

$$\begin{array}{l}(1,0,0)+(1,0,0)\,\,\,\,\,=(2,1,0),\\ (1,0,0)+(-1,0,0)\,=(0,1,0),\\ (0,1,0)+(0,1,0)\,\,\,\,\,\,=(0,0,1).\end{array}$$
(11)

Note that the phases with even winding numbers are completely decoupled from the higher order phases, i.e., we can freely stack the phases labeled by (2n, 0, 0) without affecting the higher order phases. One can check that the tuple (a, b, c) satisfying Eq. (11) is isomorphic to \(({\mathbb{Z}}\times {{\mathbb{Z}}}_{4})\). If we label the group element in \(({\mathbb{Z}}\times {{\mathbb{Z}}}_{4})\) as (g, h), the generator of \({\mathbb{Z}}\)(1, 0) corresponds to the 3d TSC with winding number 1: (1, 0, 0), and it’s inverse element (−1, 0) corresponds the phase (−1, 1, 1). The \({{\mathbb{Z}}}_{4}\) is generated by the second order phase (0, 1, 0). Due to the non-trivial stack rules, the phase labeled by (2, 0, 0) is in fact the (2, − 1) element in the group \(({\mathbb{Z}}\times {{\mathbb{Z}}}_{4})\). To better reveal the higher-order topology, it’s convenient to quotient out the subgroup generated by (2, 0, 0) and the resulting group is \({\mathfrak{C}}={{\mathbb{Z}}}_{8}\) labeled by (a, b, c) but now a is a \({{\mathbb{Z}}}_{2}\) number with the stacking rule (1, 0, 0) + (1, 0, 0) ≅ (0, 1, 0).

Fig. 4: The cell decomposition for the Topological Crystal approach.
figure 4

Cell decomposition of the unit cell for space group \(P\bar{1}\) (#2). The origin is chosen at the center of the unit cell. Colored faces are inequivalent 2-cells. Bold blue lines are inequivalent 1-cells. Red dots are inequivalent 0-cells.

Full size image
Table 2 The topological crystals and symmetry indicators
Full size table

Momentum space topological invariants

In this section, we discuss our calculation and results in Step 2, where we obtain the momentum-space topological invariants by calculating the equivariant K group \({}^{\phi }{K}_{G}^{(\tau ,c),-3}(BZ)\) using the Atiyah-Hirzebruch Spectral Sequence (AHSS)13,15,59. Please see refs. 13,15,59 for an introduction to the well-developed Equivariant K Theory and AHSS. Here, we present the essential results of our calculation. The elements \({E}_{2}^{p,-n}\) in the E2 page of the AHSS are summarized in Table 3.

Table 3 The E2 page we find from our calculation
Full size table

To show that the E2 page is the limiting page E∞, one has to calculate the third differential, which is a non-trivial task. We give a physical argument below to show that the E2 page is the limiting page.

Each element in the E2 page can be characterized by topological invariants defined in the subspaces of BZ. Topological invariants on 0-cells are classified by \({E}_{2}^{0,-3}={{\mathbb{Z}}}^{8}\). The explicit form of each \({\mathbb{Z}}\) invariant is defined in Eq. (4). By quotienting out the contribution of ASC, they give the symmetry indicators. There is a \({\mathbb{Z}}\) invariant defined on the 3-cell. This can be naturally identified with the 3D winding number N318 that characterized the 3D TSC in class DIII. It has been shown that the parity of the strong symmetry indicator \({\kappa }_{{{{\rm{sc}}}}}^{s}\) agrees with the 3D winding number N3 modulo 247,48,49:

$${\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\,\,{{\mbox{mod}}}\,\,2={N}_{3}\,\,{{\mbox{mod}}}\,\,2$$
(12)

This is consistent with the conjecture that the E2 page is the limiting page. Moreover, since there is a non-trivial relation between \({\kappa }_{{{{\rm{sc}}}}}^{s}\) and N3, the formation group extension that we use to obtain the full \({{{\mathcal{K}}}}\) group has to be non-trivial:

$$1\to {\mathbb{Z}}\to {{{\mathcal{K}}}}\to {{\mathbb{Z}}}^{8}\to 1.$$
(13)

Upon quotienting out the contribution of ASC, the \({{{\mathcal{K}}}}\) group should agree with the classification from the Topological Crystal Approach: \(({\mathbb{Z}}\times {{\mathbb{Z}}}_{4})\times {{\mathbb{Z}}}_{4}^{3}\times {{\mathbb{Z}}}_{2}^{3}\) in real space. The mixed SIs \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},x}\), \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},y}\), \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{m}}}},z}\) and weak indicators \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},xy}\), \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},yz}\), \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{w}}}},xz}\) simply correspond to the \({{\mathbb{Z}}}_{4}^{3}\times {{\mathbb{Z}}}_{2}^{3}\) factor. The \(({\mathbb{Z}}\times {{\mathbb{Z}}}_{4})\) factor that corresponds to the strong phases requires further discussion. The strong SI \({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\) itself is a \({{\mathbb{Z}}}_{8}\) number. However, strong phases are, in fact, characterized by a pair of invariants \(({N}_{3},{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}})\) with the constraint Eq. (12). Taking into account the constraint in Eq. (12), the pair \(({N}_{3},{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}})\) can be parametrized as \((2n+({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\,\,{{\mbox{mod}}}\,\,2),{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}})\), where n is an integer. One can check that the pair \((2n+({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\,\,{{\mbox{mod}}}\,\,2),{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}})\) indeed satisfies the group multiplication rule of \(({\mathbb{Z}}\times {{\mathbb{Z}}}_{4})\). If we quotient out the subgroup generated by an even winding number \(({N}_{3},{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}})=(2n,0)\), we obtain \((({\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}}\,\,{{\mbox{mod}}}\,\,2),{\kappa }_{{{{\rm{sc}}}}}^{{{{\rm{s}}}}})\in {{\mathbb{Z}}}_{8}\).

Reference Hamiltonian H ref

In this section, we review why a reference Hamiltonian Href is needed8,13 and how we make the choice of Href. In K theory, a K group can be viewed as a formal difference between two vector bundles B1 and B2. In Karoubi’s formulation, the difference is represented by a triple [B, H1, H2], where B is the vector bundle whose base space is the BZ and vector space is formed by the occupied states of the Hamiltonian Hi, and H1 and H2 are the flattened Hamiltonians.

We can further associate different equivalence classes of triples [B, H, Href] with distinct gapped phases of matter. To do so, instead of using different reference Hamiltonians H2 for different triples, it is crucial to define a ‘trivial’ Hamiltonian Href as the universal reference Hamiltonian (i.e., set H2 = Href) for all triples. For superconductors, a natural choice for the reference Hamiltonian is BdG Hamiltonian formed by a vacuum state,

$${H}_{{{{\rm{ref}}}}}={{{\rm{diag}}}}({{\mathbb{I}}}_{N},-{{\mathbb{I}}}_{N}),$$
(14)

where N is the number of normal bands. In Eq. (4), we always use Eq. (14) as the reference Hamiltonian.

Basis ambiguity of symmetry indicators

In this section, we discuss the basis ambiguity in the calculation of SIs and how to determine the basis for SIs such that the SIs have a transparent correspondence to the Majorana boundary signatures. SIs are defined as elements in the quotient group \(X={{{{\mathcal{K}}}}}^{{\prime} }/\{\,{{\mbox{AS}}}\,\}\), where \({{{{\mathcal{K}}}}}^{{\prime} }={{\mathbb{Z}}}^{8}\) is the classification of the topological invariants defined at TRIMs (the winding number is not included here since the atomic superconductors have a zero winding number), and \(\{\,{{\mbox{AS}}}\,\}={\mathbb{Z}}\times 8{\mathbb{Z}}\times {(4{\mathbb{Z}})}^{3}\times {(2{\mathbb{Z}})}^{3}\) is the classification group of the topological invariants for the atomic superconductors. More specifically, it can be written as a matrix:

$${M}_{{{{\rm{AS}}}}}=\left(\begin{array}{cccccccc}{{{{\boldsymbol{a}}}}}_{1}&{{{{\boldsymbol{a}}}}}_{2}&{{{{\boldsymbol{a}}}}}_{3}&{{{{\boldsymbol{a}}}}}_{4}&{{{{\boldsymbol{a}}}}}_{5}&{{{{\boldsymbol{a}}}}}_{6}&{{{{\boldsymbol{a}}}}}_{7}&{{{{\boldsymbol{a}}}}}_{8}\end{array}\right),$$
(15)

where each column vector contains the set of 0d invariants nk at TRIMs generated by a atomic superconductor sitting at a Wyckoff position in the real space. The explicit matrix form can be found in Ref. 6. To proceed, we compute the Smith normal form to find the linearly independent bases:

$$U{M}_{{{{\rm{AS}}}}}V=\lambda$$
(16)

where U and V are the transformation matrices for the momentum-space and real-space bases respectively, and λ is a diagonal matrix:

$$\begin{array}{r}\lambda =\left(\begin{array}{cccc}1&&&\\ &2{I}_{3}&&\\ &&2{I}_{3}&\\ &&&8\end{array}\right).\end{array}$$
(17)

From the fact that \({M}_{{f}_{0}}V={U}^{-1}\lambda\), we can now extract the linearly independent basis. Specifically, the new real-space basis vectors are given by

$$\begin{array}{r}{M}_{{{{\rm{AS}}}}}V=\left(\begin{array}{cccccccc}{{{{\boldsymbol{a}}}}}_{1}^{{\prime} }&{{{{\boldsymbol{a}}}}}_{2}^{{\prime} }&{{{{\boldsymbol{a}}}}}_{3}^{{\prime} }&{{{{\boldsymbol{a}}}}}_{4}^{{\prime} }&{{{{\boldsymbol{a}}}}}_{5}^{{\prime} }&{{{{\boldsymbol{a}}}}}_{6}^{{\prime} }&{{{{\boldsymbol{a}}}}}_{7}^{{\prime} }&{{{{\boldsymbol{a}}}}}_{8}^{{\prime} }\end{array}\right),\end{array}$$
(18)

where \(\{{{{{\boldsymbol{a}}}}}_{i}^{{\prime} }\}\) are column vectors rotated by the transformation matrix V from {ai}. The new momentum-space basis vectors are given by

$${U}^{-1}=\left(\begin{array}{cccccccc}{{{{\boldsymbol{b}}}}}_{1}^{{\prime} }&{{{{\boldsymbol{b}}}}}_{2}^{{\prime} }&{{{{\boldsymbol{b}}}}}_{3}^{{\prime} }&{{{{\boldsymbol{b}}}}}_{4}^{{\prime} }&{{{{\boldsymbol{b}}}}}_{5}^{{\prime} }&{{{{\boldsymbol{b}}}}}_{6}^{{\prime} }&{{{{\boldsymbol{b}}}}}_{7}^{{\prime} }&{{{{\boldsymbol{b}}}}}_{8}^{{\prime} }\end{array}\right),$$
(19)

where \(\{{{{{\boldsymbol{b}}}}}_{j}^{{\prime} }\}\) are column vectors rotated by U−1 from {bj} at j = TRIMs. Since the two sets of new bases are related by

$${{{{\boldsymbol{a}}}}}_{i}^{{\prime} }={{{{\boldsymbol{b}}}}}_{i}^{{\prime} }{\lambda }_{i},$$
(20)

where λi denotes the diagonal element of λ, we can span the 0D invariant group for atomic superconductors {AS} and \({{{\mathcal{K}}}}\) in the same set of linearly independent bases.

The explicit form of the symmetry indicators is given by

$${{{\boldsymbol{\kappa }}}}\equiv U\bar{{{{\bf{n}}}}},$$
(21)

where \(\bar{{{{\bf{n}}}}}\) is the set of 0d invariants. There is however a basis ambiguity in calculating the Smith normal form:

$$\lambda =U{M}_{{{{\rm{AS}}}}}V$$
(22)
$$=(U{L}^{-1})(L{M}_{{{{\rm{AS}}}}}{R}^{-1})(RV),$$
(23)
$$=\tilde{U}{\tilde{M}}_{{{{\rm{AS}}}}}\tilde{V}.$$
(24)

While the symmetry indicator group remains unchanged, the explicit form of the symmetry indicators is now given by

$$\tilde{{{{\boldsymbol{\kappa }}}}}\equiv \tilde{U}\bar{{{{\bf{n}}}}}.$$
(25)

Therefore, there is no unique explicit expression for the symmetry indicators without further input. To fix a canonical basis, we need to match with the real-space classification, and we choose a basis such that the strong, mixed, and weak phases are all separated.

Ab-initio calculations

Our DFT calculation is performed for centrosymmetric bulk β-MoTe2 (nonsymmorphic space group #11 P21/m, point group C2h) with the experimental geometry60, in the absence of spin-orbit coupling, using VASP61,62. The experimental lattice constants and angles are as follows60: a = 6.330, b = 3.469, c = 13.860 angstrom, β = 93.55∘, and α = γ = 90∘. The real space lattice vectors are \({{{{\bf{a}}}}}_{1}=a{\hat{e}}_{x}\), \({{{{\bf{a}}}}}_{2}=b{\hat{e}}_{y}\), and \({{{{\bf{a}}}}}_{3}=c\cos \beta {\hat{e}}_{x}+c\sin \beta {\hat{e}}_{z}\), where \({\hat{e}}_{x,y,z}\) are unit vectors in Cartesian coordinates. There are four Mo atoms and eight Te atoms in a monoclinic primitive unit cell. The inversion center of our atomic coordinates in the primitive unit cell is set to the origin. We use projector-augmented wave (PAW) pseudopotentials69 within the Perdew-Burke-Ernzerhof (PBE) generalized gradient approximation70. Each Mo atom has 6 valence electrons which are nominally singly occupied at the five 4d orbitals and at the 5s orbital. Each Te atom has 6 valence electrons which are nominally occupied at the 5p orbitals and 5s orbitals. The Hubbard U value of 3.0 eV is used for the Mo 4d orbitals, within the DFT+U method63 as implemented in VASP, following ref. 40. We sample k points of 10 × 20 × 5 with Γ point centered for the self-consistent calculation. The cutoff of the kinetic energy is set to 400 eV. We consider 54 bands.

In order to compute the topological indices Z8, Z4, and Z213 or the indices in Eq. (5), we calculate parity values of all bands at the eight TRIM k points, from the wave function of the self-consistent calculation using two different codes (Wang, D., https://github.com/obaica/vasp_wavecar_parity)71. The TRIM points are Γ = (0.0, 0.0, 0.0), X = (0.5, 0.0, 0.0), U = (0.5, 0.0. 0.5), Z = (0.0, 0.0, 0.5), Y = (0.0, 0.5, 0.0), S = (0.5, 0.5, 0), T = (0.0, 0.5, 0.5), and R = (0.5, 0.5, 0.5). Here the coordinates of the k points are in terms of the reciprocal lattice vectors \({{{{\bf{b}}}}}_{1}=({\hat{e}}_{x}-cot\beta {\hat{e}}_{z})2\pi /a\), \({{{{\bf{b}}}}}_{2}=2\pi /b{\hat{e}}_{y}\), and \({{{{\bf{b}}}}}_{3}=2\pi /c{\hat{e}}_{z}\). Note that the MoTe2 crystal has nonsymmorphic group. The atoms in the crystal have twofold screw symmetry along the y axis, t(b2/2)C2y (where t(b2/2) is a translation along the y or b2 axis by b2/2), mirror plane t(b2/2)σxz about the xz plane, and inversion symmetry. The twofold screw symmetry gives twofold degeneracy at k = π/b plane in the absence of spin-orbit coupling and each degenerate band consists of a band with positive parity and a band with negative parity72. Therefore, there is always twofold degeneracy at the Y, R, S, and T points with an equal number of positive-parity bands and negative-parity bands. For visualization of the Fermi surfaces, we use the c2x program73 and the XCrysDen program74.

Construction of the 44-band tight binding model

From the above VASP calculation, we then compute hopping integrals and construct a tight-binding model based on 44 Wannier functions (WFs), using WANNIER90 code version 1.264,65,66. The 44 Wannier functions consist of dxy, dyz, dxz, \({d}_{{x}^{2}-{y}^{2}}\), and \({d}_{{z}^{2}}\) orbitals of all four Mo sites and px, py, and pz orbitals of all eight Te sites in the unit cell. We use the same number of k points as the VASP calculation and set the minimum energy of disentanglement as zero. Only disentanglement65 is applied without maximum localization of the WFs. All disentagled WFs are centered at the atomic sites. We exclude the bottom most eight bands from the VASP calculation in order to generate the WF-44 tight-binding model. These eight bands have the same numbers of positive parity bands and negative parity bands. Fewer numbers of WFs than 44 orbitals would produce neither atomic-orbital-shaped WFs nor poor agreement with the VASP band structure. Figure 5a shows the comparison between the DFT band structure (solid black) and the 44-band tight-binding model (dashed red). They are in good agreement with each other. This band structure is similar to the band structure with spin-orbit coupling (Fig. 5b).

Fig. 5: Detailed band structures of 1T′-MoTe2.
figure 5

(Color online) (a) Band structure obtained using VASP (solid black) and the 44-band tight binding model (dashed red) without spin-orbit coupling where the Fermi level is 7.1515 eV. b Band structure using VASP with spin-orbit coupling where the Fermi level is 7.1544 eV. In a and b, the Fermi level is indicated as the horizontal line.

Full size image