Universality classes of topological phase transitions with higher-order band crossing

We discuss the TIs and TSCs whose low-energy effective theory is a Dirac model satisfying the following criteria:

• (1)  
  The low-energy effective Dirac Hamiltonian $H({\bf{k}})={\bf{d}}({\bf{k}})\cdot {\boldsymbol{\Gamma }}$ is constructed out of Γ^a matrices that satisfies the Clifford algebra $\left\{{{\rm{\Gamma }}}^{a},{{\rm{\Gamma }}}^{b}\right\}=2{\delta }_{{ab}}$. In addition, the low-energy dispersion takes the form

  $$\begin{eqnarray}&&{E}_{\pm }({\bf{k}})=\pm d=\pm {\left({k}^{2n}+{M}^{2}\right)}^{1/2},\end{eqnarray} \tag{ 1 }$$

  where the integer n signifies the order of band crossing at the topological phase transition M_c = 0. The mass term M is assumed to have no momentum dependence.
• (2)  
  The model is classified according to the TR, PH and chiral symmetries

  $$\begin{eqnarray}\begin{array}{rcl} & & {TH}({\bf{k}}){T}^{-1}=H(-{\bf{k}}),\\ & & {CH}({\bf{k}}){C}^{-1}=-H(-{\bf{k}}),\\ & & {SH}({\bf{k}}){S}^{-1}=-H({\bf{k}}).\end{array}\end{eqnarray} \tag{ 2 }$$

  In each spatial dimension from 1D to 3D, there are 5 symmetry classes that are topologically nontrivial. We will discuss all the 15 nontrivial dimension-symmetry classes case by case in the following sections.
• (3)  
  The points in the Brillouin zone (BZ) satisfying k₀ = −k₀ (up to a reciprocal lattice vector) are referred to as the HSPs. If the Hamiltonian has inversion, TR, or other point group symmetries, then the energy bands must have extrema at the HSPs. If the low harmonics (e.g. $\cos k$ or $\sin k$) of the lattice tight-binding Hamiltonian are dominant (i.e. their coefficients are larger than those of the higher harmonics), then the topological phase transition is accompanied by a gap-closing at a HSP. Since this is the case for most examples in the literature, we focus on models where the gap-closing occurs at a HSP.
• (4)  
  In 2D and 3D, we consider m-fold rotation symmetric Dirac Hamiltonians and choose an eigenvalue basis in which both the Dirac Hamiltonian and the rotation operator are diagonal. In this basis, the operator of the m-fold rotation symmetry can be written as

  $$\begin{eqnarray}&&{C}_{m}=\mathrm{diag}({\alpha }_{p},{\alpha }_{q}...),\ \ \ {\alpha }_{p}={{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}\left(p+\displaystyle \frac{F}{2}\right)},\end{eqnarray} \tag{ 3 }$$

  where p = 0, 1, 2 ... m − 1 are the rotational eigenvalues. The spin statistical factor F indicates whether the basis is fermionic F = 1 or bosonic F = 0. We assume that in each symmetry class, the basis can be either fermionic or bosonic. The rotational symmetry requires that

  $$\begin{eqnarray}&&{C}_{m}H({\bf{k}}){C}_{m}^{-1}=H({R}_{m}{\bf{k}}),\end{eqnarray} \tag{ 4 }$$

  where R_m rotates the momentum k about the k_z axis by 2π/m. In 1D this rotational symmetry is absent.

As will be demonstrated case by case in the following sections, the topological invariant ${ \mathcal C }$ (or a related one that judges the topological phase transition) in all the 15 dimension-symmetry classes can be cast into the form of a D-dimensional integration over the BZ

$$\begin{eqnarray}&&{ \mathcal C }={ \mathcal C }(M)={\int }_{\mathrm{BZ}}\displaystyle \frac{{{\rm{d}}}^{D}{\bf{k}}}{{\left(2\pi \right)}^{D}}\ F({\bf{k}},M),\end{eqnarray} \tag{ 5 }$$

where F(k, M) is referred to as the curvature function, the precise form of which depends on the dimension and symmetry class. Owing to the point (3) above, the curvature function is generally an even function $F({{\bf{k}}}_{0}+\delta {\bf{k}},M)=F({{\bf{k}}}_{0}-\delta {\bf{k}},M)$ around the gap-closing HSP. While the topological invariant ${ \mathcal C }$ remains constant within a topological phase in the parameter space of M, the profile of F(k, M) varies with changing M, which is key to our analysis of critical behavior.

For all the 15 dimension/symmetry classes considered, we will discuss the critical behavior of the system by means of a continuous model that respects the low energy sector of the system. In this continuous model description, the integration in equation (5) is treated as an unbounded integration over the entire D-dimensional momentum space ${\int }_{{\rm{BZ}}}{{\rm{d}}}^{D}{\bf{k}}\to {\int }_{-\infty }^{\infty }{{\rm{d}}{k}}_{1}...{\int }_{-\infty }^{\infty }{{\rm{d}}{k}}_{D}$. In such an unbounded integration treatment, even though the integration may contribute to only a fraction of the true topological invariant ${ \mathcal C }$, the jump of topological invariant ${\rm{\Delta }}{ \mathcal C }$ can be correctly captured provided the integration is properly normalized

$$\begin{eqnarray}&&{\rm{\Delta }}{ \mathcal C }={ \mathcal C }({M}_{c}^{+})-{ \mathcal C }({M}_{c}^{-})=\displaystyle \frac{1}{{ \mathcal N }}\left(\displaystyle \prod _{i=1}^{D}{\int }_{-\infty }^{\infty }{{\rm{d}}{k}}_{i}\right)\left[F({\bf{k}},{M}_{c}^{+})-F({\bf{k}},{M}_{c}^{-})\right].\end{eqnarray} \tag{ 6 }$$

The normalization factor ${ \mathcal N }$ can be fixed by the fact that ${\rm{\Delta }}{ \mathcal C }=1$ for linear Dirac models n = 1. Equation (6) is valid because near the critical point $M\to {M}_{c}$, the jump of topological invariant ${\rm{\Delta }}{ \mathcal C }$ is entirely determined by the integration of the curvature function around a small region $\sim 1/\xi$ near the HSP, so the momentum cutoff is unimportant, as will be demonstrated in section 2.2.

Another key ingredient in our analysis is the Wannier states $| {\bf{R}}n\rangle$ defined from the Bloch state by

$$\begin{eqnarray}&&| {u}_{n{\bf{k}}}\rangle =\displaystyle \sum _{{\bf{R}}}{{\rm{e}}}^{-{\rm{i}}{\bf{k}}\cdot (\hat{{\bf{r}}}-{\bf{R}})}| {\bf{R}}n\rangle ,\ \ \ | {\bf{R}}n\rangle =\displaystyle \frac{1}{N}\displaystyle \sum _{{\bf{k}}}{{\rm{e}}}^{{\rm{i}}{\bf{k}}\cdot (\hat{{\bf{r}}}-{\bf{R}})}| {u}_{n{\bf{k}}}\rangle ,\end{eqnarray} \tag{ 7 }$$

where N denotes the number of lattice sites, and $\hat{{\bf{r}}}$ is the position operator. We propose a correlation function derived through the Fourier transform of the curvature function,

$$\begin{eqnarray}&&{\lambda }_{{\bf{R}}}=\int \displaystyle \frac{{{\rm{d}}}^{D}{\bf{k}}}{{\left(2\pi \right)}^{D}}\ {{\rm{e}}}^{{\rm{i}}{\bf{k}}\cdot {\bf{R}}}\ F({\bf{k}},M).\end{eqnarray} \tag{ 8 }$$

In each of the 15 dimension-symmetry classes, we will demonstrate that the correlation function always takes the form of measuring the overlap of Wannier functions.

Through investigating the low-energy continuous models of all the 15 nontrivial dimension-symmetry classes from 1D to 3D, we found that the critical behavior of the curvature functions fall into two different scenarios.

2.2.1. Peak-divergence scenario

For linear band crossing n = 1, the curvature function F(k, M) peaks at the gap-closing HSP. The peak gradually narrows and the height increases as the system approaches the critical point $M\to {M}_{c}$, and eventually the peak diverges and flips sign across the transition. We call this critical behavior the peak-divergence scenario. In this scenario, the peak can be well fitted by an Ornstein–Zernike form

$$\begin{eqnarray}&&F({{\bf{k}}}_{0}+\delta {\bf{k}},M)=\displaystyle \frac{F({{\bf{k}}}_{0},M)}{1+{\xi }^{2}\delta {k}^{2}}.\end{eqnarray} \tag{ 9 }$$

The critical behavior described above is summarized as

$$\begin{eqnarray}\begin{array}{rcl} & & \mathop{\mathrm{lim}}\limits_{M\to {M}_{c}^{+}}F({{\bf{k}}}_{0},M)=-\mathop{\mathrm{lim}}\limits_{M\to {M}_{c}^{-1}}F({{\bf{k}}}_{0},M)=\pm \infty ,\ \ \ \\ & & \mathop{\mathrm{lim}}\limits_{M\to {M}_{c}}\xi =\infty .\end{array}\end{eqnarray} \tag{ 10 }$$

Denoting the critical exponents of F(k₀, M) and ξ as

$$\begin{eqnarray}&&| F({{\bf{k}}}_{0},M)| \propto | M-{M}_{c}{| }^{-\gamma },\xi \propto | M-{M}_{c}{| }^{-\nu },\end{eqnarray} \tag{ 11 }$$

the conservation of the topological invariant ${{ \mathcal C }}_{{\rm{div}}}$ under the diverging peak

$$\begin{eqnarray}&&{{ \mathcal C }}_{\mathrm{div}}=F({{\bf{k}}}_{0},M)\left(\displaystyle \prod _{i=1}^{D}{\int }_{-{\xi }^{-1}}^{{\xi }^{-1}}\displaystyle \frac{{\rm{d}}{{k}}_{i}}{1+{\xi }^{2}{k}_{i}^{2}}\right)=\displaystyle \frac{F({{\bf{k}}}_{0},M)}{{\xi }^{D}}\times { \mathcal O }(1)=\mathrm{const}.,\end{eqnarray} \tag{ 12 }$$

yields a scaling law

$$\begin{eqnarray}&&\gamma =D\nu ,\end{eqnarray} \tag{ 13 }$$

which is to be satisfied by any linear, isotropic Dirac model in any dimension. Because of this Ornstein–Zernike form, the Wannier state correlation function as a Fourier transform of the curvature function, as in equation (8), decays with correlation length ξ.

2.2.2. Shell-divergence scenario

For any high order band crossing n > 1, the extremum of the curvature function ${F}_{\max }({{\bf{k}}}_{\max },M)$ forms a D − 1 dimensional shell surrounding the HSP, with a radius denoted by k_max and a thickness denoted by k_wid. The critical behavior is that the shell gradually reduces it radius as $M\to {M}_{c}$, and F_max(k_max, M) gradually diverges and flips sigh across the transition. In other words,

$$\begin{eqnarray}\begin{array}{rcl} & & \mathop{{\rm{lim}}}\limits_{M\to {M}_{c}^{+}}{F}_{{\rm{\max }}}({{\bf{k}}}_{{\rm{\max }}},M)=-\mathop{{\rm{lim}}}\limits_{M\to {M}_{c}^{-1}}{F}_{{\rm{\max }}}({{\bf{k}}}_{{\rm{\max }}},M)=\pm \infty ,\ \ \ \\ & & \mathop{{\rm{lim}}}\limits_{M\to {M}_{c}}{k}_{{\rm{\max }}}=0,\ \ \ \mathop{{\rm{lim}}}\limits_{M\to {M}_{c}}{k}_{{\rm{wid}}}=0.\end{array}\end{eqnarray} \tag{ 14 }$$

For the Dirac models in section 2.1, the k_max and k_wid have the same critical exponent

$$\begin{eqnarray}&&| F({{\bf{k}}}_{0},M)| \propto | M-{M}_{c}{| }^{-\gamma },\ \ \ {k}_{\max }^{-1}\propto {k}_{{\rm{wid}}}^{-1}\propto \xi \propto | M-{M}_{c}{| }^{-\nu }.\end{eqnarray} \tag{ 15 }$$

The conservation of the portion of the topological invariant ${{ \mathcal C }}_{{\rm{div}}}$ in the diverging shell

$$\begin{eqnarray}&&{{ \mathcal C }}_{{\rm{div}}}={\rm{\Omega }}{\int }_{{k}_{{\rm{\max }}}-{k}_{{\rm{wid}}}/2}^{{k}_{{\rm{\max }}}+{k}_{{\rm{wid}}}/2}{\rm{d}}{k}\,{k}^{D-1}F({\bf{k}},M)\approx {\rm{\Omega }}{k}_{{\rm{wid}}}{k}_{{\rm{\max }}}^{D-1}{F}_{{\rm{\max }}}=\mathrm{const}.,\end{eqnarray} \tag{ 16 }$$

where Ω takes into account the angular integration, yields the same scaling law as equation (13). The two scenarios are depicted pictorially for 1D and 2D in figure 1.

Zoom In Zoom Out Reset image size

To be more specific, find that the low-energy continuous model for most of the symmetry classes has the topological invariant in equation (5) in the following scaling form

$$\begin{eqnarray}&&{ \mathcal C }={\rm{\Omega }}{\int }_{0}^{\infty }{\rm{d}}{k}\,{k}^{D-1}\displaystyle \frac{{{Mk}}^{D(n-1)}}{{\left({M}^{2}+{k}^{2n}\right)}^{(D+1)/2}}.\end{eqnarray} \tag{ 17 }$$

Expanding the curvature function in the integrand for small k yields

$$\begin{eqnarray}&&F(k,M)=\displaystyle \frac{{{Mk}}^{D(n-1)}}{{\left({M}^{2}+{k}^{2n}\right)}^{(D+1)/2}}\approx \displaystyle \frac{\mathrm{Sgn}(M)}{| M{| }^{D}}\displaystyle \frac{{k}^{D(n-1)}}{\left(1+\tfrac{D+1}{2}\tfrac{{k}^{2n}}{{M}^{2}}\right)}.\end{eqnarray} \tag{ 18 }$$

For the peak-divergence scenario at n = 1, the Ornstein–Zernike form is evident, which renders γ = D and ν = 1, satisfying equation (13). For the shell-divergence scenario at n > 1, solving for the k_max from ∂_kF(k, M) = 0 yields ${k}_{\max }\propto | M{| }^{1/n}$, and subsequently we obtain ${F}_{\max }({k}_{\max },M)\propto | M{| }^{-D/n}$. Consequently, the critical exponents are γ = D/n and ν = 1/n, in agreement with the scaling law in equation (13). This analysis indicates that the critical exponent ν for the length scale ξ are basically determined by the order of band crossing n. This statement can also be inferred by simply counting the dimensions in equation (1). From the energy spectrum, one infers that $\left[{k}^{n}\right]=\left[M\right]$ are of the same unit. Consequently the inverse of the momentum gives the characteristic length scale of the problem $\xi \sim {k}^{-1}\sim | M{| }^{1/n}$.

