2.1. Topological hierarchy superconductor from s-wave pairing SC with Rashba SO coupling and flux-superlattice

The first example of topological hierarchy superconductor is an s-wave pairing SC with Rashba SO coupling and flux-superlattice.

The s-wave pairing SC with Rashba SO coupling on the square lattice is described by^[18]

where the kinetic term , the Rashba SO coupling term , and the superconducting pairing term are given as

Here, annihilates (creates) a fermion at site j = (j_x,j_y) with spin σ = (↑,↓), μ = x or y which is a basic vector for the square lattice, t_s is the nearest neighbor (NN) hopping parameter, λ serves as the SO coupling constant, Δ serves as the s-wave SC pairing order parameter, u is the chemical potential, and h is the strength of the Zeeman field, which is illustrated in Fig. 1(a).

Firstly, we calculate the Majorana zero mode around a single vortex by the exact diagonalization numerical approach. According to the numerical results, there exists a localized state with exact zero energy around the vortex. The zero mode is localized around the vortex core within a length-scale ∼ m⁻¹ (m is the mass gap). For this zero mode, the fermionic-quasiparticle operator satisfies , which indicates that the zero mode is indeed a Majorana fermion. When there are two vortices nearby, as shown in Fig. 1(b), the intervortex quantum tunneling effect occurs and the Majorana modes on the two vortices couple. As shown in Fig. 1(c), the energy splitting δE as a function of the distance between the two vortices D oscillates and decreases exponentially. When the two vortices are well separated, the quantum tunneling effect can be ignored and we have two quantum states with exact zero energy.^[19–21]

Next, we numerically study the Majorana zero modes of a square flux-superlattice. The results are shown in Fig. 2(a). Because each vortex traps a Majorana fermion, we obtain a multi-Majorana fermion system for the system with the square flux-superlattice. The Majorana fermions on the square flux-superlattice will couple with each other due to the intervortex quantum tunneling effect.

The density of state (DOS) of the TSC with a square flux-superlattice is shown in Fig. 2(b). From Fig. 2(b), one can see that except for the energy bands of the paired electrons, there exists a mid-gap energy band in the parent TSC. In particular, the mid-gap energy band has a finite energy gap. This means that this mid-gap system shown in Fig. 2(b) may be a topological state with a non-trivial topological number.

To check the topological properties of the mid-gap system, we calculate its edge states. We consider a system on a cylinder with 12 super-unitcells along x-direction while periodic boundary along y-direction. The results are shown in Figs. 2(c) and 2(d). Figure 2(c) shows that two gapless chiral edge states are localized at the boundaries. While when the parent SC is a non-topological SC, the mid-gap energy band together with the edge states disappears (see Fig. 2(d)). So the mid-gap system induced by the flux-superlattice shown in Fig. 2(a) is really a topological state.

Then, we propose an effective tight-binding lattice model of the s-wave TSC with a square flux-superlattice, in which each vortex traps a Majorana zero mode and two Majorana zero modes couple with each other by a short range interaction. We call this effective tight-binding lattice model as the Majorana lattice model, of which the effective Hamiltonian can be written as

where t_jk is the Majorana fermion’s hopping amplitude from j to k and satisfies , and s_ij = −s_ji is a gauge factor. The pair (j,k) denotes the summation that runs over all the nearest neighbor pairs (with hopping amplitude t) and all the next-nearest neighbor (NNN) pairs (with hopping amplitude t′). From the polygon rule proposed in Ref. [14], each triangular plaquette possesses π/2 quantum flux effectively. This Hamiltonian allows Z2 gauge choice s_jk = ±1. From this result, we can derive that there always exists an energy gap in the Majorana lattice model as long as t′ ≠ 0.

To characterize the topological properties of the Majorana lattice model, we introduce the Chern number ,

where d(k) = (d_x (k),d_y (k),d_z (k)), d_x (k) = −2t sin k_x, d_y (k) = −sin2k_y (t + 2t′ cos k_x), d_z (k) = −2t sin²k_y + 4t′ cos k_x cos²k_y. In the presence of NNN hopping t′, we have . Furthermore, we calculate the edge states of the Majorana lattice model and find that a gapless mode indeed exists. Now we can conclude that the (topological) Majorana lattice model captures the key low energy physics of the parent s-wave SC with a square flux-superlattice.

