# Higher-order topological superconductors characterized by Fermi level
crossings

Hong Wang  School of Physics, MOE Key Laboratory for Non-equilibrium Synthesis
and Modulation of Condensed Matter, Xi’an Jiaotong University, Xi’an, Shaanxi
710049, China  Xiaoyu Zhu<EMAIL_ADDRESS>School of Physics, MOE Key
Laboratory for Non-equilibrium Synthesis and Modulation of Condensed Matter,
Xi’an Jiaotong University, Xi’an, Shaanxi 710049, China

###### Abstract

We demonstrate that level crossings at the Fermi energy serve as robust
indicators for higher-order topology in two-dimensional superconductors of
symmetry class D. These crossings occur when the boundary condition in one
direction is continuously varied from periodic to open, revealing the
topological distinction between opposite edges. The associated Majorana
numbers acquire nontrivial values whenever the system supports two Majorana
zero modes distributed at its corners. Owing to their immunity to
perturbations that break crystalline symmetries, Fermi level crossings are
able to characterize a wide range of higher-order topological superconductors.
By directly identifying the level-crossing points from the bulk Hamiltonian,
we establish the correspondence between gapped bulk and Majorana corner states
in higher-order phases. In the end, we illustrate this correspondence using
two toy models. Our findings suggest that Fermi level crossings offer a
possible avenue for characterizing higher-order topological superconductors in
a unifying framework.

##  I Introduction

Topological states of matter are usually endowed with a bulk-boundary
correspondence, which facilitates the identifications of topologically
protected gapless boundary modes without going into the details of the energy
spectrum at open boundaries [1, 2, 3]. Recent advancements in higher-order
topological systems have extended this correspondence to include gapped
boundaries [4, 5, 6, 7, 8, 9, 10, 11, 12, 13], with gapless corner (hinge)
modes appearing at the intersections of adjacent edges (surfaces). Tremendous
efforts have been devoted to classifying and characterizing these topological
states, mostly in crystalline-symmetry protected systems [14, 15, 16, 17, 18,
19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34]. However, it
is well known that gapless corner or hinge states persist when crystalline
symmetries are broken. This is especially evident in higher-order topological
superconductors [35, 36, 37, 38, 39, 40], where Majorana zero modes [41, 42,
43, 44, 45, 46] remain stable as mid gap states unless the bulk or boundary
gap closes. Hence it would be desirable to characterize higher-order states
regardless of whether crystalline symmetries are present.

Higher-order topology can be understood from a boundary perspective, as
different parts of the whole boundary, such as the four edges of a square
lattice, may exhibit a distinct topology in higher-order phases. For intrinsic
higher-order states, the relevant crystalline symmetry requires symmetry
related edges or surfaces to be topologically inequivalent [14, 15, 16]. A
topology change is only possible through bulk-gap closing. Consequently, bulk
invariants, such as symmetry indicators related to the crystalline symmetry,
can be defined [17, 18, 19, 20]. This stands in contrast with boundary-
obstructed topological states, which fall within an extrinsic higher-order
classification [11, 40]. Without the protections of crystalline symmetries,
the boundary topology in these states could change while the bulk gap remains
open. One may characterize the topology by Wilson loop eigenvalues of Wannier
bands that are obtained from Wilson loops of energy bands, the so-called
nested Wilson loop approach [4, 5]. However, the quantization of such
topological invariants still requires the presence of crystalline symmetries,
such as mirror symmetry [11, 38]. Establishing bulk-boundary correspondence
under broken crystalline symmetries remains an open question. Considering that
boundary topology is ultimately determined by the bulk properties for both
intrinsic and boundary-obstructed phases, it should be possible to associate a
topological invariant with it based on bulk information, which applies in both
phases.

In this paper, we focus on two-dimensional (2D) superconductors of symmetry
class D [3] and higher-order phases featuring two Majorana corner states. The
higher-order topology can be characterized by a pair of Majorana numbers,
which are intimately related to Fermi level crossings that emerge during the
continuous variation of the boundary condition along one direction, as
illustrated in Fig. 1(a). We further introduce a generic method for locating
these crossings from the bulk Hamiltonian. As a result, bulk-boundary
correspondence is established in both higher-order phases discussed earlier,
due to the robustness of the Fermi level crossings against crystalline-
symmetry-breaking perturbations.

Figure 1: (a) Higher-order topology characterized by Fermi level crossings in
the D symmetry class. Crossings occur at high-symmetry momenta $K=0,\pi$, when
the system on a torus ($\lambda_{a}=1$) is continuously deformed into a
cylinder ($\lambda_{a}=0$). (b) Topology of edges determined by Majorana
numbers. When $\mathcal{M}_{x(y)}=-1$, two opposite edges along the $y(x)$
direction exhibit a distinct topology, indicated by different colors. (c)
Schematic plots of the BdG spectrum at $K=0(\pi)$ with or without level
crossings. The fermion parity of the ground state (“$+$” for even, “$-$” for
odd) switches at each crossing. The crossing may disappear when the bulk or
boundary gap closes.

##  II General theory

To demonstrate how Fermi level crossings determine the higher-order topology
of D-class superconductors, we start from a 2D periodic lattice and modulate
its boundary condition in one direction. The resulting Bogoliubov-de Gennes
(BdG) Hamiltonian can be expressed as

$\tilde{\mathcal{H}}_{\lambda_{a}}=\sum_{k_{\bar{a}}}\mathcal{H}_{\lambda_{a}}(k_{\bar{a}})=\tilde{\mathcal{H}}_{1}-(1-\lambda_{a})\sum_{k_{\bar{a}}}\mathcal{B}_{a}(k_{\bar{a}}),$
(1)

where $\bar{a}=y(x)$ when $a=x(y)$, and the real parameter $\lambda_{a}$
controls the boundary condition in the $a$ direction, with
$\lambda_{a}=1,-1,0$ corresponding to the periodic (PBC), anti periodic
(APBC), and open boundary condition (OBC), respectively. In Eq. (1),
$\mathcal{H}_{\lambda_{a}}(k_{\bar{a}})$ represents the 1D boundary-modulated
Hamiltonian at wave vector $k_{\bar{a}}$, and $\mathcal{B}_{a}(k_{\bar{a}})$
involves all terms that cross its boundary. The lattice terminations we
consider are compatible with unit cells, thus allowing the specific form of
$\mathcal{B}$ to be directly read off from the bulk Hamiltonian
$\tilde{\mathcal{H}}_{1}$. The process of varying $\lambda_{a}$ from $1$ to
$0$ is akin to gradually cutting a torus along the $\bar{a}$ direction until
it eventually becomes a cylinder, as illustrated in Fig. 1(b) for the case of
$a=x$.

Here, we consider a gapped bulk with trivial first-order topology, which means
the cylindrical system described by $\tilde{\mathcal{H}}_{\lambda_{a}=0}$ is
fully gapped. Treating it as a quasi-1D system along the $\bar{a}$-direction,
we may characterize the higher-order topology with the Majorana number [47,
48, 49]

$\mathcal{M}_{a}=\text{sgn}\prod_{K}\text{Pf}[-iH_{\lambda_{a}=0}(K)],$ (2)

where “Pf” is shorthand for Pfaffian, $K=0,\pi$ represents the high-symmetry
momentum, and $H$ refers to the matrix representation of $\mathcal{H}$ in the
Majorana basis. In 1D, the Majorana number being $-1$ implies the presence of
a single Majorana zero mode at each end. If either $\mathcal{M}_{x}$ or
$\mathcal{M}_{y}$, or both of them, take the value of $-1$, we will instead
have two Majorana zero modes at the corners of a 2D sheet. To elaborate this
let us consider the cylindrical system in the lower left-hand panel of Fig.
1(b) with $\mathcal{M}_{x}=-1$. If we cut it along the axis, the resulting two
edges along the $x$ direction will each harbor one Majorana mode. Due to the
trivial first-order topology, these localized modes cannot propagate along the
edges and must be confined to their respective ends, i.e., the corners. If, in
addition $\mathcal{M}_{y}=-1$, the two modes would also appear at the two
edges in the $y$ direction. As a result, they can only reside at opposite
corners, as depicted in the upper right-hand panel of Fig. 1(b). If
$\mathcal{M}_{y}=1$, however, they would appear at adjacent corners along the
$y$ direction, as shown in the lower right-hand panel of Fig. 1(b).

The Majorana number defined in Eq. (2) is closely related to level crossings
at the Fermi energy $\epsilon=0$ that appear while $\lambda_{a}$ varies in the
range $[0,1]$. Notably, Eq. (2) only involves the 1D Hamiltonian at high-
symmetry momenta $K$. Therefore, we only need to consider Fermi level
crossings in these subsystems, as shown in Fig. 1(a). At each crossing, the
fermion parity of the ground state switches, indicated by the sign change of
$\text{Pf}[-iH_{\lambda_{a}}(K)]$. We can then characterize the fermion-parity
difference between PBC and OBC by the number of crossings in between, denoted
by $\eta_{a,K}$, as Fig. 1(c) demonstrates. This is formally expressed as

$(-1)^{\eta_{a,K}}=\frac{\text{sgn Pf}[-iH_{\lambda_{a}=0}(K)]}{\text{sgn
Pf}[-iH_{\lambda_{a}=1}(K)]}.$ (3)

We may also define a Majorana number for the toroidal system
$\tilde{\mathcal{H}}_{\lambda_{a}=1}$ ($\tilde{\mathcal{H}}_{1}$) similar to
Eq.(2), which due to trivial first-order topology must be positive, i.e.,

$\text{sgn}\prod_{K}\text{Pf}[-iH_{\lambda_{a}=1}(K)]=1.$ (4)

Combining Eqs. (2)-(4), we arrive at

$\mathcal{M}_{a}=\prod_{K}(-1)^{\eta_{a,K}}=(-1)^{\eta_{a}},$ (5)

where $\eta_{a}$ denotes the total number of crossings at $K=0,\pi$. An odd
value of $\eta_{x}$ or $\eta_{y}$ implies the system resides in a higher-order
phase. Fermi level crossings are protected by fermion-parity conservation and
particle-hole symmetry, making them immune to crystalline-symmetry-breaking
perturbations [50].

Intuitively, we may understand the relation between Fermi level crossings and
higher-order topology from the viewpoint of boundary topology. As shown in
Fig. 1(b), an odd value of $\eta_{a}$ ($\mathcal{M}_{a}=-1$) reveals that
opposite edges along $\bar{a}$ are topologically inequivalent (shown in
different colors). This explains the possible locations of Majorana zero
modes, which appear at the intersections of topologically distinct edges. In
some simple models, as we demonstrate later, the edge topology can be
characterized by the sign of the mass gap in the edge Hamiltonian, allowing us
to validate this argument.

To establish the bulk-boundary correspondence, we will demonstrate how the
Fermi level crossings of the 1D subsystems are identified from the bulk
Hamiltonian. For brevity, we use $\mathcal{H}_{\lambda}$ to replace
$\mathcal{H}_{\lambda_{a}}(K)$, where

$\mathcal{H}_{\lambda}=\mathcal{H}_{1}-(1-\lambda)\mathcal{B}$ (6)

represents a generic 1D Hamiltonian of D class. Following the prescription
given by Ref.[51], we first define a retarded Green’s function

$\mathcal{G}_{\lambda}(\epsilon)=(\epsilon-\mathcal{H}_{\lambda}+i\tilde{\delta})^{-1}=A^{-1}_{\lambda}(\epsilon)\mathcal{G}_{1}(\epsilon)\mathcal{,}$
(7)

where $\tilde{\delta}$ is a positive infinitesimal,
$\mathcal{G}_{1}(\epsilon)$ is the Green’s function corresponding to the bulk
Hamiltonian $\mathcal{H}_{1}$, and
$\mathcal{A}_{\lambda}(\epsilon)=1+(1-\lambda)\mathcal{G}_{1}(\epsilon)\mathcal{B}$.
Since we focus on the parameter regime in which the bulk is fully gapped, in-
gap states of $\mathcal{H}_{\lambda}$ are solely determined by the poles of
$\mathcal{A}^{-1}_{\lambda}(\epsilon)$. Consequently, level-crossing points
are identified as roots of

$\det[A_{\lambda}(\epsilon=0)]=0,$ (8)

where $A_{\lambda}$ is the matrix representation of $\mathcal{A}_{\lambda}$.
As $\mathcal{B}$ only includes intra-cell terms crossing the boundary, we then
have $[A_{\lambda}]_{IJ}=\delta_{IJ}$ if $J$ does not appear in these terms.
This enables us to calculate $\det(A)$ using a much smaller matrix $D$, which
is obtained by projecting $\mathcal{A}$ into the eigenspace of $\mathcal{B}$
and satisfies $\det(D)=\det(A)$. The entries of $D$ are given by

$[D_{\lambda}(\epsilon)]_{ij}=\delta_{ij}+(1-\lambda)\sum_{n,k}\frac{\langle
i|n,k\rangle\langle n,k|\mathcal{B}|j\rangle}{\epsilon-\epsilon_{n,k}},$ (9)

where $\epsilon_{n,k}$ denotes the energy spectrum of the bulk Hamiltonian
$\mathcal{H}_{1}$ in the Brillouin zone, with $|n,k\rangle$ being the
corresponding eigenstate, and $|i\rangle$,$|j\rangle$ represent the
eigenvectors of $\mathcal{B}$. The dimension of $D_{\lambda}$ is equal to the
rank of $\mathcal{B}$, denoted by $N_{b}$. We then obtain the characteristic
equation

$\det[D_{\lambda}(\epsilon=0)]=0,$ (10)

which has $N_{b}$ roots in total. The number of Fermi level crossings $\eta$
is half the number of real roots in the interval $[0,1]$, from which we can
readily obtain Majorana numbers according to Eq. (5).

