Emergence of topological superconductivity in doped topological Dirac semimetals under symmetry-lowering lattice distortions

Abstract

Recently, unconventional superconductivity having a zero-bias conductance peak is reported in doped topological Dirac semimetal (DSM) with lattice distortion. Motivated by the experiments, we theoretically study the possible symmetry-lowering lattice distortions and their effects on the emergence of unconventional superconductivity in doped topological DSM. We find four types of symmetry-lowering lattice distortions that reproduce the crystal symmetries relevant to experiments from the group-theoretical analysis. Considering inter-orbital and intra-orbital electron density-density interactions, we calculate superconducting phase diagrams. We find that the lattice distortions can induce unconventional superconductivity hosting gapless surface Andreev bound states (SABS). Depending on the lattice distortions and superconducting pairing interactions, the unconventional inversion-odd-parity superconductivity can be either topological nodal superconductivity hosting a flat SABS or topological crystalline superconductivity hosting a gapless SABS. Remarkably, the lattice distortions increase the superconducting critical temperature, which is consistent with the experiments. Our work opens a pathway to explore and control pressure-induced topological superconductivity in doped topological semimetals.

Introduction

Topological insulator, Dirac semimetal (DSM), Weyl semimetal, and topological superconductor are newly established quantum states of matter which are expected to have applications for dissipationless devices and quantum information technologies^{1,2,3,4,5,6,7}. Among them, topological Weyl and Dirac semimetals are characterized by relativistic quasi-particles and gapless nodes in bulk spectra^{3,6,8,9,10,11}. Because of their anomalous electromagnetic responses and topologically-protected surface Fermi arcs on the boundaries, such topological semimetals have been attracted much attention^{6,8,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26}. Moreover, due to the unique properties of Dirac and Weyl semimetals, extensive theoretical and experimental studies of their superconducting instabilities have been conducted to observe possible topological superconductivity^5,7.

Recently, the lattice-distortion induced superconductivity in DSMs of $\hbox {Cd}_3\hbox {As}_2$^27,28,29 and $\hbox {Au}_2$Pb^{30,31,32,33,34} is reported. For $\hbox {Cd}_3\hbox {As}_2$, it does not show any superconductivity at the ambient pressure until 1.8 K^27,28,29. The structural phase transition occurs near 2.6 GPa from a tetragonal lattice with $D_{4h}$ point group symmetry ($I4_1/acd$) to a monoclinic lattice with $C_{2h}$ point group symmetry ($P2_1/c$). Then, superconductivity emerges at $T_c \approx 1.8$ K under pressure higher than 8.5 GPa. When the pressure increases further, $T_c$ keeps increasing from 1.8 K to 4.0 K in the hydrostatic pressure experiment²⁸. Similarly, $\hbox {Au}_2$Pb shows superconductivity at $T_c \approx 1.2$ K after a structural phase transition from the cubic with $O_h$ symmetry (Fd3m) to the orthorhombic lattice with $D_{2h}$ symmetry (Pbnc)^30,32,34. $T_c$ increases up to 4 K at 5 GPa, then decreases with further compression³⁴. For both materials, the point-contact measurements reported that measured $T_c$ using a hard contact tip is much higher than the measured $T_c$ using a soft tip^27,29,32. The point-contact measurements for $\hbox {Cd}_3\hbox {As}_2$ showed the zero-bias conductance peak (ZBCP) and double conductance peaks symmetric around zero bias, which was interpreted as a signal of a topological Majorana surface state^27,29. Moreover, the transport data under magnetic fields reported anomalous behaviors that the conventional BCS theory cannot explain^27,29,32. At ambient pressure, the proximity-induced superconductivity in $\hbox {Cd}_3\hbox {As}_2$ is also reported³⁵.

In parallel to the experimental exploration of the superconductivity in doped DSM, several theoretical studies were conducted^36,37. In the absence of lattice distortion, the possible superconducting states in doped DSM are suggested as either fully-gapped superconductor (FGSC) or topological nodal superconductor (TNSC) hosting a flat surface Andreev bound state (SABS) on the boundary³⁷. In experiments, however, superconductivity was observed only in the presence of lattice distortion. Considering a lattice distortion (in our work, $n_1$ type lattice distortion), the topological crystalline superconductor (TCSC) hosting surface Majorana states was proposed³⁶. However, because such lattice distortion results in the orthorhombic lattice, it cannot be applied to the observed superconductivity in the monoclinic crystal structure of $\hbox {Cd}_3\hbox {As}_2$²⁸. It is, therefore, necessary to study the effect of symmetry-lowering lattice distortions on the emergence of unconventional superconductivity in doped DSM.

In this work, we systematically study possible symmetry-lowering lattice distortions and their effects on the emergence of unconventional superconductivity in doped topological DSM. As a representative model, we consider a topological DSM described by the four-band Hamiltonian having $D_{4h}$ point group symmetry in the absence of lattice distortions. While keeping time-reversal symmetry (TRS) and inversion symmetry (IS), we find four types of symmetry-lowering lattice distortions from the group-theoretical analysis, which are denoted as $n_i$ type lattice distortions ($i=1,\ldots ,4$). Two of them ($n_1$ and $n_2$ type) reduce $D_{4h}$ of the tetragonal lattice to $D_{2h}$ orthorhombic lattice, while the others ($n_3$ and $n_4$ type) transform the tetragonal lattice to $C_{2h}$ of the monoclinic lattice. They explain the structural phase transition in $\hbox {Cd}_3\hbox {As}_2$ and $\hbox {Au}_2$Pb under pressure. The symmetry-lowering lattice distortions are summarized in Table 3.

To understand the emergence of superconductivity under lattice distortions, we adopt the Bogoliubov-de Gennes (BdG) formalism and linearized gap equation, and we assume intra-orbital (U) and inter-orbital (V) electron density-density interactions which induce superconducting instabilities. From the Fermi-Dirac statistics, six possible momentum-independent superconducting pairing potentials are found³⁷. Under lattice distortions, six pairings potentials are classified according to the irreducible representation of the remaining point symmetry group. Using these pairing potentials, possible superconducting gap structures and superconducting critical temperatures ($T_c$) are calculated. By comparing critical temperatures, we obtain the superconducting phase diagram, and the dominant superconducting phases are discovered, such as fully-gapped superconductor (FGSC), topological nodal superconductor (TNSC), and topological crystalline superconductor (TCSC) depending on the lattice distortions and the ratio of U/V. Among them, FGSC is conventional superconductor, while TNSC and TCSC are unconventional.

Interestingly, the unconventional superconductors of TNSC and TCSC emerge when inter-orbital interaction V and the strength of lattice distortion are large enough while FGSC emerges in the opposite limit. Therefore, the lattice distortion and inter-orbital interaction act as physical parameters that control the phase transition between conventional and unconventional superconductivity of a topological DSM. We find that both V and lattice distortions enhance the unconventional superconducting pairings via a unique spin-orbit locking. Moreover, $T_c$ increases under the lattice distortions due to the enhancement of DOS at the Fermi surface, which is consistent with the experimentally measured $T_c$ enhancement under pressure. The unconventional superconductors host gapless SABS in mirror plane even under the lattice distortions: Under the $n_1$ or $n_2$ type lattice distortion, the superconductivity in the orthorhombic lattice with $D_{2h}$ point group symmetry hosts a gapless SABS protected by the mirror Chern number. Under the $n_3$ or $n_4$ type lattice distortion, the superconductivity in the monoclinic lattice with $C_{2h}$ point group symmetry hosts a gapless SABS protected by the unbroken mirror symmetry and a flat SABS protected by the mirror chiral winding number in specific conditions. Because there exist gapless Majorana surface states under the lattice distortions, we suggest that these states can be observed in scanning tunneling microscope (STM) or point contact Andreev reflection spectroscopy experiments.

Consequently, our theoretical work is consistent with the discovered structural phase transition and the enhancement of superconductivity in $\hbox {Cd}_3\hbox {As}_2$ and $\hbox {Au}_2$Pb under lattice distortions. Moreover, we suggest that the emergence of conventional and unconventional superconductivity in doped topological DSM can be controlled by the pressure and strength of the superconducting pairing interaction. Therefore, our woks opens a pathway to explore and control the topological superconductors in doped topological semimetals, which may have future applications in dissipationless and quantum information devices.

Results

Undistorted Dirac semimetal

Dirac semimetal (DSM) has the low energy excitations near the Fermi-level described by a massless Dirac equation. Because all bands are doubly degenerate due to the TRS and IS, a DSM is minimally described by a four-band Hamiltonian^{6,10,38,39,40}. However, TRS and IS are not enough to protect a fourfold degeneracy, so the symmetry-protected DSM is suggested, where the Dirac points are protected by TRS, IS and crystalline symmetries^{6,10,38,39,40}. DSMs are reported in many materials such as $\beta$-cristobalite $\hbox {BiO}_2$¹⁰, distorted spinels⁴¹, $\hbox {Na}_3$Bi^42,43, $\hbox {Cd}_3\hbox {As}_2$^{42,44,45,46,47,48,49}, $\hbox {Au}_2$Pb^30,50, and $\hbox {ZrTe}_5$^51,52. Among them, superconductivity is reported in $\hbox {Cd}_3\hbox {As}_2$^27,28,29 and $\hbox {Au}_2$Pb^{30,31,32,33,34}. Both materials have Dirac points protected by TRS, IS, and $C_{4}$ rotational symmetry and share the tetragonal crystal system with $D_{4h}$ point group symmetry. For this reason, we consider the undistorted topological DSM having a $D_{4h}$ point group symmetry as a representative model system.

Model and symmetry

The general $4 \times 4$ Hamiltonian representation is

$$\begin{aligned} H({\mathbf {k}}) = \sum _{i=1}^{16} a_{i} ({\mathbf {k}}) \Gamma _i. \end{aligned}$$

(1)

The coefficient function $a_{i} ({\mathbf {k}})$ are real functions and $\Gamma _i = s_j \sigma _k$ are $4 \times 4$ gamma matrices where $s_{j}$ and $\sigma _{k}$ are Pauli matrices for spin and orbital degrees of freedom in the spin $(\uparrow , \downarrow )$ and the orbital (1, 2) spaces, respectively.

The symmetry constraints can simplify the Hamiltonian’s form in Eq. (1). Due to TRS and IS, the Hamiltonian satisfies the following equations:

$$\begin{aligned} T H({\mathbf {k}}) T^{-1} = H ( - {\mathbf {k}}), ~~~~~~ P H( {\mathbf {k}}) P^{-1} = H ( - {\mathbf {k}}), \end{aligned}$$

(2)

where $T = i s_y {\hat{K}}$ is the time-reversal operator (${\hat{K}}$ is the complex conjugation operator) and P is the inversion operator. Because the inversion does not flip the spin, the inversion operator has orbital dependency only, and it can be chosen as $P= - \sigma _z$ for topological DSM without loss of generality^39,45. Then, due to TRS and IS, among sixteen $\Gamma _i$ matrices, only six $\Gamma _i$ matrices are allowed. They are $\Gamma _0 = {\mathbf {1}}_{4 \times 4}$, $\Gamma _1 = \sigma _x s_z$, $\Gamma _2 = \sigma _y s_0$, $\Gamma _3 = \sigma _x s_x$, $\Gamma _4 = \sigma _x s_y$, and $\Gamma _5 = \sigma _z s_0$. We set $a_0 ({\mathbf {k}}) = 0$ since it does not contribute to the formation of Dirac points^39,45.

Table 1 Transformation properties of gamma matrices under symmetry operations.

The $D_{4h}$ point group symmetry imposes more constraints on the Hamiltonian’s form in Eq. (1). The generators of $D_{4h}$ point group can be chosen as inversion P, fourfold rotation about the z axis $C_{4z}$, and twofold rotation about the x axis $C_{2x}$. Their matrix representations are chosen as

$$\begin{aligned} P = - \sigma _z,~~~~~ C_{4z} = \exp ( - i \frac{\pi }{2}s_z - i \frac{\pi }{4} \sigma _z s_z ),~~~~~ C_{2x} = i \sigma _z s_x, \end{aligned}$$

(3)

where we adopt the following basis set known to describe the low-energy effective Hamiltonian of $\hbox {Cd}_2\hbox {As}_3$⁴⁵.

$$\begin{aligned} \left| 1, \uparrow \right\rangle = \left| P_{J=\frac{3}{2}} , 3/2\right\rangle ,~~~~ \left| 1,\downarrow \right\rangle =\left| P_{J=\frac{3}{2}} , -3/2\right\rangle ,~~~~ \left| 2, \uparrow \right\rangle =\left| S_{J=\frac{1}{2}} , 1/2\right\rangle ,~~~~ \left| 2,\downarrow \right\rangle =\left| S_{J=\frac{1}{2}} , -1/2\right\rangle , \end{aligned}$$

(4)

where J is the total angular momentum. Other rotation and mirror symmetries are given by $C_{2z} = i \sigma _z s_z$, $M_{xy} = - i s_z$, $M_{yz} = - i s_x$, $M_{zx} = - i \sigma _z s_y$, $M_{(1 1 0)} = i ( \sigma _z s_x - s_y ) / \sqrt{2}$, and $M_{(1 {{\bar{1}}} 0)} = i ( \sigma _z s_x + s_y ) / \sqrt{2}$. The subscript in each mirror operator represents the corresponding mirror plane by using either Cartesian coordinates or Miller indices. The group elements are derived in Sec. S1 in Supplementary Information. Due to this $D_{4h}$ symmetry, the Hamiltonian in Eq. (1) satisfy

$$\begin{aligned} U H({\mathbf {k}}) U^{-1} = H ( S {\mathbf {k}} ), \end{aligned}$$

(5)

where U and S are transformation matrices for an element of $D_{4h}$ group in the spin-orbital and momentum spaces, respectively. For the group generators, the Hamiltonian in Eq. (1) satisfies

$$\begin{aligned} P H({\mathbf {k}}) P^{-1} = H( - {\mathbf {k}}),~~~~~ C_{4z} H({\mathbf {k}}) C_{4z}^{-1} = H( {\mathscr {R}}_{4z} {\mathbf {k}}),~~~~~ C_{2x} H({\mathbf {k}}) C_{2x}^{-1} = H( {\mathscr {R}}_{2x} {\mathbf {k}}), \end{aligned}$$

(6)

where ${\mathscr {R}}_{4z} {\mathbf {k}} = (-k_y, k_x, k_z)$ and ${\mathscr {R}}_{2x} {\mathbf {k}} = (k_x, -k_y, -k_z)$. Because the transformation properties of gamma matrices are given by Table 1, Eq. (5) imposes constraints to each coefficient functions $a_i ({\mathbf {k}})$, which is summarized in Table 2. Therefore, the general form of the Hamiltonian of DSM having $D_{4h}$ point group symmetry is obtained.

