Topological Crystalline Materials of $J=3/2$ Electrons: Antiperovskites, Dirac points, and High Winding Topological Superconductivity
Takuto Kawakami Tetsuya Okamura Shingo Kobayashi Masatoshi Sato


APS/123-QED

 


takuto.kawakami@yukawa.kyoto-u.ac.jp
Yukawa Institute for Theoretical Physics, Kyoto University, Kyoto 606-8502, Japan

Yukawa Institute for Theoretical Physics, Kyoto University, Kyoto 606-8502, Japan

Institute for Advanced Research, Nagoya University, Nagoya 464-8601, Japan
Department of Applied Physics, Nagoya University, Nagoya 464-8603, Japan

msato@yukawa.kyoto-u.ac.jp
Yukawa Institute for Theoretical Physics, Kyoto University, Kyoto 606-8502, Japan



 
We present a theory of the high-spin generalization of topological insulators and their doped superconducting states.
The higher-spin topological insulators involve a pair of J=3/2 bands with opposite parity, 
and are characterized  by their band inversion.
The low-energy effective theory reveals that the topological insulators host four different phases characterized by 
mirror Chern numbers, at which boundaries two different patterns of bulk Dirac points appear.
For the carrier-doped case,
it is shown that the system may host unique unconventional superconductivity 
because of its high-spin nature and additional orbital degrees of freedom intrinsic to topological insulators. 
The superconducting critical temperature is evaluated by
using density-density pairing interactions, and odd-parity Cooper pairs are shown to be naturally realized in the presence of interorbital pairing interaction. 
It is observed that even the simplest spin 0 odd-parity pairing state exhibits 
a novel class of topological superconductivity—high winding topological superconductivity. 
We also discuss the experimental signals of high winding topological superconductivity in the case of the antiperovskite superconductor Sr_3-xSnO. 










Topological crystalline materials of J=3/2 electrons:  

antiperovskites, Dirac points, and high winding topological superconductivity
    Masatoshi Sato
    December 30, 2023
======================================================================================================================================



§ INTRODUCTION

The search for topological materials is a recent trend in
condensed matter physics <cit.>. 
A promising direction to fulfil this search is a multi-orbital
material with strong spin-orbit coupling.
Such a system can be topologically non-trivial through band inversion owing to
spin-orbit coupling. 
For instance, in the case of topological
insulators (TIs), the Z_2 topological indices 
are directly related to the number of time-reversal
invariant momenta at which band inversion occurs <cit.>.
Following this insight, 
numerous TIs including Bi_2Se_3 and Bi_2Te_3 have been
discovered experimentally <cit.>.

In the case of superconductors, 
even a single orbital system may possess a 
non-trivial topological phase because the charge conjugate counterpart of the single band coexists in the Bogoliubov–de Gennes (BdG)  Hamiltonian <cit.>.  
Any spin-triplet superconductor can be topological if the Fermi surface is properly chosen <cit.>. 
However, the search for spin-triplet superconductivity itself is a challenge
because of its strongly correlated origin.
No spin-triplet superconductor has been established yet, despite the existence of several promising candidates. 

Recently, it has been recognized that 
multi-orbital systems may solve this difficulty, as they
allow another mechanism of topological superconductivity <cit.>:
Using interorbital pairing interaction and spin-orbit coupling, 
a multi-orbital system may host
odd-parity Cooper pairs even without strong correlation <cit.>,
which indicates topological superconductivity <cit.>.
In this new direction,
various types of topological superconductivity  in doped topological materials
have been discussed recently <cit.>. 


In this study, we explore another class of topological
phases that can be realized only in J=3/2 bands.
Spin-orbit coupled electrons may behave as higher-order spin states, 
owing to the mixture of spin and orbital angular momentum. 
In crystals, the higher spin state can be identified as a J=3/2 state, 
as discrete crystalline rotation allows at most fourfold degeneracy of the J=3/2 spin.
Superconductivity of J=3/2 electrons
has been discussed for Luttinger semimetals with a single J=3/2 band such as half-Heusler systems <cit.>.
While the half-Heuslers may host interesting higher-spin Cooper pairs, 
their actual realization is restricted because the systems support
only a single J=3/2 band.

Here, alternatively, we consider systems with multiple J=3/2 bands.
A class of relevant materials is antiperovskite A_3BX, where A is (Ca, Sr, La), 
B is (Pb, Sn), and X is (C, N, O). 
The antiperovskite materials support two different J=3/2
electrons near the Fermi level, i.e., d-orbital electrons of the A atom and
p-orbital ones of the B.
The first-principle calculations show that the band inversion of these
two orbitals may occur at the Γ-point, accompanied with
three-dimensional (3D) Dirac points with a tiny gap <cit.>.
Owing to the band inversion, these materials become topological
crystalline insulators <cit.>.
It was also discovered very recently that
one of the antiperovskite topological materials, Sr_3SnO, shows
superconductivity with hole doping <cit.>.
Although discussion from the
viewpoint of J=3/2 electrons was missing, 
possible topological
superconductivity in Sr_3-xSnO was suggested
theoretically <cit.>, using the analogy of
the superconducting Dirac semimetal, Cd_3As_2 <cit.>.



Herein, we develop a general theory of higher-spin topological materials. 
Considering the application to antiperovskite,
we assume cubic symmetry, which allows for four-dimensional representation corresponding to the fourfold degeneracy of J=3/2.
Based on the low-energy effective Hamiltonian, we observe that the system hosts four different topological crystalline insulating phases, where the phase boundaries are characterized by two different patterns of gapless Dirac points. 
In addition to the octahedral Dirac points discussed previously <cit.>,
cubic Dirac points appear at one of the phase boundaries.  

We thereafter study the superconducting states of doped higher-spin topological materials.   
We first classify possible momentum-independent
gap functions, based on the cubic symmetry.
The gap functions may contain spin J=2
(spin-quintet) and J=3 (spin-septet) components, and ordinary J=0 (spin-singlet) and J=1 (spin-triplet) ones,
as Cooper pairs are formed by the J=3/2 electrons.
In contrast to the single J=3/2 band case as in half-Heusler materials, 
additional orbital degrees of freedom make it possible to obtain any higher-spin gap function 
in the framework of weakly correlated constant gap functions.
Subsequently, we evaluate the superconducting critical temperature T_ c for each gap
function by using simple density-density pairing interactions. 
It is observed that the system supports odd-parity superconductivity if the
interorbital pairing interaction is dominant.
The J=0 odd-parity pairing state with
A_1u representation yields the highest T_c among odd-parity superconducting states,
but the T_1u paring state consisting of J=3 and J=1 Cooper pairs also has high T_ c.

We reveal that the simplest J=0 odd-parity pairing state shows a new class of topological superconductivity
i.e., high winding topological superconductivity. 
As a 3D time-reversal invariant superconductor,
the odd-parity pairing state has a non-zero 3D winding number similar to the ^3He-B phase. 
However, in contrast to ^3He-B, owing to the higher spin nature of J=3/2 electrons, 
it does not support the minimal non-zero value but hosts various higher winding numbers 
depending on the model parameters. 
It is predicted that the topological superconductivity with a higher winding number 
exhibits a higher value of quantized longitudinal thermal conductance, 
which is proportional to the winding number <cit.>. 
Moreover, the accompanying multiple surface Majorana fermions 
display a new physics of conformal field theory with non-Abelian gauge field <cit.>. 
Furthermore, owing to a non-trivial spin texture on the Fermi surface 
originating from the winding number, we can expect non-trivial magnetic responses.
We also identify the topological phase diagram of the superconducting state 
and determine the relevant topological numbers and the patterns of surface Majorana fermions. 
We finally discuss the experimental signals of high winding superconductivity in the case of Sr_3-xSnO. 
We observe that Sr_3-xSnO is in the vicinity of topological phase transition, 
which predicts the characteristic nodal structure of the superconducting gap 
if it realizes the J=0 odd-parity pairing state.

This paper is organized as follows. 
In Sec. <ref>, we examine the normal state of higher-spin topological materials. 
Using the general k· p Hamiltonian involving two J=3/2 bands with opposite parity, we obtain the phase diagram of the J=3/2 topological crystalline insulators.
We reproduce the results for the antiperovskite oxides <cit.>, clarifying the topological origin of Dirac points in antiperovskites. 
We also observe that a novel cubic pattern of Dirac points is possible, in
addition to the octahedral Dirac points in antiperovskites. 
In Sec. <ref>, we examine the superconductivity of doped higher-spin TIs.
We classify gap functions 
in terms of an irreducible representation of the O_h point group, 
and evaluate the superconducting critical temperature
for each gap function. 
We demonstrate that the A_1u (J=0) representation is the most stable odd-parity pairing state, 
and the T_1u (J=3, or 1) is the second one. 
In Sec. <ref>, we identify the topological crystalline superconductivity of the A_1u state. 
Sec. <ref> discusses the application of our theory 
to Sr_3-xSnO, the case without band inversion, and the property of other pairing states, 
and subsequently, the conclusion is presented in Sec. <ref>.

 


§ NORMAL STATE



 §.§ Effective Hamiltonian


Here, we formulate higher-spin TIs. 
First, we define high-spin fermions. In the presence of spin-orbit coupling, 
electrons may behave as spin J fermions with J=|l± 1/2| owing to 
a mixture of the orbital angular momentum l and the electron spin 1/2. 
For electrons in crystals, however, a more precise definition of J is necessary, 
as crystals partially break rotation symmetry. Considering that spin J in the above sense  
indicates (2J+1)-fold degeneracy at k=0, we define spin J in crystals as a
(2J+1)-dimensional representation of the point group at the Γ point.
From this definition, it is observed that 
only a J=3/2 fermion is possible as a
higher-spin fermion, as point groups allow a maximum of four-
dimensional representations at Γ.
Furthermore, cubic symmetry is necessary to obtain the four-dimensional
representation.
Therefore, we consider a cubic crystal in the following analysis.
We also assume inversion symmetry for simplicity, which specifies the
cubic symmetry as the O_h group.

For ordinary spin 1/2 electrons, a topological insulating phase 
with inversion symmetry is
obtained via band inversion between the orbitals with different parities at
the Γ point. 
We can consider a similar band inversion mechanism even for J=3/2
electrons.
The O_h group hosts two different four-dimensional, i.e., J=3/2,
representations, G_3/2g and G_3/2u, corresponding to the even- and odd-parity states, respectively.
When these bands are inverted, the gap between
them closes at Γ, and thus, we can expect a non-trivial
topological phase transition through the band inversion.

This argument specifies a system for a higher-
spin TI.
It is symmetric under the O_h group, and consists of two
different 
J=3/2 bands corresponding to 
the G_3/2g and G_3/2u representations of the O_h group.  
Moreover, time-reversal symmetry should be assumed when we consider an
insulator.  