Compared to Eq.(2), where Majorana numbers are determined by calculating the
Pfaffian of finite systems with open boundaries [52], i.e.,
$\text{Pf}[-iH_{\lambda_{a}=0}(K)]$, and the accuracy crucially depends on
system size, identifying the Fermi level crossings is computationally more
accurate and efficient for a translation-invariant system. It does not suffer
from finite-size effects, and the computational cost is similar to Wilson loop
calculations. Moreover, it provides a potential path to characterizing higher-
order topological superconductors in other symmetry classes such as the DIII
or BDI classes, where Fermi level crossings might be protected by their
topological charges. Additionally, by pinpointing the crossings directly from
the bulk Hamiltonian, we establish the correspondence between gapped bulk and
gapless corner states in higher-order phases. In the following, we shall
illustrate this in specific models.

Figure 2: (a) Geometry of the two-leg Kitaev ladder. $\lambda$ controls the
boundary condition. (b) The fermion-parity switch in a dimerized lattice
($t_{1}=0$). For the case of $|m|>|\delta m|$, the ground state switches from
the even-parity sector ($E_{e,1}$, solid lines) to the odd-parity sector
($E_{o,1}$, dashed lines) while the boundary condition varies from PBC (APBC)
to OBC. Level crossings are indicated by black arrows.

##  III toy models

First we consider a two-leg Kitaev ladder [53, 54, 55, 56] as schematically
shown in Fig. 2(a), and demonstrate how level crossings are identified from
the bulk Hamiltonian. Each unit cell contains four Majorana fermions denoted
by $\alpha_{s,j}$ and $\beta_{s,j}$, with $s=1,2$ being the chain index and
$j$ referring to the cell index. The boundary-modulated Hamiltonian with $N$
unit cells has the form $\mathcal{H}_{\lambda}=\Gamma^{T}H_{\lambda}\Gamma$ in
the Majorana basis $\Gamma=\bigoplus_{j=1}^{N}\Gamma_{j}$, where
$\Gamma_{j}=\\{\alpha_{1,j},\alpha_{2,j},\beta_{1,j},\beta_{2,j}\\}^{T}/\sqrt{2}$
and the Hamiltonian matrix is given by

$H_{\lambda}=\sum_{r=0,\pm 1}T^{r}\otimes h_{r}+\lambda(T^{N-1}\otimes
h_{1}^{\dagger}+\text{H.c.}).$ (11)

Here, $T$ denotes the translation operator that moves each cell by one site to
the left, with $T|j\rangle=|j-1\rangle$ and $T|j=1\rangle=0$ [57]. Hamiltonian
(11) includes the intra cell term $h_{0}=-t_{1}\tau_{y}-m\sigma_{y}-\delta
m\tau_{z}\sigma_{y}$, and inter cell term
$h_{1}=h_{-1}^{\dagger}=t_{2}(\tau_{y}+i\tau_{x})/2$, with $\tau$ and $\sigma$
being Pauli matrices that act in the chain and rung space separately. $t_{1}$
and $t_{2}$ represent couplings of Majorana fermions along the chain, while
$m$ and $\delta m$ are those along the rung. For brevity, we assume $t_{1}$
and $t_{2}$ to be non-negative.

In this model, $m$ and $\delta m$ determine whether level crossings occur when
$\lambda$ varies in the range $[0,1]$. This is readily seen in a perfectly
dimerized lattice ($t_{1}=0$), in which case only the boundary block shown in
Fig. 2(b) depends on $\lambda$, and its Hamiltonian has the form

$\displaystyle\mathcal{H}_{b}=$ $\displaystyle-2\lambda
t_{2}(\psi_{1}^{\dagger}\psi_{1}+\psi_{2}^{\dagger}\psi_{2}-1)$
$\displaystyle+2i(m\psi_{1}^{\dagger}\psi_{2}+\delta
m\psi_{1}^{\dagger}\psi_{2}^{\dagger}-\text{H.c.}),$ (12)

where $\psi_{s}=(\alpha_{s,1}+i\beta_{s,N})/2$ are fermionic operators. The
conservation of fermion number parity enables us to study the lowest energy
levels in the even- and odd-parity sectors separately, with
$E_{e,1}=-2\sqrt{\lambda^{2}t_{2}^{2}+\delta m^{2}}$ and $E_{o,1}=-2|m|$.
While the boundary condition goes from PBC to OBC, the two levels would cross
if $0<m^{2}-\delta m^{2}<t_{2}^{2}$, signaling a switch in the ground-state
fermion parity, as demonstrated in Fig. 2(b). This parity switch could be
observed from the zero-bias peak in an experimental setup that consists of two
quantum dots coupled by a nanowire-superconductor heterojunction [58, 59]. The
parameters $m$ and $\delta m$ are related to the electrochemical potential of
quantum dots, and $t_{2}$ or $\lambda$ is controlled by tuning the cross
Andreev reflection and elastic cotunnelling.

Figure 3: (a)-(c) Fermi level crossings in a BdG spectrum. The fermion parity
switches at each crossing. The crossing appears when $|m|>|\delta m|$. Blue
dashed lines indicate signs of the Pfaffian of a 1D Hamiltonian. (d) The
variation of crossing points $\lambda$ with $m$. $\delta m=0.5$, $t_{1}=1$,
and $t_{2}=2$.

For generic $t_{1}$, we have
$\mathcal{H}_{1}=\sum_{k}\Gamma_{-k}^{T}H_{k}\Gamma_{k}$ in $k$ space, with
the basis
$\Gamma_{k}=\\{\alpha_{1,k},\alpha_{2,k},\beta_{1,k},\beta_{2,k}\\}^{T}/\sqrt{2}$,
and the Bloch Hamiltonian

$H_{k}=(-t_{1}+t_{2}\cos k)\tau_{y}-t_{2}\sin k\tau_{x}-m\sigma_{y}-\delta
m\tau_{z}\sigma_{y}.$ (13)

The energy spectrum is given by

$\epsilon_{n,k}=\pm\left(\sqrt{t_{1}^{2}+t_{2}^{2}-2t_{1}t_{2}\cos k+\delta
m^{2}}\pm m\right),$ (14)

with $n$ being the band index. Substituting $\epsilon_{n,k}$ and
$\mathcal{B}=\Gamma_{N}^{T}h_{1}\Gamma_{1}+\text{H.c.}$ into Eqs. (9) and
(10), we obtain

$\lambda^{2}=1-\frac{2\Lambda}{\Lambda-a+2t_{2}^{2}},$ (15)

with $\Lambda=\sqrt{a^{2}-b^{2}}$, $a=t_{1}^{2}+t^{2}_{2}+\delta m^{2}-m^{2}$,
and $b=2t_{1}t_{2}$. From Eq. (15), we find that the number of crossings
$\eta=1$ when $0<m^{2}-\delta m^{2}<(t_{2}-t_{1})^{2}$ and $t_{1}<t_{2}$, as
shown in Fig. 3. This indicates that the boundary phase transition occurs at
$|m|=|\delta m|$ as in the dimerized case, which is verified by the exact
boundary spectrum (see Supplemental Material [60] and Ref. [61] therein). In
the special case where $\delta m=0$, Hamiltonian (13) is invariant under
inversion, with the corresponding operator being $\tau_{y}$, up to a gauge
factor. The inversion symmetry facilitates the direct determination of the
Fermi level crossings from the differences of the ground-state inversion
eigenvalues between PBC and APBC [60]. With the knowledge of $\eta$ in a 1D
system, we can proceed to determine the higher-order topology in a 2D system,
according to Eq. (5).

Figure 4: (a) Phase diagram of the 2D model in $(\theta,\delta m)$ space. Four
different phases are characterized by Majorana numbers
$(\mathcal{M}_{x},\mathcal{M}_{y})$. (b) Evolutions of crossing points with
$\theta$. No crossings appear at $K=\pi$ in this case and hence only those at
$K=0$ are displayed. (c) Distributions of Majorana zero modes in three
nontrivial phases that are separated by boundary phase transitions. Majorana
zero modes appear either at two adjacent corners ($\theta=0,\pi/2$), or at
opposite corners ($\theta=\pi/4$). The symbols “$+$” and “$-$” indicate the
signs of the edge gaps. $t_{1}=0.5$, $t_{2}=1$, and $m=0.4$.

The 2D Hamiltonian we consider takes the form

$\displaystyle H^{\text{2D}}_{\bm{k}}=[t_{2}(\cos k_{x}+\cos
k_{y})-t_{1}-t_{2}]\tau_{y}-\delta m\tau_{z}\sigma_{y}$ (16) $\displaystyle-
t_{2}(\sin k_{x}\tau_{x}+\sin
k_{y}\tau_{z}\sigma_{z})-m(\cos\theta\sigma_{y}+\sin\theta\tau_{y}\sigma_{x}),$

when written in the Majorana basis as in Eq.(13), and reduces to the 1D
Hamiltonian at $k_{y}=0$, $\theta=0$. This model is equivalent to the $p\pm
ip$ superconductor under an in-plane Zeeman field [62, 63]. According to Eq.
(5), Majorana numbers $(\mathcal{M}_{x},\mathcal{M}_{y})$ are determined by
Fermi level crossings of four 1D Hamiltonians, $\mathcal{H}_{\lambda_{a}}(K)$.
In Fig. 4(a), we draw the $(\theta,\delta m)$ phase diagram. Here the
crossings only occur at $K=0$ as Fig. 4(b) shows, although it is possible they
emerge at $K=\pi$ for $t_{1}$ and $t_{2}$ taking other values. Two Majorana
corner states emerge when at least one Majorana number takes $-1$, as
illustrated in Fig. 4(c).

To corroborate previous arguments concerning the relation between level
crossings and boundary topology, we obtain the mass gap for an arbitrary edge
[60], given by

$\Delta(\phi)=\delta m-m\cos(\phi-\theta),$ (17)

where $\phi$ indicates the normal direction of the edge ($\phi=0,\pi/2$ for
right and top edges respectively). The topology of the edges in D-class
systems can be characterized by the sign of the mass gap. As seen from the
three representative cases in Fig. 4(c), gaps of opposite edges along $y(x)$
indeed take different signs when $\eta_{x(y)}$ is an odd number, or
equivalently, $\mathcal{M}_{x(y)}=-1$. This can be guaranteed when inversion
symmetry is enforced, by noting that $\Delta(\phi)=-\Delta(\phi+\pi)$ in the
absence of $\delta m$. In this intrinsic higher-order phase, we always have
$\mathcal{M}_{x}=\mathcal{M}_{y}=-1$. Turning on $\delta m$ breaks inversion
symmetry and drives the system into a boundary-obstructed phase, in which
process the gap signs do not change immediately, so is the number of Fermi
level crossings. We can therefore use Fermi level crossings to characterize
the higher-order topology in both phases.

The robustness of Fermi level crossings is also reflected in their persistence
under weak disorder or boundary impurities [60]. While Eqs. (5) and (10) may
not be directly applicable due to potential broken translation symmetry, the
number of Fermi level crossings remains unchanged. This reinforces their role
as a reliable tool to characterize higher-order topological superconductors.

##  IV conclusion

In conclusion, Fermi level crossings can serve as useful indicators for
higher-order topology in the D symmetry class when the nontrivial phase
accommodates two Majorana corner states. The applicability of this approach
extends beyond the toy models introduced above, as demonstrated in the
Supplemental Material [60] for a Rashba bilayer system. The level crossings we
consider emerge while the boundary condition continuously varies from PBC to
OBC, during which two opposite edges gradually decouple. An odd number of
crossings signals a topological distinction between the two edges. From this
point of view, one may consider Fermi level crossings emerging under
variations of other twisted boundary conditions [64] when dealing with higher-
order phases with four or more Majorana corner states, where one needs to
associate the crossings with topological distinctions between neighboring
edges.

## Acknowledgements

This work was supported by National Science Foundation of China (NSFC) under
Grant No. 11704305, and the Innovation Program for Quantum Science and
Technology (2021ZD0302400).