Table 2 Symmetry constraints on $a_i({\mathbf {k}})$.

Lattice model

For concreteness, we construct an explicit lattice model that describes a class of Dirac semimetals such as $\hbox {Cd}_3\hbox {As}_2$ and $\hbox {Au}_2$Pb. The coefficient functions of Hamiltonian in Eq. (1) are given by^39,45

$$\begin{aligned} a_1 ({\mathbf {k}})&= v \sin k_x, \end{aligned}$$

(7)

$$\begin{aligned} a_2 ({\mathbf {k}})&= v \sin k_y, \end{aligned}$$

(8)

$$\begin{aligned} a_3 ({\mathbf {k}})&= ( \beta + \gamma ) \sin k_z (\cos k_y - \cos k_x), \end{aligned}$$

(9)

$$\begin{aligned} a_4 ({\mathbf {k}})&= - ( \beta - \gamma ) (\sin k_z \sin k_x \sin k_y), \end{aligned}$$

(10)

$$\begin{aligned} a_5 ({\mathbf {k}})&= M' - t_{xy}(\cos k_x + \cos k_y ) - t_z \cos k_z, \end{aligned}$$

(11)

where $M'$, $t_{xy}$, $t_z$, v, $\beta$, and $\gamma$ are material-dependent parameters. The energy eigenvalues are given by

$$\begin{aligned} E = \pm \left| a({\mathbf {k}})\right| = \pm \left( \sum _{i=1}^5 a_{i}^2 ({\mathbf {k}}) \right) ^{1/2}. \end{aligned}$$

(12)