The edge state decay length in the topologically nontrivial phase has the same critical exponent as the correlation length. This can be seen by considering the model to be defined in the half-space x > 0, and considering the edge state with zero transverse momentum, i.e. k_y = k_z = 0 in 3D and k_y = 0 in 2D. The problem can be calculated by projecting the higher-order Dirac Hamiltonian into real space and solve for the zero-energy edge state [6]

$$\begin{eqnarray}&&\left({\bf{d}}\cdot {\boldsymbol{\Gamma }}\right)\psi =\left[{k}_{x}^{n}{{\rm{\Gamma }}}^{x}+M{{\rm{\Gamma }}}^{M}\right]\psi =\left[{\left(-{\rm{i}}\right)}^{n}{\partial }_{x}^{n}{{\rm{\Gamma }}}^{x}+M{{\rm{\Gamma }}}^{M}\right]\psi =0,\end{eqnarray} \tag{ 19 }$$

where ${{\rm{\Gamma }}}^{x}$ and Γ^M represent the Γ-matrices that multiply k_xⁿ and M, respectively, whose precise form depends on the symmetry class and dimension. Multiplying the equation by Γ^x and using the ansatz ψ = χ_ηϕ(x), with χ_η an eigenstate of Γ^xΓ^Mχ_η = ηχ_η, we obtain

$$\begin{eqnarray}&&{\xi }^{n}=\displaystyle \frac{{i}^{n-2}}{\eta M}.\end{eqnarray} \tag{ 20 }$$

The eigenvalue η is chosen such that one of the roots ${\xi }_{\mathrm{loc}}$ is real and identifiable with the edge state decay length, which obviously also has critical exponent ν = 1/n. Thus, the correlation length and the edge state decay length are essentially the same length scale, although we should emphasize that the edge state only exists in the topologically nontrivial phase, whereas the Wannier state correlation function is well-defined in either the trivial or nontrivial phase.

The CRG approach has been proposed to judge topological phase transitions in a variety of systems [6, 11, 12, 21, 22]. Here we demonstrate that the method applies to all the 15 dimension-symmetry classes at any order of band crossing n and rotational symmetry m. The approach is essentially an iterative method to search for the trajectory (RG flow) in the parameter space of M along which the maximum of the curvature function is reduced. In this way, the topological phase transitions, at which the curvature function diverges, can be identified. The scaling equation demands that, at a given M, we find the new $M^{\prime}$ that satisfies

$$\begin{eqnarray}&&F({{\bf{k}}}_{0}+\delta {\bf{k}},M)=F({{\bf{k}}}_{0}+b\delta {\bf{k}},M^{\prime} ).\end{eqnarray} \tag{ 21 }$$

The mapping $M\to M^{\prime}$ yields the RG flow that identifies the topological phase transitions. Here δk is a small momentum displacement away from the HSP, and one chooses 0 < b < 1 for the peak-divergence scenario n = 1, whereas b > 1 for the shell-divergence scenario n > 1. The divergence of the curvature function is gradually reduced under the scaling procedure, as shown schematically in figure 2. We remark that a similar RG approach has also been proposed to search for the transition based on the reduced density matrix obtained from bipartition, since it shows a similar critical behavior [36].

Zoom In Zoom Out Reset image size

Since the curvature function in many symmetry classes takes the generic scaling form of equation (18), we examine here explicitly the CRG approach applied to this scaling form. As the approach concerns only a small k region near the HSP which is set to be at the origin in our continuous models k₀ = 0, we expand the scaling form in equation (18) by

$$\begin{eqnarray}&&F(\delta k,M)\approx \displaystyle \frac{\mathrm{Sgn}(M)}{| M{| }^{D}}\delta {k}^{D(n-1)}\left(1-\displaystyle \frac{D+1}{2}\displaystyle \frac{\delta {k}^{2n}}{{M}^{2}}\right).\end{eqnarray} \tag{ 22 }$$

For the peak-divergence scenario n = 1, putting equation (22) into the leading order expansion of the scaling equation, equation (21), yields the RG equation

$$\begin{eqnarray}&&\displaystyle \frac{{\rm{d}}{M}}{{\rm{d}}{\ell }}\equiv \displaystyle \frac{M^{\prime} -M}{\delta {k}^{2}}=\displaystyle \frac{1-{b}^{2}}{2}\displaystyle \frac{{\partial }_{k}^{2}F(k,M){| }_{k=0}}{{\partial }_{M}F(k,M){| }_{k=0}}=\displaystyle \frac{(1-{b}^{2})(D+1)}{2D}\displaystyle \frac{1}{M},\end{eqnarray} \tag{ 23 }$$

where we have defined the scaling parameter dℓ = δk². The corresponding RG flow flows away from the critical point M_c = 0, with a flow rate that diverges as $M\to {M}_{c}$.

For the shell-divergence scenario at n > 1, we approximate equation (22) by

$$\begin{eqnarray}\begin{array}{rcl} & & F(\delta k,M)\approx \displaystyle \frac{\mathrm{Sgn}(M)}{| M{| }^{D}}\delta {k}^{D(n-1)},\\ & & F(b\delta k,M^{\prime} )=F(\delta k+{\rm{\Delta }}k,M+\delta M)\approx \displaystyle \frac{\mathrm{Sgn}(M)}{| M+\delta M{| }^{D}}{\left(\delta k+{\rm{\Delta }}k\right)}^{D(n-1)}.\end{array}\end{eqnarray} \tag{ 24 }$$

Treating $\delta M\ll M$ and ${\rm{\Delta }}k=(b-1)\delta k\ll \delta k$ as small quantities and expanding to the leading order, and then putting equation (24) into the (21) yields the leading order RG equation

$$\begin{eqnarray}&&\displaystyle \frac{{\rm{d}}{M}}{{\rm{d}}{\ell }}\equiv \displaystyle \frac{\delta M}{({\rm{\Delta }}k/\delta k)}=(n-1)M,\end{eqnarray} \tag{ 25 }$$

where the scaling parameter is defined by ${\rm{d}}{\ell }={\rm{\Delta }}k/\delta k=b-1$. The equation corresponds to an RG flow that flows away from M_c = 0, with a flow rate that vanishes at M = M_c. In other words, the topological phase transition manifests as an unstable fixed point at n > 1.

For 3D class AIII that has S² = 1, we choose the Γ-matrices [8]

$$\begin{eqnarray}\begin{array}{rcl}{{\rm{\Gamma }}}^{a} & = & \left({\alpha }_{x},{\alpha }_{y},{\alpha }_{z},\beta ,-{\rm{i}}\beta {\gamma }^{5}\right),\\ {\alpha }_{{i}} & = & \left(\begin{array}{cc}0 & {\sigma }_{{i}}\\ {\sigma }_{{i}} & 0\end{array}\right),\ \ \ \beta =\left(\begin{array}{cc}1 & 0\\ 0 & -1\end{array}\right),\ \ \ {\gamma }^{5}=\left(\begin{array}{cc}0 & 1\\ 1 & 0\end{array}\right).\end{array}\end{eqnarray} \tag{ 26 }$$

The chiral operator is S = β, which together with the chiral symmetry in equation (2) requires that d₄(k) = 0, whereas there is no constraint on d_i(k) being even or odd in k for i = 1, 2, 3, 5. We follow the convention to choose d₅ = M as the mass term. The Hamiltonian then takes the form

$$\begin{eqnarray}H({\bf{k}})=\displaystyle \sum _{i=1,2,3,5}{d}_{i}{{\rm{\Gamma }}}^{i}=\left(\begin{array}{cccc} & & {g}^{{\prime} } & f\\ & & {f}^{* } & -{g}^{{\prime} * }\\ {g}^{{\prime} * } & f & & \\ {f}^{* } & -{g}^{{\prime} } & & \end{array}\right),\end{eqnarray} \tag{ 27 }$$

where we have denoted

$$\begin{eqnarray}\begin{array}{rcl} & & f={d}_{1}-{{\rm{i}}{d}}_{2}={k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}},\\ & & {g}^{{\prime} }={d}_{3}-{\rm{i}}{M}=g-{\rm{i}}{M}={k}_{z}^{n}-{\rm{i}}{M},\end{array}\end{eqnarray} \tag{ 28 }$$

with $n={n}_{+}+{n}_{-}$ and ${k}_{\pm }={k}_{x}\pm {{\rm{i}}{k}}_{y}$, such that the dispersion satisfies

$$\begin{eqnarray}&&{E}_{\pm }({\bf{k}})=\pm d=\pm \sqrt{| f{| }^{2}+| g{| }^{2}+{M}^{2}}=\pm \sqrt{{k}^{2n}+{M}^{2}},\end{eqnarray} \tag{ 29 }$$

with $n\in {\mathbb{Z}}$.

Using the rotational symmetry operator

$$\begin{eqnarray}&&{C}_{m}=\mathrm{diag}({\alpha }_{p},{\alpha }_{q},{\alpha }_{r},{\alpha }_{s}),\ \ \ {\alpha }_{p}={{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}\left(p+\displaystyle \frac{F}{2}\right)},\end{eqnarray} \tag{ 30 }$$

the chiral symmetry [C_m, S] = 0 does not give any constraint on C_m. The rotational symmetry in equation (4) demands

$$\begin{eqnarray}\begin{array}{rcl}{\alpha }_{p}{\alpha }_{s}^{* }f({k}_{\pm },{k}_{z}) & = & {\alpha }_{r}{\alpha }_{q}^{* }f({k}_{\pm },{k}_{z})=f({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}},{k}_{z}),\\ {\alpha }_{p}{\alpha }_{r}^{* }g^{\prime} ({k}_{\pm },{k}_{z}) & = & {\alpha }_{s}{\alpha }_{q}^{* }g^{\prime} ({k}_{\pm },{k}_{z})=g^{\prime} ({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}},{k}_{z}).\end{array}\end{eqnarray} \tag{ 31 }$$

Because $g^{\prime} =g-{\rm{i}}{M}$ contains the mass term that must be invariant under the rotational transformation, it is required that ${\alpha }_{p}{\alpha }_{r}^{* }=1$ and ${\alpha }_{s}{\alpha }_{q}^{* }=1$, and therefore r = p and q = s. The constraint on $g^{\prime}$ then gives $g({k}_{\pm },{k}_{z})=g({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}},{k}_{z})$, which is consistent with our parameterization in equation (28). The parameterization of f in equation (28), together with ${\alpha }_{p}{\alpha }_{s}^{* }f({k}_{\pm },{k}_{z})={\alpha }_{p}{\alpha }_{q}^{* }f({k}_{\pm },{k}_{z})=f({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}},{k}_{z})$ gives

$$\begin{eqnarray}&&{n}_{+}-{n}_{-}={xm}+(p-q)\in {\mathbb{Z}},\end{eqnarray} \tag{ 32 }$$

which is the relation between the parameterization and the rotational symmetry eigenvalues, and we have used ${n}_{+}-{n}_{-}\in {\mathbb{Z}}$ due to $n={n}_{+}+{n}_{-}\in {\mathbb{Z}}$. The x in equation (32) and throughout the article is an arbitrary integer. Evidently, equation (32) can be satisfied by any of the m = 2, 3, 4, 6 for either the bosonic or fermionic basis.

We define the q-matrix by (not to be confused with one of the rotational eigenvalues)

$$\begin{eqnarray}q({\bf{k}})=-\displaystyle \frac{1}{d}\left({d}_{i}{\sigma }_{i}-{\rm{i}}{M}\right)=-\displaystyle \frac{1}{d}\left(\begin{array}{cc}g-{\rm{i}}{M} & f\\ {f}^{* } & -g-{\rm{i}}{M}\end{array}\right).\end{eqnarray} \tag{ 33 }$$

Writing the parameterization in equation (28) into cylindrical coordinate $({k}_{\perp },\phi ,{k}_{z})$ by $f={k}_{\perp }^{n}{{\rm{e}}}^{{\rm{i}}({n}_{+}-{n}_{-})\phi }$ and $g={k}_{z}^{n}$, and using the derivatives

$$\begin{eqnarray}&&\displaystyle \frac{\partial }{\partial {k}_{x}}=\cos \phi \displaystyle \frac{\partial }{\partial {k}_{\perp }}-\displaystyle \frac{\sin \phi }{{k}_{\perp }}\displaystyle \frac{\partial }{\partial \phi },\ \ \ \displaystyle \frac{\partial }{\partial {k}_{y}}=\sin \phi \displaystyle \frac{\partial }{\partial {k}_{\perp }}+\displaystyle \frac{\cos \phi }{{k}_{\perp }}\displaystyle \frac{\partial }{\partial \phi },\end{eqnarray} \tag{ 34 }$$

yields the following combination in the integration of topological invariant [8]

$$\begin{eqnarray}&&{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{q}^{\dagger }{\partial }_{\mu }{{qq}}^{\dagger }{\partial }_{\nu }{{qq}}^{\dagger }{\partial }_{\rho }q\right]=-M\displaystyle \frac{12{k}_{\perp }^{2(n-1)}{k}_{z}^{n-1}({n}_{+}-{n}_{-}){n}^{2}}{{\left({k}^{2n}+{M}^{2}\right)}^{2}}.\end{eqnarray} \tag{ 35 }$$

Converting the integrand into spherical coordinates (k, θ, ϕ) yields

$$\begin{eqnarray}&&{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{q}^{\dagger }{\partial }_{\mu }{{qq}}^{\dagger }{\partial }_{\nu }{{qq}}^{\dagger }{\partial }_{\rho }q\right]=-M\displaystyle \frac{12{\left(k\sin \theta \right)}^{2(n-1)}{\left(k\cos \theta \right)}^{n-1}({n}_{+}-{n}_{-}){n}^{2}}{{\left({k}^{2n}+{M}^{2}\right)}^{2}}.\end{eqnarray} \tag{ 36 }$$

The topological invariant is given by the 3D winding number obtained from the integration of the above quantity