## References

  * Hasan and Kane [2010] M. Z. Hasan and C. L. Kane,  _Colloquium_ : Topological insulators, Rev. Mod. Phys. 82, 3045 (2010).
  * Qi and Zhang [2011] X.-L. Qi and S.-C. Zhang, Topological insulators and superconductors, Rev. Mod. Phys. 83, 1057 (2011).
  * Chiu  _et al._ [2016] C.-K. Chiu, J. C. Y. Teo, A. P. Schnyder, and S. Ryu, Classification of topological quantum matter with symmetries, Rev. Mod. Phys. 88, 035005 (2016).
  * Benalcazar  _et al._ [2017a] W. A. Benalcazar, B. A. Bernevig, and T. L. Hughes, Quantized electric multipole insulators, Science 357, 61 (2017a).
  * Benalcazar  _et al._ [2017b] W. A. Benalcazar, B. A. Bernevig, and T. L. Hughes, Electric multipole moments, topological multipole moment pumping, and chiral hinge states in crystalline insulators, Phys. Rev. B 96, 245115 (2017b).
  * Langbehn  _et al._ [2017] J. Langbehn, Y. Peng, L. Trifunovic, F. von Oppen, and P. W. Brouwer, Reflection-Symmetric Second-Order Topological Insulators and Superconductors, Phys. Rev. Lett. 119, 246401 (2017).
  * Song  _et al._ [2017] Z. Song, Z. Fang, and C. Fang, ( d - 2 ) -Dimensional Edge States of Rotation Symmetry Protected Topological States, Phys. Rev. Lett. 119, 246402 (2017).
  * Ezawa [2018] M. Ezawa, Higher-Order Topological Insulators and Semimetals on the Breathing Kagome and Pyrochlore Lattices, Phys. Rev. Lett. 120, 026801 (2018).
  * Schindler  _et al._ [2018a] F. Schindler, A. M. Cook, M. G. Vergniory, Z. Wang, S. S. P. Parkin, B. A. Bernevig, and T. Neupert, Higher-order topological insulators, Sci. Adv. 4, eaat0346 (2018a).
  * Schindler  _et al._ [2018b] F. Schindler, Z. Wang, M. G. Vergniory, A. M. Cook, A. Murani, S. Sengupta, A. Y. Kasumov, R. Deblock, S. Jeon, I. Drozdov, H. Bouchiat, S. Guéron, A. Yazdani, B. A. Bernevig, and T. Neupert, Higher-order topology in bismuth, Nature Phys 14, 918 (2018b).
  * Khalaf  _et al._ [2021] E. Khalaf, W. A. Benalcazar, T. L. Hughes, and R. Queiroz, Boundary-obstructed topological phases, Phys. Rev. Research 3, 013239 (2021).
  * Wang  _et al._ [2018] Q. Wang, C.-C. Liu, Y.-M. Lu, and F. Zhang, High-Temperature Majorana Corner States, Phys. Rev. Lett. 121, 186801 (2018).
  * Zhang  _et al._ [2019a] R.-X. Zhang, W. S. Cole, and S. Das Sarma, Helical Hinge Majorana Modes in Iron-Based Superconductors, Phys. Rev. Lett. 122, 187001 (2019a).
  * Khalaf [2018] E. Khalaf, Higher-order topological insulators and superconductors protected by inversion symmetry, Phys. Rev. B 97, 205136 (2018).
  * Geier  _et al._ [2018] M. Geier, L. Trifunovic, M. Hoskam, and P. W. Brouwer, Second-order topological insulators and superconductors with an order-two crystalline symmetry, Phys. Rev. B 97, 205135 (2018).
  * Trifunovic and Brouwer [2019] L. Trifunovic and P. W. Brouwer, Higher-Order Bulk-Boundary Correspondence for Topological Crystalline Phases, Phys. Rev. X 9, 011012 (2019).
  * Skurativska  _et al._ [2020] A. Skurativska, T. Neupert, and M. H. Fischer, Atomic limit and inversion-symmetry indicators for topological superconductors, Phys. Rev. Research 2, 013064 (2020).
  * Ono  _et al._ [2020] S. Ono, H. C. Po, and H. Watanabe, Refined symmetry indicators for topological superconductors in all space groups, Sci. Adv. 6, eaaz8367 (2020).
  * Takahashi  _et al._ [2020] R. Takahashi, Y. Tanaka, and S. Murakami, Bulk-edge and bulk-hinge correspondence in inversion-symmetric insulators, Phys. Rev. Research 2, 013300 (2020).
  * Hsu  _et al._ [2020] Y.-T. Hsu, W. S. Cole, R.-X. Zhang, and J. D. Sau, Inversion-protected Higher-order Topological Superconductivity in Monolayer WTe 2, Phys. Rev. Lett. 125, 097001 (2020).
  * Tang  _et al._ [2022] F. Tang, S. Ono, X. Wan, and H. Watanabe, High-Throughput Investigations of Topological and Nodal Superconductors, Phys. Rev. Lett. 129, 027001 (2022).
  * Yan [2019] Z. Yan, Higher-Order Topological Odd-Parity Superconductors, Phys. Rev. Lett. 123, 177001 (2019).
  * Zhang  _et al._ [2022] Z. Zhang, J. Ren, Y. Qi, and C. Fang, Topological classification of intrinsic three-dimensional superconductors using anomalous surface construction, Phys. Rev. B 106, L121108 (2022).
  * Bouhon  _et al._ [2019] A. Bouhon, A. M. Black-Schaffer, and R.-J. Slager, Wilson loop approach to fragile topology of split elementary band representations and topological crystalline insulators with time-reversal symmetry, Phys. Rev. B 100, 195135 (2019).
  * Hwang  _et al._ [2019] Y. Hwang, J. Ahn, and B.-J. Yang, Fragile topology protected by inversion symmetry: Diagnosis, bulk-boundary correspondence, and Wilson loop, Phys. Rev. B 100, 205126 (2019).
  * Kruthoff  _et al._ [2017] J. Kruthoff, J. de Boer, J. van Wezel, C. L. Kane, and R.-J. Slager, Topological Classification of Crystalline Insulators through Band Structure Combinatorics, Phys. Rev. X 7, 041069 (2017).
  * Tang  _et al._ [2019] F. Tang, H. C. Po, A. Vishwanath, and X. Wan, Comprehensive search for topological materials using symmetry indicators, Nature 566, 486 (2019).
  * Zhang  _et al._ [2019b] T. Zhang, Y. Jiang, Z. Song, H. Huang, Y. He, Z. Fang, H. Weng, and C. Fang, Catalogue of topological electronic materials, Nature 566, 475 (2019b).
  * Vergniory  _et al._ [2019] M. G. Vergniory, L. Elcoro, C. Felser, N. Regnault, B. A. Bernevig, and Z. Wang, A complete catalogue of high-quality topological materials, Nature 566, 480 (2019).
  * Zhang [2022] R.-X. Zhang, Bulk-Vortex Correspondence of Higher-Order Topological Superconductors (2022), arXiv:2208.01652 .
  * Jung  _et al._ [2021] M. Jung, Y. Yu, and G. Shvets, Exact higher-order bulk-boundary correspondence of corner-localized states, Phys. Rev. B 104, 195437 (2021).
  * Roberts  _et al._ [2020] E. Roberts, J. Behrends, and B. Béri, Second-order bulk-boundary correspondence in rotationally symmetric topological superconductors from stacked Dirac Hamiltonians, Phys. Rev. B 101, 155133 (2020).
  * Kooi  _et al._ [2021] S. Kooi, G. van Miert, and C. Ortix, The bulk-corner correspondence of time-reversal symmetric insulators, npj Quantum Mater. 6, 1 (2021).
  * Huang and Hsu [2021] S.-J. Huang and Y.-T. Hsu, Faithful derivation of symmetry indicators: A case study for topological superconductors with time-reversal and inversion symmetries, Phys. Rev. Research 3, 013243 (2021).
  * Zhu [2019] X. Zhu, Second-Order Topological Superconductors with Mixed Pairing, Phys. Rev. Lett. 122, 236401 (2019).
  * Yan  _et al._ [2018] Z. Yan, F. Song, and Z. Wang, Majorana Corner Modes in a High-Temperature Platform, Phys. Rev. Lett. 121, 096803 (2018).
  * Liu  _et al._ [2018] T. Liu, J. J. He, and F. Nori, Majorana corner states in a two-dimensional magnetic topological insulator on a high-temperature superconductor, Phys. Rev. B 98, 245413 (2018).
  * Tiwari  _et al._ [2020] A. Tiwari, A. Jahin, and Y. Wang, Chiral Dirac superconductors: Second-order and boundary-obstructed topology, Phys. Rev. Research 2, 043300 (2020).
  * Volpez  _et al._ [2019] Y. Volpez, D. Loss, and J. Klinovaja, Second-Order Topological Superconductivity in $\pi$ -Junction Rashba Layers, Phys. Rev. Lett. 122, 126402 (2019).
  * Wu  _et al._ [2020] X. Wu, W. A. Benalcazar, Y. Li, R. Thomale, C.-X. Liu, and J. Hu, Boundary-Obstructed Topological High- T c Superconductivity in Iron Pnictides, Phys. Rev. X 10, 041014 (2020).
  * Read and Green [2000] N. Read and D. Green, Paired states of fermions in two dimensions with breaking of parity and time-reversal symmetries and the fractional quantum Hall effect, Phys. Rev. B 61, 10267 (2000).
  * Wilczek [2009] F. Wilczek, Majorana returns, Nature Phys 5, 614 (2009).
  * Alicea [2012] J. Alicea, New directions in the pursuit of Majorana fermions in solid state systems, Rep. Prog. Phys. 75, 076501 (2012).
  * Stanescu and Tewari [2013] T. D. Stanescu and S. Tewari, Majorana fermions in semiconductor nanowires: Fundamentals, modeling, and experiment, J. Phys.: Condens. Matter 25, 233201 (2013).
  * Elliott and Franz [2015] S. R. Elliott and M. Franz,  _Colloquium_ : Majorana fermions in nuclear, particle, and solid-state physics, Rev. Mod. Phys. 87, 137 (2015).
  * Aguado [2017] R. Aguado, Majorana quasiparticles in condensed matter, Riv. Nuovo Cimento 40, 523 (2017).
  * Kitaev [2001] A. Y. Kitaev, Unpaired Majorana fermions in quantum wires, Phys.-Usp. 44, 131 (2001).
  * Kheirkhah  _et al._ [2021] M. Kheirkhah, Z. Yan, and F. Marsiglio, Vortex-line topology in iron-based superconductors with and without second-order topology, Phys. Rev. B 103, L140502 (2021).
  * Poduval  _et al._ [2023] P. P. Poduval, T. L. Schmidt, and A. Haller, Perfectly localized Majorana corner modes in fermionic lattices (2023), arXiv:2303.01535 .
  * Beenakker  _et al._ [2013] C. W. J. Beenakker, J. M. Edge, J. P. Dahlhaus, D. I. Pikulin, S. Mi, and M. Wimmer, Wigner-Poisson Statistics of Topological Transitions in a Josephson Junction, Phys. Rev. Lett. 111, 037001 (2013).
  * Rhim  _et al._ [2018] J.-W. Rhim, J. H. Bardarson, and R.-J. Slager, Unified bulk-boundary correspondence for band insulators, Phys. Rev. B 97, 115143 (2018).
  * Wimmer [2012] M. Wimmer, Algorithm 923: Efficient Numerical Computation of the Pfaffian for Dense and Banded Skew-Symmetric Matrices, ACM Trans. Math. Softw. 38, 1 (2012).
  * Wu [2012] N. Wu, Topological phases of the two-leg Kitaev ladder, Physics Letters A 376, 3530 (2012).
  * Chitov [2018] G. Y. Chitov, Local and nonlocal order parameters in the Kitaev chain, Phys. Rev. B 97, 085131 (2018).
  * Wakatsuki  _et al._ [2014] R. Wakatsuki, M. Ezawa, and N. Nagaosa, Majorana fermions and multiple topological phase transition in Kitaev ladder topological superconductors, Phys. Rev. B 89, 174514 (2014).
  * Yan  _et al._ [2020] Y. Yan, L. Qi, D.-Y. Wang, Y. Xing, H.-F. Wang, and S. Zhang, Topological Phase Transition and Phase Diagrams in a Two-Leg Kitaev Ladder System, Ann. Phys. 532, 1900479 (2020).
  * Alase  _et al._ [2016] A. Alase, E. Cobanera, G. Ortiz, and L. Viola, Exact Solution of Quadratic Fermionic Hamiltonians for Arbitrary Boundary Conditions, Phys. Rev. Lett. 117, 076804 (2016).
  * Dvir  _et al._ [2023] T. Dvir, G. Wang, N. van Loo, C.-X. Liu, G. P. Mazur, A. Bordin, S. L. D. ten Haaf, J.-Y. Wang, D. van Driel, F. Zatelli, X. Li, F. K. Malinowski, S. Gazibegovic, G. Badawy, E. P. A. M. Bakkers, M. Wimmer, and L. P. Kouwenhoven, Realization of a minimal Kitaev chain in coupled quantum dots, Nature 614, 445 (2023).
  * Bordin  _et al._ [2023] A. Bordin, X. Li, D. van Driel, J. C. Wolff, Q. Wang, S. L. D. ten Haaf, G. Wang, N. van Loo, L. P. Kouwenhoven, and T. Dvir, Crossed Andreev reflection and elastic co-tunneling in a three-site Kitaev chain nanowire device (2023), arXiv:2306.07696 .
  * [60] See Supplemental Material at [**] for detailed derivations of the boundary spectrum of the two-leg Kitaev ladder; the role of inversion symmetry in determining Fermi-level crossings; the effective Hamiltonian for an arbitrary edge in the 2D model; the influence of bulk disorder and boundary impurities on Fermi-level crossings; Fermi-level crossings in a Rashba bilayer superconducting system.
  * Pershoguba and Yakovenko [2012] S. S. Pershoguba and V. M. Yakovenko, Shockley model description of surface states in topological insulators, Phys. Rev. B 86, 075304 (2012).
  * Phong  _et al._ [2017] V. T. Phong, N. R. Walet, and F. Guinea, Majorana zero modes in a two-dimensional p -wave superconductor, Phys. Rev. B 96, 060505(R) (2017).
  * Zhu [2018] X. Zhu, Tunable Majorana corner states in a two-dimensional second-order topological superconductor induced by magnetic fields, Phys. Rev. B 97, 205134 (2018).
  * Song  _et al._ [2020] Z.-D. Song, L. Elcoro, and B. A. Bernevig, Twisted bulk-boundary correspondence of fragile topology, Science 367, 794 (2020).

Supplemental Material for “Higher-order topological superconductors
characterized by Fermi level crossings”

In this Supplemental Material, we provide detailed derivations of the boundary
spectrum for the two-leg Kitaev ladder, explore the role of inversion symmetry
in determining Fermi level crossings, derive the effective Hamiltonian for an
arbitrary edge in the 2D model, and investigate the robustness of Fermi level
crossings against bulk disorder and boundary impurities. We also demonstrate
in a Rashba bilayer superconducting system how Fermi level crossings
effectively identify higher-order topological phases.

##  Appendix A A. Boundary spectrum