If $t_z> (M'- 2 t_{xy}) >0$, the Hamiltonian hosts a pair of Dirac points at $(0, 0, \pm k_0)$ as shown in Fig. 1a. Here, $k_0$ is determined by $M^{\prime }{} - 2 t_{xy} - t_z \cos k_0 =0$. These Dirac points are protected by the $C_{4z}$ symmetry³⁹. Due to the $C_{4z}$, the four bands on the $k_z$ axis can have different $C_{4z}$ eigenvalues, which lead to fourfold degenerate Dirac points.

Low-energy effective Hamiltonian

Near the Dirac points $(0, 0, \pm k_0)$, the low-energy effective Hamiltonian takes the form of Dirac Hamiltonian, which is given by

$$\begin{aligned} H_{\text {Dirac}}^{(\pm )} = v k_x \Gamma _1 + v k_y \Gamma _2 \pm v_z (k_z \mp k_0) \Gamma _5. \end{aligned}$$

(13)

where $v_z = t_z k_0$. The energy spectrum shows anisotropic energy-momentum dispersion, which is given by

$$\begin{aligned} E = \pm \sqrt{ v^2 (k_x ^2 + k_y)^2 + v_z^2 (k_z \mp k_0)^2 }. \end{aligned}$$

(14)

Figure 1

Crystal systems, band structures, and Fermi surfaces of Dirac semimetal (DSM) under various lattice distortions. (a) Undistorted DSM for comparison. It has a tetragonal lattice. (b–e) Distorted crystal systems under (b) $n_1$, (c) $n_2$, (d) $n_3$, and (e) $n_4$ type lattice distortions. In (b) and (c), $n_1$ and $n_2$ type lattice distortions changes inplane lattice constants, which results in orthorhombic lattices. In (d) and (e), $n_3$ and $n_4$ type lattice distortions change the $\alpha$ and $\beta$ angles, which results in monoclinic lattices. (f–j) The corresponding 3D band structures. In (f–i) [(j)], the band structures are plotted for the $k_y$-$k_z$ ($k_x$-$k_z$) plane and the orange planes are $k_y=0$ ($k_x=0$) plane. (k–o) The corresponding Fermi surfaces. In (l–o), all Fermi surfaces are distorted according to types of lattice distortions. In (n) and (o), the Fermi surfaces are shifted as indicated by the black arrows. Each vertical orange line indicates the $k_z$ axis.

Distorted Dirac semimetal

Symmetry-lowering distortions

In the absence of lattice distortions, $\hbox {Cd}_3\hbox {As}_2$^27,28,29 and $\hbox {Au}_2$Pb^32,33,34 share the same $D_{4h}$ point group symmetry and show no superconductivity. However, both materials showed superconductivity after the structural phase transition under pressure or cooling, and the superconducting critical temperature increases with the pressure^28,34. At the high pressure, $\hbox {Cd}_3\hbox {As}_2$ becomes a monoclinic lattice having $C_{2h}$ point group symmetry²⁸ and $\hbox {Au}_2$Pb becomes an orthorhombic lattice having $D_{2h}$ point group symmetry³². Thus, IS is preserved even under lattice distortions. In addition, the superconductivity appears under the small lattice distortions in the hydrostatic experiments^28,34. Therefore, we assume that both TRS and IS are preserved under lattice distortions and the effect of the lattice distortion can be implemented as a perturbation⁵³.

We now classify the possible symmetry-lowering lattice distortions. The form of the perturbation Hamiltonian for the lattice distortions is given by

$$\begin{aligned} H_{\text {pert}} = \sum _{i=0}^{5} d_{i}({\mathbf {k}}) \Gamma _i, \end{aligned}$$

(15)

where $d_{i}({\mathbf {k}})$ is a real-valued function of momentum and $\Gamma _i$ is the gamma matrix. Because $\Gamma _{1}$, $\Gamma _{2}$, $\Gamma _{3}$, and $\Gamma _{4}$ are odd under T and P, the coefficient functions $d_{1}({\mathbf {k}})$, $d_{2}({\mathbf {k}})$, $d_{3}({\mathbf {k}})$, and $d_{4}({\mathbf {k}})$ are odd functions with respect to ${\mathbf {k}}$. Similarly, the coefficient functions $d_{0}({\mathbf {k}})$ and $d_{5}({\mathbf {k}})$ are even functions with respect to ${\mathbf {k}}$. Thus, the allowed lattice distortion terms can be either $k^{\text {odd}} ~ \Gamma _{1,2,3,4}$ or $k^{\text {even}} ~ \Gamma _{0, 5}$ types.

Because we assume TRS and IS to remain under lattice distortions, the Hamiltonians for distorted and undistorted DSM have the same form of $H = \sum _{i} a_{i}({\mathbf {k}}) \Gamma _i$. The only difference between the two Hamiltonians is the transformation properties of the coefficient function $a_{i} ({\mathbf {k}})$. In the absence of lattice distortions, $a_{i} ({\mathbf {k}})$ needs to satisfy all transformation properties under all symmetry operations of $D_{4h}$ point group in Table 2. However, in the presence of lattice distortion, $a_{i} ({\mathbf {k}})$ only needs to satisfy the transformation properties under the remaining symmetry operations, so $a_{i}({\mathbf {k}})$ is less constrained.

Lattice Hamiltonian with lattice distortions

To discuss the effect of lattice distortions explicitly, we introduce the possible symmetry-lowering lattice distortions in the lattice model in Eqs. (7-11). For weak lattice distortions, the lattice distortions are approximately proportional to $\sin k_i$ and $\cos k_i$ as only nearest neighbor hoppings are relevant. Because we are interested in the Dirac physics near the Dirac points $(0, 0, \pm k_0)$, we assume that $k_x/k_z \ll 0$ and $k_y/k_z \ll 0$, which implies that $\sin k_x$ and $\sin k_y$ are smaller than $\sin k_z$ and $\cos k_i$. Hence, $\sin k_z$ and $\cos k_i$ are dominant momentum dependent terms in the leading order, and the allowed lattice distortions are either $\sin k_z ~ \Gamma _{1,2,3,4}$ or $\cos k_i ~ \Gamma _{0, 5}$ types. Because $\cos k_i ~ \Gamma _{0, 5}$ types are included in the trivial $A_{1g}$ class of $D_{4h}$ point group, they do no break any symmetry. On the other hand, $\sin k_z ~ \Gamma _{1,2,3,4}$ types are included in $B_{1g}$, $B_{2g}$, and $E_g$, and they break the crystal symmetry properly, which are summarized in Table 3. Therefore, in the leading order, there are four types of symmetry-lowering lattice distortion, which are given by

$$\begin{aligned} H_{\text {pert}} = \sin k_z (n_1 \Gamma _3 + n_2 \Gamma _4 + n_3 \Gamma _2 + n_4 \Gamma _1), \end{aligned}$$

(16)

where $n_i$ is the strength of each lattice distortion. For convenience, each lattice distortion is denoted as $n_i$ type lattice distortions in this work. From now on, we will consider these four types of symmetry-lowering lattice distortions, and the possible higher-order terms are discussed in Sec. S2 in Supplementary Information.

Table 3 Four types of symmetry-lowering lattice distortions are classified according to the irreducible representation of $D_{4h}$ point group.

Therefore, the coefficient functions in Eq. (1) are given by

$$\begin{aligned} a_1 ({\mathbf {k}})&= v \sin k_x + n_4 \sin k_z , \nonumber \\ a_2 ({\mathbf {k}})&= v \sin k_y + n_3 \sin k_z , \nonumber \\ a_3 ({\mathbf {k}})&= ( \beta + \gamma ) (\cos k_y - \cos k_x) \sin k_z + n_1 \sin k_z , \nonumber \\ a_4 ({\mathbf {k}})&= - ( \beta - \gamma ) (\sin k_x \sin k_y \sin k_z ) + n_2 \sin k_z, \nonumber \\ a_5 ({\mathbf {k}})&= M' - t_{xy}(\cos k_x + \cos k_y ) - t_z \cos k_z. \end{aligned}$$

(17)

Under lattice distortion, the fourfold rotation symmetry is broken. Thus, the Dirac point is gapped, which can be seen from the energy eigenvalues on the $k_z$ axis, $E = \pm \sqrt{ ( n_1^2 + n_2^2 + n_3^2 + n_4^2) \sin ^2 k_z +a_{5} (0,0,k_z)^2 }$. Thus, the Dirac point is gapped unless $n_1^2+n_2^2 + n_3^2 + n_4^2 = 0$. As a result of the gap-opening, the DSM becomes a 3D topological insulator because of the band inversion at the $\Gamma$ point^36,39. Counting all the parity eigenvalues for the time-reversal-invariant momenta (TRIM) points of the bulk Brillouin zone (BZ)^1,54 gives a nontrivial ${\mathbb {Z}}_2$ invariant.

The effect of lattice distortions

The four types of symmetry-lowering lattice distortions in Eq. (16) are classified according to the irreducible representation of $D_{4h}$ group. The symmetry-lowering lattice distortions break $D_{4h}$ point group symmetry into its subgroup symmetry, which is summarized in Table 3. The $n_1$ and $n_2$ type lattice distortions are included in the one-dimensional class $B_{1g}$ and $B_{2g}$, and break $D_{4h}$ point group symmetry into $D_{2h}$ and $D_{2h}'$, respectively. The $n_3$ and $n_4$ type lattice distortions are included in the two-dimensional class $E_{u}$ and break $D_{4h}$ point group symmetry into $C_{2h}$. Note that $n_2$ type lattice distortion is related to the $n_1$ type lattice distortion via $\pi /4$ rotation, while $n_4$ type lattice distortion is related to the $n_3$ type lattice distortion via $\pi /2$ rotation.

We investigate the explicit effects of the lattice distortions on the crystal systems and the Fermi surfaces using the lattice model in Eq. (17). Figure 1 shows the crystal structures, the 3D band structures, and Fermi surfaces under various lattice distortions. Under $n_1$ type lattice distortion, the crystal system and Fermi surface are elongated along x or y direction, $C_{4z}$ symmetry is broken, the Dirac point is gapped, and the crystal system becomes orthorhombic (Fig. 1b, g). Similarly, under the $n_2$ type lattice distortion, the crystal system and Fermi surface are elongated along diagonal lines either $x=y$ or $x=-y$, $C_{4z}$ symmetry is broken, the Dirac point is gapped, and the crystal system becomes orthorhombic (Fig. 1c, h). We denote the symmetry point group of this right rhombic prism as $D_{2h}'$. Under $n_3$ type lattice distortion, the crystal structure undergoes structural phase transition from tetragonal to monoclinic (Fig. 1d). Two Dirac points in the band structure are shifted oppositely along $k_y$ direction and the centers of each Fermi surfaces are also oppositely shifted along the same $k_y$ direction (Fig. 1h). Similar effects occur under $n_4$ type lattice distortion (Fig. 1e, j) because $n_4$ type lattice distortion are related with the $n_3$ type lattice distortion via $\pi /2$ rotation. The point groups of these distorted systems under $n_3$ and $n_4$ type lattice distortions are denoted as $C_{2h(x)}$ and $C_{2h(y)}$, respectively. Therefore, the four types of symmetry-lowering lattice distortions explain the lattice distortions of $\hbox {Cd}_3\hbox {As}_2$ and $\hbox {Au}_2$Pb under pressure.

Low-energy effective Dirac Hamiltonian under lattice distortions

Near the Dirac points $(0,0, \pm k_0)$, the coefficient functions of the low-energy effective Hamiltonian can be approximated as

$$\begin{aligned} a_1 ({\mathbf {k}})&= v k_x + n_4 \sin k_0 , \\ a_2 ({\mathbf {k}})&= v k_y + n_3 \sin k_0 , \\ a_3 ({\mathbf {k}})&= ( \beta + \gamma ) \left( \frac{k_x^2 - k_y^2}{2} \right) \sin k_0+ n_1 \sin k_0 , \\ a_4 ({\mathbf {k}})&= - ( \beta - \gamma ) k_x k_y \sin k_0 + n_2 \sin k_0, \\ a_5 ({\mathbf {k}})&= \pm v_z (k_z \mp k_0) \sigma _z. \end{aligned}$$

With this low-energy effective Hamiltonian, we show that the lattice distortion acts as a Dirac mass term and increases DOS at Fermi surface. We assume that the Fermi level is slightly above the Dirac points in undistorted lattice, or near the bottom of the conduction band minima after gap-opening at the Dirac points.

For $n_1$ and $n_2$ type lattice distortions, the low-energy effective Hamiltonian is given by

$$\begin{aligned} H_{\text {Dirac}}^{(\pm )} =v k_x \Gamma _1 + v k_y \Gamma _2 \pm v_z (k_z \mp k_0) \Gamma _5 \pm n_1 \sin k_0 \Gamma _3 \pm n_2 \sin k_0 \Gamma _4. \end{aligned}$$

(18)

So, $n_1$ and $n_2$ type lattice distortion terms act as Dirac mass terms. The energy eigenvalue is given by

$$\begin{aligned} E = \pm \sqrt{ v^2 (k_x^2 + k_y^2) + v_z^2 (k_z \mp k_0)^{2} + \left| n\right| ^2 \sin ^2 k_0 }, \end{aligned}$$

(19)

where $\left| n\right| = \sqrt{n_1^2 + n_2^2}$. By the assumption of the total electron number conservation under a weak lattice distortion, the lattice distortion dependent DOS at the Fermi surface is given by

$$\begin{aligned} {\text {DOS}}(\left| n\right| ) = \frac{1}{\pi v^2 v_z}\mu _0 \sqrt{\mu _0^2 + \left| n\right| ^2 \sin ^2 k_0}, \end{aligned}$$

(20)

which indicates that DOS at the Fermi level is enhanced under the lattice distortion. Here, $\mu _0$ indicates the chemical potential of the undistorted lattice. See the detailed derivations in Sec. S2.4 in Supplementary Information.

Next, we consider the $n_3$ type lattice distortion. The $n_3$ type lattice distortion shifts the gap minima along the $k_y$ direction from $(0, 0, \pm k_0)$ to $(0, \pm k_y^{(0)}, \pm k_0)$ with $k_y^{(0)} = - n_3 \sin k_0/v$. Then, the low-energy effective Hamiltonian near the gap minima points $(0, \pm k_y^{(0)}, \pm k_0)$ is given by

$$\begin{aligned} H_{\text {Dirac}}^{(\pm )} = v k_x \Gamma _1 + v (k_y \mp k_y^{(0)}) \Gamma _2 \pm v_z (k_z \mp k_0) \Gamma _5 \pm m \Gamma _3, \end{aligned}$$

where $m = - (\beta + \gamma ) \frac{n_3^2 \sin ^3 k_0 }{2 v^2 }$ is the Dirac mass term. The energy eigenvalue is given by

$$\begin{aligned} E = \pm \sqrt{ v^2 k_x^2 + v^2 (k_y \mp k_y^{(0)} )^2 + v_z^2 (k_z \mp k_0)^2 + m^2 }. \end{aligned}$$

(21)

Similar to $n_1$ and $n_2$ type lattice distortions, DOS at the Fermi surface are given by

$$\begin{aligned} {\text {DOS}}(n_3) = \frac{1}{\pi v^2 v_z}\mu _0 \sqrt{ \mu _0^2 + m^2 }, \end{aligned}$$

(22)

which means that the DOS at the Fermi level is enhanced under $n_3$ type lattice distortion. Similarly, for $n_4$ type lattice distortion, the low-energy effective Hamiltonian and DOS are easily calculated because $n_3$ and $n_4$ type lattice distortions are related via $\pi /2$ rotation.

Multiple symmetry-lowering lattice distortions

So far, we have considered only one type of lattice distortions. However, more than two types of lattice distortions can be turned on simultaneously. In this case, the final subgroup symmetry determines the crystal system and its physical properties. When both $n_1$ and $n_3$ types lattice distortions are turned on, the remaining subgroup has P, $C_{2x}$, $M_{yz}$ symmetries. This subgroup is the same point group of the distorted Dirac semimetal under single $n_3$ type lattice distortion. In other words, under $n_3$ type lattice distortion, the addition of $n_1$ type lattice distortion is also allowed. A similar argument can be applied to $n_2$ and $n_4$ types lattice distortions. When both $n_1$ and $n_2$ type lattice distortions are turned on, the remaining symmetries are P, $C_{2z}$, $M_{xy}$ symmetries. We denote this point subgroup as $C_{2h(z)}$, and we will not consider this case seriously because there is no real material that corresponds to this case. Similarly, the other combinations such as $(n_2, n_3)$, $(n_1, n_4)$, $(n_3, n_4)$, $(n_1, n_2, n_3)$, $(n_1, n_2, n_4)$ break all crystal symmetries except the inversion, and hence these cases are not interested in this work.

Superconductivity

BdG Hamiltonian

To discuss the effects of lattice distortions on the superconductivity in doped DSM, we construct the Bogoliubov-de Gennes (BdG) Hamiltonian within mean-field approximation while keeping TRS and the crystal symmetry^55,56. The BdG Hamiltonian is given by

$$\begin{aligned} H_{\text {BdG}}&= \int d {\mathbf {k}} \Psi _{{\mathbf {k}}}^{\dagger } {\mathscr {H}} ({\mathbf {k}}) \Psi _{{\mathbf {k}}}, \end{aligned}$$

(23)

$$\begin{aligned} {\mathscr {H}}({\mathbf {k}} )&= [H ({\mathbf {k}}) - \mu ] \tau _z + \Delta ({\mathbf {k}}) \tau _x, \end{aligned}$$

(24)

where $\tau _i$ is the Pauli matrices in the Nambu space. $\Delta ({\mathbf {k}})$ and $\mu$ are a pairing potential and a chemical potential, respectively. $H ({\mathbf {k}})$ is the normal state Hamiltonian in Eq. (1). The basis is taken as

$$\begin{aligned} \Psi _{{\mathbf {k}}}^{\dagger } = ( c^{\dagger }_{1 {\mathbf {k}} \uparrow }, c^{\dagger }_{2 {\mathbf {k}} \uparrow }, c^{\dagger }_{1 {\mathbf {k}} \downarrow }, c^{\dagger }_{2 {\mathbf {k}} \downarrow }, c_{1 {\mathbf {-k}} \downarrow }, c_{2 {\mathbf {-k}} \downarrow }, -c_{1 {\mathbf {-k}} \uparrow }, -c_{2 {\mathbf {-k}} \uparrow } ). \end{aligned}$$

(25)

While the pairing mechanism of doped DSM is not known yet, we assume the following onsite density-density interaction as a superconducting pairing interaction^36,37,57,58:

$$\begin{aligned} H_{\text {int}}(x) = - U [ n_1^1(x) + n_2^2(x) ] - 2V n_1(x) n_2(x), \end{aligned}$$

(26)

where $n_{i}(x)$ is the electron density operators for ith orbital ($i = 1, 2$). U and V are intra-orbital and inter-orbital interaction strengths, respectively, and we assume that at least one of them is attractive and responsible for superconductivity. Because the pairing interaction depends on the orbital and is local in ${\mathbf {x}}$, the mean-field pairing potential is orbital dependent but momentum independent: $\Delta ({\mathbf {k}}) = \Delta$.

Symmetry of BdG Hamiltonian

The BdG Hamiltonian in Eq. (23) has time-reversal symmetry T, particle-hole symmetry C, and chiral symmetry $\Gamma$:

$$\begin{aligned} T {\mathscr {H}} ({\mathbf {k}} ) T^{-1}&= {\mathscr {H}} ( - {\mathbf {k}} ), \end{aligned}$$

(27)

$$\begin{aligned} C{\mathscr {H}} ({\mathbf {k}} ) C^{-1}&= - {\mathscr {H}} ( - {\mathbf {k}} ),\end{aligned}$$

(28)

$$\begin{aligned} \Gamma {\mathscr {H}} ({\mathbf {k}}) \Gamma ^{-1}&= {\mathscr {H}} ({\mathbf {k}}), \end{aligned}$$

(29)

where $T=is_y \sigma _0 \tau _0 {\hat{K}}$ and $C= i s_y \sigma _0 \tau _y {\hat{K}}$ are time-reversal and particle-hole symmetry operators, respectively, and $\Gamma = T C= s_0 \sigma _0 \tau _y$ is the chiral operator. ${\hat{K}}$ is the complex conjugation operator. Therefore, the BdG Hamiltonian belongs to in DIII class according to the classification table of topological insulator and superconductor⁵⁹.

If the pairing potential satisfies $P \Delta ({\mathbf {k}}) P^{-1} = \eta _P \Delta (- {\mathbf {k}})$, the BdG Hamiltonian has the inversion symmetry:

$$\begin{aligned} {{\tilde{P}}} {\mathscr {H}} ({\mathbf {k}} ) {{\tilde{P}}}^{-1}&= {\mathscr {H}} ( - {\mathbf {k}} ), ~~~ {\text {with}} ~~{{\tilde{P}}} = {\text {diag}}(P, \eta _P P), \end{aligned}$$

(30)

where P and ${{\tilde{P}}}$ are the inversion operators for the DSM and BdG Hamiltonians, respectively, and $\eta _P$ is the inversion parity. If $\eta _P = 1$ ($\eta _P = -1$), the superconducting phase is an inversion-even-parity (inversion-odd-parity) superconductor. For a single-orbital superconductor, ${{\tilde{P}}}$ is the identity operator, and an inversion-odd-parity (inversion-even-parity) pairing is equivalent to the spin-triplet (spin-singlet) pairing. However, because of the spin-orbit coupling and multi-orbital band structure, the pairings are more complex in our case.

From now on, we consider momentum independent pairing potentials, $\Delta ({\mathbf {k}}) = \Delta$, because we assume onsite pairing interaction as discussed in Eq. (26). In the absence of lattice distortions, the BdG Hamiltonian has $D_{4h}$ point group symmetry³⁷. If a pairing potential satisfies the transformation property of $U \Delta s_y U^{T} s_y = \eta _U \Delta$ under a symmetry operation of $D_{4h}$ point symmetry group, the BdG Hamiltonian satisfies the corresponding symmetry:

$$\begin{aligned} {{\tilde{U}}} {\mathscr {H}} ({\mathbf {k}}) {{\tilde{U}}}^{-1} = {\mathscr {H}}(S {\mathbf {k}}), \end{aligned}$$

(31)

where U is the symmetry operator in spin and orbital spaces, $\eta _U$ is a phase factor, and ${{\tilde{U}}} = {\text {diag}} ( U, \eta _U s_y U^* s_y)$ is the extended symmetry operator in the Nambu space.

For the generators of $D_{4h}$ point group, if the pairing potential satisfies $C_{4z} \Delta s_y C_{4z}^{T} s_y = \eta _{C_{4z}} \Delta$ with $\eta _{C_{4z}} = e^{ i \pi r /2 }$ ($r=0, \ldots , 3$) and $C_{2x} \Delta s_y C_{2x}^{T} s_y = \eta _{C_{2x}} \Delta$ with $\eta _{C_{2x}} = \pm 1$, then the BdG Hamiltonian satisfies the corresponding rotation symmetry:

$$\begin{aligned} {{\tilde{C}}}_{4z} {\mathscr {H}} ({\mathbf {k}}) {{\tilde{C}}}_{4z}^{-1}&= {\mathscr {H}} ( R_{4z} {\mathbf {k}}), \end{aligned}$$

(32)

$$\begin{aligned} {{\tilde{C}}}_{2x} {\mathscr {H}} ({\mathbf {k}}) {{\tilde{C}}}_{2x}^{-1}&= {\mathscr {H}} ( R_{2x} {\mathbf {k}}), \end{aligned}$$

(33)

where the extended symmetry operators are given by ${\tilde{C}}_{4z} = {\text {diag}} ( C_{4z}, \eta _{C_{4z}} s_y C_{4z}^* s_y )$ and ${\tilde{C}}_{2x} = {\text {diag}} ( C_{2x}, \eta _{C_{2x}} s_y C_{2x}^* s_y )$. If the pairing potential satisfies $M \Delta s_y M^{T} s_y = \eta _M \Delta$ under a mirror operator M, the BdG Hamiltonian satisfies the corresponding mirror symmetry:

$$\begin{aligned} {\tilde{M}} {\mathscr {H}} ({\mathbf {k}}_{\parallel }, {\mathbf {k}}_{\perp }) {\tilde{M}}^{-1} = {\mathscr {H}} ({\mathbf {k}}_{\parallel }, - {\mathbf {k}}_{\perp }) , \end{aligned}$$

(34)

where ${{\tilde{M}}} ={\text {diag}}~(M, \eta _M s_y M^* s_y)$ is a mirror operator for BdG Hamiltonian and ${\mathbf {k}}_{\parallel }$ (${\mathbf {k}}_{\perp }$) is the momentum vector parallel (perpendicular) to the mirror plane. The $\eta _M$ is the mirror parity of the pairing potential under the mirror operation M.

In Table 4, the transformation properties of all possible pairing potentials under the rotation and mirror operators are summarized. The details of each pairing potential will be discussed below.

Table 4 The pairing potentials are classified according to the irreducible representation of $D_{4h}$ point group.

Table 5 Pairing potentials classified according to the $D_{4h}$ point group are reclassified according to the irreducible representation of unbroken subgroup under the lattice distortions.

Pairing potentials

We now investigate the possible superconducting pairing potentials in the presence of lattice distortions. Since we are considering multi-orbital superconductivity in the basis of two spins and two orbitals, pairing potentials can be represented as a product of spin Pauli matrices and orbital Pauli matrices, which leads to sixteen matrices. Among them, only six matrices are allowed because of the fermion statistics ($\Delta s_y = s_y \Delta ^T$). We denote them as $\Delta _{1}$, $\Delta _{1}'$, $\Delta _{2}$, $\Delta _{3}$, $\Delta _{41}$, and $\Delta _{42}$, whose forms and properties are listed in Table 4. Due to Pauli’s exclusion principle, the fermion bilinear form of each pairing potential shows antisymmetric property under the particle exchange. Because the pairing potential is momentum independent, the spatial part is symmetric, while the spin-orbital part is antisymmetric under the particle exchange. Thus, if the spin part is singlet, the orbital part is triplet, and vice versa. Therefore, $\Delta _1$’s and $\Delta _{41}$ are the spin-singlet orbital-triplet pairings and $\Delta _2$, $\Delta _3$, and $\Delta _{42}$ are the spin-triplet orbital-singlet pairings as shown in the bilinear form in Table 4.

Six pairing potentials can be classified according to the irreducible representations of the unbroken point group, and the superconducting critical temperatures for the pairing potentials in the different classes are independent^{36,37,56,57,58}. In the absence of lattice distortions, the pairing potentials are classified according to the $D_{4h}$ group: $\Delta _{1}$’s, $\Delta _2$, $\Delta _3$ and $\Delta _4$’s belong to $A_{1g}$, $B_{1u}$, $B_{2u}$ and $E_u$ irreducible representations, respectively, which are summarized in Table 4.

The pairing potential belonging to a specific irreducible representation of the $D_{4h}$ group can be decomposed into a combination of different irreducible representations depending on the symmetry of the distorted lattice. Some pairing potentials in the $D_{4h}$ group’s individual representations can be included in the same representation and vice versa. As an example, in the $D_{2h}$ case, $(\Delta _{41}, \Delta _{42})$ belong to in the two-dimensional representation $E_u$ are separated into one-dimensional representations $B_{2u}$ and $B_{3u}$, respectively. Similarly, for $D_{2h}'$ case, the linear combination of $\Delta _{41}$ and $\Delta _{42}$ potential belongs to in one-dimensional representations $B_{2u}$ and $B_{3u}$. Because $D_{2h}'$ case is the $\pi /4$-rotated version of $D_{2h}$ case, $\Delta _{41}+\Delta _{42}$ ($\Delta _{42}-\Delta _{41}$) is included in $B_{3u}$ ($B_{2u}$) class when $\Delta _{41}=\Delta _{42}$ ($\Delta _{41}=-\Delta _{42}$). The reclassification of pairing potentials under various lattice distortions is summarized in Table 5.

Figure 2

Superconducting nodal structures for pairing potentials under lattice distortions. Nodal structures for (a) $D_{4h}$, (b) $D_{2h}$, (c) $D_{2h}'$, (d) $C_{2h(z)}$, and (e) $C_{2h(x)}$ cases. The orange point, line, and plane indicate nodal point and nodal line, and mirror plane ($M_{xz}$, $M_{yz}$, $M_{110}$, and $M_{1 {{\bar{1}}} 0}$), respectively. In (a–e), the $\Delta _{1}$ phases are fully gapped and the $\Delta _{1}'$ phases have two nodal rings. In (a–c, e), nodal points are located in the corresponding mirror planes. In (c), $\Delta _{42}+\Delta _{41}$ and $\Delta _{42}-\Delta _{41}$ phases are considered instead of $\Delta _{41}$ and $\Delta _{41}$ phases. In (d), the system has no mirror symmetries and hence no nodal points. These nodal structures are summarized in Table 6.

Table 6 Nodal structures of superconducting phases under lattice distortions.

Superconducting nodal structure

In this subsection, we classify the superconducting nodal structures under lattice distortions and study the symmetry and topology that guarantee the classified nodal structures.

Figure 2 shows the typical nodal structures of superconducting phases of the doped DSM under lattice distortions. There are three types of nodal structures: Full gap, point nodal, and line nodal structures, which are summarized in Table 6. For $\Delta _{1}$ and $\Delta _{1}'$ superconducting phases, $\Delta _{1}$ phase is fully gapped and $\Delta _{1}'$ phase has two nodal rings regardless of lattice distortions (Fig. 2a–e). For $\Delta _{2}$ and $\Delta _{3}$ phases, nodal points exist at the intersections between the $k_z$ axis and the Fermi surfaces in the absence of lattice distortions (Fig. 2a). These points are known to be protected by $C_{4z}$ symmetry^36,37. Even under lattice distortions, if there is an unbroken mirror symmetry, the topologically protected nodal points can exist and they are protected by the corresponding mirror symmetry (Fig. 2b–e). For $\Delta _{41}$ and $\Delta _{42}$ phases, there are accidental nodal points at the intersections between the $k_z$ axis and the Fermi surfaces in the absence of lattice distortions (Fig. 2a). However, in the presence of lattice distortions, if there is an unbroken mirror symmetry, there can exist the topologically protected nodal points in the corresponding mirror plane (Fig. 2b–e). Note that all nodal points under lattice distortions in Fig. 2(b–e) are protected by the topological mirror winding numbers, as discussed later.

We now analytically investigate the condition of nodal points in each superconducting phase. Usually, nodal points can exist where the quasi-particle energy spectrum vanishes ${\mathscr {E}}({\mathbf {k}})=0$, which gives a set of equations for the momentum variables ($k_x, k_y, k_z$). If the number of variables $N_{V}$ is greater than or equal to the number of independent equations $N_{E}$, then nodal structures can exist. That is, $N_{E} \le N_{V} = 3$ is the necessary condition for the existence of the nodes. Moreover, if there is mirror symmetry, the necessary condition changes because the number of independent variables is reduced in the corresponding mirror plane. That is, the necessary condition becomes $N_E \le N_{V} = 2$. If there is additional mirror symmetry, the necessary condition can be further reduced to $N_E \le N_{V} = 1$ on the intersection of two mirror planes.

First, we consider $\Delta _{1}$ and $\Delta _{1}'$ superconducting phases. The full gap structure of $\Delta _{1}$ phase is directly seen from the energy eigenvalues of

$$\begin{aligned} {\mathscr {E}}({{\mathbf {k}}}) = \pm \sqrt{(\left| a\right| \pm \left| \mu \right| )^2 + \left\langle \Delta _{1}\right\rangle ^2}, \end{aligned}$$

(35)

where $\left| a\right| = \sqrt{ \sum _{i=1}^5 a_i({\mathbf {k}})^2}$. Unless $\left\langle \Delta _{1}\right\rangle = 0$, $\Delta _{1}$ phase is fully gapped. For $\Delta _{1}'$ phase, the energy eigenvalues are given by

$$\begin{aligned} {\mathscr {E}}({{\mathbf {k}}}) = \pm \sqrt{ \left| a\right| ^2 + \mu ^2 + \left\langle \Delta _{1}'\right\rangle ^2 \pm 2 \sqrt{ \mu ^2 \left| a\right| ^2 + \left\langle \Delta _{1}'\right\rangle ^2 \left( \left| a\right| ^2 - a_5({\mathbf {k}})^2 \right) } }. \end{aligned}$$