$$\begin{eqnarray}\begin{array}{rcl}{ \mathcal C } & = & \displaystyle \frac{\pi }{3}\displaystyle \int \displaystyle \frac{{{\rm{d}}}^{3}{\bf{k}}}{{\left(2\pi \right)}^{3}}{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{q}^{\dagger }{\partial }_{\mu }{{qq}}^{\dagger }{\partial }_{\nu }{{qq}}^{\dagger }{\partial }_{\rho }q\right]\\ & = & \left\{\begin{array}{ll}-\tfrac{({n}_{+}-{n}_{-})}{2}\mathrm{Sgn}(M)\left[2{nB}(\tfrac{n}{2},n)\right] & {\rm{if}}\ n\in 2{\mathbb{Z}}+1,\\ 0 & {\rm{if}}\ n\in 2{\mathbb{Z}},\end{array}\right.\end{array}\end{eqnarray} \tag{ 37 }$$

where $B(a,b)={\int }_{0}^{1}{t}^{a-1}{\left(1-t\right)}^{b-1}{\rm{d}}{t}$ is the Beta function. In equation (37) and all the continuous models in the following sections, the integration of the topological invariant ${ \mathcal C }$ is over the entire D-dimensional momentum space $\int {{\rm{d}}}^{D}{\bf{k}}={\int }_{-\infty }^{\infty }{{\rm{d}}{k}}_{1}...{\int }_{-\infty }^{\infty }{{\rm{d}}{k}}_{D}$. We see that there is a topological phase transition ${\rm{\Delta }}{ \mathcal C }\ne 0$ if and only if the order of band crossing n is an odd number. Here the extra factor $2{nB}(\tfrac{n}{2},n)$ is an artifact of this continuous model, and the jump of topological invariant should be treated as ${\rm{\Delta }}{ \mathcal C }=| { \mathcal C }(M\gt 0)-{ \mathcal C }(M\lt 0)| ={n}_{+}-{n}_{-}\in 2{\mathbb{Z}}+1$, as demonstrated by the following lattice model

$$\begin{eqnarray}\begin{array}{rcl} & & f={d}_{1}-{{\rm{i}}{d}}_{2}={\left({\rm{\sin }}{k}_{x}-{\rm{i}}{\rm{\sin }}{k}_{y}\right)}^{{n}_{-}},\ \ \ g={d}_{3}={{\rm{\sin }}}^{{n}_{-}}{k}_{z},\\ & & {d}_{5}=M+{\rm{\cos }}{k}_{x}+{\rm{\cos }}{k}_{y}+{\rm{\cos }}{k}_{z}.\end{array}\end{eqnarray} \tag{ 38 }$$

The resulting topological invariant is

$$\begin{eqnarray}\begin{array}{rcl}{ \mathcal C } & = & \displaystyle \frac{\pi }{3}{\displaystyle \int }_{{\rm{BZ}}}\displaystyle \frac{{{\rm{d}}}^{3}{\bf{k}}}{{\left(2\pi \right)}^{3}}{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{q}^{\dagger }{\partial }_{\mu }{{qq}}^{\dagger }{\partial }_{\nu }{{qq}}^{\dagger }{\partial }_{\rho }q\right]\\ & = & \left\{\begin{array}{ll}{n}_{-} & {\rm{for}}\ M\gtrsim -3\ \mathrm{and}\ {n}_{-}\in 2{\mathbb{Z}}+1,\\ 0 & {\rm{for}}\ M\lesssim -3\ \mathrm{or}\ {n}_{-}\in 2{\mathbb{Z}},\end{array}\right.\end{array}\end{eqnarray} \tag{ 39 }$$

which gives ${\rm{\Delta }}{ \mathcal C }={n}_{-}\in 2{\mathbb{Z}}+1$.

The Wannier state correlation function in this class is introduced from the Chern–Simons gauge field ${{ \mathcal A }}_{\mu }^{{ab}}=\langle {u}_{a-}| {\partial }_{\mu }| {u}_{b-}\rangle$ defined from the eigenstates for the filled bands ${E}_{-}({\bf{k}})=-d$

$$\begin{eqnarray}&&| {u}_{1-}\rangle =\displaystyle \frac{1}{\sqrt{2}d}\left(\begin{array}{c}-{d}_{1}+{{\rm{i}}{d}}_{2}\\ {d}_{3}+{\rm{i}}{M}\\ 0\\ d\end{array}\right),\ \ \ | {u}_{2-}\rangle =\displaystyle \frac{1}{\sqrt{2}d}\left(\begin{array}{c}-{d}_{3}+{\rm{i}}{M}\\ -{d}_{1}-{{\rm{i}}{d}}_{2}\\ d\\ 0\end{array}\right).\end{eqnarray} \tag{ 40 }$$

We find that the integrand in equation (37) can be written as

$$\begin{eqnarray}&&{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{q}^{\dagger }{\partial }_{\mu }{{qq}}^{\dagger }{\partial }_{\nu }{{qq}}^{\dagger }{\partial }_{\rho }q\right]=8{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{{ \mathcal A }}_{\mu }{{ \mathcal A }}_{\nu }{{ \mathcal A }}_{\rho }\right]=\displaystyle \frac{12M}{{d}^{4}}{\epsilon }^{\mu \nu \rho }{\partial }_{x}{d}_{\mu }{\partial }_{y}{d}_{\nu }{\partial }_{z}{d}_{\rho }.\end{eqnarray} \tag{ 41 }$$

As a result, the topological invariant in equation (39) is equivalently

$$\begin{eqnarray}&&{ \mathcal C }=\displaystyle \frac{8\pi }{3}{\int }_{{\rm{BZ}}}\displaystyle \frac{{{\rm{d}}}^{3}{\bf{k}}}{{\left(2\pi \right)}^{3}}{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{{ \mathcal A }}_{\mu }{{ \mathcal A }}_{\nu }{{ \mathcal A }}_{\rho }\right].\end{eqnarray} \tag{ 42 }$$

In addition, the Chern–Simons gauge field ${{ \mathcal A }}_{\mu }^{{ab}}$ can be expressed as the Fourier transform of the filled-band Wannier states

$$\begin{eqnarray}&&{{ \mathcal A }}_{\mu }^{{ab}}=-{\rm{i}}\displaystyle \sum _{{\bf{R}}}{{\rm{e}}}^{{\rm{i}}{\bf{k}}\cdot {\bf{R}}}\langle {\bf{0}}a| {\hat{r}}_{\mu }| {\bf{R}}b\rangle =-{\rm{i}}\displaystyle \sum _{{\bf{R}}}{{\rm{e}}}^{-{\rm{i}}{\bf{k}}\cdot {\bf{R}}}\langle {\bf{R}}a| {\hat{r}}_{\mu }| {\bf{0}}b\rangle ,\end{eqnarray} \tag{ 43 }$$

where $| {\bf{R}}b\rangle$ is the Wannier state of filled band b localized at home cell R. We use following compact notation for the eigenstates and the corresponding Wannier states

$$\begin{eqnarray}&&| {u}_{-}\rangle =\left(\begin{array}{c}| {u}_{1-}\rangle \\ | {u}_{2-}\rangle \end{array}\right),\ \ \ | {\bf{R}}\rangle =\left(\begin{array}{c}| {\bf{R}}1\rangle \\ | {\bf{R}}2\rangle \end{array}\right).\end{eqnarray} \tag{ 44 }$$

It then follows that the topological invariant in equation (42) can be written in terms of the Wannier states entirely in real space

$$\begin{eqnarray}\begin{array}{rcl}{ \mathcal C } & = & {\rm{i}}\displaystyle \frac{8\pi }{3}\displaystyle \int {{\rm{d}}}^{3}{{\bf{R}}}_{1}\displaystyle \int {{\rm{d}}}^{3}{{\bf{R}}}_{2}{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[\langle {\bf{0}}| {\hat{r}}_{\mu }| {{\bf{R}}}_{1}\rangle \langle {{\bf{R}}}_{1}| ({\hat{r}}_{\nu }-{R}_{1\nu })| {{\bf{R}}}_{1}+{{\bf{R}}}_{2}\rangle \right.\\ & & \left.\times \langle {{\bf{R}}}_{1}+{{\bf{R}}}_{2}| ({\hat{r}}_{\rho }-{R}_{1\rho }-{R}_{2\rho })| {\bf{0}}\rangle \right],\end{array}\end{eqnarray} \tag{ 45 }$$

where we have converted the discrete sum of Bravais points into an integral. Similarly, we introduce the Wannier state correlation function from the Fourier transform of the integrand in the topological invariant

$$\begin{eqnarray}\begin{array}{rcl}{\tilde{F}}_{3{\rm{D}}}({\bf{R}}) & = & {\tilde{F}}_{3{\rm{D}}}(-{\bf{R}})={\displaystyle \int }_{{\rm{BZ}}}\displaystyle \frac{{{\rm{d}}}^{3}{\bf{k}}}{{\left(2\pi \right)}^{3}}{{\rm{e}}}^{-{\rm{i}}{\bf{k}}\cdot {\bf{R}}}{\epsilon }^{\mu \nu \rho }{\rm{Tr}}\left[{{ \mathcal A }}_{\mu }{{ \mathcal A }}_{\nu }{{ \mathcal A }}_{\rho }\right]\\ & = & {\rm{i}}\displaystyle \int {{\rm{d}}}^{3}{{\bf{R}}}_{1}\displaystyle \int {{\rm{d}}}^{3}{{\bf{R}}}_{2}{\epsilon }^{\mu \nu \rho }{\rm{Tr}}\left[\langle {\bf{0}}| {\hat{r}}_{\mu }| {{\bf{R}}}_{1}\rangle \langle {{\bf{R}}}_{1}| ({\hat{r}}_{\nu }-{R}_{1\nu })| {{\bf{R}}}_{1}+{{\bf{R}}}_{2}\rangle \right.\\ & & \left.\times \langle {{\bf{R}}}_{1}+{{\bf{R}}}_{2}| ({\hat{r}}_{\rho }-{R}_{1\rho }-{R}_{2\rho })| {\bf{R}}\rangle \right].\end{array}\end{eqnarray} \tag{ 46 }$$

These real space representations are presented diagrammatically in figure 3. In addition, equations (36) and (41) indicates that the curvature function satisfies the scaling form in equation (18) with D = 3, and hence the correlation length ξ has critical exponent ν = 1/n. Finally, the Chern–Simons invariant [8], equivalently the magnetoelectric polarizability, can as well be expressed in terms of Wannier states

$$\begin{eqnarray}\begin{array}{rcl}{{CS}}_{3} & = & -\pi \displaystyle \int \displaystyle \frac{{{\rm{d}}}^{3}{\bf{k}}}{{\left(2\pi \right)}^{3}}{\epsilon }^{\mu \nu \rho }{\rm{Tr}}\left[{A}_{\mu }{\partial }_{\nu }{A}_{\rho }+\displaystyle \frac{2}{3}{A}_{\mu }{A}_{\nu }{A}_{\rho }\right]\\ & = & {\rm{i}}\pi \displaystyle \int {{\rm{d}}}^{3}{{\bf{R}}}_{1}{\epsilon }^{\mu \nu \rho }{\rm{Tr}}\left[\langle {\bf{0}}| {\hat{r}}_{\mu }| {{\bf{R}}}_{1}\rangle \langle {{\bf{R}}}_{1}| {R}_{1\nu }({\hat{r}}_{\rho }-{R}_{1\rho })| {\bf{0}}\rangle \right]\\ & & -{\rm{i}}\displaystyle \frac{2\pi }{3}\displaystyle \int {{\rm{d}}}^{3}{{\bf{R}}}_{1}\displaystyle \int {{\rm{d}}}^{3}{{\bf{R}}}_{2}{\epsilon }^{\mu \nu \rho }{\rm{Tr}}\left[\langle {\bf{0}}| {\hat{r}}_{\mu }| {{\bf{R}}}_{1}\rangle \langle {{\bf{R}}}_{1}| ({\hat{r}}_{\nu }-{R}_{1\nu })| {{\bf{R}}}_{1}+{{\bf{R}}}_{2}\rangle \right.\\ & & \left.\times \langle {{\bf{R}}}_{1}+{{\bf{R}}}_{2}| ({\hat{r}}_{\rho }-{R}_{1\rho }-{R}_{2\rho })| {\bf{R}}\rangle \right].\end{array}\end{eqnarray} \tag{ 47 }$$

according to the same argument.

Zoom In Zoom Out Reset image size

For 3D class DIII that has T² = −1 and C² = 1, we follow [8] and use the same Γ-matrices as in equation (26). The TR and PH operators are $T={\sigma }_{y}\otimes {\tau }_{x}K$ and $C={\sigma }_{y}\otimes {\tau }_{y}K$. It follows that the TR symmetry in equation (2) requires

$$\begin{eqnarray}&&{d}_{i}({\bf{k}})=-{d}_{i}(-{\bf{k}})\ \ \ \mathrm{for}\ i=1,2,3,4,\ \ \ {d}_{5}({\bf{k}})={d}_{5}(-{\bf{k}}).\end{eqnarray} \tag{ 48 }$$

Thus only the d₅ term can be the mass term d₅ = M. On the other hand, the PH symmetry requires

$$\begin{eqnarray}&&{d}_{i}({\bf{k}})=-{d}_{i}(-{\bf{k}})\ \ \ \mathrm{for}\ i=1,2,3,\ \ \ {d}_{i}({\bf{k}})={d}_{i}(-{\bf{k}})\ \ \ \mathrm{for}\ i=4,5.\end{eqnarray} \tag{ 49 }$$

Thus the d₄ term cannot exist, and the Hamiltonian takes the same form as equation (27), with the same parameterization on f and g as in equation (28). Given the rotational symmetry operator in equation (30), demanding $\left[T,{C}_{m}\right]=0$ and $\left[C,{C}_{m}\right]=0$ requires ${\alpha }_{p}={\alpha }_{s}^{* },{\alpha }_{q}={\alpha }_{r}^{* }$.

The rotational symmetry in equation (4) requires that

$$\begin{eqnarray}\begin{array}{rcl} & & {\alpha }_{p}^{2}f({k}_{\pm },{k}_{z})={\left({\alpha }_{q}^{* }\right)}^{2}f({k}_{\pm },{k}_{z})=f({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}},{k}_{z}),\\ & & {\alpha }_{p}{\alpha }_{q}g({k}_{\pm },{k}_{z})={\alpha }_{p}^{* }{\alpha }_{q}^{* }g({k}_{\pm },{k}_{z})=g({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}},{k}_{z}),\end{array}\end{eqnarray} \tag{ 50 }$$

from which we see that ${\alpha }_{p}={\alpha }_{q}^{* }$, and hence the rotational operator is constrained to be the of form

$$\begin{eqnarray}&&{C}_{m}=\mathrm{diag}({\alpha }_{p},{\alpha }_{p}^{* },{\alpha }_{p},{\alpha }_{p}^{* }).\end{eqnarray} \tag{ 51 }$$

The constraint on f in equation (50) renders

$$\begin{eqnarray}&&{{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}(2p+F)}{k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}}={{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}({n}_{+}-{n}_{-})}{k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}},\end{eqnarray} \tag{ 52 }$$

which together with the requirement that $f={d}_{1}-{{\rm{i}}{d}}_{2}$ must be odd in momentum yields

$$\begin{eqnarray}&&{n}_{+}-{n}_{-}={xm}+2p+F\in 2{\mathbb{Z}}+1.\end{eqnarray} \tag{ 53 }$$

We see that the fermionic case F = 1 can be satisfied by any of the m = 2, 3, 4, 6, but bosonic case F = 0 can only be realized at m = 3. The topological invariant in this class is the 3D winding number discussed in section 3.1, and hence the formalism of Wannier state correlation function, correlation length, and critical exponents follows that from equation (40) to (46).

For 3D class AII that has T² = −1, we use the same representation of the Γ-matrices as in equation (26), and the TR operator T = iσ_y ⨂ IK following [8]. The TR symmetry requires that

$$\begin{eqnarray}&&{d}_{i}({\bf{k}})=-{d}_{i}(-{\bf{k}})\ \mathrm{for}\ i=1,2,3,5,\ \ \ {d}_{4}({\bf{k}})={d}_{4}(-{\bf{k}}),\end{eqnarray} \tag{ 54 }$$