We consider a semi-infinite system with boundary at $j=1$ in the two-leg
Kitaev ladder. In the Hamiltonian $H_{\lambda=0}$, $\sigma_{y}$ is a good
quantum number, and we can work in its eigenspace, where $H_{\lambda=0}$ is
block diagonal. Consequently, we can set $\sigma_{y}$ to be $\pm 1$ in the two
blocks, respectively. In this toy model, each block with $\sigma_{y}=1$
($\sigma_{y}=-1$) can be viewed as a particle (hole) version of Su-Schrieffer-
Heeger (SSH) model, with the parameter $m$ acting as a chemical potential term
that shifts the energy spectrum in corresponding block. The two blocks do not
couple due to the conservation of $\sigma_{y}$ in this simple model. It is
possible to introduce additional terms that couple the two blocks and make the
model more complicated. However, the main results do not change as we only
require particle-hole symmetry. The simplicity of this toy model allows us to
obtain analytical results in a straightforward manner, as we demonstrate in
the following.

We will now derive the condition for the appearance of gapped boundary modes
in the two-leg Kitaev ladder, as well as the boundary spectrum. Let’s first
consider the block with $\sigma_{y}=1$, and the case with $\sigma_{y}=-1$ can
be obtained by sending $m\rightarrow-m$ and $\delta m\rightarrow-\delta m$.
The Hamiltonian matrix with $\sigma_{y}=+1$ is given by

$H=\sum_{r=0,\pm 1}T^{r}\otimes\tilde{h}_{r},$ (1)

where $\tilde{h}_{0}=-t_{1}\tau_{y}-m-\delta m\tau_{z}$,
$\tilde{h}_{1}=\tilde{h}_{-1}^{\dagger}=t_{2}(\tau_{y}+i\tau_{x})/2$ and $T$
is the translation operator that moves each unit cell by one site to the left,
with $T|j\rangle=|j-1\rangle$ and $T|j=1\rangle=0$. The wavefunction
$\Psi=\bigoplus_{j=1}\phi_{j}=\bigoplus_{j=1}\\{u_{j},v_{j}\\}^{T}$ satisfies
the Schrödinger equation $H\Psi=E\Psi$, which has the form

$(\tilde{h}_{0}-E)\phi_{1}+\tilde{h}_{1}\phi_{2}=0$ (2)

and

$\tilde{h}_{1}^{\dagger}\phi_{j-1}+(\tilde{h}_{0}-E)\phi_{j}+\tilde{h}_{1}\phi_{j+1}=0$
(3)

for $j>1$. Multiplying Eq.(3) by $q^{j-1}$ and summing up all the equations,
we obtain

$[q\tilde{h}_{1}^{\dagger}+(\tilde{h}_{0}-E)+q^{-1}\tilde{h}_{1}]G(q)=(\tilde{h}_{0}-E+\tilde{h}_{1}q^{-1})\phi_{1}+\tilde{h}_{1}\phi_{2}$
(4)

where

$G(q)=\sum_{j=1}^{+\infty}q^{j-1}\phi_{j},$ (5)

and $q$ is a complex number. Utilizing Eq.(2) and (4), we can express $G(q)$
as

$G(q)=[q^{2}\tilde{h}_{1}^{\dagger}+q(\tilde{h}_{0}-E)+\tilde{h}_{1}]^{-1}\tilde{h}_{1}\phi_{1}.$
(6)

For $\Psi$ to be a localized state at $j=1$, all the poles of $G(q)$ must
satisfy $|q_{p}|>1$ [61]. Substituting the specific forms of $\tilde{h}_{0}$
and $\tilde{h}_{1}$ into Eq.(6), we have

$G(q)=\frac{u_{1}t_{2}}{D}\biggl{\\{}t_{1}-qt_{2},-i(m+\delta
m+E)\biggl{\\}}^{T},$ (7)

with the denominator

$D=(t_{1}-qt_{2})(t_{2}-qt_{1})+q[(E+m)^{2}-\delta m^{2}].$ (8)

The poles are decided from the two roots of $D=0$, which have the relation
$q_{1}q_{2}=1$, and therefore cannot both have absolute values greater than
one. So only one of the two roots can be the pole and the other has to be
eliminated from numerator. By requiring $(E+m)^{2}-\delta m^{2}=0$, we could
eliminate $t_{1}-qt_{2}$ term in the first entry of Eq.(7). However, only when
$E=-m-\delta m$ can this term be eliminated in both entries, which leads to

$G(q)=\frac{u_{1}t_{2}}{t_{2}-qt_{1}}\\{1,0\\}^{T}.$ (9)

The pole $q_{p}=t_{2}/t_{1}$, satisfies $|q_{p}|>1$ when $t_{1}<t_{2}$. A
series expansion of $G(q)$ takes the form

$G(q)=\sum_{j=1}^{+\infty}q^{j-1}\left(\frac{t_{1}}{t_{2}}\right)^{j-1}\\{u_{1},0\\}^{T}.$
(10)

Comparing this equation with the definition of $G(q)$, we immediately find
that

$\phi_{j}=\left(\frac{t_{1}}{t_{2}}\right)^{j-1}\\{u_{1},0\\}^{T},$ (11)

which is clearly localized at $j=1$ if $t_{1}<t_{2}$. For the block with
$\sigma_{y}=-1$, the energy of bound state is given by $E=m+\delta m$. Hence
the boundary mode at $j=1$ appears when $t_{1}<t_{2}$. The boundary states for
$H_{\lambda=0}$ at $j=1$ can be written in original Majorana basis, with

$\displaystyle E_{b1,+}=m+\delta m,$
$\displaystyle\Psi_{b1,+}=\frac{u_{1}}{\sqrt{2}}\bigoplus_{j=1}\left(\frac{t_{1}}{t_{2}}\right)^{j-1}\\{1,-i,0,0\\}^{T}$
(12) $\displaystyle E_{b1,-}=-(m+\delta m),$
$\displaystyle\Psi_{b1,-}=\frac{u_{1}}{\sqrt{2}}\bigoplus_{j=1}\left(\frac{t_{1}}{t_{2}}\right)^{j-1}\\{1,i,0,0\\}^{T}$
(13)

where $u_{1}$ is decided from normalization condition.

Boundary modes at the other end can be obtained in a similar way. We consider
a semi-infinite system with boundary at $j=N$. Schrödinger equation for
$\sigma_{y}=1$ block is given by

$\tilde{h}_{1}^{\dagger}\phi_{N-1}+\tilde{h}_{0}\phi_{N}=E\phi_{N}$ (14)

and

$\tilde{h}_{1}^{\dagger}\phi_{j-1}+(\tilde{h}_{0}-E)\phi_{j}+\tilde{h}_{1}\phi_{j+1}=0$
(15)

for $j\in(-\infty,N-1]$. Multiplying Eq.(15) by $q^{N-j}$ and summing up all
the equations, we have

$G(q)=\sum_{j=N}^{-\infty}q^{N-j}\phi_{j}=[\tilde{h}_{1}^{\dagger}+q(\tilde{h}_{0}-E)+q^{2}\tilde{h}_{1}]^{-1}\tilde{h}_{1}^{\dagger}\phi_{N}$
(16)

Following the same analysis in deriving boundary modes at $j=1$, we obtain the
boundary states at $j=N$, with

$\displaystyle E_{bN,+}=m-\delta m,$
$\displaystyle\Psi_{bN,+}=\frac{v_{N}}{\sqrt{2}}\bigoplus_{j=N}^{-\infty}\left(\frac{t_{1}}{t_{2}}\right)^{N-j}\\{0,0,1,-i\\}^{T}$
(17) $\displaystyle E_{bN,-}=-(m-\delta m),$
$\displaystyle\Psi_{bN,-}=\frac{v_{N}}{\sqrt{2}}\bigoplus_{j=N}^{-\infty}\left(\frac{t_{1}}{t_{2}}\right)^{N-j}\\{0,0,1,i\\}^{T}.$
(18)

So we have established that gapped boundary modes appear when $t_{1}<t_{2}$,
and the boundary spectrum is independent of $t_{1}$, being $\pm(m+\delta m)$
at the end $j=1$, and $\pm(m-\delta m)$ at $j=N$. Note that, these boundary
modes may appear in the bulk continuum if bulk gap vanishes or is smaller than
the boundary gap. For a gapped bulk, boundary phase transition occurs when
$|m|=|\delta m|$. In this context, we can treat each boundary as a zero-
dimensional gapped system that belongs to D class, which is again
characterized by a $\mathbb{Z}_{2}$ invariant [3]. One can in principle assign
a different invariant at two sides of the phase transition point based on the
number of Fermi level crossings. However, this is not related to gapless
modes.

##  Appendix B B. Inversion symmetry

Figure 1: (a) Fermion-parity switch with inversion symmetry enforced ($\delta
m=0$). Although the ground state falls in even-parity sector for both PBC and
APBC, they have distinct inversion eigenvalues (shown in the parenthesis). (b)
For $t_{1}<t_{2}$ the spectrum exhibits spectral flows with some states moving
from occupied (negative energies) to unoccupied (positive energies) bands. (c)
Due to inversion symmetry, level crossings can only disappear when bulk gap
closes and reopen.

In this section, we demonstrate that with inversion symmetry enforced, level
crossings can be inferred from ground-state difference of inversion
eigenvalues for system under PBC and APBC.

In the absence of $\delta m$, the 1D Kitaev ladder is invariant under
inversion, which transforms Majorana fermions as

$\mathcal{I}\alpha_{s,j}\mathcal{I}^{-1}=-\beta_{s,N+1-j},\
\mathcal{I}\beta_{s,j}\mathcal{I}^{-1}=\alpha_{s,N+1-j}.$ (19)

Accordingly, fermionic operators $\psi_{s}$ follow transformation
$\mathcal{I}\psi_{s}\mathcal{I}^{-1}=i\psi_{s}$. Although fermion parity is
still the same for PBC and APBC, we could discriminate the two ground states
by inversion eigenvalues. Specifically, for $t_{1}=0$, provided $m$ is finite,
the ground state would evolve from
$|11\rangle=\psi_{1}^{\dagger}\psi_{2}^{\dagger}|00\rangle$ under PBC to
$|00\rangle$ in APBC. The inversion eigenvalues of the two states are $-1$ and
$1$ respectively. This distinction leads to level crossings, as shown in Fig.
1(a).

As Fig. 1(b) shows, the level crossings persist for finite $t_{1}$, with two
states from occupied bands (negative energy) moving straight into unoccupied
bands (positive energy) as $\lambda$ varies. This suggests that ground state
under APBC should be different from that of PBC in some aspect. As we
illustrate in Fig. 1(a) for $t_{1}=0$, the difference lies in their inversion
eigenvalues. To investigate this for general $t_{1}$, we may write down the
mean-field ground state explicitly. Define fermionic operators
$c_{s,k}=(\alpha_{s,k}+i\beta_{s,k})/2$, and we could express the 1D
Hamiltonian in particle-hole basis
$\\{c_{1,k},c_{2,k},c_{1,-k}^{\dagger},c_{2,-k}^{\dagger}\\}^{T}$, which takes
the form

$\tilde{H}_{k}=(t_{1}-t_{2}\cos k)\tilde{\tau}_{z}+t_{2}\sin
k\tilde{\tau}_{y}-m\tilde{\sigma}_{y},$ (20)

with $\tilde{\tau}$ and $\tilde{\sigma}$ being Pauli matrices acting in Nambu
space and chain space. Ground state of this Bogoliubov-de Gennes (BdG)
Hamiltonian is given by [41]

$|G\rangle=\mathcal{N}\exp\left(\sum_{k>0,s}g_{k}c_{s,k}^{\dagger}c_{s,-k}^{\dagger}\right)\prod_{n\in\text{occ},K}(\phi_{n,K}^{\dagger})^{\nu_{n,K}}|0\rangle,$
(21)

where $\mathcal{N}$ is a normalization factor, $g_{k}=it_{2}\sin
k/[\epsilon_{0}+(t_{1}-t_{2}\cos k)]$ with
$\epsilon_{0}=\sqrt{t_{1}^{2}+t_{2}^{2}-2t_{1}t_{2}\cos k}$, and
$\phi_{n,K}^{\dagger}$ is the creation operator of Bogoliubov quasiparticle in
occupied bands at high symmetry momenta. Occupation number $\nu_{n,K}$ is
determined by the eigenstate of $\phi_{n,K}$, denoted by
$|\widetilde{n,K}\rangle=\\{u_{n1,K},u_{n2,K},v_{n1,K},v_{n2,K}\\}^{T}$, where
$u_{ns,K}$ and $v_{ns,K}$ are particle and hole components respectively. We
would then have
$\phi_{n,K}^{\dagger}=\sum_{s}u_{ns,K}c_{s,K}^{\dagger}+v_{ns,K}c_{s,K}$. At
$K$, there is no pairing term. Only the state that is of particle type
($v_{ns,K}=0$) will be occupied in the ground state, with $\nu_{n,K}=1$.

According to Eq.(19),
$\mathcal{I}c_{s,k}\mathcal{I}^{-1}=ic_{s,-k}e^{-ik(N+1)}$, which is differing
by a sign for PBC and APBC. Inversion symmetry requires BdG Hamiltonian (20)
to obey $U_{I}^{\dagger}(-k)\tilde{H}(k)U_{I}(-k)=\tilde{H}(-k)$, with
$U_{I}(k)=ie^{ik}\tilde{\tau}_{z}$. Each state at $K$ would be an eigenstate
of $\tilde{\tau}_{z}$, of which particle states satisfying
$\tilde{\tau}_{z}|\widetilde{n,K}\rangle=|\widetilde{n,K}\rangle$. Therefore,
$\nu_{n,K}$ is identified as the number of occupied states
$|\widetilde{n,K}\rangle$ with $\tilde{\tau}_{z}=1$. Under inversion
transformations, pairing terms in Eq.(21) remain the same while
$\phi_{n,K}^{\dagger}$ acquires a factor of $-ie^{iK}$ (under PBC). The
difference of ground-state inversion eigenvalues between PBC and APBC can then
be given by

$\mathcal{P}=\prod_{n\in\text{occ}}(-i)^{\nu_{n,0}-\nu_{n,\pi}},$ (22)