The single-particle state of the J=3/2 TI is
represented by an eight-
component spinor c_ k composed of the annihilation operator c_j_z, σ_z, k, 
where j_z=± 3/2, ± 1/2 denotes the z-component of J
corresponding to the four-dimensional representation of the O_h group,
σ_z specifies the parity of the representation, i.e., σ_z=1
(σ_z=-1) for the G_3/2g (G_3/2u) representation with even (odd) parity, 
which can be characterized as d, g, (p, f), or higher-order orbitals, 
and k is the momentum.
Generators of the O_h group 
consist of q-fold discrete crystal rotations C_q, n (q=2,3,4) and inversion P, which act on the Hamiltonian H_0( k) as

    C_q,nH_0(k) C_q,n^-1= H_0(D_C_q,n[k])  with C_q,n=e^-i2π/qJ·n,   
       PH_0(k) P^-1 = H_0(-k)  with P = σ_z.

Here, n denotes the rotation axis, D_C_q,n[k] indicates the q-fold rotation of k, J_i=x,y,z is the 4× 4 J=3/2 spin operator (see Appendix <ref>), and σ_i is the 2×2 Pauli matrix in the orbital space.
As the G_3/2u (G_3/2g) state changes (does not change) the sign under inversion, 
P is given as P=σ_z.
We also have time-reversal symmetry, 

    𝒯H_0(k) 𝒯^-1 = H_0(-k)  with 𝒯 =  C_2,y𝒦.

with the complex conjugate operator 𝒦.

In the following analysis, we use the low-energy effective Hamiltonian
for the J=3/2 electrons.  
To construct the effective Hamiltonian, we consider possible scalar matrices with respect to O_h symmetry.
First, we note that σ_z and σ_0 are scalars, whereas σ_x and σ_y are pseudo 
scalars.  
Any σ_ν (ν=0,x,y,z) is invariant under rotation, but σ_x and σ_y change their signs under inversion. 
Second, as the spin operators J_i (i=x,y,z) 
and their third-order polynomials J̃_i≡5/3∑_j≠ iJ_jJ_iJ_j-7/6J_i 
(see also Appendix <ref>)
behave as pseudo vectors for the O_h point group,
the inner product k·J and k·J̃ are pseudo scalars <cit.>.
Finally, we can construct a scalar as a product of any two of the above pseudo
scalars. 
Consequently, the effective Hamiltonian for the lower order of momentum
k is given as <cit.>

    H_0(k)	= m(k) σ_z + k· (v_1J + v_2J̃) σ_x - μσ_0,

where the σ_y term is absent owing to the time-reversal symmetry.





The effective mass m( k) in Eq. (<ref>) is given as 

    m(k) = -m_0 + α k^2 + [k·(β_1 J+β_2 J̃)] [k·(β_3 J+β_4J̃) ],

up to the second order of k, 
where 2m_0 is the energy gap between the G_3/2u and G_3/2g states at the Γ
point.
In the following analysis, we focus on the case where band inversion
occurs at the Γ point.
We consider m_0>0, α>0, and β_i=0 (i=1,2,3,4), which can
reproduce the band structure of Sr_3SnO qualitatively.
The case without band inversion will be briefly discussed in Sec. <ref>.

The second term of (<ref>) with σ_x leads to
orbital mixing between even- and odd-parity states, which is represented in 
green color in Figs. <ref>(a) and (b). 
This mixture opens the global energy gap around E=μ 
except for the parameter points displayed in Fig <ref>.
Note that the J̃ term reduces the continuous rotational symmetry to the discrete O_h one.


The σ_0 term in Eq. (<ref>) provides the Fermi level. 
In this section, to demonstrate the topological property of the band gap discussed above, 
we set μ=0, where the Fermi level lies in the band gap and the system falls into the insulating state.
The large |μ| pushes the energy gap below or above the Fermi level. 
In this case, the system is a metal rather than an insulator, 
although the topological property of the band gap remains unchanged.
This metallic state is important to induce topological superconductivity, which will be discussed later.

Note that H_0( k) cannot be an ordinary TI, 
although band inversion can occur at the Γ point.
The Z_2 number for TI is always trivial owing to the fourfold
degeneracy of the J=3/2 spin at Γ.
However, the system can be a topological crystalline
insulator, as shown below.





 §.§ Octahedral Dirac points



By diagonalizing the Hamiltonian (<ref>), we obtain the band
spectrum.
First, let us focus on the band spectrum along the primary axis,
say the k_z-axis. 
On the k_z-axis, we have [C_4, z, H_0(0,0,k_z)]=0 owing to the
fourfold rotation symmetry of the O_h group,
C_4,zH_0(-k_y,k_x,k_z)C_4,z^-1=H_0(k), and thus, 
the Hamiltonian is subdivided into the eigensector of C_4,z 
with the eigenvalue λ_C_4,z=e^-i j_z/2π, 

    H_0,j_z(k_z) = m_j_z(k_z)σ_z + b_j_zk_zσ_x,

where m_j_z(k_z)=m_-j_z(k_z) is an even function of k_z and b_j_z is the
following constant:

    b_± 3/2 = ± (3v_1-v_2)/2 ,
    
    	b_± 1/2 = ± (v_1+3v_2)/2 .

From time-reversal and inversion symmetry, 
the sectors with opposite j_z form a Kramers pair with the
same spectra.
The first term of Eq. (<ref>) provides a parabolic energy
dispersion for each orbital, and the second one mixes them.
As shown in Eq. (<ref>) (Eq. (<ref>)), the orbital mixing
in the j_z=± 3/2 (j_z=± 1/2) sector vanishes when
v_2=3v_1 (v_1=-3v_2). In such cases, Dirac points appear in
the band spectrum.

Figure <ref>(a) illustrates the band spectra
for v_2=3v_1.
The two parabolic dispersions 
with opposite orbital parity in the spin j_z=±3/2 sector
linearly cross each other at k_z=± k_0 with k_0 satisfying
m_± 3/2(k_0)=0.  
Consequently, there appear two Dirac points on the k_z axis. 
From the O_h group symmetry, their counterparts also exist on the
k_x and k_y axes, respectively, forming 
vertices of an octahedron in the momentum space. See
Fig. <ref>(c).    
Further, 
orbitals with opposite parity in the j_z=±1/2 sectors hybridize with each
other and thus, they do not form a Dirac point.
Similarly, when v_2=-v_1/3, octahedral Dirac points appear in the
j_z=±1/2 sectors.


In contrast to ordinary Dirac semimetals such as Cd_3As_2, these Dirac
points are easily gapped out by small perturbation even if 
crystalline symmetry is preserved.
Thus, they are not stable, but their existence is not accidental.

They indeed appear as gap-closing points accompanied by a topological phase
transition.
A key topological number is the mirror Chern number.
Consider one of three equivalent horizontal mirror reflections of the
O_h group, say the mirror
reflection with respect to the yz plane:

    M_xH_0(-k_x,k_y,k_z)M_x^-1=H_0(k),
 
with M_x=PC_2,x. See Fig. <ref> (c).    
On the mirror-invariant k_x=0 plane, it holds that [M_x,
H_0(0,k_y,k_z)]=0, and hence, the Hamiltonian has the block diagonal form 
H_0(0,k_y,k_z)=diag[H_0^+i(0, k_y, k_z),H_0^-i(0, k_y, k_z)] 
in the eigenbasis of M_x with the eigenvalues λ_M_x=± i.
Each sector may host its own Chern number


    ν_ Ch^± i = -i ∑_ξ_n^± i<0∫_k_x=0 d^2k ϵ^αβ⟨∂_α u_n,k^± i| ∂_βu_n,k^± i⟩,

where we use the solution of the Bloch equation H_0^± i|u_n,k^± i⟩=ξ_n^± i|u_n,k^± i⟩  
and two-dimensional Levi-Civita symbol ϵ^αβ with α,β∈{k_y,k_z}.
The mirror Chern number ν_σ_h for the horizontal mirror
reflection is defined by ν_σ_h=(ν^+i_ Ch-ν^-i_ Ch)/2.
From time-reversal symmetry, the total Chern number ν_
Ch^+i+ν_ Ch^-i vanishes, and thus, ν_σ_h is an
integer.


As summarized in Fig. <ref>(e), the mirror Chern number 
in the present system is
obtained as <cit.>

    ν_σ_h= 2sign(b_3/2b_1/2),

which jumps by four when Dirac points appear at v_2=3v_1 or
v_2=-v_1/3.
In other words, each of the four Dirac points on the k_x=0 plane change
the mirror Chern number by one. 
Therefore, these Dirac points are caused by the topological phase
transition. 
This indicates that the Dirac points must appear in the phase diagram because they separate two distinct topological phases, 
even though for any specific choice of parameters (i.e., one point in the phase diagram) they can gapped.
The small perturbation, such as mixing from other bands not considered in the low-energy effective Hamiltonian (<ref>), may open a gap in the Dirac points, but it
merely shifts the phase boundary of the topological phase transition,
which also hosts gapless Dirac points. 
Except at the topological phase transition, the system supports a
non-zero mirror Chern number, and thus,
the normal state is a topological crystalline insulator 
if the Fermi level is in the band gap.



 §.§ Cubic Dirac points


The effective Hamiltonian (<ref>) also reveals that another class
of Dirac
points is possible in higher-spin topological materials.
To demonstrate this, consider the diagonal axis, say the (111) axis. 
From the threefold rotational symmetry of the O_h group, 
the Hamiltonian along the (111) axis in momentum space 
is subdivided into eigensectors of C_3,(111) 
with the eigenvalues λ_C_3,(111)=e^-i 2/3j_111π,
where (111)=(x̂+ŷ+ẑ)/√(3) and j_111=± 1/2, ± 3/2.
Note that j_111=± 3/2 falls into the same eigensector with λ_C_3,(111)=-1.
The Hamiltonian in the subsector takes the same form as 
Eq. (<ref>), but k_z, j_z, and b_j_z in Eq. (<ref>)
are replaced with k_111=(k_x+k_y+k_z)/√(3), j_111, and b_j_111', respectively.
The interorbital coupling b_j_111' is given by

    b_± 1/2' = ±v_1-2v_2/2, 
          b_± 3/2' = [ ( 3/2v_1 + 1/3v_2 ) η_z +  5√(2)/6v_2 η_x ].

Here, η_i is the Pauli matrix in the 
j_111=± 3/2 space.

b_± 3/2' does not vanish for any v_i, whereas b_± 1/2' becomes zero 
when v_2=v_1/2.   
In the latter case, the orbitals with opposite parity do not hybridize with each
other along
the (111) axis, forming a Dirac point in
Fig. <ref>(b).   
By considering O_h symmetry, we have a total of eight Dirac
points at the vertices of a cube, as illustrated in
Fig. <ref>(d).

The cubic Dirac points are also understood as gap-closing points for
the topological phase transition.  
The relevant topological index in the present case is the mirror Chern
number for the diagonal mirror reflection of the O_h group.
For instance, consider the diagonal mirror reflection with respect to
the k_x+k_y=0 plane in Fig. <ref>(d).   
Similar to the octahedral Dirac case,
we can evaluate the mirror Chern number ν_σ_d with respect to
the diagonal
mirror plane. The result is shown in Fig. <ref>(f). 
At v_2=v_1/2, where there are four Dirac points on this plane,
the mirror Chern number jumps by four. 