(36)

From ${\mathscr {E}}({\mathbf {k}} )=0$, one can obtain the following equations:

$$\begin{aligned} \left| a\right| ^2 = \mu ^2 + \left\langle \Delta _{1}'\right\rangle ^2, ~~~~ a_5({\mathbf {k}}) = 0. \end{aligned}$$

(37)

Because the number of variable ($N_V = 3$) is larger than the number of equation ($N_E=2$), a one-dimensional solution can exist, which leads to the nodal lines. Because this argument works regardless of the lattice distortions, the nodal rings can exist for all cases in Fig. 2. On the other hand, under some lattice distortions, a mixture of $\Delta _{1}$ and $\Delta _{1}'$ phases is allowed when $\Delta _{1}$ and $\Delta _{1}'$ are in the same representation as shown in Table 5. In such case, the gap structures have full gap (nodal lines) when $\left| \left\langle \Delta _{1}\right\rangle \right| > \left| \left\langle \Delta _{1}'\right\rangle \right|$ ($\left| \left\langle \Delta _{1}\right\rangle \right| < \left| \left\langle \Delta _{1}'\right\rangle \right|$) (Fig. S2). See the detailed calculation in Sec. S3 in Supplementary Information.

Next, consider the $\Delta _{2}$ and $\Delta _{3}$ superconducting phases. In the absence of lattice distortions, the nodal points in $\Delta _{2}$ and $\Delta _{3}$ phases are protected by $C_{4z}$ symmetry^36,37. On the other hand, under lattice distortions, a mirror symmetry can protect the nodal points that appear in Fig. 2b, c, e. For $\Delta _3$ phase, the energy eigenvalues are given by

$$\begin{aligned} {\mathscr {E}}({\mathbf {k}}) = \pm \sqrt{ \left| a\right| ^2 + \mu ^2 +\left\langle \Delta _{3}\right\rangle ^2 \pm 2 \sqrt{ \mu ^2 \left| a\right| ^2 +(a_3({\mathbf {k}})^2 + a_5({\mathbf {k}})^2) \left\langle \Delta _3\right\rangle ^2} }. \end{aligned}$$

(38)

From ${\mathscr {E}}({\mathbf {k}}) =0$, we get the following equations:

$$\begin{aligned} \left| a\right| ^2 = \mu ^2 + \left\langle \Delta _{3}\right\rangle ^2,~~~ a_1({\mathbf {k}}) = a_2({\mathbf {k}}) = a_4({\mathbf {k}}) = 0. \end{aligned}$$

(39)

Because $N_E = 4$ is larger than $N_V = 3$, there seems to be no allowed nodal point. However, mirror symmetries can allow nodal points. For example, consider $D_{2h}$ point group with $M_{yz}$ and $M_{xz}$ mirror symmetries. Under the $M_{yz}$ mirror operation, $a_1({\mathbf {k}})$ and $a_4({\mathbf {k}})$ are odd according to Table 2, which gives $a_1({\mathbf {k}}) = a_4 ({\mathbf {k}})= 0$ at the mirror plane $(0, k_y, k_z)$. Similarly, $M_{xz}$ mirror symmetry gives $a_2({\mathbf {k}}) = a_4 ({\mathbf {k}})= 0$ at the mirror plane $(k_x, 0, k_z)$. Thus, along the $k_z$ axis, $a_1({\mathbf {k}}) = a_2({\mathbf {k}}) = a_4({\mathbf {k}}) = 0$ and Eq. (39) is reduced to

$$\begin{aligned} a_3^2 (k_z) + a_5^2(k_z) = \mu ^2 + \left\langle \Delta _{3}\right\rangle ^2. \end{aligned}$$

(40)

Because $N_E = 1$ is equal to $N_V = 1$, nodal points can exist as shown in Fig. 2b. However, when $M_{yz}$ and $M_{xz}$ mirror symmetries are broken, the nodal points for the $\Delta _3$ phase are not protected as shown in Fig. 2c, d.

Similarly, the nodal points in $\Delta _2$ phase can be understood using $M_{110}$ and $M_{1 {{\bar{1}}} 0}$ mirror symmetries. These mirror symmetries allow nodal points on the $k_z$ axis in Fig. 2c. On the other hand, when $M_{110}$ and $M_{1 {{\bar{1}}} 0}$ mirror symmetries are broken, the nodal points disappear as shown in Fig. 2b, d. For the $C_{2h(z)}$ case, a mixture of $\Delta _{2}$ and $\Delta _{3}$ phases is possible because $\Delta _{2}$ and $\Delta _{3}$ are included in the same $A_u$ representation. However, there is no allowed nodal point as shown in Fig. 2d because there is no helpful mirror symmetry. See the details in Sec. S3 in Supplementary Information.

Finally, consider $\Delta _{41}$ and $\Delta _{42}$ phases. Without lattice distortions, there are accidental nodal points on the $k_z$ axis (Fig. 2a). The existence of such nodal point is easily seen using four mirror symmetries $M_{xz}, M_{yz}, M_{110}$, and $M_{1 {{\bar{1}}} 0}$. These mirror symmetries force $a_i({\mathbf {k}}) = 0$ for $i=1,\cdots ,4$ on the $k_z$ axis according to Table 2. Then, the equations for nodal points are given by

$$\begin{aligned} \left| a_5(k_z)\right| ^2 = \mu ^2 + \left\langle \Delta _{41}\right\rangle ^2 + \left\langle \Delta _{42}\right\rangle ^2. \end{aligned}$$

(41)

Because $N_E=N_V=1$, the nodal points exist. Because the $\Delta _{41}$ and $\Delta _{42}$ pairing potentials included in $E_u$ representation of $D_{4h}$ point symmetry group, they break the $D_{4h}$ symmetry spontaneously to $D_{2h}$. Hence, some of non-zero $a_i({\mathbf {k}})$ ($i=1\cdots 4$) are spontaneously generated and the corresponding conditions are introduced, which makes the nodal points vanish. Thus, these nodal points are accidental. However, under lattice distortions, the nodal points can be protected by the unbroken mirror symmetry. For example, when the point group is $D_{2h}$ under the $n_1$ type lattice distortion, $\Delta _{41}$ and $\Delta _{42}$ are included in the different representations and thus we can consider each phase separately. For $\Delta _{41}$ phase, $a_1({\mathbf {k}})=a_4({\mathbf {k}})=0$ on the mirror plane $(0, k_y, k_z)$ due to $M_{yz}$ symmetry. Then, the equations for nodes are given by

$$\begin{aligned} a_2^2 (0, k_y, k_z) + a_5^2 (0, k_y, k_z) = \mu ^2+ \left\langle \Delta _{41}\right\rangle ^2, ~~~ a_3 (0, k_y, k_z) = 0. \end{aligned}$$

(42)