and hence only d₄ = M can be the mass term, and all others must be odd functions of momentum. Demanding the rotational operator in equation (30) to commute with the TR operator $\left[{C}_{m},T\right]=0$, renders ${\alpha }_{p}={\alpha }_{q}^{* }$ and ${\alpha }_{r}={\alpha }_{s}^{* }$, and consequently ${C}_{m}=\mathrm{diag}({\alpha }_{p},{\alpha }_{p}^{* },{\alpha }_{r},{\alpha }_{r}^{* })$. We may write the Hamiltonian into the form

$$\begin{eqnarray}H=\left(\begin{array}{cccc}M & & g & f\\ & M & {f}^{* } & -{g}^{* }\\ {g}^{* } & f & -M & \\ {f}^{* } & -g & & -M\end{array}\right),\end{eqnarray} \tag{ 55 }$$

where $f={d}_{1}-{{\rm{i}}{d}}_{2}$ and $g={d}_{3}-{{\rm{i}}{d}}_{5}$. The rotational symmetry in equation (4) demands that

$$\begin{eqnarray}\begin{array}{rcl}{\alpha }_{p}{\alpha }_{r}f({k}_{\pm },{k}_{z}) & = & f({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}},{k}_{z}),\\ {\alpha }_{p}{\alpha }_{r}^{* }g({k}_{\pm },{k}_{z}) & = & g({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}},{k}_{z}).\end{array}\end{eqnarray} \tag{ 56 }$$

We parameterize f and g by

$$\begin{eqnarray}&&f={k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}},\ \ \ g={k}_{z}^{n}.\end{eqnarray} \tag{ 57 }$$

With this parameterization in mind, the constraint on g in equation (56) demands that p = r, thus the rotational operator takes the form of equation (51). Consequently, we arrive at the same constraint to the rotational eigenvalues as equation (53).

The ${{\mathbb{Z}}}_{2}$ invariant ${ \mathcal C }$ for 3D lattice models in class AII can be constructed from the Pfaffian of the filled-band m-matrix ${m}_{\alpha \beta }({\bf{k}})=\langle {u}_{{\boldsymbol{\alpha }}{\rm{k}}}| T| {u}_{{\boldsymbol{\beta }}{\rm{k}}}\rangle$ at HSPs [37, 38], provided the Pfaffian is a real number [37, 38]

$$\begin{eqnarray}&&{\left(-1\right)}^{{ \mathcal C }}=\displaystyle \prod _{i}\mathrm{Sgn}\left(\mathrm{Pf}[m({{\bf{k}}}_{0,i})]\right),\end{eqnarray} \tag{ 58 }$$

where k_{0, i} is the ith HSP. The energy eigenstates for the filled bands ${E}_{-}({\bf{k}})=-d$ are

$$\begin{eqnarray}\begin{array}{rcl}| {u}_{1-}\rangle & = & \displaystyle \frac{1}{\sqrt{2d(d+M)}}\left(\begin{array}{c}-{d}_{1}+{{\rm{i}}{d}}_{2}\\ {d}_{3}+{{\rm{i}}{d}}_{5}\\ 0\\ d+M\end{array}\right),\\ | {u}_{2-}\rangle & = & \displaystyle \frac{1}{\sqrt{2d(d+M)}}\left(\begin{array}{c}-{d}_{3}+{{\rm{i}}{d}}_{5}\\ -{d}_{1}-{{\rm{i}}{d}}_{2}\\ d+M\\ 0\end{array}\right).\end{array}\end{eqnarray} \tag{ 59 }$$

The Pfaffian that enters the ${{\mathbb{Z}}}_{2}$ topological invariant in equation (58) is

$$\begin{eqnarray}&&\mathrm{Pf}\left[m({{\bf{k}}}_{0})\right]=\langle {u}_{1-}| T| {u}_{2-}\rangle =-\displaystyle \frac{M}{d}.\end{eqnarray} \tag{ 60 }$$

The critical behavior of the Pfaffian is that as $M\to {M}_{c}$, although the Pfaffian at k₀ remains at ±1, its second derivative diverges (or equivalently the Laplacian diverges since the model is isotropic around each k₀). This motivates us to consider an additional topological invariant ${{ \mathcal C }}^{{\prime} }$ that uses the Laplacian of the Pfaffian as the curvature function

$$\begin{eqnarray}&&F({\bf{k}},M)={{\rm{\nabla }}}_{{\bf{k}}}^{2}\mathrm{Pf}[m({\bf{k}},M)]=\displaystyle \sum _{i=1}^{D}{\partial }_{i}^{2}\mathrm{Pf}[m({\bf{k}},M)].\end{eqnarray} \tag{ 61 }$$

This choice of curvature function is sound because it integrates to a topological invariant, which can be seen by considering the integration over the BZ of a D-dimensional cubic lattice

$$\begin{eqnarray}&&{{ \mathcal C }}^{{\prime} }=\int \displaystyle \frac{{{\rm{d}}}^{D}{\bf{k}}}{{\left(2\pi \right)}^{D}}{{\rm{\nabla }}}_{{\bf{k}}}^{2}\mathrm{Pf}[m]=\displaystyle \sum _{j=1}^{D}{\partial }_{j}{\left.\left\{\displaystyle \prod _{i\ne j}{\int }_{-\pi }^{\pi }\displaystyle \frac{{{\rm{d}}{k}}_{i}}{2\pi }\mathrm{Pf}[m]\right\}\right|}_{{k}_{j}=-\pi }^{{k}_{j}=\pi }=0.\ \ \ \end{eqnarray} \tag{ 62 }$$

Moreover, since Pf[m] is usually an even function, its Laplacian ${{\rm{\nabla }}}_{{\bf{k}}}^{2}\mathrm{Pf}[m]$ is also an even function. Note that since ${{ \mathcal C }}^{{\prime} }$ is always zero, it is not directly related to the ${{\mathbb{Z}}}_{2}$ invariant in equation (58).

We introduce the Wannier state correlation function for TR-invariant systems by considering the Fourier transform of the curvature function in equation (61). In terms of the Wannier state $| {\bf{R}}n\rangle$ of the nth filled-band localized at home cell R, introduced in equation (7), the derivatives of the Pfaffian take the form

$$\begin{eqnarray}\begin{array}{rcl} & & {\partial }_{i}\mathrm{Pf}[m({\bf{k}})]=-{\rm{i}}\displaystyle \sum _{{\bf{R}}}{{\rm{e}}}^{-{\rm{i}}{\bf{k}}\cdot {\bf{R}}}\langle {\bf{R}}1| {R}_{i}T| {\bf{0}}2\rangle ,\\ & & {\partial }_{i}^{2}\mathrm{Pf}[m({\bf{k}})]=-\displaystyle \sum _{{\bf{R}}}{{\rm{e}}}^{-{\rm{i}}{\bf{k}}\cdot {\bf{R}}}\langle {\bf{R}}1| {R}_{i}^{2}T| {\bf{0}}2\rangle .\end{array}\end{eqnarray} \tag{ 63 }$$

The Fourier transform of the Laplacian of the Pfaffian is then

$$\begin{eqnarray}&&{\tilde{F}}_{{\rm{TR}}}({\bf{R}})=\int \displaystyle \frac{{{\rm{d}}}^{D}{\bf{k}}}{{\left(2\pi \right)}^{D}}\,{{\rm{e}}}^{{\rm{i}}{\bf{k}}\cdot {\bf{R}}}{{\rm{\nabla }}}_{{\bf{k}}}^{2}\,\mathrm{Pf}[m({\bf{k}})]=-\langle {\bf{R}}1| {R}^{2}\,T| {\bf{0}}2\rangle ,\end{eqnarray} \tag{ 64 }$$

where ${R}^{2}={\sum }_{i=1}^{D}{R}_{i}^{2}$. This corresponds to a correlation function that measures the overlap of the Wannier state of the first filled-band centered at R and that of the second filled-band centered at the origin, as a matrix element of R²T.

The Laplacian of the Pfaffian in our continuous model in spherical coordinates is given by

$$\begin{eqnarray}&&{{\rm{\nabla }}}_{{\bf{k}}}^{2}\mathrm{Pf}[m]=\displaystyle \frac{1}{{k}^{2}}\displaystyle \frac{\partial }{\partial k}\left[{k}^{2}\displaystyle \frac{\partial }{\partial k}\left(-\displaystyle \frac{M}{d}\right)\right]=\displaystyle \frac{(n-{n}^{2}){k}^{4n-2}M+(2n+1){{nk}}^{2n-2}{M}^{3}}{{\left({k}^{2n}+{M}^{2}\right)}^{5/2}},\end{eqnarray} \tag{ 65 }$$

whose Fourier transform reads

$$\begin{eqnarray}&&{\tilde{F}}_{{\rm{TR}}}({\bf{R}})=\int \displaystyle \frac{{{\rm{d}}}^{3}{\bf{k}}}{{\left(2\pi \right)}^{3}}{{\rm{e}}}^{{\rm{i}}{\bf{k}}\cdot {\bf{R}}}{{\rm{\nabla }}}_{{\bf{k}}}^{2}\mathrm{Pf}[m]=\displaystyle \frac{1}{2{\pi }^{2}R}{\int }_{0}^{\infty }{\rm{d}}{k}\,k\,{\rm{\sin }}{kR}{{\rm{\nabla }}}_{{\bf{k}}}^{2}\mathrm{Pf}[m].\end{eqnarray} \tag{ 66 }$$

At any n, the correlation function is a decaying function of argument R/ξ, with a correlation length $\xi \propto | M{| }^{-\tfrac{1}{n}}$ that has critical exponent ν = 1/n. At n = 1, the correlation function monotonically decays, whereas at n > 1, the correlation function oscillates and decays, as shown in figure 4. These features remain true for all the correlation functions discussed in the present work.

Zoom In Zoom Out Reset image size

For 3D class CII that has T² = −1 and C² = −1, the minimal model is a 8 × 8 Dirac model [7]. We choose the chiral representation, in which the seven Γ-matrices are given by [8]

$$\begin{eqnarray}\begin{array}{rcl}{{\rm{\Gamma }}}^{a} & = & {{\rm{\Gamma }}}_{4\times 4}^{a}\otimes {\eta }_{x},\ \ \ \mathrm{for}\ a=1\sim 4,\\ {{\rm{\Gamma }}}^{5} & = & {I}_{4\times 4}\otimes {\eta }_{y},\ \ \ {{\rm{\Gamma }}}^{6}={I}_{4\times 4}\otimes {\eta }_{z},\\ {{\rm{\Gamma }}}^{7} & = & {\left(-i\right)}^{3}{{\rm{\Gamma }}}^{1}{{\rm{\Gamma }}}^{2}...{{\rm{\Gamma }}}^{6},\end{array}\end{eqnarray} \tag{ 67 }$$

where ${{\rm{\Gamma }}}_{4\times 4}^{a}$ are those in equation (26). ${{\rm{\Gamma }}}^{6}={I}_{\sigma }\otimes {I}_{\tau }\otimes {\eta }_{z}=S$ is used to implement the chiral symmetry. The Hamiltonian written in the basis of the other six Γ-matrices is block-off-diagonal

$$\begin{eqnarray}H({\bf{k}})=\displaystyle \sum _{i=1,2,3,4,5,7}{d}_{i}{{\rm{\Gamma }}}^{i}=\left(\begin{array}{cc}0 & D\\ {D}^{\dagger } & 0\end{array}\right).\end{eqnarray} \tag{ 68 }$$

We will denote σ_i, τ_i, and η_i as the Pauli matrices for the spin, PH, and valley grading. The TR and PH operators may be chosen as $T={\rm{i}}{\sigma }_{y}K\otimes {I}_{\tau }\otimes {I}_{\eta }$ and $C=-{\rm{i}}{\sigma }_{y}K\otimes {I}_{\tau }\otimes {\eta }_{z}$. The block-off-diagonal part of the Hamiltonian in equation (68) expressed in terms of the d-vector yields

$$\begin{eqnarray}\begin{array}{rcl}D & = & \left(\begin{array}{cccc}h & & g & f\\ & h & {f}^{* } & -{g}^{* }\\ {g}^{* } & f & -{h}^{* } & \\ {f}^{* } & -g & & -{h}^{* }\end{array}\right),\\ f & = & {d}_{1}-{{\rm{i}}{d}}_{2},\ \ \ g={d}_{3}+{{\rm{i}}{d}}_{7},\ \ \ h={d}_{4}-{{\rm{i}}{d}}_{5}.\end{array}\end{eqnarray} \tag{ 69 }$$

The TR symmetry requires that

$$\begin{eqnarray}&&{d}_{i}({\bf{k}})=-{d}_{i}(-{\bf{k}})\ \ \ \mathrm{for}\ i=1,2,3,5,7,\ \ \ {d}_{4}({\bf{k}})={d}_{4}(-{\bf{k}}).\end{eqnarray} \tag{ 70 }$$

Consequently, the mass term is d₄ = M. The PH symmetry also gives the same condition.

We parameterize the 8 × 8 rotational operator by

$$\begin{eqnarray}{C}_{m}=\mathrm{diag}({\overline{\alpha }}_{p},{\overline{\alpha }}_{q},{\overline{\alpha }}_{r},{\overline{\alpha }}_{s}),\ \ \ {\overline{\alpha }}_{p}=\left(\begin{array}{cc}{\alpha }_{{p}_{1}} & \\ & {\alpha }_{{p}_{2}}\end{array}\right),\end{eqnarray} \tag{ 71 }$$

with ${\alpha }_{{p}_{i}}={{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}\left({p}_{i}+\tfrac{F}{2}\right)}$. The coexistence of TR and rotational symmetry $\left[{C}_{m},T\right]=0$ requires that ${\alpha }_{{p}_{1}}={\alpha }_{{p}_{2}}^{* }\equiv {\alpha }_{p}$. The coexistence of PH and rotational symmetry $\left[{C}_{m},C\right]=0$ also requires the same thing. As a result, the ${\overline{\alpha }}_{p}$ in equation (71) takes the form ${\overline{\alpha }}_{p}=\mathrm{diag}({\alpha }_{p},{\alpha }_{p}^{* })$.

We now examine the rotational symmetry. To simplify the calculation, we rewrite the Hamiltonian in equation (68) by

$$\begin{eqnarray}\begin{array}{rcl} & & H({\bf{k}})=\left(\begin{array}{cccc} & & {D}_{11} & {D}_{12}\\ & & {D}_{21} & {D}_{22}\\ {D}_{11}^{* } & {D}_{12} & & \\ {D}_{21} & {D}_{22}^{* } & & \end{array}\right),\\ & & {D}_{11}=\left(\begin{array}{cc}h & \\ & h\end{array}\right)=-{D}_{22}^{* },\ \ \ {D}_{12}=\left(\begin{array}{cc}g & f\\ {f}^{* } & -{g}^{* }\end{array}\right)={D}_{21}^{\dagger }.\end{array}\end{eqnarray} \tag{ 72 }$$

The rotational symmetry in equation (4) gives the constraint