which is valid regardless of $N$ being even or odd. $\mathcal{P}$ can take
four different values and therefore serves as a $\mathbb{Z}_{4}$ invariant.
Ground-state fermion-parity difference between PBC and APBC may also be
expressed with the occupation numbers, given by

$n_{F}=\prod_{n\in\text{occ},K}(-1)^{\nu_{n,K}}.$ (23)

For $\mathcal{P}\neq 1$, the spectrum exhibits spectral flows while $\lambda$
varies, with some states moving from occupied bands to unoccupied bands, and
therefore level crossing is inevitable. When $\mathcal{P}=\pm i$, there would
be an odd number of level crossings while $\lambda$ varies between PBC and
APBC due to the fermion-parity difference ($n_{F}=-1$), thus leaving an
unpaired Majorana zero mode at each open boundary. The gapped boundaries with
nontrivial topology is characterized by $\mathcal{P}=-1$.

##  Appendix C C. Effective Edge Hamiltonian

Figure 2: A semi-infinite 2D system with boundary (blue line) along $k_{2}$
direction. Two sets of coordinate systems are related by an in-plane rotation.
$\hat{n}$ represents the unit vector pointing along normal direction of the
edge. Figure 3: (a) The variation of lowest energy level (nonnegative) with
$\lambda$ and disorder strength $\sigma_{t}$. Level crossings are obvious at
weak disorder and become obscured when the disorder is strong enough.
Comparing energy spectrum at $\sigma_{t}=0.3$ in (b) and $\sigma_{t}=0.8$ in
(c), we find that disorder in $t_{1}$ term mainly influences the bulk gap. (d)
The variation of lowest energy level (nonnegative) with $\lambda$ and impurity
strength $V_{m}$. The impurity is added solely at site $j=1$. The level
crossings move under the variation of impurity strength and disappear when
boundary gap closes. Energy spectrum at two representative impurity strength,
$V_{m}=0.7$ and $V_{m}=1.2$, are plotted in (e) and (f) respectively. In all
the figures, $\bar{t}_{1}=1$, $t_{2}=2$, $m=0.5$, $\delta m=0.2$.

In this section, we derive the effective boundary Hamiltonian of the 2D model
for an arbitrary edge.

In the absence of $m$ and $\delta m$ terms, the bulk gap closes at $\Gamma$
point when $t_{1}=t_{2}$. Near this critical point, we write down the
continuum model by expanding the bulk Hamiltonian at $\Gamma$ point up to
second order in $\bm{k}$, which reads

$H(k_{x},k_{y})=[(t_{2}-t_{1})-\frac{t_{2}}{2}(k_{x}^{2}+k_{y}^{2})]\tau_{y}-t_{2}(k_{x}\tau_{x}+k_{y}\tau_{z}\sigma_{z}).$
(24)

To derive the effective Hamiltonian of an arbitrary edge, we consider another
coordinate system $k_{1}k_{2}$ that is obtained by rotating $k_{x}k_{y}$
system counterclockwise by an angle $\phi$, as illustrated in Fig. 2.
Coordinates in the two systems are related by

$k_{x}=k_{1}\cos\phi-k_{2}\sin\phi,\ k_{y}=k_{1}\sin\phi+k_{2}\cos\phi.$ (25)

Substituting Eq.(25) into the continuum Hamiltonian, we have

$\displaystyle
H^{\prime}(k_{1},k_{2})=[(t_{2}-t_{1})-\frac{t_{2}}{2}(k_{1}^{2}+k_{2}^{2})]\tau_{y}$
(26) $\displaystyle-t_{2}[(k_{1}\cos\phi-
k_{2}\sin\phi)\tau_{x}+(k_{1}\sin\phi+k_{2}\cos\phi)\tau_{z}\sigma_{z}].$

We further rotate the inner basis with
$U_{\phi}=e^{-i\phi/2\tau_{y}\sigma_{z}}$, and the resulting Hamiltonian for
finite $m$ and $\delta m$ takes the form

$\displaystyle\tilde{H}(k_{1},k_{2})$
$\displaystyle=UH_{1}(k_{1},k_{2})U^{\dagger}$ (27)
$\displaystyle=[(t_{2}-t_{1})-\frac{t_{2}}{2}(k_{1}^{2}+k_{2}^{2})]\tau_{y}$
$\displaystyle-t_{2}(k_{1}\tau_{x}+k_{2}\tau_{z}\sigma_{z})$
$\displaystyle-\delta
m\tau_{z}\sigma_{y}-m\sigma_{y}e^{i(\phi-\theta)\tau_{y}\sigma_{z}}$

To derive boundary Hamiltonian, we consider a semi-infinite system with edges
along $k_{2}$ direction (blue lines in Fig. 2). Hamiltonian for this semi-
infinite 2D system is obtained by making a substitution
$k_{1}\rightarrow-i\partial_{x}$ while keeping $k_{2}$ intact, which leads to

$\displaystyle\tilde{H}(-i\partial_{x},k_{2})$
$\displaystyle=[(t_{2}-t_{1})+\frac{t_{2}}{2}\partial_{x}^{2}]\tau_{y}+it_{2}\partial_{x}\tau_{x}$
$\displaystyle-\frac{t_{2}}{2}k_{2}^{2}\tau_{y}-t_{2}k_{2}\tau_{z}\sigma_{z}$
$\displaystyle-\delta
m\tau_{z}\sigma_{y}-m\sigma_{y}e^{i(\phi-\theta)\tau_{y}\sigma_{z}}$ (28)

In the absence of $m$ and $\delta m$ terms, the model supports two Majorana
zero modes at $k_{2}=0$ in topological phase, which are localized on the edge.
Eigenstates for these two modes are obtained by solving
Schr$\ddot{\text{o}}$dinger equation for $m=\delta m=0$, i.e.,

$\tilde{H}(-i\partial_{x},k_{2}=0)|\psi\rangle=0.$ (29)

The direction of an arbitrary edge is represented by a unit vector
$\hat{n}_{\phi}$ pointing along its normal direction outwards, as shown in
Fig. 2. In the following, we simply refer it to edge $\hat{n}_{\phi}$.
Eigenstates of the two zero modes at edge $\hat{n}_{\phi}$ is given by

$\displaystyle|\psi_{n}\rangle=c(e^{\eta_{1}x}-e^{\eta_{2}x})|\varphi_{n}\rangle,\
n=1,2$ (30) $\displaystyle|\varphi_{1}\rangle=(1\ 0)^{T}\otimes(1\ 0)^{T},\ \
|\varphi_{2}\rangle=(0\ 1)^{T}\otimes(0\ 1)^{T},$

where $\eta_{1/2}$ are the roots of
$\frac{t_{2}}{2}\eta^{2}-t_{2}\eta+(t_{2}-t_{1})=0$. We then obtain the
effective edge Hamiltonian by projecting the bulk Hamiltonian in Eq.(28) into
eigenspace spanned by basis $\\{|\psi_{1}\rangle,|\psi_{2}\rangle\\}^{T}$,
which reads

$H_{E}(\phi)=t_{2}k_{2}s_{z}+\delta ms_{y}-m\cos(\phi-\theta)s_{y},$ (31)

with $s$ being Pauli matrices acting in the zero mode basis. From the
effective edge Hamiltonian, we immediately obtain the edge gap
$\Delta(\phi)=\delta m-m\cos(\phi-\theta)$, whose size as well as sign depends
on edge orientation. Majorana corner states appear whenever gaps of adjacent
edges take opposite signs.

##  Appendix D D. Effects of disorder and impurity on level crossings

Figure 4: Fermi level crossings in a 2D system with disorder (upper panel) and
boundary impurities (lower panel). Boundary impurities are added uniformly on
the left edge. The system size is $40\times 40$, and $\bar{t}_{1}=1$,
$t_{2}=2$, $m=0.7$, $\delta m=0.2$, $\theta=\pi/4$. Figure 5: (a) Evolution of
energy spectrum for a finite 2D system with bulk disorder $\sigma_{t}$. (b)
and (c) Probability distributions of energy states closest to zero energy for
$\sigma_{t}=0.9$ and $\sigma_{t}=2.4$. Majorana corner states remain stable at
weak disorder, and couple with bulk states when disorder is strong enough. (d)
Evolution of energy spectrum with impurity strength $V_{m}$. The impurities
are added uniformly on the left edge. Bulk gap is not influenced. Boundary gap
of the left edge is expected to close and reopen with the variation of
$V_{m}$. Due to finite size effect, the edge gap doesn’t really close, but the
topology of left edge may change when impurity strength is strong enough.
Consequently, one of the Majorana corner states is transferred to neighboring
corner, as can be seen by comparing the cases in (e) with $V_{m}=0.7$ and (f)
with $V_{m}=1.2$. In all the figures, $\bar{t}_{1}=1$, $t_{2}=2$, $m=0.7$,
$\delta m=0.2$, $\theta=\pi/4$.

In this section, we investigate the stability of level-crossing points against
bulk disorder and boundary impurities, and study their influences on Majorana
corner states in 2D.

As we pointed out in the main text, level crossings are protected by fermion
parity conservation. A single crossing cannot disappear unless bulk or
boundary gap is closed. Therefore, if the disorder or impurity doesn’t close
the two gaps, we can expect the level crossings to persist.

First, let us consider bulk disorder of $t_{1}$ term in the 1D model, which
follows Gaussian distribution with mean value $\bar{t}_{1}$ and standard
deviation $\sigma_{t}$. We plotted the lowest energy level (nonnegative) in
$(\sigma_{t},\lambda)$ parameter space, as shown in Fig. 3(a). For weak
disorder (small $\sigma_{t}$), the level crossing points remain stable, as
verified by energy spectrum shown in Fig. 3(b). With the increase of bulk
disorder, more and more states move close to zero energy, indicating the
closure of bulk gap, as Fig. 3(c) demonstrates. In the latter case, there is
no longer any level crossing. Therefore, the disorder in $t_{1}$ term mainly
influences the bulk gap and is expected to close the gap when it is strong
enough.

In contrast to $t_{1}$ term, $m$ or $\delta m$ term would influence the
boundary gap. We consider impurities of strength $V_{m}$ at boundary site
$j=1$, which has the the same form as $m$ term, i.e.,
$V_{m}\Gamma_{1}^{T}\sigma_{y}\Gamma_{1}$. In Fig. 3(d), we find that with the
increase of $V_{m}$, the two level-crossing points at $\pm\lambda$ move
towards $\lambda=0$ and annihilate with each other where boundary gap closes.
From Fig. 3(e) and (f), we find that the bulk gap doesn’t change in this
process, but the boundary gap closes and reopens, accompanied by the
disappearance of level crossings.

Turning to 2D system, the topological invariant introduced in the main text is
determined from level crossings at high symmetry momenta, and hence relies on
translation symmetry. When disorder or impurities break translation symmetry,
we can no longer say that the crossing appears at 1D subsystem with $K=0$ or
$\pi$, but have to look at the 2D spectrum instead. We should emphasize that
$\eta_{a}$ is the sum of the crossings at $K=0$ and $K=\pi$, and doesn’t
necessarily equal the number of crossings that appear while a toroidal system
is deformed into a cylindrical one. This is because $K=\pi$ is not allowed
when a periodic system has an odd number of unit cells but is allowed under
anti-periodic boundary condition. When the numbers of unit cells along both
directions are even, we could safely say that the two numbers are equal to
each other. In Fig. 4, we show the energy spectrum of a $40\times 40$ lattice,
and Fermi level crossings indeed survive under weak bulk disorder and boundary
impurities that are added uniformly on the left edge.

Considering that the higher-order topology is intimately related to Fermi
level crossings, we could expect Majorana corner states to be robust under
these perturbations. In Fig. 5(a) and (d) we plotted the variations of energy
spectrum at open boundaries with bulk disorder of $t_{1}$ term, as well as
boundary impurities $V_{m}$ that are uniformly distributed on the left edge.
Indeed, Majorana corner states survive under weak disorder, and disappear when
disorder becomes so strong that the bulk gap closes, as shown in Fig. 5(b) and
(c). Adding impurities on one edge only influences the boundary spectrum. The
topology of this particular edge would change when the impurities are strong
enough, so is the topological difference between it and neighboring edges. As
a result, Majorana corner states do not disappear but may hop from one corner
to an adjacent one, as shown in Fig. 5(e) and (f).

Figure 6: Fermi Level crossings in trivial, Weyl superconductor (SC) and
higher-order topological superconductor (HOTS) phases. (a) $\phi=\pi$. Both
Majorana numbers take $-1$ in higher-order phase, suggesting the two Majorana
corner states (shown in the inset) sit at opposite corners. (b) $\phi=3\pi/4$.
Only $\mathcal{M}_{x}=-1$. The two Majorana corner states sit at adjacent
corners. Note that level-crossings only occur at $K=0$ with given parameters.
In all the figures, $t=\alpha=\Delta=2\Delta_{Z}=1$, $\mu=0$, $\theta=\pi/3$.

##  Appendix E E. Rashba Bilayer

In this section, we apply the topological invariant proposed in the main text
to a Rashba bilayer superconducting system that is known to support higher-
order phases [39]. There are three key ingredients that make it a higher-order
superconductor: Rashba spin-orbit coupling, in-plane Zeeman field and phase
difference of $s$-wave pairing between the two layers. The model Hamiltonian
in $k$-space is given by

$\displaystyle H(\bm{k})=$ $\displaystyle[\mu+4t-2t(\cos k_{x}+\cos
k_{y})]\tau_{z}$ (32) $\displaystyle+$ $\displaystyle\alpha(\sin
k_{y}s_{x}-\sin k_{x}\tau_{z}s_{y})\sigma_{z}$ $\displaystyle+$
$\displaystyle\frac{\Delta}{2}\tau_{y}s_{y}[(\sigma_{0}+\sigma_{z})+e^{i\phi\tau_{z}}(\sigma_{0}-\sigma_{z})]$
$\displaystyle+$
$\displaystyle\Gamma\tau_{z}\sigma_{x}+\Delta_{Z}(\cos\theta\tau_{z}s_{x}+\sin\theta
s_{y}),$