Notably, two of the octahedral Dirac points appear on
the diagonal mirror plane when v_2=3v_1 or v_2=-v_1/3. 
Consequently, the diagonal mirror Chern number in Fig. <ref>(f)
jumps by two at v_2=3v_1
and v_2=-v_1/3.  





§ SUPERCONDUCTIVITY





 §.§ O_h classification of the gap function


With carrier doping, TIs can be superconductors at
low temperature. 
Here, we consider the superconducting states of doped higher-spin TIs.

To describe the superconducting states,
we use the Nambu space spanned by the basis 
(c_k, c̅_-k) 
with the eight-component spinor c_k of the annihilation operator 
and its time-reversal hole partner c̅_- k with the component 
c̅_j_z, σ_z, -k=∑_j_z' (C_2,y)_j_z,j_z'c_j_z', σ_z, -k^†.
In this basis, the BdG Hamiltonian is written as 

    H(k) = H_0( k)τ_z + Re(Δ( k))τ_x + Im(Δ( k))τ_y

where τ_μ represents the 2× 2 matrices acting on the particle and time-reversal hole space, 
and we have used time-reversal symmetry for the one-particle Hamiltonian 
𝒯H_0( k)𝒯^-1=H_0(- k). 

Note that, in this basis, the point group and time-reversal act on the gap
function Δ( k) in the same manner as Eqs. (<ref>), (<ref>), and (<ref>) on the normal Hamiltonian
H_0( k)
(see Appendix <ref> for more details).

We classify the multi-component gap function Δ( k).
For simplicity, we assume that Δ( k) is independent of the
momentum k as in the conventional BCS  theory. 
The higher spin and multi-orbital natures of the system enable us to consider
various unconventional Cooper pairs even in the weakly correlated
case.

The gap function is expanded as

    Δ = ∑_i, νΔ_i,νΦ_i,ν , 
      Φ_i,ν = φ_i(J) ⊗σ_ν/N,

where φ_i(J) is a set of 4× 4 matrices spanning the
spin space,
σ_ν (ν=0,x,y,z) is
the Pauli matrix in the orbital space, and 
N=√(tr[(φ_i ⊗σ_ν)^2]/8) is the
normalization constant.
A convenient basis of φ_i( J) is the 
spherical harmonic Y_l,m( J), which is defined by Y_l,m±1=∓[J_±,Y_l,m] with
J_±=J_x± iJ_y.  
In the present case, the spin of a Cooper pair is given by 3/2⊗ 3/2=0⊕
1⊕ 2⊕ 3, and hence, the azimuthal and magnetic
quantum numbers l and m in Y_l,m assume the values l=0,1,2,3 and
m=-l,…, l.
Reconstructing the basis Y_l,m for each l in terms of
the representations of the O 
group, we obtain φ_i(J) in
Table <ref>. 
Here, we use the O group rather than O_h, as inversion acts on the spin space trivially as the identity operator.






Notably, a possible combination of φ_i( J) and σ_ν in 
the right-hand side of Eq. (<ref>) is restricted by the
constraint 

    (Δ( k) C_2,y)^T =-Δ(- k) C_2,y,

from the Fermi statistics of c_ k.
As φ_i( J) for even (odd) l satisfies 
[φ_i(J)C_2,y]^T= -φ_i(J)C_2,y
([φ_i(J)C_2,y]^T= +φ_i(J) C_2,y),
it is combined only with the symmetric Pauli matrices   
σ_0, σ_x, and σ_z (the anti-symmetric Pauli matrix
σ_y). 
Consequently, we have 28 different Φ_i,ν in
Eq. (<ref>).     
In the group theory, they are classified into eight irreducible
representations of the O_h group, as shown in Table. <ref>. 


Here, we compare the possible gap functions in Table <ref> with those in the superconductor in the single J=3/2 band.
In single J=3/2 systems such as half-Heuslers, 
only even-parity superconductivity corresponding to A_1g, E_g, and T_2g in the Table <ref> is  
possible within the k independent gap function <cit.>,
because they do not have degrees of freedom of σ_ν associated with orbital parity.
In addition, we have odd-parity superconductivity in the Table <ref>.
These odd-parity pairing states are expected to be topologically non-trivial <cit.>.
In the present system, these superconductivities 
are realized with a constant gap function
in contrast to the k-dependent ones in single J=3/2 states <cit.>.














 §.§ Critical temperatures

Now, we evaluate the superconducting transition temperature T_c for each representation in Table <ref>.  
Here, we assume the following density-density pairing interactions, 

    H_int = -∫ d^3 r	[ U(n_+^2+n_-^2) + 2Vn_+n_- ],

where n_σ_z(r)=∑_j_zψ^†_j_z,σ_z(r)ψ_j_z,σ_z(r) is the density of the electron with orbital σ_z
defined by the field operator ψ_j_z,σ_z(r)=∫ d^3 k e^-ik·rc_j_z,σ_z,k.


The transition temperature T_c is obtained by solving the linearized gap equation.
For each representation in Table <ref>, 
we have <cit.> (see also Appendix <ref>),

    detX(T_c) = 0,

where the matrix element of X(T_c) is 

    X_αα'(T_c)       ∼δ_αα' +  
    	∑_nI(T_c)/2∫_k_F d^2k' D_n(k') V_α
          ×∑_n'=n,n̅⟨ u_n',k'|Φ_α|u_n,k'⟩⟨ u_n,k'|Φ_α' |u_n',k'⟩.