$$\begin{eqnarray}\begin{array}{rcl} & & {\overline{\alpha }}_{p}{D}_{12}({\bf{k}}){\overline{\alpha }}_{s}^{* }={\overline{\alpha }}_{r}{D}_{12}{\overline{\alpha }}_{q}^{* }={D}_{12}({R}_{m}{\bf{k}}),\\ & & {\overline{\alpha }}_{q}{D}_{21}({\bf{k}}){\overline{\alpha }}_{r}^{* }={\overline{\alpha }}_{s}{D}_{21}{\overline{\alpha }}_{p}^{* }={D}_{21}({R}_{m}{\bf{k}}),\\ & & {\overline{\alpha }}_{p}{D}_{11}{\overline{\alpha }}_{r}^{* }={D}_{11}({R}_{m}{\bf{k}}),\\ & & {\overline{\alpha }}_{q}{D}_{22}{\overline{\alpha }}_{s}^{* }={D}_{22}({R}_{m}{\bf{k}}).\end{array}\end{eqnarray} \tag{ 73 }$$

The first two equations in equation (73) yields ${\overline{\alpha }}_{p}={\overline{\alpha }}_{q}={\overline{\alpha }}_{r}={\overline{\alpha }}_{s}$, so the final form of the rotational operator is

$$\begin{eqnarray}{C}_{m}=\mathrm{diag}({\overline{\alpha }}_{p},{\overline{\alpha }}_{p},{\overline{\alpha }}_{p},{\overline{\alpha }}_{p}),\ \ \ {\overline{\alpha }}_{p}=\left(\begin{array}{cc}{\alpha }_{p} & \\ & {\alpha }_{p}^{* }\end{array}\right).\end{eqnarray} \tag{ 74 }$$

Using the parameterization

$$\begin{eqnarray}&&f={d}_{1}-{{\rm{i}}{d}}_{2}={k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}},\ \ \ g={d}_{3}+{{\rm{i}}{d}}_{7}={d}_{3}={k}_{z}^{n},\ \ \ h={d}_{4}-{{\rm{i}}{d}}_{5}={d}_{4}=M,\end{eqnarray} \tag{ 75 }$$

such that the dispersion satisfies equation (29), the constraint on f yields

$$\begin{eqnarray}&&{n}_{+}-{n}_{-}={xm}+2p+F\in 2{\mathbb{Z}}+1,\end{eqnarray} \tag{ 76 }$$

because $n={n}_{+}+{n}_{-}\in 2{\mathbb{Z}}+1$ due to equation (70). For fermions F = 1, the condition in equation (76) can be by any of the m = 2, 3, 4, 6, where as for bosons F = 0 it can only be satisfied at m = 3.

We now introduce a winding number and the corresponding Wannier state correlation function. Defining the 4 × 4 q-matrix from the off-diagonal part of the Hamiltonian as [7]

$$\begin{eqnarray}&&q({\bf{k}})=-\displaystyle \frac{1}{d}\left({d}_{\mu }{\alpha }_{\mu }+M\beta \right)={q}^{-1}({\bf{k}})={q}^{\dagger }({\bf{k}}),\end{eqnarray} \tag{ 77 }$$

we found that ${\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{q}^{\dagger }{\partial }_{\mu }{{qq}}^{\dagger }{\partial }_{\nu }{{qq}}^{\dagger }{\partial }_{\rho }q\right]=0$, and hence the simple version of the winding number in equation (37) vanishes. However, if we add a ${{\rm{\Gamma }}}_{4\times 4}^{5}$ (defined in equation (26)) into the definition

$$\begin{eqnarray}&&{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{{\rm{\Gamma }}}_{4\times 4}^{5}{q}^{\dagger }{\partial }_{\mu }{{qq}}^{\dagger }{\partial }_{\nu }{{qq}}^{\dagger }{\partial }_{\rho }q\right]=24M\displaystyle \frac{{k}_{\perp }^{2(n-1)}{k}_{z}^{n-1}({n}_{+}-{n}_{-}){n}^{2}}{{\left({k}^{2n}+{M}^{2}\right)}^{2}},\end{eqnarray} \tag{ 78 }$$

then it yields a non-vanishing winding number ${ \mathcal C }^{\prime}$ (see equation (37))

$$\begin{eqnarray}\begin{array}{rcl}{ \mathcal C }^{\prime} & = & \displaystyle \frac{\pi }{3}\displaystyle \int \displaystyle \frac{{{\rm{d}}}^{3}{\bf{k}}}{{\left(2\pi \right)}^{3}}{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{{\rm{\Gamma }}}_{4\times 4}^{5}{q}^{\dagger }{\partial }_{\mu }{{qq}}^{\dagger }{\partial }_{\nu }{{qq}}^{\dagger }{\partial }_{\rho }q\right]\\ & = & \left\{\begin{array}{ll}({n}_{+}-{n}_{-})\mathrm{Sgn}(M)\left[2{nB}\left(\tfrac{n}{2},n\right)\right] & \mathrm{if}\ n\in 2{\mathbb{Z}}+1\\ 0 & \mathrm{if}\ n\in 2{\mathbb{Z}}\end{array}\right.,\end{array}\end{eqnarray} \tag{ 79 }$$

where B(a, b) is the Beta function. Although this winding number does not correspond to the true ${{\mathbb{Z}}}_{2}$ invariant in this class, the integrand in equation (79) satisfies all the critical properties of the curvature function discussed in section 2.2. In addition, from the eigenstates of the filled bands

$$\begin{eqnarray}\begin{array}{rcl}| {u}_{1-}\rangle & = & \displaystyle \frac{1}{\sqrt{2}d}\left(\begin{array}{c}-{d}_{1}+{{\rm{i}}{d}}_{2}\\ {d}_{3}\\ 0\\ M\\ 0\\ 0\\ 0\\ d\end{array}\right),\ \ \ | {u}_{2-}\rangle =\displaystyle \frac{1}{\sqrt{2}d}\left(\begin{array}{c}-{d}_{3}\\ -{d}_{1}-{{\rm{i}}{d}}_{2}\\ M\\ 0\\ 0\\ 0\\ d\\ 0\end{array}\right),\\ | {u}_{3-}\rangle & = & \displaystyle \frac{1}{\sqrt{2}d}\left(\begin{array}{c}0\\ -M\\ -{d}_{1}+{{\rm{i}}{d}}_{2}\\ {d}_{3}\\ 0\\ d\\ 0\\ 0\end{array}\right),\ \ \ | {u}_{4-}\rangle =\displaystyle \frac{1}{\sqrt{2}d}\left(\begin{array}{c}-M\\ 0\\ -{d}_{3}\\ -{d}_{1}-{{\rm{i}}{d}}_{2}\\ d\\ 0\\ 0\\ 0\end{array}\right).\end{array}\end{eqnarray} \tag{ 80 }$$

The Berry connection ${{ \mathcal A }}_{\mu }^{{ab}}=\langle {u}_{a-}| {\partial }_{\mu }| {u}_{b-}\rangle$ renders

$$\begin{eqnarray}&&{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{{\rm{\Gamma }}}_{4\times 4}^{5}{q}^{\dagger }{\partial }_{\mu }{{qq}}^{\dagger }{\partial }_{\nu }{{qq}}^{\dagger }{\partial }_{\rho }q\right]=-8{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{{\rm{\Gamma }}}_{4\times 4}^{5}{{ \mathcal A }}_{\mu }{{ \mathcal A }}_{\nu }{{ \mathcal A }}_{\rho }\right]=\displaystyle \frac{24M}{{d}^{4}}{\epsilon }^{\mu \nu \rho }{\partial }_{x}{d}_{\mu }{\partial }_{y}{d}_{\nu }{\partial }_{z}{d}_{\rho },\end{eqnarray} \tag{ 81 }$$

similar to equation (41). Therefore the formalism of constructing correlation function in section 3.1 can be directly applied to 3D class CII, with the inclusion of ${{\rm{\Gamma }}}_{4\times 4}^{5}$ in the argument of all the traces $\mathrm{Tr}[...]\to \mathrm{Tr}[{{\rm{\Gamma }}}_{4\times 4}^{5}...]$. We denote the resulting correlation function as ${\tilde{F}}_{3D}^{{\prime} }({\bf{R}})$ in table 1.

For 3D class CI that has T² = 1, C² = −1, we use the minimal 8 × 8 model that writes the Hamiltonian in the chiral basis [7]

$$\begin{eqnarray}H=\left(\begin{array}{cc} & D\\ {D}^{\dagger } & \end{array}\right),\ \ \ D({\bf{k}})={\rm{i}}{\sigma }_{y}\beta \left({d}_{i}{\alpha }_{i}-{{\rm{i}}{d}}_{4}{\gamma }^{5}\right)={D}^{T}(-{\bf{k}}),\end{eqnarray} \tag{ 82 }$$

where the 4 × 4 matrices $\left\{{\alpha }_{i},\beta ,{\gamma }^{5}\right\}$ are those in equation (26). The D matrix expressed in terms of the d = (d₁, d₂, d₃, d₄) vector is

$$\begin{eqnarray}D=\left(\begin{array}{cc} & {D}_{12}\\ {D}_{21} & \end{array}\right)=\left(\begin{array}{cccc} & & {f}^{* } & -g\\ & & -{g}^{* } & -f\\ -{f}^{* } & g & & \\ {g}^{* } & f & & \end{array}\right),\end{eqnarray} \tag{ 83 }$$

with $f={d}_{1}-{{\rm{i}}{d}}_{2}$ and $g={d}_{3}+{{\rm{i}}{d}}_{4}$. The 8 × 8 TR and PH operators are $T=I\otimes I\otimes {\sigma }_{x}K$ and $C=I\otimes I\otimes (-{\rm{i}}{\sigma }_{y})K$. The TR and PH symmetry require

$$\begin{eqnarray}&&{d}_{i}({\bf{k}})=-{d}_{i}(-{\bf{k}})\ \ \ \mathrm{for}\ i=1,2,3,\ \ \ {d}_{4}({\bf{k}})={d}_{4}(-{\bf{k}}).\end{eqnarray} \tag{ 84 }$$

Consequently, only d₄ = M can be the mass term.

Now we consider the 8 × 8 rotational operator in equation (71). The coexistence of TR and rotational symmetry, and the coexistence of PH and rotational symmetry imply ${\overline{\alpha }}_{p}={\overline{\alpha }}_{r}^{* }$ and ${\overline{\alpha }}_{q}={\overline{\alpha }}_{s}^{* }$, so the rotational operator is constrained to take the form

$$\begin{eqnarray}&&{C}_{m}=\mathrm{diag}({\overline{\alpha }}_{p},{\overline{\alpha }}_{q},{\overline{\alpha }}_{p}^{* },{\overline{\alpha }}_{q}^{* }).\end{eqnarray} \tag{ 85 }$$

We use equations (82) and (83) to write the Hamiltonian into the form

$$\begin{eqnarray}H({\bf{k}})=\left(\begin{array}{cccc} & & & {D}_{12}({\bf{k}})\\ & & -{D}_{12}({\bf{k}}) & \\ & -{D}_{12}^{\dagger }({\bf{k}}) & & \\ {D}_{12}^{\dagger }({\bf{k}}) & & & \end{array}\right).\end{eqnarray} \tag{ 86 }$$

Writing down each component in the 8 × 8 matrix explicitly and using the definition in equation (71), we see that rotational symmetry in equation (4) requires

$$\begin{eqnarray}\begin{array}{rcl} & & {\alpha }_{{p}_{1}}{f}^{* }({\bf{k}}){\alpha }_{{q}_{1}}={\alpha }_{{p}_{2}}^{* }{f}^{* }({\bf{k}}){\alpha }_{{q}_{2}}^{* }={f}^{* }({R}_{m}{\bf{k}}),\\ & & {\alpha }_{{p}_{2}}g({\bf{k}}){\alpha }_{{q}_{1}}={\alpha }_{{p}_{1}}^{* }g({\bf{k}}){\alpha }_{{q}_{2}}^{* }={\alpha }_{{p}_{1}}g({\bf{k}}){\alpha }_{{q}_{2}}=g({R}_{m}{\bf{k}}).\end{array}\end{eqnarray} \tag{ 87 }$$

To satisfy the Dirac form of the dispersion, we proceed to parameterize

$$\begin{eqnarray}&&f={d}_{1}-{{\rm{i}}{d}}_{2}={k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}},\ \ \ g={d}_{3}+{\rm{i}}{M}={k}_{z}^{n}+{\rm{i}}{M},\end{eqnarray} \tag{ 88 }$$

with $n={n}_{+}+{n}_{-}\in 2{\mathbb{Z}}+1$ due to equation (84). Since g is invariant under rotation g(k) = g(R_mk), it requires that ${\alpha }_{{p}_{2}}{\alpha }_{{q}_{1}}={\alpha }_{{p}_{1}}^{* }{\alpha }_{{q}_{2}}^{* }={\alpha }_{{p}_{1}}{\alpha }_{{q}_{2}}=1$. A detailed analysis shows that these conditions cannot be satisfied at any nonzero m for either fermions F = 1 or bosons F = 0, hence 3D class CI is not compatible with the m-fold rotational symmetry.

We proceed to discuss higher-order Dirac model n ≥ 1 without considering the rotational symmetry. Defining the q-matrix

$$\begin{eqnarray}&&q({\bf{k}})=-\displaystyle \frac{D({\bf{k}})}{d},\ \ \ {q}^{\dagger }({\bf{k}})={q}^{-1}({\bf{k}})=-\displaystyle \frac{{D}^{\dagger }({\bf{k}})}{d},\end{eqnarray} \tag{ 89 }$$

the winding number is calculated from the integration of

$$\begin{eqnarray}&&{\epsilon }^{\mu \nu \rho }\mathrm{Tr}\left[{q}^{\dagger }{\partial }_{\mu }{{qq}}^{\dagger }{\partial }_{\nu }{{qq}}^{\dagger }{\partial }_{\rho }q\right]=-M\displaystyle \frac{24{k}_{\perp }^{2(n-1)}{k}_{z}^{n-1}({n}_{+}-{n}_{-}){n}^{2}}{{\left({k}^{2n}+{M}^{2}\right)}^{2}},\end{eqnarray} \tag{ 90 }$$

which is twice of that in equation (35) for 3D class AIII. As a result, the winding number defined in equation (37) is always an even number ${ \mathcal C }\in 2{\mathbb{Z}}$. The notion of correlation function and critical exponents then follows the discussion in section 3.1.

The minimal model for a 2D class A model is described by the Hamiltonian

$$\begin{eqnarray}&&H({\bf{k}})=f({\bf{k}}){\sigma }_{+}+{f}^{* }({\bf{k}}){\sigma }_{-}+M{\sigma }_{z}.\end{eqnarray} \tag{ 91 }$$

The rotational operation reads C_m = diag(α_p, α_q), and hence the rotational symmetry in equation (4) requires

$$\begin{eqnarray}&&{\alpha }_{p}{\alpha }_{q}^{* }f({k}_{\pm })=f\left({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}}\right).\end{eqnarray} \tag{ 92 }$$

The parameterization

$$\begin{eqnarray}&&f({k}_{\pm })={k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}},\end{eqnarray} \tag{ 93 }$$

yields