where $\tau$, $s$ and $\sigma$ are Pauli matrices that act in Nambu, spin and
layer space respectively. In this model, $\mu$ represents chemical potential,
$t$ is amplitude of nearest-neighboring hopping, $\alpha$ represent the
strength of Rashba spin-orbit coupling, $\Delta$ is the pairing of upper
layer, $\phi$ is phase difference between the two layers, $\Gamma$ denotes the
coupling between two layers, and $\Delta_{Z}$ is in-plane Zeeman field with
$\theta$ being its direction.

For $\phi=\pi$ and $\Delta_{Z}<\Delta$, the model realizes higher-order phase
with two Majorana corner states at opposite corner when
$\Gamma>|\Delta+\Delta_{Z}|$ and is in trivial phase when
$\Gamma<|\Delta-\Delta_{Z}|$. In between the two phases, the system becomes a
Weyl superconductor. Indeed, Fermi level crossings (real root
$\lambda\in(0,1)$) appear in the higher-order phase, as can be seen in Fig.
6(a). In the Weyl superconductor phase, level crossings appear exactly at
$\lambda=0$. It should be noted that the higher-order phase persists for a
wide range of $\phi$, in which case Majorana corner states may reside at
neighboring corners instead of opposite corners, as shown in the inset of Fig.
6(b).


# Symmetry dictated universal helicity redistribution of Dirac fermions in
transport

Jun-Yin Huang  Lanzhou Center for Theoretical Physics, Key Laboratory of
Theoretical Physics of Gansu Province, and Key Laboratory of Quantum Theory
and Applications of MoE, Lanzhou University, Lanzhou, Gansu 730000, China
Rui-Hua Ni  Lanzhou Center for Theoretical Physics, Key Laboratory of
Theoretical Physics of Gansu Province, and Key Laboratory of Quantum Theory
and Applications of MoE, Lanzhou University, Lanzhou, Gansu 730000, China
Hong-Ya Xu  Corresponding author<EMAIL_ADDRESS>Lanzhou Center for
Theoretical Physics, Key Laboratory of Theoretical Physics of Gansu Province,
and Key Laboratory of Quantum Theory and Applications of MoE, Lanzhou
University, Lanzhou, Gansu 730000, China  Liang Huang  Corresponding author:
<EMAIL_ADDRESS>Lanzhou Center for Theoretical Physics, Key Laboratory of
Theoretical Physics of Gansu Province, and Key Laboratory of Quantum Theory
and Applications of MoE, Lanzhou University, Lanzhou, Gansu 730000, China

###### Abstract

Helicity is a fundamental property of Dirac fermions. Yet, the general rule of
how it changes in transport is still lacking. We uncover, theoretically, the
universal spinor state transformation and consequently helicity redistribution
rule in two cases of transport through potentials of electrostatic and mass
types, respectively. The former is dictated by Lorentz boost and its complex
counterpart in Klein tunneling regime, which establishes miraculously a
unified yet latent connection between helicity, Klein tunneling, and Lorentz
boost. The latter is governed by an abstract rotation group we construct,
which reduces to SO(2) when acting on the plane of effective mass and
momentum. They generate invariant submanifolds, i.e., leaves, that foliate the
Hilbert space of Dirac spinors. Our results provide a basis for unified
understanding of helicity transport, and may open a new window for exotic
helicity-based physics and applications in mesoscopic systems.

Introduction.—Helicity, the projection of the spin onto the direction of
momentum, is known as an intrinsic and measurable property of Dirac fermions
in relativistic quantum mechanics, and plays a crucial role in understanding
the nature of fundamental particles. For instance, the helicity nature of
massless Dirac particles turns out to be responsible for a highly anisotropic
tunneling, i.e., chiral tunneling, through an electrostatic potential barrier
[1, 2, 3, 4, 5]. Flip and conservation of helicity has been recognized as a
key issue in addressing scattering of massive Dirac particles by a magnetic
monopole [6, 7] and an Aharonov-Bohm potential [8, 9, 10]. As an independent
and exploitable degree of freedom, it attracts attention for electronic
applications in mesoscopic systems recently [11, 12]. However, before stepping
forward further, it is crucial to understand the redistribution rules of
helicity in typical scattering, tunneling, and transport processes.

Generally, the transport process can be abstracted as an operation on the
Hilbert space of Dirac spinors mapping the initial to the final spinor states
characterized by energy and momentum. Yet, the helicity degree of freedom is
undetermined as the coefficients of the two helicity components can still take
different values, which typically depends on the transport details. One
fundamental question is then under commonly encountered transport situations,
how does the helicity change?

By utilizing the (3+1)-D Dirac equation formalism in the single-particle
framework, we consider two typical transport processes, i.e., through
piecewise-constant electrostatic [13, 14, 15, 16, 17] or mass potentials [18,
19, 20, 21, 22, 23, 24]. For each process, as the potential of the concerned
region changes, the final spinor evolves out a curve started from the given
initial state. We find, strikingly, that the final spinor depends only on the
potential height of this region, i.e., process-independent. Thus the
corresponding abstract operation forms a one-parameter transformation group,
whose exertion on any initial state generates a one-dimensional invariant
submanifold, which forms leaves foliating the Hilbert space of Dirac spinors.

Figure 1:  Schematics of the two transport processes and the rule of spinor
state transformations. (a) Transport through a piecewise-constant
electrostatic potential, with momentum $\mathbf{p}_{j}$ and relative helicity
components $|t_{j}^{(\pm)}|^{2}$ (red for positive, blue for negative) in each
region $j$ demonstrated. $\ket{\psi_{i}}$ is the initial (incident) spinor
state with completely positive helicity, and $\ket{\psi_{f}}$ is the final
state in region $2$. (b,c) Visualizing the Lorentz boost $\Lambda(w)$ dictated
transformation rule we uncovered from $\ket{\psi_{i}}$ (hollow circle) to
$\ket{\psi_{f}}$s. Projection of the leaf consisted of all $\ket{\psi_{f}}$s
onto subspaces parameterized by effective energy $E-V_{2}$, momentum $p_{x}$
and $p_{{}_{\bot}}$ (b), and by $p_{x}/p$ and helicity polarization $P$ (c).
Solid (dashed) curves are for transmitted (reflected) parts. Red and blue
indicate positive and negative effective energy branches, respectively, which
are connected by $\Lambda(\mathrm{i}\pi)$. (d-f) Mass potential case. The
transformation rule is now governed by an abstract rotation group
$\Gamma(\mu)$ we identified. (e) Projection of the leaf in the ($m+U_{2}$,
$p_{x}$, $p_{{}_{\bot}}$) subspace, and (f) on the ($p_{x}/p$, $P$) plane. Red
and yellow curves are for positive and negative effective mass $m+U_{2}$,
respectively.

Our main analytical results are shown in Fig. 1. The left panels plot the
schematics of the two processes, where the potentials are along the $x$-axis,
and the perpendicular momentum component $p_{{}_{\bot}}$ is conserved. The
final state $\ket{\psi_{f}}$ is chosen in the middle region to simply include
the reflected wavefunction. The middle and right panels show the foliation
structures. For transport through an electrostatic potential with a given
initial energy $E$ [Fig. 1(a)], the Hilbert space of spinors can be fully
parameterized by the effective energy $E-V_{2}$, momentum $p_{x}$ and
$p_{{}_{\bot}}$, and helicity polarization $P$ [25]. As the potential $V_{2}$
changes, from any initial state $\ket{\psi_{i}}$, a one-dimensional leaf grows
out in the Hilbert space, i.e., a curve in the four dimensional parameter
space. Particularly, the projection of this leaf into the ($E-V_{2}$, $p_{x}$,
$p_{{}_{\bot}}$) subspace is a hyperboloid characterized by $E$ [Fig. 1(b)].
Different $p_{{}_{\bot}}$, as determined by the incident angle, corresponds to
different sections, i.e., hyperbolic curves. The projection onto the
($p_{x}/p$, $P$) plane is shown in Fig. 1(c), where
$p=\sqrt{p_{x}^{2}+p_{{}_{\bot}}^{2}}$. Surprisingly, the one-parameter
transformation group for the spinors and consequently the redistribution of
helicity is just the spinor representation of Lorentz boost $\Lambda(w)$ along
the $x$-axis. When $V_{2}$ is large enough, the branch for negative energy
states ($E-V_{2}<0$) associated with Klein tunneling emerges, which can be
mapped from the original branch by a complex Lorentz boost
$\Lambda(\mathrm{i}\pi)$. This builds miraculously a latent, unified
connection between different essential properties, i.e., helicity, Klein
tunneling, and Lorentz boost, of Dirac fermions. The universal redistribution
rule of the helicity is corroborated by extensive numerics in more general
cases.

For the mass potential case with initial energy $E$ [Fig. 1(d)], the effective
mass $m+U_{2}$ plays a similar role of the effective energy. The leaf can then
be projected onto the ($p_{x}$, $m+U_{2}$) plane, which, in the natural units
($c=1$), is a circle. All these circles for possible $p_{{}_{\bot}}$s
construct a sphere in the ($m+U_{2}$, $p_{x}$, $p_{{}_{\bot}}$) space [Fig.
1(e)]. The projection on the ($p_{x}$, $P$) plane is shown in Fig. 1(f). As
$U_{2}$ varies, the one-parameter transformation group $\Gamma(\mu)$ acting on
the ($p_{x}$, $m+U_{2}$) plane can be conveniently described by an SO(2)
rotation. Its spinor representation for the abstract spinor state
transformation has also been constructed.

These findings uncover universal redistribution rules of helicity, dictated by
symmetries described by Lorentz boost and an abstract rotation group,
respectively, in transport processes through electrostatic or mass potentials.
Our results provide conveniently a modulation strategy of the helicity by
simply configuring the electrostatic or mass potential profiles, holding great
promise to exploit the new degree of freedom for novel applications.

Model.—The transport problem of the (3+1)-D massive Dirac fermions can be
described by the Dirac Hamiltonian in the position representation [26]

$\hat{H}=\bm{\alpha}\cdot\hat{\mathbf{p}}+\beta[m+U(x)]+V(x),$ (1)

where $\bm{\alpha}=\begin{pmatrix}0&\bm{\sigma}\\\ \bm{\sigma}&0\end{pmatrix}$
with Pauli matrices $\bm{\sigma}=(\sigma_{x},\sigma_{y},\sigma_{z})$,
$\hat{\mathbf{p}}$ is the 3-momentum operator,
$\beta=\begin{pmatrix}\mathds{1}_{2}&0\\\ 0&-\mathds{1}_{2}\end{pmatrix}$, $m$
is the Dirac mass, $U$ is the mass potential, and $V$ is the electrostatic
potential. Natural units $c=\hbar=e=1$ are employed. For the sake of
simplicity, we consider two typical cases, i) $U(x)=0$, and ii) $V(x)=0$. We
exploit a piecewise-constant potential with $n$ regions of width $\\{d_{j}\\}$
along the $x$-axis, where $d_{1}=d_{n}=\infty$, and $V_{1}=U_{1}=0$, which, in
the large $n$ limit, can approach arbitrary continuous potential profile.

Note that the whole Hamiltonian $\hat{H}$ does not commute with $\hat{h}$.
Nevertheless, in each region $j$, the corresponding Hamiltonian $\hat{H}_{j}$
commutes with the helicity operator [26]

$\hat{h}=\begin{pmatrix}\bm{\sigma}&0\\\
0&\bm{\sigma}\end{pmatrix}\cdot\hat{\mathbf{p}}/p_{j},$ (2)

where $p_{j}=\sqrt{(E-V_{j})^{2}-(m+U_{j})^{2}}$, and the incident energy
$E>m$. The transmitted wavefunctions (chosen as the common eigenstates of
$\hat{H}$, $\hat{\mathbf{p}}$, and $\hat{h}$) are plane wave solutions [27]

$\displaystyle\psi_{j}^{(h)}=\begin{pmatrix}\chi^{(h)}_{j}\\\
hk_{j}\chi^{(h)}_{j}\end{pmatrix}e^{\mathrm{i}(\mathbf{p}_{j}\cdot\mathbf{x}-Et)},$
(3)

where $k_{j}=p_{j}/(E-V_{j}+m+U_{j})$, $h=\pm 1$ is the eigenvalue of
$\hat{h}$, and
$\chi^{(+)}_{j}=[\cos\frac{\theta_{j}}{2}e^{-\mathrm{i}\varphi_{j}/2},\sin\frac{\theta_{j}}{2}e^{\mathrm{i}\varphi_{j}/2}]^{T}$,
$\chi^{(-)}_{j}=[-\sin\frac{\theta_{j}}{2}e^{-\mathrm{i}\varphi_{j}/2},\cos\frac{\theta_{j}}{2}e^{\mathrm{i}\varphi_{j}/2}]^{T}$
[27], with
$\mathbf{p}_{j}=(p_{j,x},p_{y},p_{z})=p_{j}(\sin\theta_{j}\cos\varphi_{j},\sin\theta_{j}\sin\varphi_{j},\cos\theta_{j})$
as $p_{y}$ and $p_{z}$ are conserved in different regions [13, 28],
$\theta_{j}=\cos^{-1}(p_{1}\cos\theta_{1}/p_{j})$,
$\varphi_{j}=\cos^{-1}(\lambda_{j}\sqrt{p_{j}^{2}-p_{y}^{2}-p_{z}^{2}}\big{/}p_{j}\sin\theta_{j})$,
and $\lambda_{j}={\rm sgn}(E-V_{j})$ is for positive or negative energy
states. Due to the rotational symmetry with respect to the $x$-axis, the angle
$(\theta_{j},\varphi_{j})$ can be characterized by a single angle
$\varrho_{j}=\cos^{-1}(\sin\theta_{j}\cos\varphi_{j})$ between
$\mathbf{p}_{j}$ and the positive $x$-axis. For the reflected plane wave
$\overline{\psi}_{j}^{(h)}$, one only needs to change $p_{j,x}$ to $-p_{j,x}$
and $\varphi_{j}$ to $\pi-\varphi_{j}$. Without loss of generality, we assume
that the incident Dirac fermions are completely positively helicity-polarized
plane waves.