Here, |u_n,k⟩ is the normalized single-particle state for the Fermi surface of n-th band 
defined by H_0(k)|u_n k⟩=0. 
(Hereafter, all the ket vectors |u_n,k⟩ indicate the wave function of the normal state unless otherwise noted.)
The band index n̅ represents the state 
|u_n̅k⟩ = P𝒯 |u_n, k⟩, 
which is degenerated with |u_n k⟩.
D_n(k')=dk_F/dξ_n is the density of states on the Fermi surface. 
The interaction is V_α=U (V) for intraorbital (interorbital) pairing. 
We also use I(T_c)=∫_-ħω_0^ħω_0 dξ1/2ξtanh(ξ/2k_BT_c)∼ln (2e^γħω_0/π k_BT_c), Euler's constant γ,
cutoff frequency ω_0, and Boltzmann constant k_B. 
Here, α and α' specify a different basis φ_i( J)⊗σ_ν/N≡Φ_α' (if existing) in the representation.  
For instance, the A_1g gap has two different bases, A_1^(0)⊗σ_0 and A_1^(0)⊗σ_z, and hence, α' runs α'=1,2.

We show the numerical solutions of Eq. (<ref>) in Fig. <ref>. 
The results for even-parity states and those for odd ones are shown
separately as 
T_c's of even- (odd-) parity states are independent of the
interorbital (intraorbital) pairing interaction V (U).
It is observed that the spin-singlet A_1g (A_1u) representation
exhibits the highest T_c among even- (odd-) parity states.




These results can be understood using the spin and orbital textures of
the Fermi surfaces in the normal state.
In Fig. <ref>, we show the expectation values of ⟨ J⟩ and ⟨σ⟩ on the Fermi
surfaces with
respect to the single-particle state |u_n k⟩. 
The hole-doped system with μ<-m_0
supports 
two Fermi surfaces around the Γ point, and 
each Fermi surface has twofold degeneracy with different spin and
orbital textures, owing to the coexistence of inversion and time-reversal symmetries (see also Appendix <ref>).  
The orbital texture ⟨σ⟩ indicates mixing between orbitals with even and odd parity.
When ⟨σ⟩∥ẑ, the
corresponding single-particle state satisfies σ_z=± 1, and hence, it
should be an orbital with either even or odd parity  and there is no orbital
mixing.
In contrast, when ⟨σ⟩∥x̂, 
orbitals with even and odd parity are fully mixed in equal weight.





First, consider the A_1g state.
As it is a spin-singlet with the angular momentum J=0, the Cooper pair has an anti-parallel
spin configuration. Thus, it is formed by an electron and its
time-reversal partner, which have opposite spins.
In this case, they can form an intraorbital pairing without obstruction. 	
Indeed, as shown in Fig. <ref>, as time-reversal does not
change the orbital, 
the electron and the
partner have the same ⟨σ⟩. Consequently, 
the A_1g state is naturally realized in the presence of
the intraorbital pairing interaction. 
It corresponds to the usual s-wave superconductivity, where the bulk is topologically trivial [see also Sec. <ref>]. 


Similarly, the A_1u state is naturally realized in the
presence of the interorbital pairing interaction.
The A_1u state is also spin-singlet with J=0, and hence, its Cooper pair is also
formed between an electron and its time-reversal partner. However,
for the interorbital pairing state, orbital mixing is
necessary, which results in an interesting nodal structure of 
the gap depending on the model parameters.
For instance, when v_2= 3v_1, the outer Fermi surface 
does
not show orbital mixing at the intersections with the primary
axis (i.e., the k_i=x,y,z axis).
See Fig. <ref>.
Correspondingly, there are gap nodes at the intersections in
Fig. <ref> (b).
We will show later that this nodal structure is not accidental 
but has a topological origin.
Note that, except for the special point 
v_2= 3v_1, v_1/2, and -v_1/3, the A_1u state is full-gap as shown in Fig. <ref> (c).
In Fig. <ref>, we show the phase diagram obtained through our calculation.
It is observed that the A_1u phase is realized when V is much larger than U.








§ TOPOLOGICAL SUPERCONDUCTIVITY IN THE A_1U STATE

In this section, 
we study the topological properties of the A_1u (J=0) superconducting state,
in terms of the topological indices and numerically obtained surface states. 




 §.§ 3D winding number
 
 
We first examine the topological properties independent of the crystal symmetry.
For the time-reversal symmetric superconductivity generally described by the BdG Hamiltonian (<ref>),
we have chiral symmetry, 

    Γ H(k) Γ=-H(k)  with Γ=τ_y,

which corresponds to the combination of time-reversal and particle-hole symmetries [see Appendix <ref> for more details].
As these symmetries do not depend on any particular crystal structure, 
the chiral symmetry specifies the most general symmetry-protected topological number <cit.>.
In the present case, the general topological number is the 3D winding number,

    w_3D = ∫d^3k/48π^2ϵ^αβγtr[Γ
    H^-1∂_α H 
    H^-1∂_β H 
    H^-1∂_γ H 
    ].

Here, we calculate w_3D for the A_1u state with Δ=Δ_0 σ_x. 
In particular, we consider the weak pairing limit, Δ_0/μ≪ 1, similar to ordinary superconductors. 

First, we numerically evaluate w_ 3D.
It is observed that, when the absolute value of the chemical potential |μ| is less than the critical value μ_c,  

    μ_c=√(m_0^2+Δ_0^2)

w_ 3D is zero i.e.,

    w_3D(|μ|<μ_c)= 0.

At |μ|=μ_c, the BdG Hamiltonian in Eq. (<ref>) becomes gapless at k=0; thus 
topological phase transition takes place at μ=μ_c, then 
for |μ|>μ_c, we have

    w_3D(|μ|>μ_c)= 
    	{[ -2 sign(v_1)             ;  4 sign(v_1)             ; -4 sign(v_1)             ;  2 sign(v_1)              ]. .

For |μ|>μ_c, the system supports two Fermi surfaces (inner and outer Fermi surfaces) around the Γ point, as illustrated in Fig. <ref>.  
When w_ 3D changes at v_2/v_1=-1/3, 3 (v_2/v_1=1/2), six (eight) point nodes appear on the outer Fermi surface, 
each of which contributes to the change of w_ 3D by 1 (-1) with the increase in v_2. 
These results indicate that the A_1u state hosts nontrivial topological superconductivity when |μ|>μ_c. 
In particular, it supports a characteristic higher winding number |w_ 3D|≥ 1.
As discussed immediately below, the higher winding number of the topological superconductivity is a direct consequence of the higher spin of the present system. 

Here, we consider the origin of the higher winding number.
Usually, the pairing interaction is too weak to mix the energetically separated bands. 
In this weak pairing case, the topological number of the superconductor 
is attributed to that of each band forming Fermi surfaces,
where the electron and hole states are degenerated.
The contribution from each band remains unchanged unless 
the gap of the system closes or the Fermi surfaces contact each 
other.[In the strong pairing case, the large gap function may mix the contributions from each band, 
but their summation is conserved.]
In the present case, we have <cit.>

    w_ 3D=1/2∑_n sign[v^n_ F] sign[d_n] ν_ Ch^n,

where n is summed for all band indices forming disconnected Fermi surfaces, 
v_ F^n is the Fermi velocity on the n-th Fermi surface, and
d_n is the expectation value of the gap function 

    d_n=⟨ u_n, k|Δ|u_n, k⟩,

the ket vector |u_n, k⟩ is the solution of Bloch equation H_0(k)|u_n,k⟩=0 on the n-th Fermi surface,
and ν_ Ch is the first Chern number of the Fermi surface,

    ν^n_ Ch= -i∫_k=k_F^n d^2k' ϵ^αβ⟨∂_k_α' u_n,k| ∂_k_β' u_n,k⟩

where k_α' represents the two-dimensional momenta on the Fermi surface. 
Here, we consider Kramers partners separately in the summation of n.
See Appendix <ref> for more details.

The relation between the J=3/2 spin and the higher winding number becomes evident 
in the spherical symmetric case at v_2=0.
In this case, the solution of the Bloch equation H_0(k)|u_j_z,±,k⟩ = ξ_j_z,±|u_j_z,±,k⟩ for the normal state is
given by

    |u_j_z,±,k⟩ = R(k̂,J)|j_z> ⊗|Ψ_±(ρ_k)>_σ.

Here, R(k̂, j)=e^-iJ_zϕ_ke^-iJ_yθ_k 
with the polar and azimuthal angles (θ_k, ϕ_k) of k is the rotation matrix that diagonalizes the spin-dependent part of H_0( k) as 
R^†(k̂,J) k· JR(k̂,J)=k J_z , 
and |j_z>∈{|± 3/2>, |± 1/2>} is the eigenstate of J_z. 
The orbital-dependent part of the wave functions is given by

    |Ψ_±(ρ_k)>_σ = e^-iσ_y ρ_k|±>,

with ρ_k=1/2arctan[v_1kj_z /m(k)] and the eigenstate |±> of σ_z.
It is observed that the states |u_j_z,+,k⟩ (|u_j_z,-,k⟩) 
yield the electron (hole) branches of the normal spectra, and hence, they define the Fermi surfaces 
for μ> μ_c (μ< -μ_c).
Therefore, Eq. (<ref>) leads to   

    w_ 3D=∑_j_z w_ 3D^j_z,

with

    w_3D^j_z = 
    ±1/2 sign(d_j_z,± ) 
    ν^j_z_Ch,

where the double sign ± corresponds to the case with μ≷±μ_c, d_j_z,± = ⟨ u_j_z,±,k|Δ_0σ_x |u_j_z,±,k⟩, 
and
ν^j_z_Ch is the first Chern number of |u^±_j_z,  k⟩ on each Fermi surface. 
Equation (<ref>) yields 

    sign(d_j_z,±,k) = ± sign(j_z v_1).

Using the polar coordinates of the Fermi surface, we can also evaluate ν_ Ch^j_z analytically,

    ν^j_z_Ch   =   1/2π∫ dθ_k dϕ_k 
    	(-i ⟨∂ u^±_j_z,k/∂θ_k|∂ u^±_j_z,k/∂ϕ_k⟩+ h.c.)
    	
       =   2 j_z.

Thus, we have

    w_ 3D^j_z= sign(v_1)|j_z|,

both for μ>μ_c and μ<μ_c.  


The last expression indicates that a higher spin provides a higher winding number.
Using this expression, we can also evaluate the winding number for each Fermi surface.
For v_2=0, the outer (inner) Fermi surface consists of j_z=± 1/2 (j_z=± 3/2) components, and hence, we obtain

    w_out = w_3D^1/2+w_3D^-1/2=  sign(v_1), 
    
           w_in = w_3D^3/2+w_3D^-3/2= 3  sign(v_1),

where w_ out (w_ in) denotes the winding number of the outer (inner) Fermi surface.
These winding numbers retain the same values for -1/3<v_2/v_1<1/2 as either a gap-closing or a contact of Fermi surfaces does not occur until v_2/v_1 reaches the boundary of the region at v_2/v_1=1/2 or -1/3.
Here, we note that the total winding number w_ out+w_ in reproduces the numerical result of Eq. (<ref>) in the same region of -1/3<v_2/v_1<1/2.


Using the above result, we can also identify the winding numbers (w_ out, w_ in) for the other regions of v_2/v_1.
First, to evaluate them in the region 1/2<v_2/v_1<3, we use a “duality” relation. 
As explained in Appendix <ref>, 
a unitary transformation maps the parameter (v_1, v_2) to (-3v_1/5-4v_2/5, -4v_1/5+3v_2/5).
This transformation exchanges the region of -1/3<v_2/v_1<1/2 with that of 1/2<v_2/v_1<3, retaining the energy spectra and reversing the winding numbers. 
Thus, combining the duality relation with Eq. (<ref>), we obtain
 
    w_ out=- sign(v_1),
      
    w_ in=-3 sign(v_1),

for 1/2<v_2/v_1<3.
Subsequently, to determine the winding numbers for the remaining two regions v_2/v_1<-1/3 and v_2/v_1>3, we use the properties of the topological phase transitions at v_2/v_1=-1/3 and v_2/v_1=3.  
As mentioned above, 
at the topological phase transitions, six point nodes appear on the outer Fermi surface, 
which change |w_ out| by 6 [We observe that |w_ out| changes only by 2 at v_2/v_1=1/2 although there appear eight point nodes on the outer Fermi surface. This is because the inner and outer Fermi surfaces contact each other at v_2/v_1=1/2, and hence, they exchange their winding numbers.].
Provided  also that the total winding number w_ out+w_ in reproduces Eq. (<ref>), we can uniquely determine the winding numbers as

    w_ out=-5 sign(v_1),
      
    w_ in=3 sign(v_1),
 
for v_2/v_1<-1/3, and 

    w_ out=5 sign(v_1),
      
    w_ in=-3 sign(v_1),
 
for v_2/v_1>3. 
In both cases, the Fermi surfaces have characteristic higher winding numbers.







 §.§ Mirror Chern numbers

The present system has the O_h group crystalline symmetry, which provides more
detailed topological information <cit.>.
Specifically, we first consider
the horizontal (or equivalently vertical) mirror reflection M_x with
respect to the k_x-axis.




The A_1u gap function is odd under the mirror reflection.
In this case, the mirror reflection operator for the BdG Hamiltonian is
given by Mτ_z, where M=PC_2,n is the mirror operator for the normal
Hamiltonian H_0( k), and n is the normal vector of the reflection plane <cit.>.
On the mirror-invariant plane k·n=0, the BdG Hamiltonian 
can be block diagonal in the diagonal basis of Mτ_z. Thus, the mirror
Chern number ν̃_σ_h,d for the BdG Hamiltonian is defined 
in the same manner as ν_σ_h,d in Sec.<ref>.

First, let us evaluate ν̃_σ_h in the limit μ=0 and Δ_0 = 0 
for the mirror reflection with respect to three equivalent horizontal planes: x=0, y=0, or z=0.  
In this limit, the BdG Hamiltonian reduces to H_0( k)τ_z, and hence,
ν̃_σ_h can be evaluated as
ν̃_σ_h=2ν_σ_h.  
Therefore, from Eq. (<ref>), we have

    ν̃_σ_h(μ=0)= 4 sign(b_3/2b_1/2).

For small Δ_0 and μ, ν̃_σ_h retains the same
value unless the gap of the system closes. 


When |μ|=μ_c, the gap of the system closes at the Γ point, and
ν̃_σ_h changes.
Using the fourfold rotation symmetry around the k_x axis, 
we obtain the following useful formula to calculate
ν̃_σ_h <cit.>.

    ν̃_σ_h= ∑_j_xN_j_x,+(Γ) j_x  4,

where N_j_x,+(Γ) is the number of negative energy states with
the spin j_x at the Γ point in the eigensector of M_xτ_z 
with the eigenvalue λ_M_xτ_z=± i.
In the weak coupling limit Δ_0→ 0, 
the particle and hole sectors are decoupled, so that λ_M_xτ_z=λ_M_xλ_τ_z. 
In this case, a negative energy
state in the λ_M_xτ_z=i sector reduces to either an electron state (λ_τ_z=+1) in the
λ_M_x=i sector or a hole state (λ_τ_z=-1) in the λ_M_x=-i
sector below the Fermi level.
Therefore, N_j_x,+(Γ) can be evaluated as the number of
these electron and hole states. 
Thus, as summarized in Fig. <ref>, 
when |μ| exceeds m_0, the j_x=3/2, -1/2 electron (j_x=1/2,
-3/2 hole) bands in the λ_M_x=i (λ_M_x=-i) sector move above (below) the Fermi level at Γ. 
Therefore, the mirror Chern number jumps by 

    Δν̃_σ_h   =    -(3/2-1/2) + (1/2-3/2)  4 
       =    -2  4

We can also evaluate ν̃_σ_h numerically as

    ν̃_σ_h(|μ|>μ_c)
    = 2 sign(b_3/2b_1/2),

which is consistent with Eq. (<ref>).

According to the bulk-edge correspondence, the non-zero mirror Chern
number indicates the existence of surface states. 
Here, we calculate the energy spectrum in the slab geometry
in Fig. <ref> with a finite size L along the z axis. 
This system has horizontal (diagonal) mirror reflection symmetry with respect to
the x (x+y) direction even in the presence of the surface.
Along the high symmetric line k_x=0 (k_x+k_y=0), 
we numerically solve the BdG equation for each mirror sector by replacing 
the momentum k_z→ -i∂_z, and using the Gauss–Lobatto expansion method <cit.>.








As shown in left column of Fig. <ref>, when |μ|<μ_c, the
spectrum along the k_x=0 line exhibits four branches of surface states. 
The number of branches coincides with |ν̃_σ_h|.
Furthermore, the surface states change the chirality (or the
connectivity with the bulk states) when
ν̃_σ_h in the bulk changes the sign. 
For -1/3<v_1/v_2<3 [see Fig. <ref>(b) and (c)], each surface state connects the bulk bands upward
with the increase in k_y,
but for v_1/v_2<-1/3 or for v_1/v_2>3 [see Fig. <ref>(a) or (d)],
it connects them downward.
A similar bulk-edge correspondence holds for the higher-doped phase |μ|>μ_c as shown in Fig. <ref>.
In this case, when the surface state cuts the zero energy upward (downward) with the increase in k_y,
it accumulates the mirror Chern number by 1 (-1).
The total mirror Chern number counted from Fig. <ref>
is consistent with Eq. (<ref>).



Similarly, the diagonal mirror Chern number ν̃_σ_d
with respect to the k_x+k_y direction can be defined 
for the BdG Hamiltonian.
For μ=0, it is evaluated as twice the value of that in
Fig. <ref>(f), which persists as long as μ reaches
the critical value μ_c,

    ν̃_σ_d(|μ|<μ_c)=
    	{[ +4  for  -1/3<v_2/v_1<1/2,;     -4  for 1/2<v_2/v_1<3,;             0  otherwise .;                            ].

In addition, for |μ|>μ_c, we numerically determine that 

    ν̃_σ_d(|μ|>μ_c)=
    	{[ +2  for  -1/3<v_2/v_1<1/2,;     -2  for 1/2<v_2/v_1<3,;             0  otherwise .;                            ].

The quasi-particle spectra along the mirror-invariant line k_x+k_y=0 for the slab geometry
are shown in right columns of Figs. <ref> and <ref>, 
all of which are consistent with the above values of the diagonal mirror Chern number.
In particular, for μ=-0.75 m_0, 
the surface dispersions of Figs. <ref>(a) and (d) cut 
the zero energy upward as many times as downward, from right to left, 
and hence, there is no topological protection.

In Fig <ref>, all the surface states for
|μ|> μ_c pass through the high symmetry point k_y=0
(k_11̅0=0) at E=0.
The origin of this property will be explained in the next section.













 §.§ 1D winding numbers

The system also has fourfold rotation symmetry, from which we can
define another topological invariant. 
The A_1u gap function is invariant under the fourfold
rotation along the primary axis, and hence, the BdG Hamiltonian trivially realizes
the fourfold symmetry as 

    C̃_4, z
    H(k_y, -k_x, k_z)C̃_4, z^-1=H( k),
      C̃_4, z=C_4, zτ_0.

On the primary axis k_x=k_y=0, the BdG Hamiltonian is block
diagonal in the eigenbasis of C̃_4, z with the eigenvalue λ_C̃_4,z=e^-ij_z/2π,

    H_j_z(k_z) = H_0,j_z(k_z)τ_z + Δ_0σ_x τ_x,

where H_0,j_z(k_z) is given in Eq. (<ref>). We also observe that each sector H_j_z(k_z) retains chiral symmetry Γ=i T C in Eq. (<ref>), as it holds that
[Γ, C̃_4,z]=0. 
Thus, on the basis where Γ is diagonal, 
H_j_z(k_z) takes the form 

    H_j_z = 
    	([            0   h_j_z(k_z); h_j_z^†(k_z)            0 ]).

Thus, we can introduce the 1D winding number 

    w^j_z_ 1D = 1/2πIm{∫ dk_z  ∂_k_zln[ (h_j_z)]},

with

    (h_j_z)=μ^2-m_j_z^2(k_z)+(b_j_z k_z)^2-Δ^2_0+2i b_j_zk_zΔ_0,

where b_j_z is given by Eqs. (<ref>) and (<ref>).
It can be calculated analytically as 

    w_1D^j_z=
    	{[              0 for  |μ| < μ_c;    sign(b_j_z) for  |μ| > μ_c ]. .

It is also observed that the total winding number vanishes, ∑_j_z
w_1D^j_z=0, because of b_-j_z=- b_j_z, and 
the 1D winding number is related to the mirror Chern
number as 

    2w_ 1D^3/2w_ 1D^1/2 = ν̃_σ_h 4

When |μ|>μ_c, 
|w_ 1D^j_z|=1 for each j_z. This indicates
that there must be four zero modes at k_x=k_y=0 in each mirror
subsector, which explains why the surface states in one mirror sector in
Figs. <ref> pass the high symmetric points k_y=0 or k_11̅0=0 when
μ=-1.25m_0.    





 §.§ Topological phase diagram

We summarize the obtained topological phase diagram in Fig. <ref>.
On the bottom of Fig. <ref> with |μ|=0, 
 the topological indices and the corresponding surface states originate from those in the normal state. 
As already discussed in Sec. <ref>, octahedral or cubic Dirac points appear
at topological phase transition points in this regime.
In contrast, for |μ|>μ_c, as shown in Fig. <ref>(f) [(e) and (g)],
eight (six) point nodes of the superconducting gap appear on the outer Fermi surface at v_2=v_1/2 (3v_1 and -v_1/3).
These point nodes are also caused by the topological phase transitions where the bulk topological number w_3D changes by eight (six). 
These kinds of nodal structures have been overlooked in the group theoretical classification of gap nodes <cit.>,
but our topological analysis clarifies that the presence of the point nodes is not accidental
at the same level as the bulk Dirac points in the normal states.
The presence of point dips near the topological phase transitions provides a chance to confirm the A_1u state through low-temperature measurements of
the specific heat, nuclear magnetic resonance spectrum, etc. 

In Fig.<ref>, we also illustrate patterns of helical Majorana fermions composed of electrons of the outer and inner Fermi surfaces.
The plus and minus signs in Fig.<ref> denote the helicities of the Majorana fermions, where 
the total helicities for the outer (inner) Fermi surface should be the same as the winding number w_ out (w_ in).
As explained below, these configurations are determined by comparing the surface state spectra in Fig.<ref> with the winding numbers and considering the requirement from symmetry.  

Let us first consider the region -1/3<v_2/v_1<1/2, which includes the spherical symmetric point at v_2=0. 
See Fig. <ref>(b) and Fig.<ref>(b). 
In this region, the outer Fermi surface with w_out=1 
provides a single helical Majorana fermion at the center of the surface Brillouin zone, 
which we denote as the blue symbols and cones in Figs. <ref>(b) and <ref>(b), respectively, at k=0 (k=k_x or k_11̅1).  
The inner Fermi surface has an unusual higher winding number 
w_in = 3.
In the spherical symmetric limit v_2=0, we have a single surface state at the center of the surface Brillouin zone, but the 
surface spectrum shows a higher-order k dependence near k=0.  
Moreover, when v_2≠ 0, the surface state breaks into several helical Majorana fermions because the spherical symmetry reduces to the point group O_h. 
For consistency with C_4v symmetry on the surface,
there should be at least a single helical Majorana  fermion at the center of the surface Brillouin zone and four others at other points, as illustrated in Fig.<ref> (b). 
In terms of the winding number, this splitting corresponds to w_in = - 1+4. 
From the actual surface spectrum in Fig.<ref>(b), 
it is observed that the four Majorana fermions are located on the diagonal mirror lines in the surface Brillouin zone.
We also observe that the two center helical Majorana fermions on the inner and outer Fermi surfaces do not mix, although they are located at the same position with opposite helicities.
This is because they also have different topological numbers, w_ 1D^1/2=w_ 1D^3/2=1.

To generate the center and satellite Majorana fermions simultaneously,
the surface state has non-monotonic dispersion. We can observe this behavior in Fig. <ref>(b) for example, 
where a surface state connecting the bulk bands from k_11̅0∼ -0.7k_0 to 0.7k_0 cuts 
the zero energy upward in k_11̅0<0, downward at k_11̅0=0, and upward again in k_11̅0>0. 
This non-monotonic behavior enhances the density of states around E=0,
providing a characteristic zero bias conductance peak in 
the tunnel conductance spectra <cit.>.

We can also determine the configurations of surface Majorana fermions in other regions.
Similar to the argument for the bulk winding numbers in Sec. <ref>, we use the duality relation to determine the configuration in the region 1/2<v_2/v_1<3. 
As the duality retains the spectrum but reverses the winding numbers, the configuration in Fig.<ref> (c) is obtained from Fig. <ref>(b).
Finally, the configurations in Figs.<ref>(a) and (d) are determined by considering the topological phase transitions at v_2/v_1=-1/3 and 3.
At these transition points, six point nodes appear on the outer Fermi surface.  
Upon being projected on the surface Brillouin zone, these point nodes change the helicity of the center Majorana fermions by two, and create four additional satellite Majorana fermions in the outer Fermi surface excitations, as shown in Figs.<ref> (a) and (d).
The obtained Majorana configurations in Figs.<ref> (a), (c), and (d) are consistent with the spectra in Figs.<ref> (a), (c) and (d), respectively.




§ DISCUSSION




 §.§ Application to Sr_3-xSnO

Notably, a class of antiperovskite materials A_3BX with 
A=(Ca, Sr, La), B=(Pb. Sn), and X=(C, N, O) realizes a higher-spin
system. Near the Fermi level, they host a d-orbital band of the A atom and
a p-orbital band of the B, which belong to the G_3/2g and
G_3/2u representations of the O_h group, respectively.  
Therefore, our theory is directly applicable to antiperovskite
materials. [In antiperovskites, the opposite parities of two J=3/2 states
originate from those of the atomic orbitals of A and B. 
However, in general, not only the atomic orbital degrees of freedom 
but also sublattice ones can form the two states with different parities discussed here <cit.>. 
In fact, in the spin 1/2 TIs of the Bi_2Se_3 family, 
orbital states with odd and even parity both originate from p-orbitals forming 
bonding and anti-bonding states between different sublattices.
The two J=3/2 states with opposite parity can also be obtained through this mechanism in general.]

The first-principle calculations for Sr_3SnO 
were carried out based on different packages, VASP <cit.>, WIEN2k, and AkaiKKR <cit.>. 
These results have scheme dependence in the detailed band structure, such as behavior near the Brillouin zone boundary. 
However, at least, the presence of band inversion at the Γ point is consistent in all the schemes.

The first-principle calculations show that the antiperovskites support
the octahedral Dirac points with a tiny gap in the j_z=3/2 sector <cit.>. 
Therefore, they are near the topological phase transition at
v_2=3v_1 <cit.>.
In addition, in the antiperovskite with valency (A^+2)_3(B^-4)(X^-2), 
p-orbitals of B^-4 ions are almost filled but d-orbitals of A^+2 are almost empty <cit.>.
Therefore, in the antiperovskites with stoichiometric ratio, the Dirac point is pinned near the Fermi level, 
corresponding to μ=0 in our model Hamiltonian in Eq. (<ref>).


Recently, the superconductivity of antiperovskites has been experimentally observed in Sr_3-xSnO with the deficiency x <cit.>.
This deficiency yields hole doping. The Fermi level is 
expected to be much lower than the fourfold degenerate point at Γ.
Therefore, this system can realize the high winding superconductivity in |μ|>μ_c discussed in the present work 
if the interorbital attractive interaction is dominant.

Finally, we discuss possible experimental signals of high winding topological superconductivity in the case of Sr_3-xSnO. 
As mentioned before, Sr_3-xSnO corresponds to our model with v_2/v_1∼ 3 and μ<-μ_c. 
This indicates that 
if the system realizes the J=0 odd-parity superconductivity, 
Sr_3-xSnO is in the vicinity of the topological phase transition 
between the w_out=-1 phase and the w_out=5 one 
in our topological phase diagram in Fig. <ref>. 
Therefore, as shown in Figs <ref> (a)-(c), there should be point dips on the outer Fermi surface. 
Note that there is no gap node on the inner Fermi surface, 
as the topological phase transition retains the value of w_in.
The point node structure can be detected via low-temperature 
experiments such as heat capacity measurement.
Notably, the difference in behaviors of the superconducting gap 
between the outer and inner Fermi surfaces might explain the double superconducting transition 
reported in Sr_3-xSnO <cit.>.




 §.§ Systems without band inversion

In the present work, for discussing the topological property of normal and A_1u superconducting states, 
we have assumed m_0α>0 where the system exhibits a band inversion of G_3/2g and the G_3/2u representation at the Γ point.
In contrast, when m_0α<0, without band inversion, the normal state is topologically trivial with the mirror Chern number 
ν_h=ν_d=0 in any parameter region of v_1 and v_2.

For the A_1u superconducting state in this case, 
the low-doped phase |μ|<μ_c is also topologically trivial as it can be derived in a manner similar to Sec <ref>.
However, in the case without band inversion, the topological phase transition 
by doping, which shrinks the Fermi surface to one point, occurs at |μ|=μ_c.
When the system changes from the low-doped phase |μ|<μ_c to the high-doped one |μ|>μ_c, 
the mirror Chern number jumps with the opposite sign to the case with m_0α>0 
as the angular momenta of valence and conduction bands at the Γ point (see Fig. <ref> for m_0α>0)
exhibit opposite signs for the cases with and without band inversion.
In addition, applying the analysis in Sec. <ref> and <ref>
we can derive the 3D and 1D winding numbers with the same value as m_0α>0. 
Consequently, all the topological indices for the high-doped phase of m_0α<0 
are the same as those in m_0α>0, as summarized in Fig. <ref>.





In contrast, when m_0α<0, without band inversion, 
the orbital mixing in both the Fermi surfaces is strongly suppressed.
This reduces the transition temperature of the interorbital pairing state A_1u. 
Thus, although the topological property of the high-doped phase is unchanged regardless of the sign of m_0α, 
the band inversion plays a key role in obtaining the new topological superconductivity 
discussed in this work. 



 §.§ Other pairing states

Depending on the pairing interaction, we cannot avoid the possibility of the s-wave superconductivity 
belonging to the A_1g representation in Table <ref>. 
As discussed in Sec. <ref>, it exhibits the highest transition temperature in even-parity superconductivity.
In this case, not only the 3D winding number, but also the symmetry-protected 1D winding numbers and mirror Chern numbers
in the superconducting state are trivial.


However, the s-wave superconducting state may host surface topological superconductivity. 
Owing to the non-zero mirror Chern numbers in the normal state, 
the system supports surface Dirac fermions in the normal state. 
In the s-wave superconducting state, the surface Dirac fermion 
exhibits surface topological superconductivity <cit.>.

The A_1u and A_1g states with J=0 mainly discussed so far 
are caused by the Cooper pair between time-reversal partners.
However, this is not the only possible partner.  
An electron may also 
form a Cooper pair with its inversion partner. 
Such a Cooper pair inevitably has a parallel spin configuration with
J=1, 3 as inversion does not flip the spin.
As shown in Table <ref>, 
the A_2u, T_1u, and T_2u states contain this class of Cooper pairs. 

As shown in Fig. <ref>, T_1u exhibits a relatively higher T_c. 
Notably, T_c of T_1u is comparable to that of
A_1u, which suggests that T_1u could win by improving the spin and orbital 
textures slightly. 
Indeed, a similar exchange of the pairing symmetry has been reported in the case of doped TIs with spin 1/2. 
The higher-order term with respect to k enhances the T_c of 
the spin-triplet E_u representation of D_3d symmetry more than that of the spin-singlet A_1u one <cit.>.
It supports the experimental observation of the anisotropic gap <cit.>.
Similarly, the present calculation does not exclude the possibility of T_1u.

Let us consider the gap structure and spin property of the T_1u state. 
It hosts J=1 and J=3 components, 
and hence, the gap function is given as 

    Δ = 2Δ_0/√(5)(χ_1 J_z + χ_3 J̃_z)σ_y

with the normalization condition of the mixing ratio χ_1^2 + χ_3^2 = 1.
Projecting Eq. (<ref>) onto the band basis of the normal state,
we obtain the energy gap at the Fermi surface of n-th band as [see Appendix <ref>]

    d'_n(k)       = -⟨ u_n, k| Δ | u_n̅, k⟩
          = χ_1 d^'(1)_n(k) + χ_3 d^'(3)_n(k),

where we use the fact that the Cooper pair of the T_1u state is formed by the electron and its inversion partner. 

Along the main axis k_z, the eigenvalue λ_C_4=e^-i j_z/2π of C_4 rotation is a good quantum number of the normal state. 
The band n and n̅ in Eq. (<ref>) belong to a different eigensector with λ_C_4 and -λ_C_4. 
In addition, as J_z and J̃_z preserve λ_C_4, 
there is a point node d'_n(k_z)=0 at the intersection of the k_z axis and the two Fermi surfaces. 
The presence of the point node is consistent with the group theoretical classification of the gap nodes <cit.>.
Along the other momentum direction, adjusting 
the ratio of these two components χ_1 and χ_2, the 
T_1u state can optimize the gap function so as to be consistent with the spin and
orbital textures on the Fermi surfaces, as much as possible.
We illustrate the optimized superconducting gap of the T_1u state for v_2=3 v_1 (v_2=0)
in Figs. <ref> (c) and (d) (Fig. <ref> (a) and (b)).

In Fig. <ref>(e), we show the parameter dependence of the optimum mixing ratio
of the J=1 and J=3 components at T_c in our numerical solution of the gap
equation. It is observed that the higher-spin J=3 component cannot be
neglected in the whole region of v_2/v_1, 
and it can even dominate the gap function in some regions.
Therefore, the T_1u state realizes a higher-spin pairing state.
The higher-spin nature of the pairing state could be tested via spin-
sensitive experiments.

In particular, at the rotational symmetric point v_2=0, 
the amplitude of the total momentum is a good quantum number.
Thus, the J=1 state and J=3 state do not mix with each other. 
Let us compare the stability of these two superconductivities in terms of the gap structure.
From the analytical form of wave function (<ref>), the
J=1 component is given as 

    d_j_z,±^'(1)(k̂) = ±2Δ_0/√(5)cos(2ρ_k) 
    	⟨ j_z| J·k̂ | -j_z ⟩,

for two Fermi surfaces with the angular momentum j_z=± 3/2 and ± 1/2 along k.
Here, as the spin operators J_x and J_y shift j_z only by ± 1 and J_z preserves j_z, 
the J=1 component cannot couple the j_z=±3/2 states, i.e., 

    d_±3/2,±^'(1)(k̂)=0

for all the directions of k. 
In contrast, the J=3 component opens the energy gap in this Fermi surface 
because J̃_x and J̃_y have finite matrix elements between j_z= ±3/2 states.
[see Appendix <ref>] 
Therefore, it is observed that the J=3 state is more stable than J=1 at this point.





Notably, the T_1u state realizes a nematic
superconducting state, which spontaneously breaks the O_h group symmetry. 
As illustrated in Fig. <ref>(a-d), T_1u supports an
anisotropic superconducting gap breaking the O_h group.
The nematic feature can also be used to identify the T_1u state
through thermal transport measurements <cit.>.



§ CONCLUSIONS

We have investigated topological superconductivity in
doped TIs of J=3/2 electrons. 
Antiperovskites realize this class of topological materials, whose normal state hosts Dirac points at the topological quantum phase
transition.
In contrast to half-Heuslers with a single J=3/2 band, the antiperovskites support
multiple J=3/2 bands, which allows us to realize
various unconventional superconductivities within the BCS
type of constant gap functions.  
Using the k· p Hamiltonian, 
we have demonstrated that the J=0 odd-parity pairing state
with the A_1u representation of O_h group exhibits the highest transition temperature T_c if the
interorbital pairing interaction is dominant.
The A_1u state shows a new class of topological superconductivity
with a higher winding number w_3D=±2 or ±4.
Moreover, it displays topological crystalline superconductivity
with respect to the mirror reflections and fourfold rotation of the
O_h group. 

We have also revealed that the nematic T_1u
state with J=3 Cooper pairs has a high T_c, which is 
comparable to that of A_1u. 
This result suggests an interesting possibility of spin-septet
nematic superconductivity in antiperovskite materials.  
Similar to the case of doped Bi_2Se_3, 
the T_1u phase might be stabilized against the A_1u phase
by considering the higher-order terms of the k· p Hamiltonian. 





This work was supported by JSPS KAKENHI Grant numbers JP16K17755, JP16H06861, JP17H02922, JP17J08855, and
the JSPS Core-to-Core program. 
S.K. is also supported by the Building of Consortia for the Development of Human Resources in Science and Technology.




§ SPIN MATRICES IN J=3/2 SPACE

Here, we provide the explicit forms of the matrices useful to describe J=3/2 systems.
First, the spin matrices in the J=3/2 space are given by 

    J_x=
    1/2([    0 √(3)    0    0; √(3)    0    2    0;    0    2    0 √(3);    0    0 √(3)    0 ]), 
    
    J_y=
    i/2([     0 -√(3)     0     0;  √(3)     0    -2     0;     0     2     0 -√(3);     0     0  √(3)     0 ]), 
    
    J_z=
    1/2([  3  0  0  0;  0  1  0  0;  0  0 -1  0;  0  0  0 -3 ]).


In general, the 4× 4 Hermitian matrix in the J=3/2 space is expanded by 
the basis given in Table <ref>. 
Among them , the T_1u representation including spin matrices itself {J_x, J_y, J_z} 
behaves as a vector under the spin rotation C_q n defined in Eq. (<ref>). 
Here, the other set of T_1u basis is given by the third-order polynomials of spin matrix, 
which has an explicit form as follows:

    J̃_x    ≡   5/3(J_yJ_xJ_y+J_zJ_xJ_z) - 7/6J_x 
       =   1/4([    0 √(3)    0   -5; √(3)    0   -3    0;    0   -3    0 √(3);   -5    0 √(3)    0 ]), 
    J̃_y    ≡   5/3(J_zJ_yJ_z+J_xJ_yJ_x) - 7/6J_y 
       =   i/4([     0 -√(3)     0    -5;  √(3)     0     3     0;     0    -3     0 -√(3);     5     0  √(3)     0 ]), 
    J̃_z    ≡   5/3(J_xJ_zJ_x+J_yJ_zJ_y) - 7/6J_z 
       =   1/2([ -1  0  0  0;  0  3  0  0;  0  0 -3  0;  0  0  0  1 ]).

By this definition, the two sets of T_1u basis satisfy the orthogonal condition 

    J·J = J̃·J̃ = 15/4,    J·J̃=0.




§ SUPERCONDUCTIVITY WITH TIME-REVERSAL AND INVERSION SYMMETRIES



 §.§ Nambu space and symmetry

In this section, we summarize the general properties of the BdG Hamiltonian
for the systems with time-reversal symmetry.
Here, we use the Nambu space spanned by the basis 
(c_k, c̅_-k), 
where the spinor c_k consists of the annihilation operator c_a, k
with the spin and orbital indices a= (σ_z, j_z), and 
c̅_-k is its time-reversal hole partner with the components
c̅_a,-k=∑_a'(C_2,y)_aa'c_a',-k^†.
In this basis, the BdG Hamiltonian is written as 

    H(k) = 
    	([  H_0(k)    Δ(k);  Δ^†(k) -H_0(k);         ]).
 
where we have used the time-reversal symmetry of the one-particle Hamiltonian in the normal state (<ref>). 

The point group operation for the particle space is given as 
c_a k→∑_a G_aa'c_a', k
with G=C_q,n, P, or PC_q,n 
corresponding to discrete rotation, inversion, and their combined operation 
acting on the internal degrees of freedom a=(σ_z,j_z).
 The transformation G for the time-reversal hole is 
c̅_a,-k→∑_a',b'(C_2,y)_aa' (𝒦 G_a'b'𝒦^-1) c^†_b',-k  
=∑_a',b' G_a, a'(C_2,y)_a'b'c^†_b',-k = ∑_a' G_aa'c̅_a',-k,
because the time-reversal operator T=C_2,y𝒦 commutes with any point group operation G.
The transformations for the particle and time-reversal hole components are similar. 
𝒯 and G act on the gap function as 

    𝒯Δ(-k)𝒯^-1 and GΔ(D_G^-1[k])G,

where D_G[k] indicates the G operation of k. When 𝒯Δ(-k)𝒯^-1=Δ(k), 
the BdG Hamiltonian is time reversal invariant,

    𝒯̃ H(k) 𝒯̃^-1 = H(-k)   with  𝒯̃ = 𝒯τ_0









The Fermi statistics of c_ k result in a restriction on the gap function  

    (Δ( k) C_2,y)^T =-Δ(- k) C_2,y.

To satisfy both (<ref>) and (<ref>), the gap function is Hermitian:

    Δ^†(k) = Δ(k).

Consequently, the BdG Hamiltonian (<ref>) for the time-reversal symmetric superconductivity is reduced to

    H(k) = H_0(k) τ_z + Δ(k) τ_x,

where τ_i represents the Pauli matrices acting on the particle and time-reversal hole space. 
 

We also consider the particle-hole symmetry of the BdG Hamiltonian 

    𝒞H(k)𝒞^-1=-H(-k)   with  𝒞=
    	([       0 C_2,y^†;   C_2,y       0;         ])𝒦,

inherent to superconductors. 
By combining (<ref>) and (<ref>), chiral symmetry is given as

    Γ H(k) Γ=-H(k)   with  Γ=-i𝒯̃𝒞=τ_y.




 §.§ Band representation

In this section, we consider the BdG Hamiltonian in terms of the band basis of the normal state.
This basis is given by the Bloch equation of the normal state

    H_0(k) |u_n,k⟩ = ξ_n(k) |u_n,k⟩.

with the band index n and eigenenergy ξ_n(k). 
From orthonormality and completeness of basis, 

    ⟨u_n, k | u_n', k⟩ = δ_n n',   ∑_n |u_n, k⟩⟨u_n, k| = 1̂,

the one-particle Hamiltonian is written as a diagonal form on this basis

    H_0(k) = ∑_nξ_n(k) |u_n,k⟩⟨u_n,k|.


For later convenience, we consider the time-reversal and inversion symmetry of the basis. 
From (<ref>) for the one-particle Hamiltonian of the normal state, 
the time-reversal partner of |u_n,k⟩ is 

    |u_n̅,-k⟩ = 𝒯 |u_n,k⟩,

with n̅≠ n because of 𝒯^2=-1. 
From (<ref>), 
the inversion partner is 

    |u_n,-k⟩ = P |u_n,k⟩.

In the presence of both symmetries, the energy of the state

    |u_n̅,k⟩ = P𝒯 |u_n,k⟩

is degenerated with that of |u_n,k⟩, i.e., ξ_n̅=ξ_n. 

Let us consider the band representation of the gap function. 
As both |u_n, k⟩ and |u_n̅, k⟩ have completeness (<ref>), 
we can expand the gap function as

    Δ(k) 
    	= -∑_n, n'd_n,n'(k) |u_n,k⟩⟨u_n̅',k|

with the matrix element
d_n,n'=-⟨u_n,k| Δ(k) |u_n̅',k⟩ = -⟨u_n,k| Δ𝒯|u_n',-k⟩. 
The d_n,n' provides the superconducting correlation
for Cooper pairs consisting of |u_n,k⟩ and |u_n',-k⟩. 
In the weak coupling limit where the gap function Δ is much smaller than the level 
spacing among the different bands, 
the matrix element d_n,n' is non-zero 
only for the degenerated bands n'=n or n'=n̅ [see Eq.(<ref>)]. 
This indicates that the time-reversal or inversion partners can form the Cooper pair.
These pairing states are described by the matrix elements 

    d_n,n'=
    	{[ d_n δ_n,n̅' for time-reversal pair;     d_n' δ_n,n' for inversion pair ]. .


In particular, when the time-reversal partners form the Cooper pair,
not only the normal Hamiltonian (<ref>)  
but also the gap function is diagonal in the band representation [see Eq. (<ref>)]. 
Therefore, the BdG Hamiltonian is also diagonal

    H(k) =∑_n(ξ_n(k)τ_z - d_n(k)τ_x ) |u_n,k⟩⟨ u_n,k|,

where d_n=⟨ u_n,k| Δ |u_n,k⟩ is real owing to hermiticity (<ref>). 
Note that the time-reversal and inversion symmetries for the gap function are

    d_n̅(-k) = d_n(k) 
    
    	d_n(-k) = ± d_n(k),

where the double sign indicates even- and odd-parity superconductivities. 



 §.§ 3D winding number for time-reversal pairs

The time-reversal and particle-hole symmetries defined by Eqs. (<ref>) and (<ref>) 
do not depend on any particular crystal structure. 
They thus specify the most general symmetry-protected topological number <cit.>.
In the present system with 𝒯^2=-1 and 𝒞^2=+1,
the general topological number is the 3D winding number

    w_3D=∫d^3k/48π^2ϵ^αβγTr[Γ
    	(H^-1∂_α H)
    	(H^-1∂_β H)
    	(H^-1∂_γ H)
    	].

with the 3D Levi-Civita symbol ϵ^αβγ. 
For the time-reversal pairing state, the diagonal form of the Hamiltonian (<ref>) provides a simple analytic form of w_3D. 
Here, we derive it by applying the procedure in Ref. <cit.> 
to the system with the additional inversion symmetry. 

First, we consider the unitary transformation (τ_x,τ_y,τ_z)→ (τ_y,τ_z,τ_x). 
It yields the eigenbasis of the chiral operator with Γ=τ_z, convenient to analyze w_3D.
In this basis, the BdG Hamiltonian (<ref>) is written as

    H(k)    =   ∑_n(ξ_n(k)τ_x - d_n(k)τ_y ) 
    	|u_n,k⟩⟨ u_n,k|
       =   ([      0   h(k); h^†(k)      0;        ])

where the off-diagonal element is 

    h(k)=∑_n(ξ_n(k) + id_n(k))|u_n,k⟩⟨ u_n,k|.





Second, we use the so-called spectral flattening technique <cit.>.
The eigenenergy of (<ref>) is given by the absolute values of the
off-diagonal elements E=±|ξ_n+id_n|.
Therefore, the replacement 

    h(k)→ q(k)   =   ∑_nξ_n(k)+id_n(k)/|ξ_n(k)+id_n(k)||u_n,k⟩⟨u_n,k| 
       =   ∑_n e^iϑ_n(k)|u_n,k⟩⟨u_n,k|,

provides the adiabatic transformation of the original Hamiltonian to that with flat bands with E=±1 without changing w_3D.
The phase factor ϑ(k) is defined as in Fig. <ref>.
In this technique, the Hamiltonian is replaced as

    H→ Q (k) = 
    	([       0    q(k); q^† (k)       0 ]),

and its inverse is replaced as 

    H^-1→ Q^-1 (k) =  Q(k),

where we have used q^-1(k)=q^†(k).
Through this replacement, we can evaluate the 3D winding number (<ref>) as

    w_3D   =   ∫d^3k/48π^2ϵ^αβγTr[Γ
    (Q ∂_α Q)
    (Q ∂_β Q)
    (Q ∂_γ Q)
    ] 
       =   ∫d^3k/24π^2ϵ^αβγTr[(q∂_α q^†) (q∂_β q^†) (q∂_γ q^†)]. 
       =   ∫d^3k/4π^2ϵ_αβγ∑_n,n'[  (i∂_γϑ_n) 
    	⟨∂_α u_n,k|∂_β u_n,k⟩δ_n,n'
          -  (i∂_γϑ_n) cosϑ_n'n⟨∂_αu_n,k| u_n',k⟩⟨ u_n',k|∂_β u_n,k⟩
          -  isinϑ_n'n⟨∂_α u_n',k| ∂_β u_n,k⟩⟨∂_γ u_n,k| u_n',k⟩]

with ϑ_n'n=ϑ_n'-ϑ_n. 
To obtain the last equality in Eq. (<ref>), 
we have used the time-reversal symmetry for the phase
ϑ_n̅(-k)=ϑ_n(k) given by (<ref>) 
and for the eigenstates 
⟨∂_αu_n̅,-k|∂_βu_n̅',-k⟩=
⟨∂_βu_n',k|∂_αu_ n,k⟩ and 
⟨∂_αu_n̅,-k|u_n̅',-k⟩=
⟨ u_n',k|∂_αu_ n,k⟩
by (<ref>) in addition to the usual algebraic expansion. 

Subsequently, we use the weak coupling limit with d_n→ 0, 
where the phase factor is given by the step function

    ϑ_n(k) =
    	{[                 0  for  k>k_F;                π for  k < k_F;  π/2  for  k = k_F and  d_n>0; -π/2  for  k = k_F and  d_n<0 ]. .

The contribution from the third term in Eq. (<ref>) is zero in this limit. 
In addition, away from the Fermi surface, the derivative of the phase factor ∂_γϑ_n=0. 

Let us consider the non-trivial contribution for the 3D winding number in (<ref>) from momenta close to the Fermi surface. 
Here, we consider the coordinate in momentum space (k_α,k_β,k_γ) 
as (k,k'_1,k'_2) with momenta k and k'_α=1,2 perpendicular and parallel to 
the Fermi surface, respectively.  The weak coupling limit (<ref>) yields 
∂_k'_αϑ = 0. In addition, the derivative with respect to k is 

    ∂_kϑ_n  = -sign[v_F^nΔ_n] πδ(k-k_F),

where v_F^n is the Fermi velocity of band n. 
Consequently, the 3D winding number (<ref>) is reduced to 

    w_3D   =   -isign[v_F^nΔ_n] ∑_n∫d^2k'/4πϵ_αβ[ 
    	⟨∂_k_α' u_n,k|∂_k_β' u_n,k⟩
          
    	-∑_n'cosϑ_n'n(k_F^n) 
    		⟨∂_k_α'u_n,k| u_n',k⟩⟨ u_n',k|∂_k_β' u_n,k⟩]

Equation (<ref>) and the time-reversal and inversion symmetries ϑ_n̅(k)=±ϑ_n(k) 
yield the phase factor at the Fermi surface as 

    cos (ϑ_n(k_F^n)-ϑ_n'(k_F^n)) = 
    	{[    1  for n'= n; ± 1  for  n'=n̅;    0  otherwise;                 ].

where the double sign indicates even- and odd-parity superconductivities.
When n'=n, the second term of (<ref>) is zero owing to ϵ_αβ⟨∂_k_α'u_n,k| u_n,k⟩⟨ u_n,k|∂_k_β' u_n,k⟩ =0. 
In addition, even for n'=n̅, the second term of (<ref>) is 

    f_n = ϵ_αβ⟨∂_k_α'u_n,k| u_n̅,k⟩⟨ u_n̅,k|∂_k_β' u_n,k⟩
          
    	=
    	-ϵ_αβ⟨∂_k_α'u_n̅,k| u_n,k⟩⟨ u_n,k|∂_k_β' u_n̅,k⟩
    	= -f_n̅

where we have used ϵ_αβ=-ϵ_βα and ⟨∂_αu_n,k|u_n',k⟩ = -⟨ u_n,k|∂_αu_n',k⟩. 
Therefore, the contribution from f_n̅ cancels that from f_n when performing summation over n and n̅. 
Consequently, only the first term in (<ref>) yields the non-trivial contribution to w_3D. 

Consequently, the 3D winding number is given by

    w_3D = 1/2∑_nsign[v_Fd_n] ν_Ch^j_z

with the first Chern number

    ν_Ch^j_z = -i/2π∫_k_F d^2k' ϵ^αβ⟨∂_k'_αu_n,k| ∂_k'_βu_n,k⟩.




 §.§ Linearized gap equation

In this section, we derive the linearized gap equation for the system with time-reversal and inversion symmetries.
The mean field Hamiltonian of the superconducting system is given as 

    Ĥ(k) = 
    	([   c_k^† c̅_-k^† ])
    	H(k)
    	([   c_k; c̅_-k ]),

where the spinor c_k consists of the annihilation operator c_j_z, σ_z, k 
with the indices j_z and σ_z for spin and orbital, respectively, 
its time-reversal hole partner c̅_-k is given as c̅_j_z,σ_z,-k=∑_j_z'(C_2,ŷ)_j_zj_z'c_j_z',σ_z,-k^†, 
and H(k) is the BdG Hamiltonian given in Eq. (<ref>).
Here, we define the thermal and anomalous Green's functions

    𝒢_k(τ,τ') = -<T_τ[c_k(τ)c_k^†(τ')]>, 
    ℱ_k(τ,τ') = -<T_τ[c̅_-k(τ)c_k^†(τ')]>,

where c_ k(τ)=e^H_0τ/ħc_ke^-H_0τ/ħ is the Heisenberg operator for imaginary time τ. 
From the Heisenberg equation of c_ k(τ), we obtain the Gor'kov equation of the Green's function <cit.>, 

    [  [iħω_N- H_0(k)] 𝒢_k(ω_N) -Δℱ_k(ω_N) =ħ,; [iħω_N+ H_0(k)] ℱ_k(ω_N) -Δ𝒢_k(ω_N) = 0, ]

where we have used the time-reversal symmetry of the normal state (<ref>) and of the superconducting gap function (<ref>).
The Fourier component of thermal and anomalous Green's functions is given by 

    𝒢_k(τ,τ') = 1/βħ∑_Ne^iω_N(τ-τ')𝒢_k(ω_N),

with fermionic Matsubara frequency ω_N = (2N+1)π /ħβ and β= 1/k_BT.
In terms of the Matsubara anomalous Green's function, the gap equation is written as

    Δ = -lim_η→ 01/ħβ∑_N∫ d^3k Ṽ e^iω_Nηℱ_k(ω_N),

where the interaction is Ṽ=U (V) when the gap function Δ is intended for intraorbital (interorbital) coupling [see also Eq. (<ref>)].

Near the T_c where Δ→ 0, we can linearize the gap equation. 
By solving (<ref>) for 𝒢_k and ℱ_k and 
substituting them successively, we approximate the anomalous Green's function in the lowest order of Δ as  

    ℱ_ k(ω_N) = -1/ħ𝒢^0_ k(-ω_N) Δ𝒢^0_ k(ω_N) + 𝒪(Δ^2),

with the single-particle Green's function 𝒢^0_ k(ω_N) = ħ[ iħω_N - H_0(k) ]^-1. 
It is written in terms of the band representation (<ref>) as 

    𝒢^0_ k(ω_N) = ∑_nħ |u_n,k⟩⟨u_n,k|/iħω_N - ξ_n(k).

By using Eq. (<ref>) and (<ref>), we can linearize the gap equation (<ref>) as 

    Δ   =    -lim_η→ 01/β∑_N,n,n'∫ d^3k Ṽ e^iω_Nη
          ×|u_n,k⟩⟨ u_n,k|Δ |u_n',k⟩⟨ u_n',k|/(-iħω_N-ξ_n)(iħω_N-ξ_n'),

In addition, the matrix element d_n,n'=⟨ u_n,k|Δ |u_n',k⟩ is non-zero for n'=n or n̅ [see also Eq. (<ref>)]
as the time-reversal or inversion partners form the Cooper pair. Hence, we reduce Eq. (<ref>) as 

    Δ      = -lim_η→ 01/β∑_N,n∫ d^3k Ṽ e^iω_Nη
          ×∑_n'=n,n̅|u_n,k⟩⟨ u_n,k|Δ |u_n',k⟩⟨ u_n',k|/ħ^2 ω_N^2-ξ_n^2, 
          = - ∑_n∫ dξ_n d^2k' D_n(k) Ṽ1/4ξ_ntanhβξ_n/2
          ×∑_n'=n,n̅ |u_n,k⟩⟨ u_n,k|Δ |u_n',k⟩⟨ u_n',k|,

where we have used ξ_n̅=ξ_n, the expansion formula z^-1tanh(z/2)=4∑_n=0^∞[(2n+1)^2π^2+z^2]^-1, 
and the momentum coordinates (k,k_1',k_2') with momenta k and k'_1,2 perpendicular and parallel to the Fermi surface, respectively.
D_n(k)=d k/d ξ_n is the density of states. 

In the weak coupling limit where β_c^-1=k_BT_c is much smaller than any other energy scale, 
the factor ξ_n^-1tanh(β_cξ_n/2) has a sharp peak around the Fermi level ξ_n=0.
Hence, we can approximate that the other factors are represented by their values at the Fermi surface ξ_n=0. 
Therefore, Eq. (<ref>) is rewritten as  

    Δ      = - ∑_nI(β)/2∫_k_F d^2k' D_n(k') Ṽ
          ×∑_n'=n,n̅ |u_n,k'⟩⟨ u_n,k'|Δ |u_n',k'⟩⟨ u_n',k'|,

where the energy integral is evaluated as 
I(β)=∫_-ħω_0^ħω_0 dξ1/2ξtanh(βξ/2)= ln (2e^γβħω_0/π) with Euler's constant γ 
and the cutoff frequency ω_0 <cit.>.

In the system discussed in the main text, 
the gap function is expanded as Δ=∑_αΔ_αΦ_α [see Eq. (<ref>)].
By using the orthonormality of the basis tr[Φ_αΦ_α']=δ_αα' in Eq. (<ref>), 
we obtain the gap equation for Δ_α as 

    ∑_β X_αα'(β) Δ_α' = 0,

with the coefficient

    X_αα'(β)       = δ_αα' +  
    	∑_nI(β)/2∫_k_F d^2k' D_n(k') Ṽ
          ×∑_n'=n,n̅⟨ u_n',k'|Φ_α|u_n,k'⟩⟨ u_n,k'|Φ_α' |u_n',k'⟩.

The transition temperature is given by det X(β_c)=0. 









§ DUALITY RELATION FOR K· P HAMILTONIAN

It can be easily verified that the identical vectors J and J̃ in the spin space hold that

    D_J J_i D_J^† = -3/5 J_i - 4/5J̃_i, 
     
    	D_JJ̃_i D_J^† = -4/5 J_i + 3/5J̃_i,

with the unitary matrix acting on the spin J=3/2 space

    D_J = ([  0  0  1  0;  0  0  0 -1; -1  0  0  0;  0  1  0  0 ]).

As the A_1u gap function is invariant under D_J, 
this unitary transformation maps the system with a set of parameters v_1 and v_2 to that with another set 

    ([ v_1; v_2 ])
    	([ v_1'; v_2' ]) = 
    	([ - [3v_1 + 4 v_2]/5; - [4v_1 - 3 v_2]/5 ]),

which yields the “duality” of the BdG Hamiltonian.
In addition, the operator D_J commutes with the chiral operator 

    [ Γ, D_J ] = 0.

Therefore, the systems at the two sets of parameters in (<ref>) yield 
the same winding number w_3D.

Similarly, the inversion operator for the superconducting state
maps the parameters as 

    ([ v_1; v_2 ]) 
    	([ v_1'; v_2' ])= 
    	([ - v_1; - v_2 ]) 
    	 with P̃  = σ_z τ_z.

Note that this mapping does not change v_2/v_1.
As the operator P̃ anti-commutes with the chiral operator Γ, 

    {Γ, P̃} = 0,

the systems at the two sets of parameters in (<ref>) yield 
the winding number w_3D with the opposite sign.

Finally, in Fig. <ref>, 
we summarize these relations in the parameter space spanned by v_1 and v_2/v_1 used in the main text.
For v_2/v_1>-3/4, the unitary transformation D_J maps v_2/v_1 to the opposite side of v_2/v_1=1/2 
denoted by the dashed line in Fig. <ref>. Simultaneously, it changes the sign of v_1. 
Furthermore, we can map a system with v_1 to another with -v_1 using P̃. 
Therefore, at the opposite side with respect to v_2/v_1=1/2 and with the same sign as v_1 
the winding number w_3D exhibits the opposite sign. 
In contrast, for v_2/v_1<-3/4, the unitary transformation D_J connects the 
parameters with the same sign as v_1 (see Fig. <ref>), 
which have the same winding number w_3D. 