$$\begin{eqnarray}&&{n}_{+}-{n}_{-}={xm}+(p-q)\in {\mathbb{Z}},\end{eqnarray} \tag{ 94 }$$

since there is no nonspatial symmetry in this class to constrain the order of band crossing, one has $n={n}_{+}+{n}_{-}\in {\mathbb{Z}}$. This condition can be satisfied by any m for either bosons or fermions.

The topological invariant ${ \mathcal C }$ in 2D class A is given by the integration of the Berry curvature ${\rm{\Omega }}({k}_{x},{k}_{y})$. The valance band Berry connection as ${\bf{A}}=\langle {u}_{{\bf{k}}}| {\rm{i}}{{\boldsymbol{\nabla }}}_{{\bf{k}}}| {u}_{{\bf{k}}}\rangle$, which yields the Berry curvature

$$\begin{eqnarray}&&{\rm{\Omega }}({k}_{x},{k}_{y})={\partial }_{{k}_{k}}{A}_{y}-{\partial }_{{k}_{y}}{A}_{x}=\displaystyle \frac{1}{2{d}^{3}}\left\{{\bf{d}}\cdot {\partial }_{{k}_{x}}{\bf{d}}\times {\partial }_{{k}_{y}}{\bf{d}}\right\}.\end{eqnarray} \tag{ 95 }$$

The derivatives may be converted into the polar coordinates by

$$\begin{eqnarray}&&{\partial }_{{k}_{x}}{d}_{i}=\cos \phi {\partial }_{k}{d}_{i}-\displaystyle \frac{\sin \phi }{k}{\partial }_{\phi }{d}_{i},\ \ \ {\partial }_{{k}_{y}}{d}_{i}=\sin \phi {\partial }_{k}{d}_{i}+\displaystyle \frac{\cos \phi }{k}{\partial }_{\phi }{d}_{i},\end{eqnarray} \tag{ 96 }$$

which yields the Berry curvature

$$\begin{eqnarray}&&{\rm{\Omega }}({k}_{x},{k}_{y})=\displaystyle \frac{({n}_{+}^{2}-{n}_{-}^{2}){{Mk}}^{2n-2}}{2{\left[{k}^{2n}+{M}^{2}\right]}^{3/2}}.\end{eqnarray} \tag{ 97 }$$

The topological invariant is calculated by converting the integration over the BZ into an integration over the entire polar plane

$$\begin{eqnarray}&&{ \mathcal C }={\int }_{{\rm{BZ}}}\displaystyle \frac{{{\rm{d}}}^{2}{\bf{k}}}{2\pi }{\rm{\Omega }}({k}_{x},{k}_{y})\approx \displaystyle \frac{1}{2\pi }{\int }_{0}^{2\pi }{\rm{d}}\phi {\int }_{0}^{\infty }{\rm{d}}{k}\,k\,{\rm{\Omega }}({k}_{x},{k}_{y})=\displaystyle \frac{{n}_{+}-{n}_{-}}{2}\mathrm{Sgn}(M).\end{eqnarray} \tag{ 98 }$$

The above invariant is normalized such that the change of topological invariant across the topological phase transition ${\rm{\Delta }}{ \mathcal C }={ \mathcal C }(M\gt 0)-{ \mathcal C }(M\lt 0)={n}_{+}-{n}_{-}$ is consistent with ${\rm{\Delta }}{ \mathcal C }=1$ for the usual linear Dirac model (n₊, n₋) = (1, 0).

The Fourier transform of the Berry curvature can be expressed in terms of the Wannier states by [6, 39–41]

$$\begin{eqnarray}\begin{array}{rcl}{\tilde{F}}_{2{\rm{D}}}({\bf{R}}) & = & {\displaystyle \int }_{{\rm{BZ}}}\displaystyle \frac{{{\rm{d}}}^{2}{\bf{k}}}{{\left(2\pi \right)}^{2}}{{\rm{e}}}^{{\rm{i}}{\bf{k}}\cdot {\bf{R}}}{\rm{\Omega }}({k}_{x},{k}_{y})=-{\rm{i}}\langle {\bf{R}}| {\left({\bf{R}}\times \hat{r}\right)}_{z}| {\bf{0}}\rangle \\ & = & -{\rm{i}}\displaystyle \int {{\rm{d}}}^{2}{\bf{r}}({R}_{x}{r}_{y}-{R}_{y}{r}_{x}){W}^{* }({\bf{r}}-{\bf{R}})W({\bf{r}}),\end{array}\end{eqnarray} \tag{ 99 }$$

which measures the overlap of the Wannier functions centering at ${\bf{R}}$ and at the origin ${\bf{0}}$, sandwiched by the factor ${R}_{x}{r}_{y}-{R}_{y}{r}_{x}$. Within our continuous model formalism

$$\begin{eqnarray}&&{\tilde{F}}_{2{\rm{D}}}({\bf{R}})\approx \displaystyle \frac{1}{2\pi }{\int }_{0}^{\infty }{\rm{d}}{k}\,k\,{J}_{0}({Rk}){\rm{\Omega }}({k}_{x},{k}_{y}),\end{eqnarray} \tag{ 100 }$$

where J₀(Rk) is the 0th-order Bessel function. Using the Berry curvature in equation (97), we find that the Fourier transform at any $\left\{{n}_{+},{n}_{-}\right\}$ gives a decaying function λ(R) with a correlation length proportional to $\xi \propto {\left|M\right|}^{-1/n}\propto | M{| }^{-\nu }$. Moreover, because $\left\{{n}_{+},{n}_{-}\right\}$ are positive numbers, and hence $n\geqslant | {n}_{+}-{n}_{-}|$, we have the inequality

$$\begin{eqnarray}&&\displaystyle \frac{1}{\nu }\geqslant | {\rm{\Delta }}{ \mathcal C }| ,\end{eqnarray} \tag{ 101 }$$

i.e. the order of band crossing can only be greater than or equal to the jump of topological invariant. Finally, we remark that out of all the correlation functions introduced in the present work, only ${\tilde{F}}_{2{\rm{D}}}({\bf{R}})$ is gauge invariant (because Berry curvature is gauge invariant) and hence measurable. It can be measured by, for instance, performing a Fourier transform of the Berry curvature measured in the cold atom experiments [42–44], from which the critical behavior discussed above may be verified.

We now discuss 2D class C that has only PH symmetry C² = −1. The minimal model is a 2 × 2 Dirac model with the PH operator C = σ_yK. Using the rotational operator C_m = diag(α_p, α_q), the requirement [C, C_m] = 0 yields ${\alpha }_{q}={\alpha }_{p}^{* }$, and hence the rotational operator has to take the form ${C}_{m}=\mathrm{diag}({\alpha }_{p},{\alpha }_{p}^{* })$. On the other hand, the PH symmetry requires ${d}_{i}({\bf{k}})={d}_{i}(-{\bf{k}})$, i.e. all the three components are even in momentum. Without loss of generality we choose d₃ to be the mass term, and parameterize the Hamiltonian by

$$\begin{eqnarray}H({\bf{k}})=\left(\begin{array}{cc}M & f\\ {f}^{* } & -M\end{array}\right),\ \ \ f={d}_{1}-{{\rm{i}}{d}}_{2}={k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}}.\end{eqnarray} \tag{ 102 }$$

The requirement that d_i(k) is even in momentum dictates $n\in 2{\mathbb{Z}}$, indicating that the order of band crossing can only be even.

The rotational symmetry in equation (4) renders

$$\begin{eqnarray}&&{\left({\alpha }_{p}\right)}^{2}f({k}_{\pm })=f({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}}),\end{eqnarray} \tag{ 103 }$$

and consequently

$$\begin{eqnarray}&&{n}_{+}-{n}_{-}={xm}+2p+F\in 2{\mathbb{Z}},\end{eqnarray} \tag{ 104 }$$

must be satisfied. As a result, when the system is fermionic F = 1, it can only have m = 3 fold rotational symmetry, whereas the bosonic case F = 0 can be satisfied with any m.

The topological invariant is calculated from the same Chern number as in 2D class A in section 4.1. Consequently, the Wannier state correlation function and the critical exponent of the correlation length are introduced exactly in the same way, using the formalism from equation (95) to (101). However, it should be emphasized that in 2D class C here the jump of topological invariant ${\rm{\Delta }}{ \mathcal C }$, order of band crossing, and inverse of the critical exponent are all even numbers, in contrast to 2D class A in which they can be either even or odd.

The calculation of 2D class D, where the PH symmetry satisfies C² = 1, is similar to that in 2D class C presented in section 4.2. We follow the formalism in [8] and use the PH operator $C={\sigma }_{x}K$. The PH symmetry requires

$$\begin{eqnarray}&&{d}_{1}({\bf{k}})=-{d}_{1}(-{\bf{k}}),\ \ \ {d}_{2}({\bf{k}})=-{d}_{2}(-{\bf{k}}),\ \ \ {d}_{3}({\bf{k}})={d}_{3}(-{\bf{k}}).\end{eqnarray} \tag{ 105 }$$

We thereby choose d₃ = M as the mass term. Using the Hamiltonian and the lowest order expansion of f in equation (102), we see that the order of band crossing is odd $n\in 2{\mathbb{Z}}+1$.

The commutation of rotational operator ${C}_{m}=\mathrm{diag}({\alpha }_{p},{\alpha }_{q})$ with the PH operator [C, C_m] = 0 yields ${\alpha }_{q}={\alpha }_{p}^{* }$. As a result, the rotational eigenvalues take the form ${C}_{m}=\mathrm{diag}({\alpha }_{p},{\alpha }_{p}^{* })$. Equation (4) renders the following condition on f

$$\begin{eqnarray}&&{\left({\alpha }_{p}\right)}^{2}f({k}_{\pm })={{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}(2p+F)}f({k}_{\pm })=f({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}})={k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}}{{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}({n}_{+}-{n}_{-})}.\end{eqnarray} \tag{ 106 }$$

Therefore,

$$\begin{eqnarray}&&{n}_{+}-{n}_{-}={xm}+2p+F\in 2{\mathbb{Z}}+1,\end{eqnarray} \tag{ 107 }$$

must be satisfied. If the system is fermionic F = 1, this condition can be satisfied in any m-fold rotational symmetry, whereas if the system is bosonic F = 0, then only m = 3 is allowed.

The topological invariant is calculated from the Chern number in section 4.1, and the Wannier state correlation function and the critical exponent of the correlation length are introduced by the formalism from equation (95) to (101). However, it should be emphasized that in 2D class D here the jump of topological invariant ${\rm{\Delta }}{ \mathcal C }$, order of band crossing, and inverse of the critical exponent are all odd numbers.

In this section we discuss 2D class AII that has T² = −1, with the consideration of m-fold rotational symmetry [24]. The minimal model is a 4 × 4 Dirac Hamiltonian. We use the representation for the Bernevig–Hughes–Zhang model [45, 46]

$$\begin{eqnarray}&&{{\rm{\Gamma }}}_{a}=\left\{{s}_{x}\otimes {\sigma }_{z},{s}_{y}\otimes I,{s}_{z}\otimes I,{s}_{x}\otimes {\sigma }_{x},{s}_{x}\otimes {\sigma }_{y}\right\}.\end{eqnarray} \tag{ 108 }$$

The TR operator in this basis reads $T=-{\rm{i}}{I}\otimes {\sigma }_{y}K$. The coexistence of TR and rotational symmetry ${C}_{m}=\mathrm{diag}({\alpha }_{p},{\alpha }_{q},{\alpha }_{r},{\alpha }_{s})$ requires ${\alpha }_{r}={\alpha }_{p}^{* }$ and ${\alpha }_{s}={\alpha }_{q}^{* }$. Consequently, ${C}_{m}=\mathrm{diag}({\alpha }_{p},{\alpha }_{q},{\alpha }_{p}^{* },{\alpha }_{q}^{* })$, in agreement with [24].

Now the representation in equation (108) indicates that the Hamiltonian takes the form

$$\begin{eqnarray}H=\left(\begin{array}{cccc}M & f & & g\\ {f}^{* } & -M & {g}^{* } & \\ & g & M & -{f}^{* }\\ {g}^{* } & & -f & -M\end{array}\right),\end{eqnarray} \tag{ 109 }$$

where we have defined f and g, and we assume their leading order expansions are

$$\begin{eqnarray}&&f={d}_{1}-{{\rm{i}}{d}}_{2}={k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}},\ \ \ g={d}_{4}-{{\rm{i}}{d}}_{5}={k}_{+}^{{{\ell }}_{+}}{k}_{-}^{{{\ell }}_{-}},\ \ \ {d}_{3}=M.\end{eqnarray} \tag{ 110 }$$

Note that TR invariance requires

$$\begin{eqnarray}&&{d}_{i}(-{\bf{k}})=-{d}_{i}({\bf{k}})\ \ \ \mathrm{for}\ i=1,2,4,5,\ \ \ {d}_{3}({\bf{k}})={d}_{3}(-{\bf{k}}),\end{eqnarray} \tag{ 111 }$$

so we use d₃ = M as the mass term. The oddness in equation (111) requires the expansion in equation (110) to satisfy $n\in 2{\mathbb{Z}}+1$ and ${{\ell }}_{+}+{{\ell }}_{-}=2{\mathbb{Z}}+1$. Equation (4) yields the following constraints on f and g

$$\begin{eqnarray}\begin{array}{rcl} & & {\alpha }_{p}{\alpha }_{q}^{* }f({k}_{\pm })={{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}(p-q)}f({k}_{\pm })=f({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}2\pi /m})={k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}}{{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}({n}_{+}-{n}_{-})},\\ & & {\alpha }_{p}{\alpha }_{q}g({k}_{\pm })={{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}(p+q+1)}g({k}_{\pm })=g({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}2\pi /m})={k}_{+}^{{{\ell }}_{+}}{k}_{-}^{{{\ell }}_{-}}{{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}({{\ell }}_{+}-{{\ell }}_{-})},\\ & & {\alpha }_{p}{\alpha }_{q}{g}^{* }({k}_{\pm })={{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}(p+q+1)}{g}^{* }({k}_{\pm })=g({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}2\pi /m})={k}_{-}^{{{\ell }}_{+}}{k}_{+}^{{{\ell }}_{-}}{{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}({{\ell }}_{-}-{{\ell }}_{+})}.\end{array}\end{eqnarray} \tag{ 112 }$$

The last two equations of the above equation can never be satisfied at the same time, so we conclude that $g={d}_{4}-{{\rm{i}}{d}}_{5}=0$. Thus rotational symmetry requires that Γ₄ and Γ₅ terms cannot exist. The constaint on f then gives

$$\begin{eqnarray}&&{n}_{+}-{n}_{-}={xm}+(p-q)\in 2{\mathbb{Z}}+1,\end{eqnarray} \tag{ 113 }$$

which can be satisfied by any m for either fermionic or bosonic cases.

Denoting $d=\sqrt{{d}_{1}^{2}+{d}_{2}^{2}+{d}_{3}^{2}}$, the four eigenstates in a gauge choice that makes the Pfaffian real and equal to ±1 at HSPs are