The wavefunction in each region can be decomposed into the eigenmodes in terms
of Eq. (3) and its reflected counterpart, i.e.,
$\psi_{j}=t_{j}^{(h)}\psi_{j}^{(h)}+r_{j}^{(h)}\overline{\psi}_{j}^{(h)}$. For
transmitted flow, the helicity polarization is defined as [29, 30]

$\displaystyle
P_{j}=(J_{j}^{(+)}-J_{j}^{(-)})\big{/}(J_{j}^{(+)}+J_{j}^{(-)}),$ (4)

where
$J_{j}^{(h)}=|t_{j}^{(h)}|^{2}|\psi_{j}^{(h){\dagger}}\bm{\alpha}\psi_{j}^{(h)}|$
is the magnitude of the probability current. For reflected flow,
$\overline{J}_{j}^{(h)}$ and $\overline{P}_{j}$ can be defined similarly. The
transmission and reflection coefficients ($t_{j}^{(h)},r_{j}^{(h)}$) are then
numerically calculated via the transfer matrix method by matching them at the
boundary of each region [31, 32, 33]. The overall transmission probability is
then

$T^{(h)}=|t_{n}^{(h)}|^{2}k_{n}\sin\theta_{n}\cos\varphi_{n}\big{/}k_{1}\sin\theta_{1}\cos\varphi_{1}.$

As a test, by approximating a smooth linear potential with large $n$, our
numerical calculation matches perfectly with the analytical result [34]
[Supplemental Material (SM) Sec. I].

Theory.—We derive the redistribution rules of helicity in transport as
follows.

Case i), only consider the electrostatic potential, i.e., for any $j\in[1,n]$,
$U_{j}=0$. We have

$\displaystyle\hat{H}_{j}\psi_{j}^{(h)}=\lambda_{j}|E-V_{j}|\psi_{j}^{(h)},\
\hat{h}\psi_{j}^{(h)}=h\psi_{j}^{(h)}.$ (5)

For $V_{j}\in(V_{-},V_{+})$ with
$V_{\pm}=E\pm\sqrt{p_{y}^{2}+p_{z}^{2}+m^{2}}$, the transmitted part has an
imaginary $p_{j,x}$, which is evanescent and no longer has a well-defined
helicity. Therefore we only consider $V_{j}<V_{-}$ scattering or $V_{j}>V_{+}$
Klein tunneling.

For $n=3$, the normalized transmission and reflection coefficients
$t_{j}^{(h)}$ ($j=2,3$) and $r_{j}^{(h)}$ ($j=1,2$) can be obtained
analytically (SM Sec. II), e.g.,

$\displaystyle\begin{pmatrix}t_{j}^{(+)}\\\ t_{j}^{(-)}\\\
\end{pmatrix}=|\mathcal{M}|^{-1/2}\begin{pmatrix}\mathcal{M}_{(+,+)}&\mathcal{M}_{(+,-)}\\\
\mathcal{M}_{(-,+)}&\mathcal{M}_{(-,-)}\\\
\end{pmatrix}\cdot\begin{pmatrix}t_{1}^{(+)}\\\ t_{1}^{(-)}\\\ \end{pmatrix},$
(6)

where
$\mathcal{M}_{(h^{\prime},h)}\equiv\bra{\psi_{j}^{(h^{\prime})}}\alpha_{x}\ket{\psi_{1}^{(h)}}$
with normalized ket $\ket{\psi_{j}^{(h)}}=\psi_{j}^{(h)}/|\psi_{j}^{(h)}|$.
For positively helicity-polarized initial state, $t_{1}^{(+)}=1$,
$t_{1}^{(-)}=0$. Plugging into Eq. (4), the exact expression of the helicity
polarization for the transmitted flow is

$\displaystyle P_{j}=(\eta_{j}-1)\big{/}(\eta_{j}+1),$ (7)

where

$\displaystyle\eta_{j}=\left|t_{j}^{(+)}\big{/}t_{j}^{(-)}\right|^{2}=\left(\frac{k_{1}+k_{j}}{k_{1}-k_{j}}\right)^{2}\cdot\frac{1+\cos(\varrho_{1}+\varrho_{j})}{1-\cos(\varrho_{1}+\varrho_{j})}.$
(8)

For reflected flow, one only needs to change $\varrho_{j}$ to
$\pi-\varrho_{j}$ to get $\overline{P}_{j}$.

Remarkably, the spinor state $\psi_{j}^{(h)}$ in region $j$ is determined
completely by $V_{j}$ and $(E,\mathbf{p}_{1})$ [Eq. (3)], and is irrelevant to
the potential in other regions. So do the matrix $\mathcal{M}$ and the
helicity polarization $P_{j}$. Regarding the transformation from
$\ket{\psi_{1}}$ to $\ket{\psi_{j}}$ as an abstract operation, it is additive
and forms a one-parameter ($V_{j}$) transformation group (SM Sec. III-A),
generating leaves in the spinor Hilbert space.

Alternatively, disregard the actual transport process, only consider the two
sets of energy and momentum of the initial and final states. Now imagine that
they are for the same state but observed in two different inertial reference
frames, i.e., $O$ and $O^{\prime}$ in Fig. 2, and are thus connected by a
Lorentz boost $\Lambda(w)$ with rapidity $w$ (defined by $\cosh
w=(1-v^{2})^{-1/2}$ with $v$ being the relative velocity of the two frames),

$\displaystyle\begin{pmatrix}E-V_{j}\\\
p_{j,x}\end{pmatrix}=\begin{pmatrix}\cosh w&-\sinh w\\\ -\sinh w&\cosh w\\\
\end{pmatrix}\cdot\begin{pmatrix}E\\\ p_{1,x}\end{pmatrix}.$ (9)

In this regard, the spinor state will be transformed to
$\hat{S}\ket{\psi_{1}}$ with

$\displaystyle\hat{S}[\Lambda(w)]=\cosh(w/2)\mathds{1}_{4}-\sinh(w/2)\alpha_{x},$
(10)

being the spinor representation of the Lorentz boost [26]. Surprisingly, we
find $\hat{S}\ket{\psi_{1}}$ from the imagined inertial frame transformation
is the same as the final state $\ket{\psi_{j}}$ in the actual transport
problem. As such, the one-parameter transformation group characterizing the
transport process is just the Lorentz boost. Retrospectively, as the
electrostatic potential only shifts the energy of the Dirac fermion, the
equation structure is preserved, and the corresponding abstract operation is
equivalent to the inertial frame transformation (Lorentz boost). This simple,
universal behavior unveils yet another concealed attribute of the miraculous
Dirac equation.

Figure 2:  The change of effective energy, momentum, and helicity from
$\ket{\psi_{1}}$ to $\ket{\psi_{j}}$ in the actual transport problem can be
described by a Lorentz boost $\Lambda(w)$ with real rapidity $w$ of the
imagined inertial frame transformation from $O$ to $O^{\prime}$. The Klein
tunneling can be described by the combined operation
$\Lambda(\mathrm{i}\pi)\Lambda(w)$, i.e.,
$O\xrightarrow{\Lambda(w)}O^{\prime}\xrightarrow{\Lambda(\mathrm{i}\pi)}\widetilde{O}$.

For Klein tunneling ($V_{j}>V_{+}$), we will take the convention $p_{j,x}<0$,
so that the $x$-component of the group velocity
${v}_{x}=\lambda_{j}p_{j,x}/|E-V_{j}|$ is always positive for transmitted part
[35, 36, 37]. Conventional Lorentz boost fails to link the positive and
negative energy states. By replacing $w$ with $\mathrm{i}\pi$, the operation
$\Lambda(\mathrm{i}\pi)$ flips the sign of effective energy and $p_{x}$, as
from $O^{\prime}$ to $\widetilde{O}$ in Fig. 2.
$\hat{S}[\Lambda(\mathrm{i}\pi)]=-\mathrm{i}\alpha_{x}$ flips the helicity
polarization, exactly the same as that in the actual tunneling process. The
combined operation from $O$ to $\widetilde{O}$ is then
$\Lambda(\mathrm{i}\pi)\Lambda(w)$ (SM Sec. III-B), which is traceable [38,
39, 40, 41, 42]. The complex Lorentz boost $\Lambda(\mathrm{i}\pi)$ enforces
an antisymmetry between the helicity polarization and the effective energy
$E-V$ [Fig. 3(a)].

Indeed, for Dirac spinors, introducing $\hat{S}[\Lambda(\mathrm{i}\pi)]$ is
required by the $\mathcal{PCT}$ symmetry. The $\mathcal{PCT}$ symmetry
interchanges the positive and negative energy states, with
$\mathcal{P}=\mathcal{P}_{x}\mathcal{P}_{y}\mathcal{P}_{z}$, $\mathcal{C}$,
$\mathcal{T}$ being the parity, charge conjugate, and time-reversal
operations. Due to the rotational symmetry with respect to the $x$-axis,
$\mathcal{P}_{y}\mathcal{P}_{z}=\mathds{1}$, the $\mathcal{PCT}$ symmetry
reduces to $\mathcal{P}_{x}\mathcal{CT}$, which equals to
$\hat{S}[\Lambda(\mathrm{i}\pi)]$ disregarding an unobservable phase (SM Sec.
III-B).

The above results do not explicitly depend on the fact of $n=3$, and are valid
for arbitrary $n$ regions.

Case ii), the mass potential. Requiring $p_{j,x}$ to be real leads to
$U_{j}\in(U_{-},U_{+})$ with $U_{\pm}=-m\pm\sqrt{E^{2}-p_{y}^{2}-p_{z}^{2}}$.
Note that the negative mass potential has been widely used in models of
nuclear physics and condensed matter physics [21, 18, 19, 20, 22, 23].
Similarly, Eqs. (6-8) have the same form except that now $\psi_{j}^{(h)}$,
$k_{j}$, and $\varrho_{j}$ are determined by $U_{j}$. The corresponding map
from $\ket{\psi_{1}}$ to $\ket{\psi_{j}}$ also forms a one-parameter
transformation group, denoted as $\Gamma(\mu)$ (SM Sec. IV). In the plane of
effective mass $m+U_{j}$ and $x$-momentum, the action of $\Gamma(\mu)$ reduces
to SO(2) rotation

$\displaystyle\begin{pmatrix}m+U_{j}\\\
p_{j,x}\end{pmatrix}=\begin{pmatrix}\cos\mu&-\sin\mu\\\ \sin\mu&\cos\mu\\\
\end{pmatrix}\cdot\begin{pmatrix}m\\\ p_{1,x}\end{pmatrix},$ (11)

where $\mu$ is an angle parameter determined by $U_{j}$. The transformation of
the spinor states from $\ket{\psi_{1}}$ to $\ket{\psi_{j}}$ in the actual
transport process can be constructed as
$\ket{\psi_{j}}=\hat{S}_{m}\ket{\psi_{1}}$, with

$\displaystyle\hat{S}_{m}[\Gamma(\mu)]=\cos(\mu/2)\mathds{1}_{4}-\sin(\mu/2)\beta\alpha_{x}$
(12)

being the spinor representation of $\Gamma(\mu)$. Thus the redistribution rule
of the helicity is derived, and is completely different from that due to
electrostatic potentials. To the best of our knowledge, this group
$\Gamma(\mu)$ and its spinor representation are unknown before.

Figure 3:  Plots of helicity polarization for transmitted ($P_{j}$, solid
curves) and reflected ($\overline{P}_{j}$, dashed) waves. (a) vs effective
energy $E-V_{j}$ in units of $E$ for electrostatic potential, and (b) vs
effective mass $m+U_{j}$ in units of $m$ for mass potential. The incident
energy and angle: $E=1.08m$, $\varrho_{1}=2\pi/5$. The curves are from Eqs.
(7-8). The circles are numerical results from the random piecewise-constant
potential model with $n=10$ and $j$ a randomly chosen region, e.g., $j=5$.

Numerical verification.—To be as general as possible, we consider a piecewise-
constant potential model with $n=10$ and random widths and potential heights
for each region: $d_{j}\in(0,1/m)$, $V_{j}$ or $U_{j}\in(-10m,10m)$. We then
systematically vary the potential value in a randomly chosen region, say,
$j=5$, and plot the helicity polarization of the transmitted and reflected
waves in Fig. 3 as symbols. They agree with the theory well for both of the
electrostatic and mass potentials. The numerical simulation corroborates
clearly that the helicity redistribution only depends on the potential in this
region, but is independent of the transport processes before or after this
region.

Figure 4:  Density maps for helicity polarization in the output region
($j=n=2$) through step potentials. (a) vs incident energy $E$ and
electrostatic potential $V_{2}$, and (b) vs $E$ and mass potential $U_{2}$,
all in units of $m$. The incident angle $\varrho_{1}=\pi/3$. The helicity
polarization of the incident waves is $P=1$ for (a) and $P=0$ ($t_{1}^{(+)}=1$
and $t_{1}^{(-)}=\mathrm{i}$) for (b).

Potential applications.—Figure 4 demonstrates the variation of the helicity
polarization $P$ vs the key parameters. For electrostatic potential, $P$ can
be reversed from $1$ (incident) to $-1$ (output). For mass potential, with a
non-polarized incident wave of $P=0$, the output can be modulated to either
$1$ or $-1$. That is, complete polarization of helicity can be generated from
completely non-polarized incident waves by applying a properly adjusted mass
potential. Due to the inherent connection to spin polarization, the helicity,
with the unveiled redistribution rule, may lead to appealing applications for
spin-based electronics.