2.2. Topological hierarchy superconductor from topological p_x + ip_y SC with flux-superlattice

The second example of topological hierarchy superconductor is a topological p_x + ip_y SC with flux-superlattice.

Topological p_x + ip_y superconductivity is encoded in the following lattice Hamiltonian H (Fig. 3(a)):

where u is the chemical potential, is the electron pairing function, and t is the hopping strength. In the following parts, we set t to be the energy unit. The operator creates/destroys an electron on lattice site r and satisfies the anti-commutation statistics .

The topological properties of the topological p_x + ip_y SC and those of the topological s-wave pairing SC with Rashba SO coupling are similar. If we do the unitary transformation

on the BdG Hamiltonian of the s-wave pairing SC, it will map into a p_x + ip_y SC.^[22] The Rashba SO coupling corresponds to the pairing terms for the p_x + ip_y SC through the unitary transformation. Thus, the quantized vortex in p_x + ip_y SC also traps a Majorana zero mode. The energy splitting between two energy levels versus the flux-distance d oscillates and decreases exponentially, which is shown in Fig. 3(b).

We then study the p_x + ip_y topological SC with a flux-superlattice. We choose 8a × 4a sites to be a unit cell (8a along x-direction and 4a along y-direction), which indicates that the flux-superlattice is anisotropy. To show the topological properties of the mid-gap spectra, we put the system on a cylinder (open boundary condition along x-direction and periodic boundary condition along y-direction). The result of this case is plotted in Fig. 3(c), which is similar to Fig. 2(c). The energy bands from the bulk consist of two kinds of fermion spectra: the Bogoliubov quasi-particles spectrum and the vortex core fermion spectrum induced by the flux-superlattice.^[23] Both energy bands have energy gaps. Except for the energy bands from the bulk, there exist gapless edge states. There are two types of edge spectra: one comes from the parent topological superconductor due to its nontrivial topological properties, and the other comes from the flux-superlattice. This means that the mid-gap spectra induced by the flux-superlattice indeed have nontrivial topological properties.

Due to the unit cell of the vortex-superlattice is a rectangle, the effective tight-binding Majorana lattice model is given by

where γ_i is the operator of Majorana fermion in the i-th vortex core obeying the self-conjugate condition and the canonical commutate relation {γ_i,γ_j} = 2δ_ij. The indices l = 1,2 denote the nearest neighbor couplings and the next-next nearest neighbor couplings, respectively, t_lx and t_ly are the corresponding coupling strengths along x and y directions, and is the next nearest neighbor coupling strength. Each s_ij connecting bond ⟨ij⟩ has Z2 gauge degree of freedom which does not affect any physical conclusions. In Fig. 3(d), we give the illustration of the gauge choice of the square Majorana lattice. Motivated by the decompose rules that two Majorana fermions fuse into a complex fermion, and no net flux through the primitive cell, we adopt the condition that there is a π/2 flux in each triangular plaquette, and the detail argument about the flux is shown in Ref. [14].

We obtain the tunneling parameters by fitting the energy dispersion of the mid-gap spectra induced by the flux-superlattice. As shown in the top panel of Fig. 4, the mid-gap spectra are obtained by a numerical approach.^[23,24] The fitting dispersion of the effective tight-binding Majorana lattice model is shown in the bottom panel. We calculate the topological invariant of the effective tight-binding Majorana lattice model which gives the Chern number . So the effective tight-binding Majorana lattice model is a “topological SC state”. This hierarchical effect has been employed to understand the topological quantum transition induced by vortex excitation in the context of Kitaev’s honeycomb model supporting non-Abelian anyon excitation.^{[15,16,25–27]} In Refs. [28] and [29], the flux-superlattice in the p_x + ip_y topological SC has been studied, in which the effective tight-binding Majorana lattice model has no nontrivial topological properties.

In the above section, we illustrate the s-wave pairing SC with Rashba SO coupling and the topological p_x + ip_y SC with a square flux-superlattice. Due to the inter-vortex tunneling, there exist mid-gap energy bands induced by the flux-superlattice. We find that such mid-gap energy bands always have nontrivial topological properties including the gapless edge states and non-zero winding number. We can write down an effective tight-binding Majorana lattice model to characterize the mid-gap states and obtain the tunneling parameters by fitting the low energy spectrum with numerical calculations. Therefore, we call this type of topological superconductors with topological MG states topological hierarchy superconductors (THSC). In summary, a THSC is given by the following equation:

3.1. Topological hierarchy insulators induced by flux-superlattice

Based on the Haldane model on the honeycomb lattice, we review the topological hierarchy insulator from a topological insulator with a flux-superlattice.

The Hamiltonian of the Haldane model is given by^[3]

Here, the fermionic operator ĉ_i annihilates a fermion on lattice site i. T is the NN hopping amplitude and T′ is the NNN hopping amplitude, which are illustrated in Fig. 5(a). ⟨i,j⟩ and ⟪i,j⟫ denote the NN and the NNN links, respectively. e^iϕ_ij is a complex phase along the NNN link, and we set the direction of the positive phase |ϕ_ij| = π/2 clockwise. μ is the chemical potential, which is set to be zero in this paper. In the following, we take the lattice constant a ≡ 1. The Haldane model is an integer quantum Hall insulator without Landau levels. It breaks the time reversal symmetry without any net magnetic flux through the unit cell of a periodic two-dimensional honeycomb lattice. There exists a topological invariant for the Haldane model, the Thouless–Kohmoto–Nightingale–den Nijs (TKNN) number (or the Chern number).^[31] Thus, the Haldane model in Eq. (7) is a typical topological band insulator for the case of T′ ≠ 0 at half-filling.

By using the exact diagonalization numerical approach on a 72 × 72 lattice, we find that there exists a fermionic zero mode around each π-flux. The particle density is mainly localized around the π-flux. The length-scale of the wave-function of the zero mode is ξ ∼ v_F/Δf, where v_F = 3aT/2 is the Fermi velocity. When there are two fluxes nearby, the inter-flux quantum tunneling effect occurs and the two zero modes couple. The energy splitting δE between two energy levels versus the flux-distance L oscillates and decreases exponentially, which is shown in Fig. 5(b).

Next we superpose a triangular flux-superlattice on the TI. By using the numerical calculations, we obtain the DOS for the Haldane model with a triangular flux-superlattice. See numerical results in Fig. 6(a) with flux-superlattice constant for the case of T = 1.0, T′ = 0.1. From Fig. 6(a), one may see that the mid-gap bands appear. In Fig. 6(b), the mid-gap bands in Fig. 6(a) are zoomed in. Now, we find that there exist four points with van Hove singularity, and the mid-gap bands also have an energy gap.

To illustrate the topological properties of the mid-gap states induced by the flux-superlattice of the parent TI, we study the edge states of the TI with flux-superlattice for the case of T = 1.0, T′ = 0.1. We put the Haldane model with a finite flux-superlattice along x-direction on a torus. The parent topological insulator has the periodic boundary condition along both x-direction and y-direction. While the flux-superlattice has the periodic boundary condition along y-direction and the open boundary condition along x-direction. See the numerical results in Fig. 6(c). There exist edge states along the boundaries of the flux-superlattice. The existence of the gapless zero modes on the boundary of the flux-superlattice indicates that the mid-gap system is indeed an induced “topological insulator” on the parent TI.

Finally, we also write down an effective tight-binding model to describe the mid-gap states. The quantum states of the fermionic zero mode around a π-flux can be formally described in terms of the fermion Fock states {|0⟩, |1⟩}. Here, |0⟩ and |1⟩ denote the unoccupied state and the occupied state, respectively. These quantum states are localized around the flux within a length-scale ξ ∼ v_F/Δ_f. For the case of ξ < L, we can consider each flux as an isolated “atom” with localized electronic states and use the effective tight-binding model to describe these quantum states on the fluxes.

Now, we superpose the localized states to obtain the sets of Wannier wave functions w(R) with R denoting the position of the flux. Due to the formation of the triangular flux-superlattice, both the quantum-tunneling-strength δE between NN fluxes and the quantum-tunneling-strength δE′ between NNN fluxes exist. Owing to the inter-flux quantum tunneling effect, we obtain energy splitting |δE_RR′|, which is just the particle’s hopping amplitude |t_RR′| between two localized states on two π-fluxes at R and R′. The effective tight-binding model of the two fermionic zero modes then takes the form of with |t_RR′| = |δ E_RR′|. With considering the NN and NNN hoppings, the effective tight-binding model of the localized states around the triangular flux-superlattice is given by