$$\begin{eqnarray}\begin{array}{rcl}| {u}_{\mathrm{1,3}}\rangle & = & \displaystyle \frac{1}{\sqrt{2d(d\mp {d}_{3})}}\left(\begin{array}{c}0\\ 0\\ -{d}_{3}\pm d\\ {d}_{1}-{{\rm{i}}{d}}_{2}\end{array}\right),\\ | {u}_{\mathrm{2,4}}\rangle & = & \displaystyle \frac{1}{\sqrt{2d(d\mp {d}_{3})}}\left(\begin{array}{c}{d}_{3}\mp d\\ {d}_{1}+{{\rm{i}}{d}}_{2}\\ 0\\ 0\end{array}\right),\end{array}\end{eqnarray} \tag{ 114 }$$

where the upper signs are for the occupied states $\left\{1,2\right\}$ which have eigenenergies ${E}_{1}={E}_{2}=-d$, and the lower signs are for the unoccupied states $\left\{3,4\right\}$ which have eigenenergies E₃ = E₄ = d. The Pfaffian of the m-matrix of the two occupied states is

$$\begin{eqnarray}&&\mathrm{Pf}(m)={m}_{12}=\langle {u}_{1}| T| {u}_{2}\rangle =\displaystyle \frac{{d}_{3}}{d}=\displaystyle \frac{M}{{\left({k}^{2n}+{M}^{2}\right)}^{1/2}}.\end{eqnarray} \tag{ 115 }$$

Thus the sign of ${d}_{3}/d$ at HSPs determines the topology of the system according to equation (58).

Using the Pfaffian in equation (115), its Laplacian is

$$\begin{eqnarray}&&{{\rm{\nabla }}}_{{\bf{k}}}^{2}\mathrm{Pf}[m({\bf{k}})]=\displaystyle \frac{1}{k}\displaystyle \frac{\partial }{\partial k}\left(k\displaystyle \frac{\partial }{\partial k}\displaystyle \frac{{d}_{3}}{d}\right)+\displaystyle \frac{1}{{k}^{2}}\displaystyle \frac{{\partial }^{2}}{\partial {\phi }^{2}}\displaystyle \frac{{d}_{3}}{d}=\displaystyle \frac{{n}^{2}{{Mk}}^{4n-2}-2{n}^{2}{M}^{3}{k}^{2n-2}}{{\left[{k}^{2n}+{M}^{2}\right]}^{5/2}}.\end{eqnarray} \tag{ 116 }$$

The correlation function in our continuous model renders

$$\begin{eqnarray}&&{\tilde{F}}_{{\rm{TR}}}({\bf{R}})\approx \displaystyle \frac{1}{2\pi }{\int }_{0}^{\infty }{\rm{d}}{k}\,{{kJ}}_{0}({Rk}){{\rm{\nabla }}}_{{\bf{k}}}^{2}\mathrm{Pf}[m({\bf{k}})],\end{eqnarray} \tag{ 117 }$$

which displays the same behavior as other correlation functions.

The discussion of 2D class DIII follows that in 3D class DIII with a dimensional reduction k_z = 0, using the same Γ-matrices and TR and PH operators described in section 3.2. The analysis from equation (48) to (51) remains true in 2D. We now translate the 3D Hamiltonian in equation (27) into 2D and parameterize f and $g^{\prime}$ by

$$\begin{eqnarray}&&f={k}_{+}^{{n}_{+}}{k}_{-}^{{n}_{-}}={d}_{1}-{{\rm{i}}{d}}_{2},\ \ \ g^{\prime} ={k}_{+}^{{{\ell }}_{+}}{k}_{-}^{{{\ell }}_{-}}-{\rm{i}}{M}={d}_{3}-{\rm{i}}{M}.\end{eqnarray} \tag{ 118 }$$

First let us look at the parameterization of $g^{\prime}$. The constraint on $g^{\prime}$ is $g^{\prime} ({k}_{\pm })=g^{\prime} ({k}_{\pm }{{\rm{e}}}^{\pm {\rm{i}}\tfrac{2\pi }{m}})$. If one parameterizes $g^{\prime} ={k}_{+}^{{{\ell }}_{+}}{k}_{-}^{{{\ell }}_{-}}-{\rm{i}}{M}$, this means ${k}_{+}^{{{\ell }}_{+}}{k}_{-}^{{{\ell }}_{-}}={{\rm{e}}}^{{\rm{i}}\tfrac{2\pi }{m}({{\ell }}_{+}-{{\ell }}_{-})}{k}_{+}^{{{\ell }}_{+}}{k}_{-}^{{{\ell }}_{-}}$, and hence ${{\ell }}_{+}={{\ell }}_{-}={\ell }/2$, which implies $g={k}^{{\ell }}-{\rm{i}}{M}$ only depends on the module of momentum, which contradicts the requirement that the ${d}_{3}({\bf{k}})$ in $g^{\prime} ={d}_{3}-{\rm{i}}{M}$ must be odd in momentum. Thus we conclude that the d₃ term cannot exist, and we have eventually $g^{\prime} =-{\rm{i}}{M}$. The oddness of d₁ and d₂ again yields the same constraint between the ${n}_{+}-{n}_{-}$ and the rotational eigenvalues as in equation (53), which can be satisfied by any m for fermions, and m = 3 for bosons.

Using the eigenstates in equation (40), with setting d₃ = 0 as discussed above, the ${{\mathbb{Z}}}_{2}$ topological invariant can be constructed from the Pfaffian of the m-matrix using the formalism in section 4.4. The TR operator in section 3.2 acting on the filled-band eigenstates yields

$$\begin{eqnarray}&&T| {u}_{1-}\rangle =\displaystyle \frac{1}{\sqrt{2}d}\left(\begin{array}{c}-{\rm{i}}{d}\\ 0\\ -M\\ -{{\rm{i}}{d}}_{1}+{d}_{2}\end{array}\right),\ \ \ T| {u}_{2-}\rangle =\displaystyle \frac{1}{\sqrt{2}d}\left(\begin{array}{c}0\\ {\rm{i}}{d}\\ {{\rm{i}}{d}}_{1}+{d}_{2}\\ M\end{array}\right).\end{eqnarray} \tag{ 119 }$$

The Pfaffian of the m-matrix is then

$$\begin{eqnarray}&&\mathrm{Pf}\left[m({\bf{k}})\right]=\langle {u}_{1-}| T| {u}_{2-}\rangle =-\langle {u}_{2-}| T| {u}_{1-}\rangle =\displaystyle \frac{M}{d}.\end{eqnarray} \tag{ 120 }$$

So we see that the sign of M determines the ${{\mathbb{Z}}}_{2}$ index in our continuous model. Obviously, because of the ${{\mathbb{Z}}}_{2}$ invariant, we always have ${\rm{\Delta }}{ \mathcal C }=1$. The rest of the correlation function, critical exponents, and universality class follows that of section 4.4 for 2D class AII.

Following [8], we model the 1D class BDI systems, which have T² = 1 and C² = 1, by $T={\sigma }_{z}K$ and $C={\sigma }_{x}K$. The TR symmetry requires

$$\begin{eqnarray}&&{d}_{1}(k)=-{d}_{1}(-k),\ \ \ {d}_{2}(k)={d}_{2}(-k),\ \ \ {d}_{3}(k)={d}_{3}(-k).\end{eqnarray} \tag{ 121 }$$

On the other hand, the PH symmetry requires

$$\begin{eqnarray}&&{d}_{1}(k)=-{d}_{1}(-k),\ \ \ {d}_{2}(k)=-{d}_{2}(-k),\ \ \ {d}_{3}(k)={d}_{3}(-k).\end{eqnarray} \tag{ 122 }$$

Combining the TR and PH symmetry, we see that the d₂ term cannot exist, and we may use d₃ = M as the mass term. The oddness of d₁ implies that we can parameterize the Hamiltonian by

$$\begin{eqnarray}&&H(k)={k}^{n}{\sigma }_{x}+M{\sigma }_{z},\ \ \ {E}_{\pm }(k)=\pm \sqrt{{k}^{2n}+{M}^{2}},\end{eqnarray} \tag{ 123 }$$

with an odd order of band crossing $n\in 2{\mathbb{Z}}+1$.

Denoting the valance band eigenstate of the Hamiltonian in equation (123) by $| {u}_{-}\rangle$, the corresponding Berry connection actually vanishes $A(k)=\langle {u}_{-}| {\rm{i}}{\partial }_{k}| {u}_{-}\rangle =0$, but the eigenstates do not return to themselves after a cycle. If we, however, choose a basis in which the Hamiltonian is off-diagonal, then the topological invariant can be computed solely from the integral of the Berry connection. For this purpose, we rotate the Hamiltonian by

$$\begin{eqnarray}\begin{array}{rcl} & & R={{\rm{e}}}^{{\rm{i}}{\sigma }_{x}\pi /4}=\displaystyle \frac{1}{\sqrt{2}}\left(\begin{array}{cc}1 & {\rm{i}}\\ {\rm{i}} & 1\end{array}\right),\\ & & \tilde{H}(k)={RH}(k){R}^{-1}={d}_{1}{\sigma }_{x}+{d}_{3}{\sigma }_{y}={k}^{n}{\sigma }_{x}+M{\sigma }_{y}.\end{array}\end{eqnarray} \tag{ 124 }$$

That is, we rotate the basis such that the mass term resides in the σ_y component. The eigenstate in this basis $\tilde{H}(k)| {\tilde{u}}_{-}\rangle =E(k)| {\tilde{u}}_{-}\rangle$ is related to the original basis by a rotation $| {\tilde{u}}_{-}\rangle ={V}_{k}| {u}_{-}\rangle$, and gives a Berry connection

$$\begin{eqnarray}&&\tilde{A}(k)=\langle {\tilde{u}}_{-}| {\rm{i}}{\partial }_{k}| {\tilde{u}}_{-}\rangle =\langle {u}_{-}| {\rm{i}}({V}_{k}{\partial }_{k}{V}_{k})| {u}_{-}\rangle =\displaystyle \frac{{d}_{3}{\partial }_{k}{d}_{1}-{d}_{1}{\partial }_{k}{d}_{3}}{2{d}^{2}}.\end{eqnarray} \tag{ 125 }$$

Using the parameterization in equation (123), we obtain the Berry connection

$$\begin{eqnarray}&&\tilde{A}(k)=\displaystyle \frac{{{nMk}}^{n-1}}{2({k}^{2n}+{M}^{2})}.\end{eqnarray} \tag{ 126 }$$

The topological invariant is the Zak phase calculated from the integration of Berry connection, which in our continuous model reads

$$\begin{eqnarray}&&{ \mathcal C }={\int }_{-\infty }^{\infty }\displaystyle \frac{{\rm{d}}{k}}{\pi }\tilde{A}(k)=\displaystyle \frac{1}{2}\mathrm{Sgn}(M),\end{eqnarray} \tag{ 127 }$$

normalized such that jump of topological invariant across the critical point M_c = 0 is ${\rm{\Delta }}{ \mathcal C }={ \mathcal C }(M\gt 0)-{ \mathcal C }(M\lt 0)=1$, consistent with the result of the usual linear Dirac model. Following the recipe in section 2.1, we consider our model as an low-energy effective theory of a lattice model, and introduce the Wannier states $| R\rangle$ constructed from $| {\tilde{u}}_{-}\rangle$. The topological invariant constructed from the Wannier state reads, according to the theory of charge polarization [13, 14],

$$\begin{eqnarray}&&\displaystyle \frac{{ \mathcal C }}{2}=\langle 0| \hat{r}| 0\rangle .\end{eqnarray} \tag{ 128 }$$

The correlation function is given by the Fourier transform of the Berry connection

$$\begin{eqnarray}&&{\tilde{F}}_{1{\rm{D}}}(R)={\int }_{{\rm{BZ}}}\displaystyle \frac{{\rm{d}}{k}}{2\pi }{{\rm{e}}}^{{\rm{i}}{k}{R}}\tilde{A}(k)=\langle 0| \hat{r}| R\rangle ,\end{eqnarray} \tag{ 129 }$$

which measures the overlap of the Wannier state centering at the origin and that centering at R, sandwiched by a position operator $\hat{r}$. These Wannier state representations are shown diagrammatically in figure 5.

Zoom In Zoom Out Reset image size

For n = 1, the Berry connection $\tilde{A}(k)$ displays an extremum ${\tilde{A}}_{\max }$ at the HSP k = 0, whereas for any n > 1, the extremum is located at a finite momentum k_max, i.e. showing a double peak structure in momentum space. It is straight forward to obtain

$$\begin{eqnarray}&&{k}_{\max }={\left(\displaystyle \frac{n-1}{n+1}{M}^{2}\right)}^{1/2n},\ \ \ {\tilde{A}}_{\max }=\displaystyle \frac{n}{2}{\left(\displaystyle \frac{n-1}{n+1}\right)}^{\tfrac{n-1}{2n}}{M}^{-1/n}\end{eqnarray} \tag{ 130 }$$

indicating a critical exponent α = 1/n. Moreover, the half-width-at-half-maximum scales as ${k}_{{\rm{HWHM}}}\propto | M{| }^{1/n}$, implying the correlation length $\xi \propto | M{| }^{-1/n}$ and a critical exponent ν = 1/n.

For 1D class AIII that satisfies the chiral symmetry S² = 1, we follow [8]. The chiral operator is S = σ_y, and hence the chiral symmetry demands that the d₂ term must vanish. We choose d₃ to be the mass term and write

$$\begin{eqnarray}&&H(k)={k}^{n}{\sigma }_{x}+M{\sigma }_{z}.\end{eqnarray} \tag{ 131 }$$

The order of band crossing n can be either even or odd without violating the symmetry. Since the Hamiltonian takes exactly the same form as that for 1D class BDI discussed in section 5.1, the formalism therein also applies. We make the same rotation of basis by equation (124), and introduce the Berry connection and Zak phase. The resulting Wannier state correlation function takes the same form as in equation (129). The critical exponent of the correlation length is again $\nu =1/n$. However, it is obvious from equation (126) that even values of n will give a Berry connection that is odd in k, and thereby integrates to zero. In other words, only odd values of n can be topologically nontrivial.

We close this section by revisiting an assumption made in section 2.1, namely throughout the present work, we only consider systems in which the gap-closing takes place at an HSP. In appendix, we use a 1D class AIII model to demonstrate that, even for systems that have gap-closing at a non-HSP, the critical behavior may still be entirely determined by the low-energy Dirac model of the system. In other words, the critical exponents and the jump of topological invariant may still be determined by the order of band crossing and the symmetry class, and possibly other symmetries compatible with the non-HSP. The clarification of this conjecture would require a systematic study of this type of systems, which is left for future investigations.

The 1D class DIII formalism can be obtained from 3D class DIII in section 3.2 by a dimensional reduction k_y = k_z = 0. Since rotational symmetry is absent in 1D, the issue now is reduced to how to parameterize $f={d}_{1}-{{\rm{i}}{d}}_{2}$ and $g^{\prime} ={d}_{3}-{\rm{i}}{M}$ in equation (27) using 1D momentum k_x = k, provided ${d}_{i}(k)=-{d}_{i}(-k)$ for i = 1, 2, 3. Given the oddness of d_i, we should parameterize