Discussion and conclusion.— Another transport setup is through magnetic
barriers. In this case, the helicity is conserved [43, 44, 45, 46] as it
commutes with the (3+1)-D Dirac Hamiltonian in a static magnetic field [47, 8,
48, 49]. Thus the redistribution rule is trivial.

Helicity is one of the most fundamental and intriguing properties of Dirac
fermions. Due to the complexity of its combined nature of spin and momentum, a
unified understanding of how it changes in transport is still lacking. We
tackle this issue and have uncovered the universal helicity redistribution
rule, from which process-independent helicity polarizations arise. The
underlying principle of the universality is unravelled by abstracting the
physical transport of helicity into imagined transforming operations on the
spinor states. In such way, we identify one-parameter transformation groups as
well as their spinor representations, which are described by the Lorentz boost
$\Lambda(w)$ for electrostatic potentials, and the abstract rotation group
$\Gamma(\mu)$ for mass potentials. This explicitly reveals an unknown
connection between the physical rule of helicity in transport that is highly
nontrivial and the underlying symmetries of the Dirac equation. Our results
not only are fundamental to the understanding of relativistic quantum
scattering, tunneling, and transport, but also offer promises of generating
desired helicity polarization for exotic applications in mesoscopic systems
[12, 50, 51, 52].

###### Acknowledgements.

Acknowledgments.— We thank Profs. Yong-Shi Wu, Zhong-Zhou Ren, Li-Sheng Geng
for insightful discussions. This work was supported by NSFC under Grants No.
12175090, No. 12105125, No. 11775101, and No. 12247101, and by the 111 Project
under Grant No. B20063.

## References

  * Katsnelson  _et al._ [2006] M. Katsnelson, K. Novoselov, and A. Geim, Chiral tunnelling and the Klein paradox in graphene, Nat. Phys. 2, 620 (2006).
  * Zeb  _et al._ [2008] M. A. Zeb, K. Sabeeh, and M. Tahir, Chiral tunneling through a time-periodic potential in monolayer graphene, Phys. Rev. B 78, 165420 (2008).
  * Tudorovskiy  _et al._ [2012] T. Tudorovskiy, K. J. A. Reijnders, and M. I. Katsnelson, Chiral tunneling in single-layer and bilayer graphene, Phys. Scr. T146, 014010 (2012).
  * He  _et al._ [2013] W.-Y. He, Z.-D. Chu, and L. He, Chiral tunneling in a twisted graphene bilayer, Phys. Rev. Lett. 111, 066803 (2013).
  * Habib  _et al._ [2015] K. M. M. Habib, R. N. Sajjad, and A. W. Ghosh, Chiral tunneling of topological states: Towards the efficient generation of spin current using spin-momentum locking, Phys. Rev. Lett. 114, 176801 (2015).
  * Zwanziger [1968] D. Zwanziger, Exactly soluble nonrelativistic model of particles with both electric and magnetic charges, Phys. Rev. 176, 1480 (1968).
  * Kazama  _et al._ [1977] Y. Kazama, C. N. Yang, and A. S. Goldhaber, Scattering of a Dirac particle with charge Ze by a fixed magnetic monopole, Phys. Rev. D 15, 2287 (1977).
  * Vera and Schmidt [1990] F. Vera and I. Schmidt, Helicity and the Aharonov-Bohm effect, Phys. Rev. D 42, 3591 (1990).
  * Coutinho and Perez [1994] F. A. B. Coutinho and J. F. Perez, Helicity conservation in the Aharonov-Bohm scattering of Dirac particles, Phys. Rev. D 49, 2092 (1994).
  * Araujo  _et al._ [2001] V. S. Araujo, F. A. B. Coutinho, and J. F. Perez, On the most general boundary conditions for the Aharonov-Bohm scattering of a Dirac particle: helicity and Aharonov-Bohm symmetry conservation, J. Phys. A 34, 8859 (2001).
  * Gu  _et al._ [2011] N. Gu, M. Rudner, and L. Levitov, Chirality-assisted electronic cloaking of confined states in bilayer graphene, Phys. Rev. Lett. 107, 156603 (2011).
  * Zhao  _et al._ [2015] L. Zhao, J. Wang, J. Liu, Y. Xu, B.-L. Gu, Q.-K. Xue, and W. Duan, Electronic analog of chiral metamaterial: Helicity-resolved filtering and focusing of Dirac fermions in thin films of topological materials, Phys. Rev. B 92, 041408(R) (2015).
  * De Leo and Rotelli [2009] S. De Leo and P. Rotelli, Planar Dirac diffusion, Eur. Phys. J. C 63, 157 (2009).
  * Setare and Jahani [2010] M. Setare and D. Jahani, Klein tunneling of massive Dirac fermions in single-layer graphene, Phys. B: Condens. Matter 405, 1433 (2010).
  * De Leo and Rotelli [2012] S. De Leo and P. Rotelli, Relative helicity phases in planar Dirac scattering, Phys. Rev. A 86, 032113 (2012).
  * Bittencourt  _et al._ [2015] V. A. S. V. Bittencourt, S. S. Mizrahi, and A. E. Bernardini, SU (2) $\otimes$ SU (2) bi-spinor structure entanglement induced by a step potential barrier scattering in two-dimensions, Ann. Phys. 355, 35 (2015).
  * Navarro-Giraldo and Quimbay [2020] J. Navarro-Giraldo and C. Quimbay, Two-dimensional Klein tunneling for massive Dirac fermions with a defined helicity, Ann. Phys. 412, 168022 (2020).
  * Jackiw [2012] R. Jackiw, Fractional and Majorana fermions: the physics of zero-energy modes, Phys. Scripta 2012, 014005 (2012).
  * Hunt  _et al._ [2013] B. Hunt, J. D. Sanchez-Yamagishi, A. F. Young, M. Yankowitz, B. J. LeRoy, K. Watanabe, T. Taniguchi, P. Moon, M. Koshino, P. Jarillo-Herrero, and R. C. Ashoori, Massive Dirac fermions and Hofstadter butterfly in a van der Waals heterostructure, Science 340, 1427 (2013).
  * Shen [2017] S.-Q. Shen,  _Topological Insulators: Dirac Equation in Condensed Matter_ , Springer Series in Solid-State Sciences (Springer Singapore, 2017).
  * Zhao  _et al._ [2022] Q. Zhao, Z. Ren, P. Zhao, and J. Meng, Covariant density functional theory with localized exchange terms, Phys. Rev. C 106, 034315 (2022).
  * Cheng  _et al._ [2022] H. Cheng, J. Yang, Z. Wang, and L. Lu, Monopole topological resonators (2022), arXiv:2210.09056 [cond-mat.mes-hall] .
  * Wang  _et al._ [2022] M. Wang, Q. Ma, S. Liu, R.-Y. Zhang, L. Zhang, M. Ke, Z. Liu, and C. T. Chan, Observation of boundary induced chiral anomaly bulk states and their transport properties, Nat. Commun. 13, 5916 (2022).
  * Chen  _et al._ [2023] K. Chen, F. Komissarenko, D. Smirnova, A. Vakulenko, S. Kiriushechkina, I. Volkovskaya, S. Guddala, V. Menon, A. Alù, and A. B. Khanikaev, Photonic Dirac cavities with spatially varying mass term, Sci. Adv. 9, eabq4243 (2023).
  * [25] The relative phase of the two spinor components is irrelevant for the helicity polarization.
  * Greiner [2000] W. Greiner,  _Relativistic Quantum Mechanics: Wave Equations_ (Springer, Berlin, 2000).
  * Strange [1998] P. Strange,  _Relativistic Quantum Mechanics: With Applications in Condensed Matter and Atomic Physics_ (Cambridge University Press, Cambridge, 1998).
  * [28] Note that $(\theta,\varphi)$ are spherical coordinate parameters in general sense [37], different from the well addressed incident angles $(\alpha,\beta)$ in graphene [1, 53] and Weyl semimetal [54]. There is a simple relation between them: for $\lambda=1$, $(\theta,\varphi)=(\beta,\alpha)$; and for $\lambda=-1$, $(\theta,\varphi)=(\pi-\beta,\alpha-\pi)$, i.e., $\mathbf{p}=\lambda p(\sin\beta\cos\alpha,\sin\beta\sin\alpha,\cos\beta)$ [55].
  * Bandyopadhyay and Cahay [2008] S. Bandyopadhyay and M. Cahay,  _Introduction to Spintronics_ (CRC press, Boca Raton, FL, 2008).
  * Eschrig [2015] M. Eschrig, Spin-polarized supercurrents for spintronics: a review of current progress, Rep. Prog. Phys. 78, 104501 (2015).
  * Davies [1998] J. H. Davies,  _The Physics of Low-Dimensional Semiconductors: An Introduction_ (Cambridge University Press, Cambridge, 1998).
  * Markos and Soukoulis [2008] P. Markos and C. M. Soukoulis,  _Wave Propagation_ (Princeton University Press, NJ, 2008).
  * Cotăescu  _et al._ [2007] I. I. Cotăescu, P. Gravila, and M. Paulescu, Applying the Dirac equation to derive the transfer matrix for piecewise constant potentials, Phys. Lett. A 366, 363 (2007).
  * Sauter [1931] F. Sauter, Über das Verhalten eines Elektrons im homogenen elektrischen Feld nach der relativistischen Theorie Diracs, Z. Phys. 69, 742 (1931).
  * Klein and Nishina [1929] O. Klein and Y. Nishina, Über die Streuung von Strahlung durch freie Elektronen nach der neuen relativistischen Quantendynamik von Dirac, Z. Phys. 52, 853 (1929).
  * Castro Neto  _et al._ [2009] A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S. Novoselov, and A. K. Geim, The electronic properties of graphene, Rev. Mod. Phys. 81, 109 (2009).
  * Nguyen and Nguyen [2009] H. C. Nguyen and V. L. Nguyen, Tunneling of Dirac electrons through one-dimensional potentials in graphene: a T-matrix approach, J. Phys.: Condens. Matter 21, 045305 (2009).
  * Jost [1957] R. Jost, Eine Bemerkung zum CTP theorem, Helv. Phys. Acta. 30, 153 (1957).
  * Bergknoff and Thacker [1979] H. Bergknoff and H. B. Thacker, Structure and solution of the massive Thirring model, Phys. Rev. D 19, 3666 (1979).
  * Thacker [1981] H. B. Thacker, Exact integrability in quantum field theory and statistical systems, Rev. Mod. Phys. 53, 253 (1981).
  * Froggatt and Nielsen [1991] C. D. Froggatt and H. B. Nielsen,  _Origin of Symmetries_ (World Scientific, Singapore, 1991) pp. 86–91.
  * Lehnert [2016] R. Lehnert, CPT symmetry and its violation, Symmetry 8, 114 (2016).
  * De Martino  _et al._ [2007] A. De Martino, L. Dell’Anna, and R. Egger, Magnetic confinement of massless Dirac fermions in graphene, Phys. Rev. Lett. 98, 066802 (2007).
  * Wu  _et al._ [2010] Z. Wu, F. M. Peeters, and K. Chang, Electron tunneling through double magnetic barriers on the surface of a topological insulator, Phys. Rev. B 82, 115211 (2010).
  * Johnson  _et al._ [1997] M. Johnson, B. R. Bennett, M. J. Yang, M. M. Miller, and B. V. Shanabrook, Hybrid Hall effect device, Appl. Phys. Lett. 71, 974 (1997).
  * Nogaret  _et al._ [2000] A. Nogaret, S. J. Bending, and M. Henini, Resistance resonance effects through magnetic edge states, Phys. Rev. Lett. 84, 2231 (2000).
  * Goldhaber [1977] A. S. Goldhaber, Dirac particle in a magnetic field: symmetries and their breaking by monopole singularities, Phys. Rev. D 16, 1815 (1977).
  * Andreev and Spivak [2018] A. V. Andreev and B. Z. Spivak, Longitudinal negative magnetoresistance and magnetotransport phenomena in conventional and topological conductors, Phys. Rev. Lett. 120, 026601 (2018).
  * Wang  _et al._ [2021] H.-W. Wang, B. Fu, and S.-Q. Shen, Helical symmetry breaking and quantum anomaly in massive Dirac fermions, Phys. Rev. B 104, L241111 (2021).
  * Wang  _et al._ [2016] J. Wang, N. Wang, H. Huang, and W. Duan, Electronic properties of SnTe-class topological crystalline insulator materials, Chin. Phys. B 25, 117313 (2016).
  * Gogin  _et al._ [2022a] L. Gogin, L. Rossi, F. Rossi, and F. Dolcini, The Dirac paradox in 1 + 1 dimensions and its realization with spin–orbit coupled nanowires, New J. Phys. 24, 053045 (2022a).
  * Gogin  _et al._ [2022b] L. Gogin, F. Rossi, and F. Dolcini, Electron transport in quantum channels with spin–orbit interaction: effects of the sign of the rashba coupling and applications to nanowires, New J. Phys. 24, 093025 (2022b).
  * Cheianov and Fal’ko [2006] V. V. Cheianov and V. I. Fal’ko, Selective transmission of Dirac electrons and ballistic magnetoresistance of n-p junctions in graphene, Phys. Rev. B 74, 041403(R) (2006).
  * Hills  _et al._ [2017] R. D. Y. Hills, A. Kusmartseva, and F. V. Kusmartsev, Current-voltage characteristics of Weyl semimetal semiconducting devices, Veselago lenses, and hyperbolic Dirac phase, Phys. Rev. B 95, 214103 (2017).
  * Calogeracos and Dombey [1999] A. Calogeracos and N. Dombey, History and physics of the Klein paradox, Contemp. Phys. 40, 313 (1999).