where is the fermionic annihilation operator of a localized state on a flux R, and t_RR′ is the hopping parameter between NN (NNN) sites R and R′. The NN (NNN) hopping term is denoted as . In particular, owing to the polygon rule (each fermion gains an accumulated phase shift |ϕ| = (n − 2)π/2 encircling around a smallest n-polygon^[14]), the total phase around each plaquette of the triangular flux-lattice is ±π/2 for fermions. The ratio between the NN hopping strength t and the NNN hopping strength t′ is indeed the same as the ratio between |δE| and |δE′|. Because people can tune |δE′|/|δE| (t′/t) by changing the flux-superlattice constant L, we may regard the flux-superlattice model of the mid-gap states as an emergent controllable topological system.

Now, we study the topological properties of the effective flux-superlattice model. According to the gauge shown in Fig. 7, we choose if fermions hop along the directions of the arrows. Now we label the fermionic annihilation operators of the localized states on the two sub-flux-superlattices by Â_R, B̂_R. By the Fourier transformation, we obtain the effective Hamiltonian of the tight-binding model of the TI with triangular flux-superlattice in momentum space as

where and with h = (h_x, h_y, h_z) and τ being the Pauli matrix. The three components of h are

where η₁ = L(1, 0), , . In the following, for simplicity, we set the flux-superlattice constant to be L ≡ 1. The energy spectrums of the effective tight-binding model of the TI with triangular flux-superlattice are then obtained as

The topological quantum phase transition occurs when the band gap closes E_VL(k) = 0. See Fig. 8(a). We find that there exists a quantum critical point at t′/t = 1 that separates two quantum phases, 0 < t′/t < 1, t′/t > 1. To characterize the two quantum phases, we introduce the Chern number^[31]

with n = h/|h|. According to the Chern number, we find that in the region of 0 < t′/t < 1, the Chern number is 1, and in the region of t′/t > 1, the Chern number is −3. The tight-binding model of the TI with triangular flux-superlattice always has nontrivial topological properties. As a result, we call it topological flux-superlattice model that can describe the mid-gap states of the Haldane model with triangular flux-superlattice.

In Fig. 8(b), we show the DOS of the flux-superlattice model for the case of t/T = 0.04, t′/T = 0.014. From the DOS, we can see that the effective flux-superlattice model has a finite energy gap and there exist four points with van Hove singularity. Figures 8(c) and 8(d) show the edge states for the cases of t′ = 0.1t (the Chern number is C = 1) and t′ = 3.1t (the Chern number is C = −3), respectively.

Hence, the low energy physics of the Haldane model with flux-superlattice can be described by an effective flux-superlattice model. The effective flux-superlattice model shows nontrivial topological properties, including a nontrivial topological invariant and gapless edge states. In this sense, the effective flux-superlattice model is really an emergent “topological insulator” on the parent TIs. Our results are more exotic than those in Ref. [32], in which a (non-topological) free-fermion model on a honeycomb lattice with flux-superlattice is studied.

We have studied the Haldane model with triangular flux-superlattice. Using a similar approach, we studied the Haldane model with other types of flux-superlattices such as the square flux-superlattice and honeycomb flux-superlattice, and found similar topological properties. In this sense, the topological mid-gap states always exist in a TI with flux-superlattice. In addition, we need to point out that the Kane–Mele model^[5] or the spinful Haldane model^[33] with flux-superlattice exhibit similar topological features. Due to the spin degree of freedom, a π-flux on these models traps two zero modes. The overlap of the zero modes also gives rise to a topological mid-gap system inside the band gap of the parent TIs. When one considers the interaction between the two-component fermions, the topological mid-gap system may lead to quite different physics consequences.

As a result, we call this type of topological insulators with topological MG states topological hierarchy insulators (THIs). In summary, THI from flux-superlattice is given by the following equation:

3.2. Topological hierarchy insulators induced by vacancy-superlattice

Based on the Haldane model on the square lattice, we review the topological hierarchy insulator from a topological insulator with a vacancy-superlattice.

The Hamiltonian of the (spinless) Haldane model on the square lattice is given by^[3,34]