$$\begin{eqnarray}&&f={a}_{1}{k}^{n}-{{\rm{i}}{a}}_{2}{k}^{n},\ \ \ g^{\prime} ={a}_{3}{k}^{n}-{\rm{i}}{M},\ \ \ n\in 2{\mathbb{Z}}+1,\end{eqnarray} \tag{ 132 }$$

such that the dispersion is

$$\begin{eqnarray}&&{E}_{\pm }(k)=\pm \sqrt{| f{| }^{2}+| g^{\prime} {| }^{2}}=\pm {\left[({a}_{1}^{2}+{a}_{2}^{2}+{a}_{3}^{2}){k}^{2n}+{M}^{2}\right]}^{1/2}.\end{eqnarray} \tag{ 133 }$$

Therefore, the most generalized form is that all (d₁, d₂, d₃) are nonzero, and we denote $d=\sqrt{{d}_{1}^{2}+{d}_{2}^{2}+{d}_{3}^{2}+{M}^{2}}$. The filled band eigenstates are that in equation (40), since we use the same representation of the Γ-matrices. Applying the TR operator in section 3.2 to the filled band eigenstates gives the same result as equation (119), and the Pfaffian is that in equation (120). So once again the sign of the mass term M determines the ${{\mathbb{Z}}}_{2}$ topological invariant in this continuous model.

The Wannier state correlation function in this class is the 1D version of that in section 4.4. The correlation function is again calculated from the Fourier transform of the second derivative of the Pfaffian on the m-matrix

$$\begin{eqnarray}&&{\partial }_{k}^{2}\mathrm{Pf}\left[m\right]=\displaystyle \frac{(n+{n}^{2}){{Mk}}^{4n-2}+(n-2{n}^{2}){M}^{3}{k}^{2n-2}}{{\left({k}^{2n}+{M}^{2}\right)}^{5/2}}.\end{eqnarray} \tag{ 134 }$$

The Fourier transform of which gives

$$\begin{eqnarray}&&{\tilde{F}}_{{\rm{TR}}}(R)={\int }_{{\rm{BZ}}}\displaystyle \frac{{\rm{d}}{k}}{2\pi }{{\rm{e}}}^{{\rm{i}}{k}{R}}{\partial }_{k}^{2}\mathrm{Pf}[m(k)]=-\langle R1| {R}^{2}T| 02\rangle \approx \displaystyle \frac{1}{\pi }{\int }_{0}^{\infty }{\rm{d}}{k}\,{{\rm{e}}}^{{\rm{i}}{k}{R}}{\partial }_{k}^{2}\mathrm{Pf}[m(k)],\end{eqnarray} \tag{ 135 }$$

where $| R1\rangle$ denotes the Wannier state of band 1 centered at R, and $| 02\rangle$ denotes the Wannier state of band 2 centered at the origin.

For 1D class D with ${C}^{2}=1$, we follow [10] to use the PH operator $C={\sigma }_{x}K$, which requires that

$$\begin{eqnarray}&&{d}_{i}(k)=-{d}_{i}(-k)\ \ \ \mathrm{for}\ i=1,2,\ \ \ {d}_{3}(k)={d}_{3}(-k),\end{eqnarray} \tag{ 136 }$$

and hence only d₃ = M can be the mass term, and $\left\{{d}_{1},{d}_{2}\right\}$ must be of odd powers of k. In this class, the basis is the Bogoliubov quasiparticles

$$\begin{eqnarray}{ \mathcal H }=\displaystyle \sum _{k}\left({c}_{k}^{\dagger }\ {c}_{-k}\right)\left(\begin{array}{cc}M & {d}_{1}-{{\rm{i}}{d}}_{2}\\ {d}_{1}+{{\rm{i}}{d}}_{2} & -M\end{array}\right)\left(\begin{array}{c}{c}_{k}\\ {c}_{-k}^{\dagger }\end{array}\right).\end{eqnarray} \tag{ 137 }$$

To have a dispersion of the form ${E}_{\pm }(k)=\pm \sqrt{{k}^{2n}+{M}^{2}}$, we parameterize d_i by

$$\begin{eqnarray}&&{d}_{1}={a}_{1}{k}^{n},\ \ \ {d}_{2}={a}_{2}{k}^{n},\end{eqnarray} \tag{ 138 }$$

where $| a| =\sqrt{{a}_{1}^{2}+{a}_{2}^{2}}=1$ is rescaled to unity. We may further eliminate the d₁ component by a gauge transformation such that the Hamiltonian takes the form of the well-studied 1D Majorana chain [47], and the correlation function can be introduced [48]. We rewrite the Hamiltonian by

$$\begin{eqnarray}&&\begin{array}{l}({d}_{1}-{{\rm{i}}{d}}_{2}){c}_{k}^{\dagger }{c}_{-k}^{\dagger }=({a}_{1}-{{\rm{i}}{a}}_{2}){k}^{n}{c}_{k}^{\dagger }{c}_{-k}^{\dagger }=| a| {{\rm{e}}}^{{\rm{i}}\varphi }{k}^{n}{c}_{k}^{\dagger }{c}_{-k}^{\dagger }\\ \quad =\,-{{\rm{i}}{k}}^{n}{{\rm{e}}}^{{\rm{i}}(\varphi +\pi )}{c}_{k}^{\dagger }{c}_{-k}^{\dagger }=-{\rm{i}}{\tilde{d}}_{2}{\tilde{c}}_{k}^{\dagger }{\tilde{c}}_{-k}^{\dagger }\to -{{\rm{i}}{d}}_{2}{c}_{k}^{\dagger }{c}_{-k}^{\dagger },\end{array}\end{eqnarray} \tag{ 139 }$$

such that the Hamiltonian now reads

$$\begin{eqnarray}{ \mathcal H }=\displaystyle \sum _{k}\left({c}_{k}^{\dagger }\ {c}_{-k}\right)\left(\begin{array}{cc}M & -{{\rm{i}}{d}}_{2}\\ {{\rm{i}}{d}}_{2} & -M\end{array}\right)\left(\begin{array}{c}{c}_{k}\\ {c}_{-k}^{\dagger }\end{array}\right),\end{eqnarray} \tag{ 140 }$$

with

$$\begin{eqnarray}&&{d}_{2}={k}^{n},\ \ \ n\in 2{\mathbb{Z}}+1,\end{eqnarray} \tag{ 141 }$$

which is now in a more convenient form.

For the case that the basis is fermionic, we introduce the Majorana basis by assuming that our continuous Hamiltonian is the low-energy effective theory of a lattice model, and introduce the Majorana fermions from the real space lattice electron creation and annihilation operators at site i [10]

$$\begin{eqnarray}&&{\lambda }_{i}={c}_{i}+{c}_{i}^{\dagger },\ \ \ {\lambda }_{i}^{{\prime} }=\displaystyle \frac{1}{i}({c}_{i}-{c}_{i}^{\dagger }),\ \ \ {{\rm{\Lambda }}}_{i}=\left(\begin{array}{c}{\lambda }_{i}\\ {\lambda }_{i}^{{\prime} }\end{array}\right).\end{eqnarray} \tag{ 142 }$$

Then the Hamiltonian in equation (140) can be written as

$$\begin{eqnarray}\begin{array}{rcl}{ \mathcal H } & = & \displaystyle \frac{{\rm{i}}}{4}\displaystyle \sum _{k}{{\rm{\Lambda }}}^{T}(k)X(k){\rm{\Lambda }}(-k)\\ & = & \displaystyle \sum _{k}{{\rm{\Lambda }}}^{T}(k)\left(\begin{array}{cc}0 & {d}_{2/4}-{\rm{i}}{M}/4\\ {d}_{2}/4+{\rm{i}}{M}/4 & 0\end{array}\right){\rm{\Lambda }}(-k),\\ X(k) & = & -{\rm{i}}{M}{\sigma }_{y}-{{\rm{i}}{d}}_{2}{\sigma }_{x}.\end{array}\end{eqnarray} \tag{ 143 }$$

That is, d₂/4 is effectively the σ_x component and M/4 the σ_y component in this Majorana basis.

At HSPs ${k}_{0}=-{k}_{0}$, the PH symmetry ${CH}({k}_{0}){C}^{-1}=-H(-{k}_{0})=-H({k}_{0})$ implies that d₁ and d₂ in equation (137) must vanish, as it does for our continuous model that sets k₀ = 0 at the origin. The X(k₀) in equation (143) is

$$\begin{eqnarray}X({k}_{0}=0)=-{\rm{i}}{M}{\sigma }_{y}=\left(\begin{array}{cc} & -M\\ M & \end{array}\right),\end{eqnarray} \tag{ 144 }$$

which is antisymmetric. The ${{\mathbb{Z}}}_{2}$ topological invariant is calculated from the sign of the Pfaffian of X(k₀) at the two HSPs k₀ = 0 and k₀ = π in the original lattice model. In our continuous model

$$\begin{eqnarray}&&\mathrm{Sgn}\left\{\mathrm{Pf}\left[X({k}_{0}=0)\right]\right\}=-\mathrm{Sgn}(M),\end{eqnarray} \tag{ 145 }$$

so there is a topological phase transition at M_c = 0.

Alternatively, the topological invariant can be calculated from the Berry connection in the Majorana basis as in section 5.1. The eigenstates in the Majorana basis calculated from equation (143) is

$$\begin{eqnarray}&&| {u}_{\pm }\rangle =\displaystyle \frac{1}{\sqrt{2}d}\left(\begin{array}{c}\pm d\\ {d}_{2/4}+{\rm{i}}{M}/4\end{array}\right),\ \ \ d=\displaystyle \frac{1}{4}\sqrt{{d}_{2}^{2}+{M}^{2}}.\end{eqnarray} \tag{ 146 }$$

The Berry connection of the filled band in this gauge is

$$\begin{eqnarray}&&A(k)=\langle {u}_{-}| {\rm{i}}{\partial }_{k}| {u}_{-}\rangle =\displaystyle \frac{M{\partial }_{k}{d}_{2}}{2\left({d}_{2}^{2}+{M}^{2}\right)}.\end{eqnarray} \tag{ 147 }$$

After using d₂ = kⁿ, we obtain the same form of A(k) as in equation (126), and thus the same Zak phase ${ \mathcal C }$ as in equation (127). We see that ${\rm{\Delta }}{ \mathcal C }=1$ regardless of the order of band crossing n. The correlation function follows from equation (129) which measures the overlap of Wannier function sandwiched by a position operator. Because the basis of this Wannier state is the Majorana fermions, the correlation function has been referred to as the Majorana-Wannier state correlation function, as previously calculated for the 1D Kitaev chain which is an explicit example for n = 1 in this class [48]. It should be reminded that although we formulate this section in the fermionic language, our formalism also applies to bosonic models.

For 1D class CII with T² = −1 and C² = 1, we use the Γ-matices [49]

$$\begin{eqnarray}&&{{\rm{\Gamma }}}^{a}=\left\{{\sigma }_{x}\otimes {\tau }_{z},{\sigma }_{y}\otimes {\tau }_{z},I\otimes {\sigma }_{x},I\otimes {\tau }_{y},{\sigma }_{z}\otimes {\tau }_{z}\right\}.\end{eqnarray} \tag{ 148 }$$

The TR and PH operators are $T={\sigma }_{y}\otimes {IK}$ and $C=I\otimes {\tau }_{y}K$. The TR and PH symmetries require the mass term to be d₃ = M, and the only allowed kinetic term is ${d}_{2}({\bf{k}})=-{d}_{2}(-{\bf{k}})$. Thus a general higher-order Dirac Hamiltonian in this class takes the form

$$\begin{eqnarray}&&H={d}_{2}{{\rm{\Gamma }}}^{2}+{d}_{3}{{\rm{\Gamma }}}^{3}={k}^{n}{{\rm{\Gamma }}}^{2}+M{{\rm{\Gamma }}}^{3},\end{eqnarray} \tag{ 149 }$$

with $n\in 2{\mathbb{Z}}+1$. The topological invariant is given by the winding number [49]

$$\begin{eqnarray}&&{ \mathcal C }=\displaystyle \frac{1}{2{\rm{i}}}{\int }_{0}^{2\pi }\displaystyle \frac{{\rm{d}}{k}}{2\pi }\mathrm{Tr}\left[{\sigma }_{y}\otimes {\tau }_{y}{H}^{-1}(k){\partial }_{k}H(k)\right],\end{eqnarray} \tag{ 150 }$$

which in our continuous model of equation (149), with the replacement in the integral ${\int }_{0}^{2\pi }\to {\int }_{-\infty }^{\infty }$, yields ${ \mathcal C }=-\mathrm{Sgn}(M)$, and hence the jump of topological invariant across the transition is ${\rm{\Delta }}{ \mathcal C }=2$ for any order of band crossing n.

The filled-band eigenstates read

$$\begin{eqnarray}&&| {u}_{1-}\rangle =\displaystyle \frac{1}{\sqrt{2}d}\left(\begin{array}{c}-{{\rm{i}}{d}}_{2}\\ -d\\ 0\\ {d}_{3}\end{array}\right),\ \ \ | {u}_{2-}\rangle =\displaystyle \frac{1}{\sqrt{2}d}\left(\begin{array}{c}-d\\ {{\rm{i}}{d}}_{2}\\ {d}_{3}\\ 0\end{array}\right).\end{eqnarray} \tag{ 151 }$$

To construct the correlation function, we consider a modified Berry connection

$$\begin{eqnarray}&&{\overline{A}}^{{ab}}=\langle {u}_{a-}| {\sigma }_{y}\otimes {\tau }_{y}{\rm{i}}{\partial }_{k}| {u}_{b-}\rangle ,\end{eqnarray} \tag{ 152 }$$

such that the integrand in equation (150) can be written as

$$\begin{eqnarray}&&\mathrm{Tr}\left[{\sigma }_{y}\otimes {\tau }_{y}{H}^{-1}(k){\partial }_{k}H(k)\right]=4{\rm{i}}\mathrm{Tr}\left[{\overline{A}}^{{ab}}\right],\end{eqnarray} \tag{ 153 }$$

meaning that the topological invariant can as well be expressed as an integration of $\mathrm{Tr}\left[{\overline{A}}^{{ab}}\right]$. The Fourier transform of $\mathrm{Tr}\left[{\overline{A}}^{{ab}}\right]$ yields a correlation function

$$\begin{eqnarray}&&{\tilde{F}}_{1{\rm{D}}}^{{\prime} }(R)={\int }_{0}^{2\pi }\displaystyle \frac{{\rm{d}}{k}}{2\pi }{{\rm{e}}}^{{\rm{i}}{k}{R}}{\rm{Tr}}\left[{\overline{A}}^{{ab}}\right]={\rm{Tr}}\left[\langle 0| {\sigma }_{y}\otimes {\tau }_{y}\hat{r}| R\rangle \right],\end{eqnarray} \tag{ 154 }$$

that measures the overlap of Wannier functions weighted by ${\sigma }_{y}\otimes {\tau }_{y}\hat{r}$.