Because $N_E=N_V=2$, there can exist nodal points (Fig. 2b). For $\Delta _{42}$ phase, nodal points also can exist due to $M_{xz}$ mirror symmetry (Fig. 2b). When the point group is $D_{2h'}$ under the $n_2$ type lattice distortion, nodal points can exist due to $M_{110}$ or $M_{1 {{\bar{1}}} 0}$ mirror symmetries (Fig. 2c). For $C_{2h(z)}$, a mixture of $\Delta _{41}$ and $\Delta _{42}$ phases is possible. However, there is no allowed nodal point due to the lack of mirror symmetry (Fig. 2d). When the point group is $C_{2h(x)}$ under the $n_3$ type lattice distortion, nodal points can exist due to $M_{yz}$ mirror symmetry (Fig. 2e). See the detailed calculations in Sec. S3 in Supplementary Information.

Figure 3

Topologically protected nodal structures and chiral winding numbers. The orange ring, point, plane, and vertical line indicate nodal ring, nodal point, mirror plane, and $k_z$ axis, respectively. Each winding number is defined along each blue loop. (a) The chiral winding numbers ($W=\pm 2$) protect nodal rings. (b, c) The mirror chiral winding numbers ($W_M=\pm 2$) protect nodal points on the mirror planes. (d) Evolution of nodal points in $\Delta _{42}$ phases and the corresponding mirror chiral winding number under the $n_1$ type lattice distortion. For clarity, the blue winding loops are omitted. For $n_1=0$, nodal points with $W_M=0$ are located on $k_z$ axis. These are fine-tuned accidental nodal points because $D_{4h}$ is spontaneously broken into $D_{2h}$ due to $\Delta _{42}$ pairing [see the main text below Eq.  (41)]. As $n_1$ increases, the nodal points split into two nodal points with $W_M=\pm 2$. The bottom plot shows the evolution of the energy dispersion along $k_x$ axis. As $n_1$ increases, the blue (orange) band moves downward (upward), which results in two Dirac points.

Stability of nodal structures

There are two types of nodes in Table 6, which are symmetry-protected node and topologically-protected node. In this subsection, we investigate the stability of them.

Chiral winding number

Because of the chiral symmetry of the BdG Hamiltonian, the nodal lines can be protected by a chiral winding number^4,59,60. The chiral winding number is defined along a path ${\mathscr {C}}$ enclosing a singular point in the Brillouin zone as shown in Fig. 3a:

$$\begin{aligned} W = \frac{1}{4 \pi i} \oint _{{\mathscr {C}}} {\text {Tr}} \left[ \Gamma {\mathscr {H}}^{-1} ({\mathbf {k}}) d {\mathscr {H}}({\mathbf {k}}) \right] , \end{aligned}$$

(43)

where $\Gamma$ is the chiral operator. As shown in Sec. S4 in Supplementary Information, the transformation property of the winding number under PT symmetry is given by

$$\begin{aligned} W = - \eta _{\Gamma , {{\tilde{P}}} T} W, \end{aligned}$$

(44)

where the parity $\eta _{A,B} = \pm 1$ is determined by the relation $AB = \eta _{A,B} BA$. For the inversion-even-parity (inversion-odd-parity) pairing potential, $\eta _{\Gamma , {{\tilde{P}}} T}$ is $-1$ ($+1$). Thus, the chiral winding number is zero for the inversion-odd-parity superconductor and only the inversion-even-parity superconducting phases having $\Delta _{1}$ and $\Delta _{1}'$ pairing potentials can have a nontrivial chiral winding number.

$\Delta _{1}$ and $\Delta _{1}'$ phases

Because $\Delta _{1}$ phase is fully gapped, the chiral winding number is zero. On the other hand, two nodal rings in $\Delta _{1}'$ phase are topologically protected by the chiral winding numbers. The calculated chiral winding numbers around the nodal rings are $W=\pm 2$ (Fig. 3). These chiral winding numbers do not change even under the lattice distortions because chiral winding number depends only on T, C, P, and $\Gamma$ symmetries. Thus, the topologically-protected nodal rings in $\Delta _{1}'$ phase exist regardless of the lattice distortion (Fig. 2).

Mirror chiral winding number

If there is mirror symmetry, the BdG Hamiltonian commutes with the mirror symmetry operator in the mirror plane:

$$\begin{aligned}{}[ {\tilde{M}}, {\mathscr {H}} ({\mathbf {k}}_{M})] = 0. \end{aligned}$$

(45)

where ${\tilde{M}}$ is a mirror operator and ${\mathbf {k}}_{M}$ is the momentum vector located in the mirror plane. Then, the BdG Hamiltonian can be block-diagonalized according to the mirror eigenvalues $\lambda = \pm i$. Besides, if the mirror operator commutes with the chiral operator,

$$\begin{aligned}{}[ \Gamma , {\tilde{M}}] = 0 , \end{aligned}$$

(46)

the chiral operator also can be block diagonalized according to the same mirror eigenvalue. Then, the winding number $W_{\lambda }$ in each mirror eigenvalue sector can be defined. The condition in Eq. (46) is satisfied only when the pairing potential is mirror even. The reason is as follows: In our convention, the mirror operator for BdG Hamiltonian is defined as ${\tilde{M}} = {\text {diag}}[ M, \eta _M s_y M^* s_y ]$ where M and $s_y M^* s_y$ are mirror operators for electron part and hole part, respectively. $\eta _M=\pm 1$ is the mirror parity of a pairing potential, which is given in Table 4. Because the mirror operator commutes with the time-reversal operator $[ T , M] = 0$, all the mirror operator satisfies $s_y M^* s_y = M$. Then, ${\tilde{M}} = M \tau _0$ (${\tilde{M}} = M \tau _z$) for the mirror-even-parity (mirror-odd-parity) pairing potential. Thus, only the mirror-even-parity superconducting phase satisfies the condition of Eq. (46).

Furthermore, the mirror chiral winding number can be defined as $W_{M} = W_{i} - W_{-i}$, where $W_{\lambda }$ is the chiral winding number for each block having a mirror eigenvalue $\lambda$. The mirror chiral winding number $W_{M}$ can also be defined for a path ${\mathscr {C}}$ that encloses the Dirac point in the mirror plane as shown in Fig. 3b. When the path ${\mathscr {C}}$ is parametrized by $\theta \in [0, 2 \pi )$, the mirror chiral winding number is given by^37,61

$$\begin{aligned} W_{M} = \frac{-1}{4 \pi } \oint _{{\mathscr {C}}} d \theta {\text {Tr}} \left[ {\tilde{M}} \Gamma {\mathscr {H}}^{-1} ({\mathbf {k}}(\theta )) d {\mathscr {H}}({\mathbf {k}}(\theta )) \right] . \end{aligned}$$

(47)

$\Delta _2$ and $\Delta _3$ phases

In the absence of lattice distortions, the $C_{4z}$ symmetry protects the nodal points by assigning different eigenvalues^36,37. The same nodal points are also topologically protected by the mirror chiral winding number in Eq.(47) because the $\Delta _2$ and $\Delta _3$ pairing potentials are mirror-even. For $\Delta _{3}$ pairing potential, which is mirror-even under $M_{xz}$ and $M_{yz}$, the calculated mirror chiral winding numbers around the nodal points are $\pm 2$ (Fig. 3c). Similarly, the nodal points in the $\Delta _{2}$ phase are topologically protected by $M_{110}$ and $M_{1 {{\bar{1}}} 0 }$ mirror chiral winding numbers.

Even though $C_{4z}$ symmetry is broken under lattice distortions, the mirror chiral winding number topologically protects the nodal points if the corresponding mirror symmetry is unbroken. For example, consider $D_{2h}$ point group which has $M_{xz}$ and $M_{yz}$ mirror symmetries. Among $\Delta _{2}$ and $\Delta _{3}$ pairings, $\Delta _{3}$ pairing is mirror even under $M_{xz}$ and $M_{yz}$. Thus, the nodal points in the $\Delta _{3}$ phase are topologically protected by the corresponding mirror chiral winding numbers (Fig. 3c).

Furthermore, the nodal points are positioned on the $k_z$ axis because $C_{2z}$ symmetry gives an additional constraint as follows: Let $W_M( {\mathbf {k}})$ denote a mirror chiral winding number at ${\mathbf {k}}$. Then, the mirror chiral winding number at $C_{2z} {\mathbf {k}}$ is related with that at ${\mathbf {k}}$ by

$$\begin{aligned} W_M( C_{2z} {\mathbf {k}}) = \eta _{C_{2z}}W_M( {\mathbf {k}}), \end{aligned}$$

(48)

where $\eta _{C_{2z}}$ is the parity of the pairing potential under $C_{2z}$ transformation. The detail derivation is in Sec. S4 in Supplementary Information. Since $\eta _{C_{2z}}=1$ for $\Delta _{2}$ and $\Delta _{3}$, $W_M( {\mathbf {k}}) =W_M( C_{2z} {\mathbf {k}})$, which means that the mirror chiral winding numbers are the same for the two nodal points that are related by $C_{2z}$ rotation. Now, let us assume that a nodal point on the $k_z$ axis in the absence of lattice distortions deviates from the $k_z$ axis under the $n_2$ type lattice distortion. Due to the $C_{2z}$ symmetry, there exists another nodal point having the same mirror chiral winding number. Thus, the total mirror winding number under lattice distortion becomes twice the original winding number, which is a contraction with the topological charge conservation. Therefore, the nodal points should be located on the $k_z$ axis under the $n_2$ type lattice distortion.

A similar argument can be applied to the $D_{2h}'$ case having $M_{110}$ and $M_{1{{\bar{1}}}0}$ mirror symmetries. The nodal points in the $\Delta _{2}$ phase is topologically protected by the $M_{110}$ and $M_{1{{\bar{1}}} 0}$ mirror chiral winding numbers and the nodal points are located on the $k_z$ axis due to the $C_{2z}$ symmetry (Figs. 2c and 3b). For $C_{2h(x)}$ case, $M_{yz}$ is unbroken while $C_{2z}$ is broken. Thus, nodal points on $M_{yz}$ plane in $\Delta _3$ phase are protected by the $M_{yz}$ mirror chiral winding number and can be deviated from $k_z$ axis due to the $C_{2z}$ symmetry breaking (Fig. 2e).

$\Delta _{41}$ and $\Delta _{42}$ phases

In the absence of lattice distortions, the nodal points in each $\Delta _{41}$ and $\Delta _{42}$ phases (Fig. 2a) are accidental nodal points because a single phase, either $\Delta _{41}$ or $\Delta _{42}$ phase, would break the $D_{4h}$ point group symmetry spontaneously. Only if we neglect such lattice symmetry breaking, the accidental nodal points can be understood to be protected by the different eigenvalues of $C_{2z}$ and $s_z$ symmetry operators (see the details in Sec. S3 in Supplementary Information). Note that the existence of the accidental point nodes also can be verified via $C_{4z}$ symmetry^36,37. In the viewpoint of topological winding numbers, the mirror chiral winding numbers are zero in the absence of lattice distortions (Fig. 3d) . Due to the $C_{2z}$ symmetry, Eq. (48) gives

$$\begin{aligned} W_M( 0, 0, k_z) = - W_M( 0, 0, k_z), \end{aligned}$$

(49)

which implies that $W_M=0$ on the $k_z$ axis. Here, $\eta _{C_{2z}} = -1$ is used for $\Delta _{41}$ and $\Delta _{42}$. Thus, the nodal points are not topologically protected for $D_{4h}$ case.

Figure 4

Surface band structures of superconducting phases under distortions. Surface band structures on the (010) surface for (a–d) $D_{4h}$, (e–h) $D_{2h}$, (i–l) $C_{2h(z)}$ and (m–p) $C_{2h(x)}$. In each panel, the upper figure indicates the close-up view of the band structure near $E=0$ corresponding to the red box in the lower figure. The red vertical arrows indicate the nodal points of the bulk superconducting states. In the insets of (e, h, i, l), the bulk states are gapped. The cyan vertical arrows indicate the gapped surface states. In (b, f, j, k), red horizontal lines show the surface flat bands. The nature of gapless surface state (GSS) is distinguished by the colored circle: Red ones in (a, b, e, f), green ones in (a, b, d–f, h, i, l), and black ones in (d, h, i, l) indicate GSS’s protected by mirror Chern numbers, zero-dimensional topological numbers, and mirror eigenvalues, respectively. In (j), GSS’s are accidental. The details are in Table 7 and in the main text. Region I, II, and III are $(0,k_1)$− $(0,k_2)$, $(k_2, 0)$− (0, 0), (0, 0) − $(\pi /2,0)$, respectively, where $k_1$ and $k_2$ ($k_1>k_2>0$) indicate two intersecting points between the upper Fermi surface and the $k_z$ axis.

However, under lattice distortions, nodal points can be topologically protected by the mirror chiral winding number. Let us consider the $D_{2h}$ point group under the $n_1$ type lattice distortion. Since $\Delta _{41}$ and $\Delta _{42}$ pairings are mirror-even under $M_{yz}$ and $M_{xz}$ operations, the corresponding mirror chiral winding number protects nodal points in each mirror plane (Fig. 2b). The calculated mirror chiral winding numbers are $W_M=\pm 2$ (Fig. 3d). Note that the nodal points are off the $k_z$ axis and the calculated mirror chiral winding numbers satisfy Eq. (49). For the $D_{2h}'$ case, $M_{110}$ and $M_{1 {{\bar{1}}} 0}$ mirror chiral winding numbers ($W_M = \pm 2$) protect nodal points in corresponding mirror planes for the superconducting phases having $\Delta _{41} \pm \Delta _{42}$ pairing potentials (Fig. 2c). For the $C_{2h(z)}$ case, all the relevant mirror symmetries are broken and hence there are no topologically-protected nodal points (Fig. 2d). For the $C_{2h(x)}$ case, there are unbroken $M_{yz}$ and $C_{2x}$. Thus, $M_{yz}$ mirror chiral winding numbers ($W_M=\pm 2$) protect the nodal points but the nodal points need not to be located symmetrically with respect to the $k_z$ axis (Fig. 2e). These nodal points in the $C_{2h(x)}$ case can be understood from the nodal points in the $D_{2h}$ case: Among four nodal points in the $D_{2h}$ case, two nodal points are pair-annihilated, and only two nodal points survive in the $C_{2h(x)}$ case.

Finally, we discuss a gap structure change of $\Delta _{42}$ phase under $n_1$ type lattice distortion (Fig. 3d). When $n_1=0$, each nodal points has $W_M=0$ and a quadratic energy-momentum dispersion relation along the $k_x$. With the increasing lattice distortion, nodal points with $W_M=\pm 2$ are created pairwise from a nodal point with $W_M=0$, and linear energy-momentum dispersion relation for all three momentum directions appears. Similar gap structure changes occur under the other lattice distortions.

Table 7 Gapless surface Andreev bound state (SABS) on (010) surface.

Surface spectrum

Surface Andreev bound state (SABS) in superconducting phases of the topological DSM have been studied in the absence of lattice distortion³⁷. In this subsection, we systematically investigate SABS in superconducting phases under lattice distortions. There are four types of gapless surface Majorana states under lattice distortions. Three types are topologically protected by mirror chiral winding, mirror Chern, and zero-dimensional winding numbers. The fourth type is protected by mirror symmetry and corresponding eigenvalues.

Using the Möbius transformation based method⁶², we calculate the surface band structures. Figure 4 shows the numerically obtained surface spectra for (010) surface in various superconducting phases under lattice distortions. For $\Delta _{1}$ and $\Delta _{1}'$ phases, there is no SABS; $\Delta _{1}$ phase is fully gapped and topologically trivial, and $\Delta _{1}'$ phase has two nodal lines having opposite chiral winding numbers as shown in Fig. 3a, which does not have protected SABS because of the positions and shapes of two nodes in momentum space. On the other hand, $\Delta _{2}$, $\Delta _{3}$, $\Delta _{41}$, and $\Delta _{42}$ have various types of SABS (Fig. 4), which are summarized in Table 7.

Without loss of generality, we will focus on the (010) surface and the surface Brillouin zone $(k_{x}, k_{z})$. A similar analysis for the (010) surface can be easily applied to the other surfaces such as (100), (110) planes, because the results for the other plane only depend on the mirror symmetries and the transformation properties of the pairing potentials under the unbroken symmetries. For convenience, we consider the surface states in the three regions: Region I, II, and III, which are $(0,k_1)$-$(0,k_2)$, $(k_2, 0)$-(0, 0), (0, 0)-$(\pi /2,0)$, respectively. Here, $k_1$ and $k_2$ ($k_1>k_2>0$) indicate two intersecting points between the upper Fermi surface and the $k_z$ axis.

First, we consider the flat SABS in the Region I, which is topologically protected by the nontrivial mirror chiral winding number in Eq. (47). For example, let us consider $\Delta _{3}$ phase and $M_{yz}$ mirror symmetry. For $D_{4h}$, $D_{2h}$, and $C_{2h(x)}$ cases, $M_{yz}$ mirror is unbroken and $\Delta _{3}$ has odd parity under $M_{yz}$, which leads to the opposite mirror chiral winding numbers ($W_M=\pm 2$) for two nodal points near the upper Fermi sphere as shown in Fig. 3c. Then, there exists a flat SABS on (010) surface as shown in Fig. 4b, f, j. To understand such SABS on (010) surface, the mirror winding number $W_{M}(k_z)$ along the mirror invariant $k_z$ axis is defined as³⁷

$$\begin{aligned} W_{M}(k_z) = \frac{-1}{4 \pi i } \int _{-\pi }^{\pi } d k_y {\text {Tr}} \left[ {\tilde{M}} \Gamma {\mathscr {H}}^{-1} d_{k_y} {\mathscr {H}} \right] , \end{aligned}$$

(50)

which is nontrivial between nodal points. Therefore, between the nodal points, there exists a flat SABS. Similarly, for $\Delta _{41}$ phase, $M_{yz}$ mirror symmetry gives nontrivial mirror chiral winding numbers, which guarantees the existence of the zero-energy flat SABS in the Region I on (010) surface (Fig. 4k). Note that, under the $n_3$ type lattice distortion, the mixture of $\Delta _{3}$ and $\Delta _{41}$ phases are allowed. But the flat SABS is still present due to the $M_{yz}$ mirror chiral winding number.

Second, we consider the gapless SABS protected by the mirror Chern number $C_M$. The topological mirror superconducting phases^36,63 are allowed for $\Delta _{2}$ and $\Delta _{3}$ phases because $\Delta _{2}$ and $\Delta _{3}$ pairing potentials are $M_{xy}$ mirror-odd and the corresponding mirror Chern numbers for each mirror eigenvalue block are nontrivial. Under the $n_1$ ($n_2$) type lattice distortion, $\Delta _{2}$ ($\Delta _{3}$) phase is fully gapped, and the mirror Chern number defined in $M_{xy}$ plane is nontrivial ($C_{M} = \pm 2$), which leads to a topologically-protected Majorana states on $M_{xy}$ plane. For example, see the surface spectra in the Region III in Fig. 4e.

Third, we consider the gapless SABS protected by the zero-dimensional topological number. Since $\Delta _{2}$ and $\Delta _{42}$ pairings are odd under $M_{yz}$, a zero-dimensional topological number $\rho (k_x)$ can be defined using $\Gamma M_{yz}$^36,37. Then, the zero-dimensional topological number protects the gapless state in the Region III. See the surface spectra at the Region III in Fig. 4d, h, i, l and Table 7. Similarly, $\Delta _{2}$ and $\Delta _{3}$ pairings are odd under $M_{yz}$, a zero-dimensional topological number $\rho (k_z)$ is defined using $C M_{yz}$^36,37, which protects the gapless states in the Region II for $D_{4h}$ and $D_{2h}$ cases. See the surface spectra at the Region II in Fig. 4a, b, e, f and Table 7.

Fourth, we consider the gapless SABS protected by mirror eigenvalues. If the pairing potential has an odd parity under the mirror operation, the mirror eigenvalues for the electron and hole bands are different, which protects the band crossing of surface states^36,37. For example, consider $\Delta _{2}$ phase and $M_{yz}$ symmetry. Because $(k_x, k_y, k_z) \rightarrow (-k_x, k_y, k_z)$ under $M_{yz}$, the mirror eigenvalues are properly defined on the $k_z$ axis. Moreover, $\Delta _{2}$ pairing has odd parity under $M_{yz}$ symmetry. Hence, the different mirror eigenvalues protect the gapless states in the Region I. See Fig. 4a, e, i. Similarly, $\Delta _{42}$ phases has odd parity under $M_{yz}$, which protects the gapless states in the Region I and II. See Fig. 4d, h, l.

In summary, we find the various types of surface states depending on the pairing potentials and lattice distortions. Even under the lattice distortions, most of the inversion-odd-parity superconducting phases have gapless SABS, which may be observed as zero bias conductance peak (ZBCP) in experiments.

Superconducting critical temperature and phase diagram

In this subsection, we study superconducting critical temperatures and their enhancements under lattice distortions. We also investigate the phase diagram for the various superconducting phases under lattice distortion.

In the weak-coupling limit, the superconducting critical temperature $T_c$ can be calculated by solving the linearized gap equation and a phase diagram for various pairing potentials is obtained by comparing the critical temperatures^{37,55,56,57,58}. The linearized gap equation can be expressed using the pairing susceptibility^{37,55,56,57,58}. The pairing susceptibility $\chi _i$ for each pairing potential $\Delta _i$ is given by

$$\begin{aligned} \chi _i(T) = - \frac{1}{\beta } \sum _{\omega _n} \sum _{{\mathbf {k}}} {\text {Tr}} [ (\Delta _i \tau _x) G_0({\mathbf {k}}) (\Delta _i \tau _x) G_0({\mathbf {k}})]. \end{aligned}$$

(51)

Here, $\beta = 1/(k_B T)$ is the inverse temperature, $k_B$ is the Boltzmann constant, $\omega _n$ is the Matsubara frequency, and $\Delta _i$ is the matrix representation of a pairing potential listed in Table 4. $G_0 ({\mathbf {k}}) = \frac{ {\mathscr {P}}_{{\mathbf {k}}}}{i \omega _n - \varepsilon _{{\mathbf {k}}}}$ is the single-particle Green’s function of the normal state and ${\mathscr {P}}_{{\mathbf {k}}} \equiv \sum _{m=1,2} \left| \phi _{m,{\mathbf {k}}} \right\rangle \left\langle \phi _{m,{\mathbf {k}} }\right|$ is the projection operator onto the two degenerate Bloch states in the conduction bands. Here, $\varepsilon _{{\mathbf {k}}} =\left| a({\mathbf {k}} )\right| - \mu$. Then, the superconducting susceptibility has the following generic form:

$$\begin{aligned} \chi _i(T) = \int \frac{ d ^3 {\mathbf {k}} }{(2\pi )^3} f_i({\mathbf {k}}) \frac{ \tanh ( \beta \varepsilon _{{\mathbf {k}}}/2) }{ 2 \varepsilon _{{\mathbf {k}}} }, \end{aligned}$$

(52)

where $f_{i} ({\mathbf {k}})$ is momentum dependent form factor. The explicit expressions for form factors $f_i ( {\mathbf {k}})$ are given in Sec. S5 in Supplementary Information.

With these susceptibilities, we now solve the linearized gap equation. The linearized gap equations are obtained by minimizing the mean-field free energy in the weak coupling limit. Since superconducting critical temperatures with pairing potentials in the same classes are not independent, the $\Delta _{i}$’s in the same class can appear in the same linearized gap equation.

First, consider the gap equation in the absence of lattice distortions. According to the irreducible representation of $D_{4h}$, $\Delta _{1}$’s, $\Delta _2$, $\Delta _3$ and $\Delta _4$’s belong to $A_{1g}$, $B_{1u}$, $B_{2u}$ and $E_u$ irreducible representations (see Table 5). Then, the gap equations are given by

$$\begin{aligned}&\left| \begin{array}{cc} U \chi _{1} (T_c) - 1 &{} U \chi _{1,1'} (T_c) \\ U \chi _{1,1'} (T_c) &{} U \chi _{1'} (T_c) -1 \end{array} \right| = 0, ~~~{\text {for}}~{ \Delta _{1}}~ {\text {and}}~{ \Delta _{1}'}~ {\text {phases}}, \end{aligned}$$

(53)

$$\begin{aligned}&V \chi _2 (T_c) - 1 =0,~~~ V \chi _3 (T_c) - 1 = 0, ~~~~~{\text {for}}~ {\Delta _{2}} ~{\text {and}}~{ \Delta _{3}}~{\text { phases}}, \end{aligned}$$

(54)

$$\begin{aligned}&\left| \begin{array}{cc} V \chi _{41} (T_c) - 1 &{} V \chi _{41,42} (T_c) \\ V \chi _{41,42} (T_c) &{} V \chi _{42} (T_c) -1 \end{array} \right| = 0, ~~~~~{\text {for}}~{ \Delta _{41}}~{\text { and}}~{ \Delta _{42}}~{\text { phases}}, \end{aligned}$$

(55)

where $\chi _{i, j}$ is the generalized superconducting susceptibility for mixed pairings $\Delta _{i}$ and $\Delta _{j}$ by replacing the second $\Delta _{i}$ with $\Delta _{j}$ in Eq. (51). Using the low-energy effective Hamiltonian in Eq. (13), the superconducting susceptibility can be further simplified and hence one can solve the gap equation analytically. Using an ellipsoidal coordinate, the superconducting susceptibility can be represented as a product of two independent integrals (see more details in Sec. S5 in Supplementary Information):

$$\begin{aligned} \chi _i = {\mathscr {R}} (\beta _c) \Omega _i(\mu ). \end{aligned}$$

(56)

Here, the radial integral part ${\mathscr {R}} (\beta _c)$ is given by

$$\begin{aligned} {\mathscr {R}} (\beta _c) = \int _{-\omega _D}^{\omega _D} d E ~ \frac{ \tanh ( \beta _c E /2) }{ E }, \end{aligned}$$

(57)

where E is an integration variable and $\omega _D$ is the energy cutoff of the pairing potential. The angular integral part $\Omega _i(\mu )$ is given by

$$\begin{aligned} \Omega _i(\mu ) = \int _{0}^{ \pi } \int _{0}^{ 2 \pi } \frac{d\theta d\phi }{(2\pi )^3} \left| \frac{ 2 \mu ^2 \sin \theta }{ v^2 v_z } \right| f_i( r = \mu ,\theta , \phi ), \end{aligned}$$

(58)

where the form factor $f_i ({\mathbf {k}})$ is represented as a function of $r, \theta$ and $\phi$ in the ellipsoidal coordinates. After the integration over $\theta$ and $\phi$, the susceptibilities can be obtained as follows:

$$\begin{aligned} \chi _{1} = 4 \pi C_0 {\mathscr {R}}(T_c), ~ \chi _{1'} = \frac{4 \pi }{3} C_0 {\mathscr {R}}(T_c), ~ \chi _{2} = \chi _{3} = \frac{8 \pi }{3} C_0 {\mathscr {R}}(T_c), ~ \chi _{41} = \chi _{42} =\frac{4 \pi }{3} C_0 {\mathscr {R}}(T_c), ~ \chi _{1,1'} = \chi _{41,42}=0, \end{aligned}$$

(59)

where $C_0 = \frac{2}{(2\pi )^3}\frac{\mu ^2}{v^2 v_z }$. Then, the linearized gap equations are given by

$$\begin{aligned}&\chi _{1} (T_c) = \chi _{1'}(T_c) = 1/U, \end{aligned}$$

(60)

$$\begin{aligned}&\chi _2 (T_c) = \chi _3 (T_c) = \chi _{41} (T_c) = \chi _{42} (T_c) = 1/V. \end{aligned}$$

(61)

If we denote the critical temperature $T_c^{(i)}$ for a pairing potential $\Delta _i$, then the gap equations are given by

$$\begin{aligned}&{\mathscr {R}}(T^{(1)}_c) = \frac{1}{3} {\mathscr {R}}(T^{(1')}_c) = \frac{1}{4 \pi U C_0}, \end{aligned}$$

(62)

$$\begin{aligned}&{\mathscr {R}}(T^{(2)}_c) = {\mathscr {R}}(T^{(3)}_c) = \frac{1}{2} {\mathscr {R}}(T^{(41)}_c) = \frac{1}{2} {\mathscr {R}}(T^{(42)}_c) = \frac{3}{ 8\pi V C_0}. \end{aligned}$$

(63)

Because ${\mathscr {R}}(x)$ is a monotonically decreasing function with respect to x, $T_{c}^{(1)} > T_{c}^{(1')}$ and $T_{c}^{(2)} = T_{c}^{(3)} > T_{c}^{(41)} = T_{c}^{(42)}$. Thus, the highest $T_c$ is determined among $T_{c}^{(1)}, T_{c}^{(2)}$ and $T_{c}^{(3)}$. Because the critical temperatures are same at the phase boundary, the phase boundary in Fig. 5a is determined by the equation ${\mathscr {R}}(T^{(1)}_c) = {\mathscr {R}}(T^{(2)}_c) = {\mathscr {R}}(T^{(3)}_c)$, which gives the critical value of $U/V = 2/3$.

Figure 5

Phase diagrams for the tetragonal and orthorhombic crystal systems. (a) Superconducting phase diagram in the U and V plane in the absence of lattice distortions when $\mu /t_z=0.1$. In the orange (blue) region, $\Delta _2$ or $\Delta _3$ ($\Delta _1$) phase is dominant. The slope of the phase boundary is approximately $U/V=2/3$. The white region indicates a non-superconducting phase. (b) The numerically calculated critical value of U/V ratio as a function of the chemical potential in the absence of lattice distortions. Since $\mu = 0.75 t_z$ is the band inversion point, there is a local maximum due to Van Hove singularity near $\mu = 0.75 t_z$. (c, e) Phase diagrams with respect to (c) $n_1$ and (e) $n_2$ type lattice distortions when $U=0.045 t_z$. The corresponding point groups are (c) $D_{2h}$ and (e) $D'_{2h}$. Each black arrow indicates the possible phase transition from an inversion-even-parity to inversion-odd-parity superconducting phases. (d, f) The normalized critical temperature $T_c / T_0$ for various pairing potentials with respect to (d) $n_1$ and (f) $n_2$ type lattice distortions. In both figures, $U/V=0.75$ and $U=0.045 t_z$, which corresponds to the black arrows in (c, e). $T_0$ is the critical temperature of the $\Delta _{1}$ phase in the absence of the lattice distortions.

When the chemical doping is low, the superconducting phase diagram for undistorted Dirac semimetal is shown in Fig. 5a. When the intra-orbital interaction U is strong, the conventional s-wave superconductivity with pairing potential $\Delta _{1}$ is the dominant phase. However, with the increasing inter-orbital interaction V, the unconventional superconducting phase with inter-orbital pairing potential $\Delta _{2}$ or $\Delta _{3}$ can emerge. Figure 5b shows the numerically obtained critical value of U/V ratio using the lattice Hamiltonian. Thus, by controlling the U/V ratio, both conventional and unconventional superconductivity can emerge for for the large range of chemical doping. The calculated value of U/V ratio is similar with 2/3 using the low-energy effective Hamiltonian, which means that $\Delta _{2}$ or $\Delta _{3}$ phase can emerge for the large range of chemical doping.

Next, consider the effect of $n_1$ and $n_2$ types of lattice distortions on the superconducting temperatures and the phase diagrams. When $n_1$ type lattice distortion is turned on, the point group becomes $D_{2h}$. In this case, only $\Delta _{1}$ and $\Delta _{1}'$ belong to the same $A_{g}$ class, and the others are belong to different classes (see Table 5). So the linearized gap equation is given by

$$\begin{aligned}&\left| \begin{array}{cc} U \chi _{1} (T_c) - 1 &{} U \chi _{1,1'} (T_c) \\ U \chi _{1,1'} (T_c) &{} U \chi _{1'} (T_c) -1 \end{array} \right| = 0, ~~~{\text {for}}~{ \Delta _{1}}~{\text { and}}~{ \Delta _{1}'}~{\text { phases}}, \end{aligned}$$

(64)

$$\begin{aligned}&V \chi _2 (T_c) = V \chi _3 (T_c) = V \chi _{41} (T_c) = V \chi _{42} (T_c) = 1,~~~{\text {for}}~ \Delta _{2}, \Delta _{3}, \Delta _{41}, {\text {and}}~ \Delta _{42}~{\text {phases}}. \end{aligned}$$

(65)

Similar to $D_{4h}$ case, the susceptibility can be analytically calculated when the chemical doping level is small. The relevant gap equations that determine the phase map are given by

$$\begin{aligned} {\mathscr {R}}(T^{(1)}_c) = \frac{1}{4 \pi U C_0},~~~~~~ {\mathscr {R}}(T^{(2)}_c) = \frac{1}{\frac{4\pi }{3} (2 + \frac{n_1^2 \sin ^2 k_0}{\mu ^2}) V C_0}. \end{aligned}$$

(66)

Thus, the phase boundary is given by

$$\begin{aligned} \frac{U}{V} = \frac{2 + (n_1^2 \sin ^2 k_0)/{\mu ^2} }{3}. \end{aligned}$$

(67)

Similarly, the other cases can be calculated. See the details in Supplementary Information.

Figure 5a shows the numerically calculated phase map in the absence of lattice distortions, which is consistent with the previous work³⁷. The superconducting phases are separated by the $U/V=2/3$ line as shown in Fig. 5a and the U/V ratio of the phase boundary depends on the chemical potential as shown in Fig. 5b. However, such superconducting phase is not reported in the real materials of $\hbox {Au}_2$Pb and $\hbox {Cd}_2\hbox {As}_3$ in the absence of lattice distortion. Such discrepancy might happen because either Tc is very low or the interaction strength is repulsive in real materials at ambient pressure.

Figure 5c–f shows the numerically calculated phase maps under the $n_1$ and $n_2$ types of lattice distortions using the low-energy effective Hamiltonian. The phase diagrams are plotted in the plane of the U/V ratio versus strength of $n_1$ or $n_2$ type lattice distortion. In each diagram, the dominant phases are conventional spin-singlet $\Delta _{1}$ phase and unconventional spin-triplet $\Delta _2$ or $\Delta _3$ phase depending on the parameters. When U/V is small (large) enough, $\Delta _2$ or $\Delta _3$ ($\Delta _{1}$) phase emerges. Remarkably, the unconventional superconductivity can emerge with increasing lattice distortions. As an example, near the phase boundary of $U/V \approx 0.7$, there is a phase transition between conventional superconducting $\Delta _{1}$ and unconventional superconducting $\Delta _{2}$ phases when $n_1$ increases (see the black arrow in Fig. 5c). To see this phase transition more clearly, we plot the normalized superconducting critical temperatures along the black arrow (Fig. 5d). When $n_1=0$, the $\Delta _{1}$ phase is dominant. With increasing $n_1$, the superconducting critical temperatures for $\Delta _{2}$ are increasing, which leads to the $\Delta _{2}$ superconducting phase under enough lattice distortion. Note that $T_c$’s for $\Delta _{1}$, $\Delta _{2}$, $\Delta _{41}$, and $\Delta _{42}$ increase while $T_c$ for $\Delta _{3}$ decreases with the increasing $n_1$ (Fig. 5d). This can be explained by the expectation values of the Cooper pairings and spin-orbital texture at the Fermi surface, which will be discussed later. Because $n_1$ and $n_2$ type lattice distortions are related with $\pi /4$ rotation, similar features are observed except for the exchange of $\Delta _2$ and $\Delta _3$ phases (Fig. 5e, f).

Figure 6

Phase diagrams for the monoclinic crystal system. (a–c) Phase diagrams with respect to U/V ratio and $n_3$ type lattice distortions for (a) $n_1=0.0$, (b) $n_1=0.05$, and (c) $n_1=0.1$. (d–f) The normalized critical temperature $T_c / T_0$ along the black arrows in (a–c). Here, $U/V=0.7$ and $T_0$ is the critical temperature of the $\Delta _{1}$ phase in the absence of the lattice distortions. In (d), the red and orange lines for $\Delta _{2}$ and $\Delta _{3}$ overlap.

For the $n_3$ type lattice distortion, similar features can be observed in Fig. 6. Under the $n_3$ type lattice distortion, $n_1$ type lattice distortions also can be involved as discussed before. Thus, we plot three representative phase diagrams for $n_1=0.0$, 0.05, and 0.1. Surprisingly, when $n_1=0$, $\Delta _{2}$ and $\Delta _{3}$ phases are degenerate, and they are dominant unconventional phases as shown in Fig. 6a,d. With increasing $n_1$, the region of the unconventional phase $\Delta _{2}$ increases (Fig. 6a–c) and the degenerate $\Delta _{2}$ and $\Delta _{3}$ phases become distinguishable.

Under $n_1$, $n_2$, and $n_3$ lattice distortions, the $T_c$’s of $\Delta _{2}$, $\Delta _{3}$, $\Delta _{41}$, and $\Delta _{42}$ increases much more than that of $\Delta _{1}$ (Figs. 5d, f, 6d–f), and hence the unconventional superconducting phases emerge. The mechanism of this will be discussed below.

Mechanism for $T_c$ enhancement of unconventional superconductivity

The $T_c$ enhancement of unconventional superconductivity under lattice distortions can be understood by the enhancement of DOS at Fermi surface and the enhancement of the expectation values of unconventional pairings at Fermi surfaces due to the unique spin-orbital texture.

First, we consider the increment of DOS at the Fermi surface. Under the lattice distortions, the DOS’s at the Fermi surface increase as shown in Eqs. (20) and (22). Then, the superconducting critical temperature increases under lattice distortions because $T_c \propto e^{- \frac{1}{g N(0)} }$. Here, g is the strength of the pairing potential in the standard BCS theory and N(0) is the DOS at Fermi surface. Due to this enhancement of DOS, most of the superconducting temperatures increase under the lattice distortions (see Figs. 5d, f and 6d–f). However, some unconventional superconducting temperatures decrease while some unconventional superconducting temperatures increase under lattice distortions. To understand this, we investigate the pairing expectation values for each superconducting pairing potentials.

As a representative example, we calculate the normalized expectation values for the $\Delta _{1}$, $\Delta _2$, and $\Delta _3$ pairings at the Fermi surface with and without the $n_1$ type lattice distortion (Fig. 7a, b). For a clear comparison, the differences $\Delta _i^{\text {diff}} \equiv \left\langle \Delta _i\right\rangle _{n_1 \ne 0} - \left\langle \Delta _i\right\rangle _{n_1=0}$ are calculated (Fig. 7c). Without lattice distortions, $\left\langle {\Delta }_{1}\right\rangle$ is uniform while $\left\langle {\Delta }_{2}\right\rangle$ and $\left\langle {\Delta }_{3}\right\rangle$ show zeros on the $k_z$ axis. With the $n_1$ type lattice distortion, $\left\langle {\Delta }_{2}\right\rangle$ increases while $\left\langle {\Delta }_{3}\right\rangle$ decreases, which leads to $\Delta _{2}^{\text {diff}} > 0$ and $\Delta _{3}^{\text {diff}} <0$ (Fig. 7c, d). On the other hand, $\Delta _{1}^{\text {diff}}= 0$. These behaviors of the expectation values of $\left\langle {\Delta }_{{i}}\right\rangle$ explains that the tendency of $T_c$ under lattice distortions. $T_c$ of $\Delta _{2}$ phase increase greater than that of $\Delta _{1}$ phase while $T_c$ of $\Delta _{3}$ phase decreases under $n_1$ type lattice distortion (Fig. 5a). Similarly, the effect of the other types of lattice distortions on $T_c$ can be understood by the expectation value change of the pairing potentials.

Figure 7

Expectation values of pairing potentials at the upper Fermi surface under the $n_1$ type lattice distortion. (a, b) The normalized expectation values $\left\langle {\Delta }_{i}\right\rangle$ of $\Delta _{1}$, $\Delta _{2}$, and $\Delta _{3}$ (a) without and (b) with the $n_1$ type lattice distortion are plotted at the upper Fermi surface of DSM in the $k_x$-$k_z$ plane. (c) The differences $\Delta _i^{\text {diff}} \equiv \left\langle {\Delta }_{i}\right\rangle _{n_1 \ne 0} -\left\langle {\Delta }_{i}\right\rangle _{n_1=0}$ are plotted. In (a–c), the black arrows indicate the points having zero expectation values. (d) The normalized integrated expectation values of each pairing potentials, $\int _{\text {FS}} d^2 k \Delta _i^{\text {diff}} / \int _{\text {FS}} d^2 k \left\langle { \Delta }_{i}\right\rangle _{n_1=0}$, are plotted with respect to $n_1$. Note that the upper Fermi surfaces encloses the Dirac point $( 0, 0, k_0$) as shown in Fig. 1(k–o).

Microscopically, we can understand the emergence of unconventional superconducting phases under lattice distortions as a result of the enhancement of inter-orbital pairing at the Fermi surface. Even though our argument can be applied to all distortions, we discuss the effect of $n_1$ type lattice distortion for convenience. We consider two Fermi surfaces encapsulating Dirac points $( 0, 0, \pm k_0)$ which are related by time-reversal and inversion. On the upper Fermi surface near the Dirac point $( 0, 0, + k_0)$, the Dirac Hamiltonian in Eq.(18) in the $k_x$-$k_z$ plane is given by

$$\begin{aligned} H_{\text {Dirac}} = ( n_1 \sin k_0 s_x + v k_x s_z ) \sigma _x + v_z (k_z-k_0) \sigma _z. \end{aligned}$$

(68)

The spin and orbital parts can be diagonalized separately and the total wavefunction can be represented by the product of spin and orbital wavefunctions³⁷:

$$\begin{aligned} \left| \Psi \right\rangle = \left| \phi \right\rangle _{\text {orbital}} \otimes \left| \psi \right\rangle _{\text {spin}}. \end{aligned}$$

(69)

Let us diagonalize the spin part. The spin part of Hamiltonian is given by

$$\begin{aligned} H_{\text {spin}} = ( n_1 \sin k_0 s_x + v k_x s_z ) = {\mathbf {h}} \cdot {\mathbf {s}}, \end{aligned}$$

(70)

where ${\mathbf {h}} = (n_1 \sin k_0 , 0, v k_x)$. Since this Hamiltonian is a product of momentum and spin operators, the spin wavefunction can be represented in the helicity basis $\left| \lambda \right\rangle _{\text {spin}}$ with $\lambda =\pm 1$:

$$\begin{aligned} H_{\text {spin}} \left| \lambda \right\rangle _{\text {spin}}&= \lambda \left| {\mathbf {h}} \right| \left| \lambda \right\rangle _{\text {spin}}. \end{aligned}$$

(71)

Next, we diagonalize the remaining orbital part. Depending on the spin helicity $\lambda$, the Hamiltonian in Eq. (68) can be written as follows:

$$\begin{aligned} H_{\text {orbital},~\lambda } = \lambda \left| {\mathbf {h}}\right| \sigma _x + v_z (k_z-k_0) \sigma _z = {\mathbf {d}}_{\lambda }\cdot \varvec{\sigma }, \end{aligned}$$

(72)

where ${\mathbf {d}}_{\lambda } = ( \lambda \left| {\mathbf {h}}\right| , 0, d_z)$ and $d_z = v_z (k_z-k_0)$. The orbital wavefunction can be represented by the pseudo-spin along ${{\hat{\mathbf{d}}}_{\lambda }}$. For each spin helicity $\lambda$, there are two orbital wavefunctions $\left| \kappa {{\hat{\mathbf{d}}}_{\lambda }}\right\rangle _{\text {orbital}}$ with $\kappa =\pm 1$ that satisfy the following equations:

$$\begin{aligned} H_{\text {orbital},~\lambda } \left| \kappa {{\hat{\mathbf{d}}}_{\lambda }} \right\rangle _{\text {orbital}} = \kappa d \left| \kappa {{\hat{\mathbf{d}}}_{\lambda }} \right\rangle _{\text {orbital}}, \end{aligned}$$

(73)

where $d = \sqrt{ \left| {\mathbf {h}}\right| ^2 + d_z^2 }$. When the chemical potential is positive, two degenerate wavefunctions located in conduction bands participate in the superconducting pairing. These wavefunctions are given by

$$\begin{aligned} \left| \Psi _1\right\rangle = \left| {{\hat{\mathbf{d}}}} _{+}\right\rangle _{\text {orbital}} \otimes \left| + \right\rangle _{\text {spin}},~~~~~~~~ \left| \Psi _2\right\rangle = \left| {{\hat{\mathbf{d}}}} _{-}\right\rangle _{\text {orbital}} \otimes \left| -\right\rangle _{\text {spin}}, \end{aligned}$$

(74)

which form a Kramer’s pair due to the PT symmetry regardless of lattice distortions: PT operation conserves the momentum while it flips helicity and the x-component of the orbital because $T= i s_y {\hat{K}}$ and $P=-\sigma _z$.

Since we have obtained the spin and orbital texture in one Fermi surface, we can obtain the spin and orbital texture of the other Fermi surface by applying either time-reversal or inversion operator. Let $\Psi ({\mathbf {k}})$ be a wavefunction on the Fermi surface. Because there is no $\sigma _y$ in the Hamiltonian Eq. (68), the time-reversal partner $T \Psi ({\mathbf {k}})$ has the same orbital direction and the opposite spin direction regardless of lattice distortions comparing with $\Psi ({\mathbf {k}})$. On the other hand, since $P=-\sigma _z$, the inversion partner $P \Psi ({\mathbf {k}})$ has the opposite $d_x$ while keeping $d_z$ and spin direction comparing with $\Psi ({\mathbf {k}})$. Figure 8 shows the numerically calculated spin and orbital textures using the lattice model. The P and T symmetry operators connects spin and orbital wavefunctions in Fig. 8. The red and green arrows indicate time-reversal and inversion pairs, respectively.

Figure 8

Spin and orbital textures without and with lattice distortion. (a, d) [(c,d)] Numerically calculated spin (orbital) textures at two Fermi surface surfaces. The $n_1$ type lattice distortion is absent in (a, c) and present in (b, d). In (a, b) [(c,d)], the spin (orbital) textures are represented by the small black (blue) arrows. In (a–d), the textures in left and right panels correspond to the spin helicity up and down wavefunctions, respectively. In (a, b), the red and green arrows indicate time-reversal and inversion pairs, respectively. In (c, d), the orange and blue arrows indicate the possible Cooper pairing between two electrons with opposite momenta. Note that the orbital pseudo-spin vectors connected by orange arrows are parallel regardless of the lattice distortion. On the other hand, the orbital pseudo-spin vectors connected by cyan arrows are parallel in (c) while non-parallel in (d).

Using these spin and orbital textures, let us investigate how the lattice distortions promote the unconventional pairings. The conventional $\Delta _{1}$ pairing is not affected by the lattice distortion. The expectation value of $\Delta _{1}$ is constant over the entire Fermi surface regardless of lattice distortion as shown in Fig. 7a, b. Because $\Delta _{1}= c^{\dag }_{1\uparrow }c^{\dag }_{1\downarrow }+c^{\dag }_{2\uparrow }c^{\dag }_{2\downarrow } + H.c.$ connects two wavefunctions that are related by time-reversal, the expectation value of $\Delta _{1}$ is constant due to TRS. In other words, because $\Delta _{1}$ is represented by the identity matrix $1_{4 \times 4}$, the expectation value of the $\Delta _{1}$ over the Fermi surface is constant even under the lattice distortions.

On the other hand, $n_1$ type lattice distortion can increase the expectation values of the inter-orbital pairing $\Delta _{2}$. For example, let us consider two wavefunctions located at the south pole of the upper Fermi surface ($k_z = + k_s$ with $k_s<k_0$) and the north pole of the lower Fermi surface ($k_z=-k_s$). Two wavefunctions are indicated by the orange and cyan arrows in Fig. 8c, d. At $k_z = \pm k_s$, the Dirac Hamiltonian in Eq.(18) is given by

$$\begin{aligned} H_{\text {Dirac}} ^{(\pm )} = \pm n_1 \sin k_0 s_x \sigma _x - v_z (k_0 - k_s) \sigma _z, \end{aligned}$$

(75)

where $H_{\text {Dirac}} ^{(+)}$ and $H_{\text {Dirac}} ^{(-)}$ correspond $k_z=k_s$ and $k_z=-k_s$, respectively. When $n_1=0$, wave functions on the conduction bands at $k_z=\pm k_s$ are given by

$$\begin{aligned}&\left| 2 \right\rangle _{\text {orbital}} \otimes \left| \uparrow _{x} \right\rangle _{\text {spin}}, ~~~ \left| 2 \right\rangle _{\text {orbital}} \otimes \left| \downarrow _{x} \right\rangle _{\text {spin}}, ~~~ {\text {for}} ~ k_z= + k_s, \end{aligned}$$

(76)

$$\begin{aligned}&\left| 2 \right\rangle _{\text {orbital}} \otimes \left| \uparrow _{x} \right\rangle _{\text {spin}}, ~~~ \left| 2 \right\rangle _{\text {orbital}} \otimes \left| \downarrow _{x} \right\rangle _{\text {spin}}, ~~~ {\text {for}} ~ k_z=-k_s, \end{aligned}$$

(77)

where $\left| 1 \right\rangle _{\text {orbital}}$ and $\left| 2 \right\rangle _{\text {orbital}}$ indicate the orbital basis for $\sigma$ matrix as defined before. $\left| \uparrow _{x} \right\rangle _{\text {spin}}$ and $\left| \downarrow _{x} \right\rangle _{\text {spin}}$ indicate the spin up and down along x-direction. Thus, the expectation value of inter-orbital pairing is zero for these wavefunctions because the orbital states of the wavefunctions in Eqs. (76) and, (77) are same. On the other hand, when $n_1 \ne 0$, the x-component of the orbital pseudo-spin is generated (indicated in the large cyan arrows in Fig. 8d). The wave functions at $k_z=\pm k_s$ are given by

$$\begin{aligned}&\left| \hat{{\mathbf {d}}}_{+}\right\rangle _{\text {orbital}}\otimes \left| \uparrow _{x}\right\rangle _{\text {spin}},~~~ \left| \hat{{\mathbf {d}}}_-\right\rangle _{\text {orbital}}\otimes \left| \downarrow _{x}\right\rangle _{\text {spin}},~~~ {\text {for}} ~ k_z=+k_s, \end{aligned}$$

(78)

$$\begin{aligned}&\left| \hat{{\mathbf {d}}}_-\right\rangle _{\text {orbital}}\otimes \left| \uparrow _{x}\right\rangle _{\text {spin}},~~~ \left| \hat{{\mathbf {d}}}_+\right\rangle _{\text {orbital}}\otimes \left| \downarrow _{x}\right\rangle _{\text {spin}},~~~ {\text {for}} ~ k_z=-k_s, \end{aligned}$$

(79)

where ${\mathbf {d}}_{\pm } = ( \pm n_1 \sin k_0 , 0, - v_z(k_0 - k_s))$. Therefore, under the lattice distortion, the expectation value of the inter-orbital pairing is allowed and $\Delta _{2}$ pairing is enhanced. This mechanism for the enhancement of unconventional pairings can be applied to the other cases. In summary, the emergence of unconventional superconductivity under lattice distortion can be understood due to the enhancement of inter-orbital pairings and DOS at Fermi surfaces.

Topological superconductivity of doped Dirac semimetal

Table 8 Possible topological superconductivity in doped DSM under lattice distortions.

As summarized in Table 8, we characterize possible superconducting states in doped Dirac semimetal by the gap structures, topological winding numbers, and surface spectra.

First, the conventional superconducting phase having $\Delta _{1}$ pairing potential can emerge. Because $T_c$ of the $\Delta _{1}$ phase increases under lattice distortions as shown in Figs. 5 and 6, conventional fully-gapped s-wave superconductivity can emerge.

Second, we consider the inversion-odd-parity superconductor. The BdG Hamiltonian in Eq. (23) are included in the DIII class according to 10-fold Altland-Zirnbauer classes^4,59 because $T^2 = -1, C^2 = +1$, and $\Gamma ^2 = +1$. With the additional inversion symmetry, the DIII class superconductor can be an inversion-odd-parity topological superconductor⁵⁷ classified by $Z_2$ invariants $(-1) ^ { w_{ \text {DIII} } }$, where

$$\begin{aligned} w_{\text {DIII}} = - \int \frac{d^3 k}{48 \pi ^3} \varepsilon _{\mu \nu \rho } {\text {Tr}}[ \Gamma (Q \partial _{\mu } Q) (Q \partial _{\nu } Q) (Q \partial _{\rho } Q) ]. \end{aligned}$$

(80)

Here, $\Gamma$ is the chiral operator, and Q is the so-called Q-matrix^4,59 (or projection matrix). The sufficient condition for realizing the inversion-odd-parity topological superconductor is that it has an inversion-odd-parity pairing with a full gap and its Fermi surface encloses an odd number of time-reversal-invariant momenta. In the absence of lattice distortions, the inversion-odd-parity pairings, $\Delta _{2}$, $\Delta _{3}$, $\Delta _{41}$, and $\Delta _{42}$, are not fully gapped (Fig. 2a) and cannot be such a topological superconductor. However, under the lattice distortions, these inversion-odd-parity phases can be fully gapped, and the sufficient condition above can be satisfied for the large chemical potential ($\mu > M_0$) because the Fermi surface can enclose only (0, 0, 0) in BZ. However, when the chemical potential is large with a lattice distortion, the band structure near the Fermi energy is far from that of DSM. Because we are discussing the Dirac physics, we do not consider such a superconducting phase in this work.

Third, topological mirror superconducting phases^36,63 can exist under lattice distortions. Topological DSM has a nontrivial mirror Chern number defined in the $M_{xy}$ plane and the corresponding surface states on the mirror-symmetric boundary^36,39. Similarly, topological mirror superconductivity for $\Delta _{2}$ and $\Delta _{3}$ phases can exist under lattice distortions. Under the $n_1$ ($n_2$) type lattice distortion, $\Delta _{2}$ ($\Delta _{3}$) phase is fully gapped, the $\Delta _{2}$ ($\Delta _{3}$) pairing potential is mirror-odd under the $M_{xy}$ symmetry, and the mirror Chern number defined in $M_{xy}$ plane is nontrivial ($C_{M} = \pm 2$), which leads to topological mirror superconductivity with a topologically-protected Majorana states on the mirror symmetric boundary. For example, see the gapless surface spectra of $\Delta _{2}$ and $\Delta _{3}$ phases in Region III in Fig. 4a, b, e, f. Due to the TRS and IS, this topological mirror superconductor is classified as $2{\mathbb {Z}}$.

Fourth, topological line nodal superconducting phases can exist under lattice distortions. As discussed in Fig. 3b, c, the inversion-even-parity $\Delta _{1}'$ pairing allows a topologically-protected nodal lines protected by the chiral winding number in Eq. (43). According to this chiral winding number, in general, the topological line nodal superconductor in doped topological DSM is classified as $2{\mathbb {Z}}$. The reason is as follows. Since there are PT and PC, the nodal points are fourfold degenerate, which means that there are even number of winding source at the same points. Therefore, our generic model has a topological winding number of even integers. Note that the topological class of a line node in 3D DIII superconductor using Clifford algebra⁶¹ is $2{\mathbb {Z}}$, which is consistent with our result. However, there is no surface state because $\Delta _{1}'$ phase has two nodal lines having opposite chiral winding numbers (Fig. 3a).

Fourth, topological point nodal superconducting phases can exist under lattice distortions. For an inversion-odd-parity and mirror-even-parity pairing potential, we have a topological point nodal superconductor of which nodal points are protected by the mirror chiral winding number in Eq. (47). Because the chiral winding number is zero for inversion-odd-parity superconductor ($W = W_{\lambda =i} + W_{\lambda =-i} = 0$), the mirror chiral winding number is given by $W_{M} = W_{\lambda =i} - W_{\lambda =-i} = 2 W_{\lambda =i}$. From this mirror chiral winding number, this topological point nodal superconductor is classified as $2{\mathbb {Z}}$. Note that the classification of a point node using Clifford algebra^4,61 is $M{\mathbb {Z}}$ considering one mirror sector, which is consistent with our results.

Discussion

Now, we compare our results with experimental works in doped DSM of $\hbox {Au}_2$Pb^{30,31,32,33,34} and $\hbox {Cd}_3\hbox {As}_2$^27,28,29. $\hbox {Au}_2$Pb shows superconductivity at $T_c \approx 1.2$K with $D_{2h}$ symmetry at the ambient pressure^30,32,34. This structural transition corresponds to the $n_1$ or $n_2$ type lattice distortion. $T_c$ increases to 4 K until 5 GPa under compression³⁴. The point-contact measurements also reported that $T_c \approx$ 2.1 K using a hard contact tip is higher than the measured $T_c \approx$ 1.13 K using a soft tip. Assuming that the hard tip induces higher pressure than the soft tip, the experimental results are consistent with our result that $T_c$ is enhanced with increasing $n_1$ or $n_2$ lattice distortion (Fig. 5). The experiments reported that the superconductivity is either conventional^33,34 or unconventional³² depending on the physical situations. From our analysis, the superconducting phase of $\hbox {Au}_2$Pb is expected to be either a conventional fully gapped or unconventional topological mirror superconductor with a gapless SABS depending on physical parameters.

Similarly, in $\hbox {Cd}_3\hbox {As}_2$, the structural phase transition occurs near 2.6 GPa, resulting in a monoclinic lattice $C_{2h}$. Then, a superconductivity emerges at $T_c \approx 1.8$ K under pressure higher than 8.5 GPa. This structural transition corresponds to $n_3$ or $n_4$ type lattice distortion. When the pressure increases further, $T_c$ keeps increasing from 1.8 K (8.5 GPa) to 4.0K (21.3 GPa), which is consistent with the enhancement of $T_c$ under lattice distortions (Fig. 6). In this case, $n_1$ or $n_2$ can also be added without breaking the symmetry further. From our analysis, the superconducting phases of $\hbox {Cd}_3\hbox {As}_2$ are expected to be either a conventional or topological mirror superconductor with a gapless SABS. We emphasize that the topological nodal superconductor having a flat SABS can appear only if either $n_3$ or $n_4$ lattice distortion is turned on. The point-contact measurements for $\hbox {Cd}_3\hbox {As}_2$ showed the zero-bias conductance peak (ZBCP) and double conductance peaks symmetric around zero bias, which was interpreted as a signal of a Majorana surface states^27,28,29. Even though our result cannot directly explain the result of the point-contact measurement, the unconventional superconductivity having gapless Majorana fermion can emerge regardless of the lattice distortions according to the surface spectra (see Fig. 4), which seems to support the measured conductance peaks. Further experimental studies that reveal the nature of superconductivity are necessary, and our theoretical results will be a helpful guideline to interpret the experimental result and search for the possible topological superconductivity in DSM.

The other way to induce superconductivity in DSM is to use the proximity effect. Recently, an 1D proximity-effect-induced superconductivity in $\hbox {Cd}_2\hbox {As}_3$ Dirac semimetal nanowire-based Josephson junctions is reported, where the superconductivity is induced by the proximity effect from conventional s-wave superconductor^64,65. On the other hand, the lattice-distortion-induced superconductivity in this work is intrinsic and three-dimensional. Interestingly, a strain-induced topological superconductivity with Majorana bound states was also reported in 2D Dirac semimetals⁶⁶. We expect that combining these approaches would be helpful to generate and manipulate Majorana bound states, which is compatibly applicable to future quantum information technologies.

Summary

In this work, we have studied the possible symmetry-lowering lattice distortions and their effects on the emergence of unconventional superconductivity in doped topological DSM. From the group theoretical analysis, four types of symmetry-lowering lattice distortions that reproduce the crystal systems present in experiments are identified. We investigated the possible superconductivity under such symmetry-lowering lattice distortions considering inter-orbital and intra-orbital electron density-density interactions. We found that both conventional and unconventional superconductivity can emerge depending on the lattice distortion and electron density-density interaction. Remarkably, the unconventional inversion-odd-parity superconductivity hosts gapless surface Andreev bound states (SABS) even under lattice distortions. We found that the lattice distortion enhances the superconducting critical temperature. Therefore, our work is consistent with the observed structural phase transition and the enhancement of superconductivity in $\hbox {Cd}_3\hbox {As}_2$ and $\hbox {Au}_2$Pb under pressure. We also suggest that enhanced conventional and unconventional superconductivity in doped topological DSM can be controlled by physical parameters such as the pressure and strength of the superconducting pairing interaction. Thus, our work will provide a valuable tool to explore and control the superconductivity in topological materials.

Methods

To study the effects of symmetry-lowering lattice distortions, we assume a minimal $4 \times 4$ Hamiltonian that describes representative topological Dirac semimetals^39,45, where the lattice distortions are implemented as a perturbation⁵³. To study the superconductivity, we construct the Bogoliubov-de Gennes (BdG) Hamiltonian within the mean-field approximation while keeping TRS and the crystal symmetry^55,56. The momentum independent pairing potentials are classified using irreducible representations of the unbroken point group^{36,37,56,57,58}. The nodal structures, chiral winding number in Eq. (43), and chiral mirror winding number in Eq. (47) are calculated using the BdG Hamiltonian. The surface Green’s functions are calculated using a Möbius transformation-based method⁶². The superconducting critical temperature $T_c$ is calculated by solving the linearized gap equation in the weak-coupling limit^{37,55,56,57,58}. All the details are provided in the main text and Supplementary Information.