where ĉ_i is the annihilation operator of the fermions at site i. A and B label the sublattices. t and t′ are the nearest-neighbor and the next-nearest-neighbor hopping parameters, respectively. For the Hamiltonian in Eq. (14), there exists a π-flux in each square plaquette and a π/2-flux in each triangular lattice. The illustration of the Haldane model on the square lattice is shown in Fig. 9(a). The Haldane model on the square lattice is a TI with quantum anomalous Hall (QAH) effect. The Chern number C_parent of the (parent) Haldane model is 1 and the quantized Hall conductivity is e²/h.

Vacancy is a non-topological defect and cannot induce a fermion zero mode when we consider a vacancy in TI. However, we have pointed out that if the Hamiltonian has the particle–hole (PH) symmetry, then a vacancy can induce a zero mode. Under PH transformation, we have

where is the PH transformation operator.^[33,35] Here, 𝓡 is an operator that leads to and is the complex conjugate operator. As a result, each energy level with positive energy E is paired with an energy level with negative energy −E. If the Hamiltonian is on a bipartite lattice, the quantum levels of the system with a single vacancy become an odd number. As a result, there must exist an unpaired electronic state when we remove a lattice site to create a vacancy. Because of the PH symmetry, the corresponding unpaired electronic state must have exactly zero energy.^[36]

The Hamiltonian in Eq. (14) satisfies the particle–hole (PH) symmetry. So, each vacancy induces a zero mode. We have calculated the fermionic zero mode around a single vacancy by the exact diagonalization numerical approach on a 60 × 60 square lattice and obtain the particle density distribution around the vacancy in Fig. 9(b). The particle density is localized around the vacancy center within a length-scale ∼ (Δ_f)⁻¹, where Δ_f is the fermion energy gap. The localization effect is similar to that of a flux on TI/TSC. We denote the wave-function of the zero mode by ψ₀(r_i − R), where R is the position of the vacancy. The quantum states of the fermionic zero mode around a vacancy can be formally described in terms of the fermion Fock states {|0⟩,|1⟩}. Here, |0⟩ and |1⟩ denote the empty state and the occupied state, respectively.

We next study the Haldane model with vacancy-superlattice. When there are two vacancies nearby, the inter-vacancy quantum tunneling effect occurs and the fermionic zero modes on the two vacancies couple. The lattice constant of the vacancy-superlattice (VSL) is denoted by L, which is the distance between two nearby vacancies. When the vacancies are located on different sublattices, the energy splitting ΔE is finite (L/a is an odd number); when the vacancies are located on the same sublattice, the energy splitting vanishes, ΔE = 0 (L/a is an even number). The energy splitting ΔE between two vacancies as a function of L is shown in Fig. 9(c). When two vacancies are well separated (or L → ∞), the quantum tunneling effect can be ignored and we have two quantum states with exact zero energy. In contrast, for the smaller L (odd number times of a), the coupling between two zero modes becomes stronger and the energy splitting cannot be neglected.

For the case of ξ < L, we can consider each vacancy as an isolated “atom” with localized electronic states and use the effective tight-binding model to describe these quantum states induced by the VSL. We now construct the superpositions of the localized states to obtain the sets of Wannier wave functions ψ₀(r_i − R). In general, the effective tight-binding model of the localized states induced by the VSL becomes

where d̂_R is the fermionic annihilation operator of a localized state on vacancy R. is the hopping parameter between NN (NNN, next next nearest neighbor (NNNN)) sites R and R′ and is determined by the energy splitting ΔE, . In general, because , |t_RR′|, we only consider the NN and NNN hopping terms. Because the parent Hamiltonian in Eq. (14) has the PH symmetry, the effective tight-binding model of the localized states induced by the VSL Ĥ_1−VL in Eq. (16) also has the PH symmetry,

Because the energy splitting ΔE only determines the particle’s hopping amplitude |t_RR′|, we need to settle down the phase of T_RR′. In particular, to preserve the PH symmetry, the total phase around a plaquette with four NN vacancies can be either π or 0. One approach is to count the total flux number inside the plaquette. For the mid-gap states induced by VSL, every four NN vacancies form an L × L square, of which the total flux number is just L²/2 and the total phase around a plaquette is πL². After considering the compact condition, we have 0/π-flux inside each plaquette if L² is an even/odd number. In addition, we also use a numerical approach to check the above prediction of the flux number inside the plaquette.

We then consider the Haldane model with a square-VSL and the lattice constant L set to 3a. See the illustration in Fig. 10(a). From the above discussion it can be seen that there exists a π-flux in each square plaquette of VSL and a π/2-flux in each triangular plaquette of VSL. From the configuration of zero modes around vacancies in Fig. 10(a), one can see that the quantum states induced by VSL can be regarded as a new generation of lattice model, of which the vacancies play the role of the “atoms”.

As a result, we obtain an effective Hamiltonian that describes the vacancy-induced square lattice, written in the form

where T and T′ are the effective NN and NNN hopping parameters, respectively. Ā and B̄ label the sublattices of VSL. By comparing Eq. (14) with Eq. (18), one can see the corresponding relationship, ĉ_i ⟷ d̂_R, t ⟷ T, and t′ ⟷ ′T′. For the defect-induced quantum states, the energy spectrum and DOS are obviously similar to those of the parent Hamiltonian in Eq. (14).

From the DOS of the Haldane model with square-VSL for t′/t = 0.2 and L = 3a in Fig. 10(b), the mid-gap energy bands appear near the chemical potential. It is obvious that the MG states are induced by the VSL. In particular, we find that the mid-gap states also have an energy gap and exist two points with van Hove singularity. The DOS of the MG states is very similar to that of the parent TI. By fitting a curve to the DOS of the vacancy-induced quantum states, we obtain the effective hopping parameters as T = 0.22t and T′ = 0.07t.

To illustrate the topological properties of the Haldane model with square-VSL, we study its edge states. As a result, when the system has periodic boundary condition along y-direction but open boundary condition along x-direction, there exist the gapless edge states. On the other hand, we consider the half filling case, of which the system has zero Chern number. When the parent topological insulator has periodic boundary condition along both x-direction and y-direction, while the VSL has periodic boundary condition along y-direction but open boundary condition along x-direction, there may exist gapless edge states. See the numerical results in Fig. 10(c) for the case of t′/t = 0.2 and L = 3a, there indeed exist gapless edge states on the boundaries of VSL.

In addition, we also calculate the Haldane model on the honeycomb lattice with the honeycomb-VSL. The Haldane model on the honeycomb lattice also has the PH symmetry, and honeycomb-VSL can also induce a topological hierarchy insulator.

In general, the Haldane model with a square-VSL of odd number L² has nontrivial topological properties and is different from the traditional TIs. As a result, we call this type of topological insulator with topological MG states topological hierarchy insulators (THIs). In summary, THI from vacancy-superlattice is given by the following equation:

There also exist a particular type of topological hierarchy insulator, topological fractal insulator (TFI), which is a topological hierarchy insulator with infinite generations of self-similar L × L-square-VSLs. To construct a TFI, we first consider a topological hierarchy insulator with self-similar MG states. Because we can tune t′/t (or ΔE′/ΔE) by changing the NNN hopping parameter t′ of the parent TI or the VSL constant L, the MG states can be self-similar to the parent topological states for the case of T/T′ = t/t′. From numerical calculations, we obtain the self-similarity case, where smaller/bigger t/t′ leads to bigger/smaller T/T′. Now, under the self-similar condition T/T′ = t/t′, the effective Hamiltonian that describes the MG states becomes

where α = T/t (α < 1) is an energy-scaling ratio.

Next, we regard the vacancy-induced MG states as a parent TI and introduce an additional VSL on it. Using the recursive relation, we can see that there exist additional MG states inside the energy gap of the first generation MG states (which we refer to as generation-1 MG states). Consequently, we obtain a new generation of topological MG states (we call such MG states in the energy gap of the MG states as generation-2 MG states). Using the same procedure, we can construct a THI with n-generation MG states and eventually a THI with infinite generations of self-similar MG states that is just a TFI. Then, under the self-similar condition T/T′ = t/t′, the effective Hamiltonian of the generation-n MG states is given by

From the DOS in Fig. 10(d), we can see that the DOS of the MG states is always self-similar.

It is known that a self-similar object looks “roughly” the same on different scales and fractal is a particularly self-similar object that exhibits a repeating pattern displaying at every scale. The topological fractal insulator provides a unique example of topological matters with self-similarity, in which the topological properties always look the same on different energy scales (α²ⁿ) or different length scales (L²ⁿ).