Twisted bands and topological invariants

In an insulator, the energy gap at the Fermi level separates the electronic states into occupied bands and unoccupied bands. At zero temperature, there are no states at the Fermi level, and therefore no electronic conductivity. From this point of view, all insulators should behave in the same way. However, Chern insulators and TIs inform us that insulators can be further classified into trivial and non-trivial ones according to their topological properties of occupied electronic states in momentum space. Topologically non-trivial insulators have twisted bands and can be characterized by the non-zero topological invariants, such as Z or Z ₂. Insulators with infinitely large size (i.e., without a boundary) still look the same; however, in the presence of a boundary (say to vacuum), the topologically non-trivial insulators always support gapless boundary states protected by their topological properties, as distinguished from the trivial insulators. The existence of gapless (conducting) boundary states is responsible for the dissipationless quantum transport properties, such as the QAHE^{Reference Yu, Zhang, Zhang, Zhang, Dai and Fang23,Reference Chang, Zhang, Feng, Shen, Zhang, Guo, Li, Ou, Wei, Wang, Ji, Feng, Ji, Chen, Jia, Dai, Fang, Zhang, He, Wang, Lu, Ma and Xue24} and the QSHE,^{Reference Bernevig, Hughes and Zhang28,Reference König, Wiedmann, Brüne, Roth, Buhmann, Molenkamp, Qi and Zhang29} which are useful for emerging technologies with ultralow power consumption. There are several in-depth review articles on TIs^{Reference Hasan and Kane1,Reference Qi and Zhang2} in which readers can learn the fundamental mathematics and physics associated with the topological nature of insulators. See the Introductory article in this issue.

In a periodic solid crystal, the topology of electron bands and its twist in momentum space are due to the non-vanishing Berry phase^{Reference Berry46} and the gauge field,^{Reference Wilczek and Shapere47} and can be also related to the Wannier functions^{Reference Wannier48,Reference Wannier49} (a set of localized orbitals that represent the Bloch wave function in a solid) and their centers.^{Reference Qi and Zhang2,Reference Fu and Kane50} Let us consider a 2D system with a periodic boundary condition, whose occupied eigen wave functions are given as the Bloch states $\psi _{nk} \left( r \right) = e^{ik \cdot r} u_{nk} \left( r \right)e^{i\phi _n \left( k \right)} $, with n being the band index and u _nk(r) the cell-periodic part. We have explicitly written out the arbitrary phase ϕ_n(k), which is not assigned by the Schrödinger equation and depends on the explicit gauge choice in each calculation. The momentum k = (k _x, k _y) is defined in the 2D BZ, which is a torus (a closed manifold without boundary) due to the periodic boundary condition. Considering the Berry connection A _n(k) = i〈u _nk ⎢∇_k⎟ u _nk〉 (which is gauge dependent) and the Berry curvature Ω_n(k) = ∇_k × A _n(k) (which is gauge independent), the following equation holds as claimed by the Chern theorem:^{Reference Chern18}

(1)$$\sum\limits_n {\int_{BZ} {\Omega _n \left( k \right) \cdot dS} } = 2\pi Z.$$

Namely, the surface integration of Berry curvature over a torus must be an integer multiple of 2π. The integer number Z is a topological invariant, the Chern number. An insulator with a non-zero Z number is a Chern insulator, which supports chiral edge states^{Reference Haldane21} and the QAHE.^{Reference Yu, Zhang, Zhang, Zhang, Dai and Fang23,Reference Chang, Zhang, Feng, Shen, Zhang, Guo, Li, Ou, Wei, Wang, Ji, Feng, Ji, Chen, Jia, Dai, Fang, Zhang, He, Wang, Lu, Ma and Xue24} In such a case, the Z number is physically meaningful and is the number of chiral edge states, similar to that in the QHE.^{Reference Prange and Girvin17,Reference Stone20}

The band topology discussed previously can be understood intuitively as a band twist. Here we consider the simplest case: a 2D lattice with only one single occupied band (the band index n can therefore be neglected). Let us chose a special gauge that is smooth and periodic along the k _x direction, but not necessarily along the k _y direction (this choice is always possible^{Reference Soluyanov and Vanderbilt51}). Thus, we can do the integration in Equation 1, explicitly for the k _x direction with fixed k _y and obtain

(2)$$\eqalign{ & 2\pi Z = - \int_0^{2\pi } {\int_0^{2\pi } {dk_x dk_y \left( {\partial _x A_y - \partial _y A_x } \right)} } \cr & = \int_0^{2\pi } {dk_y \partial _y \left( {\int_0^{2\pi } {dk_x A_x \left( {k_x ,k_y } \right)} } \right)} \cr & = \int_0^{2\pi } {d\theta \left( {k_y } \right)} . \cr} $$

Here, $\theta \left( {k_y } \right) = \int_0^{2\pi } {dk_x A_x } $ is an angle calculated from the 1D integration of A _x(k _x, k _y) along the k _x axis for each fixed k _y. We can plot θ (k _y) over a cylinder surface (in the cylinder coordinates), as shown in Figure 1. k _y is plotted along the longitudinal direction, and the azimuth is the angle θ (k _y) for each fixed k _y. Moving from k _y = 0 to k _y = 2π, we can see the difference between the trivial insulator (Z = 0) and the Chern insulator (Z ≠ 0). The winding number of θ (k _y) over the cylinder surface is zero for the former (Figure 1a), and non-zero for the latter (Figure 1b). In this way, the Chern number can be related to the winding number of eigen wave functions, and it can be seen that the Chern insulator has “twisted” energy bands.

Figure 1.

Schematic illustration of band twist in a 2D Brillouin zone. The evolution of θ(k _y) is plotted over a cylinder surface, whose azimuth is θ and longitudinal axis is k _y. (a) Trivial insulator with Z = 0; (b) non-trivial Chern insulator with Z = 1; (c) trivial 2D insulator with Z ₂ = 0; (d) non-trivial 2D TI with Z ₂ = 1; and (e) trivial 2D insulator with |Z| = 2, but Z ₂ = 0 (see text for detailed explanation.) These figures are re-plotted based on Reference 52.

The above is for Chern insulators, which break the time-reversal symmetry (TRS). In the presence of TRS, however, the situation is different, because the total Z number (after summation over the band index n) should always be zero due to Kramers degeneracy—the degeneracy of spin up and down eigenstates in a time-reversal invariant system. In this case, the total number of occupied electronic states must be even, and we can, in principle, divide them into two classes, ψ^I and ψ^II, which are TRS related and are called Kramers pairs. Evaluating the Chern number for each class independently, if one class has Chern number Z, the other one must have Chern number –Z (due to the TRS), and the total Chern number of the system is zero. The question is whether |Z| can be used as the definition of the topological properties. The answer is no. This is because there are many different ways to divide the occupied wave functions into two classes in the presence of band degeneracies. A different way of dividing may lead to different |Z| numbers, which are different by an even number.^{Reference Roy32,Reference Fu and Kane50} Therefore, although |Z| is not a well-defined number, |Z| mod 2 is a topological invariant, which is called the Z ₂ number and can be used for the definition of TIs with TRS.^{Reference Kane and Mele26,Reference Fu, Kane and Mele30–Reference Roy32}

Let us again consider the simplest example of a 2D insulator with only two occupied eigenstates, ψ_1k and ψ_2k, which are related by TRS. The two states form a Kramers pair and must be degenerate at the time-reversal-invariant momentum (TRIM) of the BZ. We can evaluate the Berry connection and check the winding number of each occupied eigenstate (as we have done previously for the Chern insulator). In general, there are three situations (as shown in Figure 1c–e). First, if the winding numbers of two states are zero (i.e., |Z| = 0), this is the trivial situation with Z ₂ index v = 0 (Figure 1c). Second, if the winding numbers of two states are ±1 (i.e., |Z| = 1), this is the non-trivial topological insulator with a Z ₂ index v = 1, where the crossing of two θ(k _y) curves at k _y = π is protected by TRS (Figure 1d). Third, if the winding numbers of two states are ±2 (i.e., |Z| = 2), the two θ(k _y) curves must cross at some k _y other than the TRIM. Such crossings are not protected by TRS and can be removed by small perturbations, which drive the crossings into anti-crossings. As a result, the original winding number |Z| = 2 is reduced to zero, and the system becomes equivalent to the trivial case with a Z ₂ index v = 0 (Figure 1e).

The 2D Z ₂ topological insulator with TRS is different from a Chern insulator in the sense that it supports helical (rather than chiral) edge states and the QSHE.^{Reference Kane and Mele25–Reference Bernevig, Hughes and Zhang28} Here, the helicity of edge states means the orientation of its spin and momentum is locked, and the scattering between counter-propagating channels is forbidden within TRS. The Z ₂ number can be physically related to the number of helical edge state pairs mod 2. The Chern number is defined only for 2D insulators, while the idea of the Z ₂ number can be extended from 2D to 3D compounds^{Reference Fu, Kane and Mele30–Reference Roy32} by using the 2D Z ₂ numbers as building blocks. For 3D insulators, there are four independent indices for Z ₂ numbers, v ₀(v _xv _yv _z), as discussed in Reference 30. Here v ₀ is for the total Z ₂ number, and v _x,y,z are for the Z ₂ numbers in the x, y, z directions, respectively. Using our scheme presented prievously for the 2D Z ₂ number, we can obtain the three weak Z ₂ indices v _x,y,z directly from the three 2D Z ₂ numbers evaluated for the k _x = π, k _y = π, and k _z = π planes, respectively. For the total Z ₂ index v ₀, however, we have to consider two planes together. We can choose the k _x = 0 and π planes, or k _y = 0 and π planes, or k _z = 0 and π planes. If the Z ₂ numbers of two planes are different, we have v ₀ = 1, otherwise, if they are the same, we have v ₀ = 0.

Wilson loop method for the evaluation of topological invariants

In real materials, the band structures are usually very complicated, and the determination of topological numbers becomes essential and computationally demanding. The band degeneracy, either accidentally or due to symmetry, makes the numerical determination of phases of wave functions a tough task. Here the Wilson loop method is presented for determining topological indices efficiently.^{Reference Soluyanov and Vanderbilt51–Reference Ringel and Kraus53} This method is computationally easy, and is equally applicable for Chern insulators, Z ₂ TIs (2D and 3D), and the topological crystalline insulators,^{Reference Fu54,Reference Hsieh, Lin, Liu, Duan, Bansi and Fu55} which have attracted a lot of recent interest.

For band insulators with inversion symmetry, the Z ₂ indices can be easily computed as the product of half of the parity (Kramers pairs have identical parities) numbers for all the occupied states at the TRIM points.^{Reference Fu and Kane56} In such a case, parity is a good quantum number for the TRIM points and can be computed from the eigen wave functions obtained from first-principles calculations. Up until now, this method has been frequently used and is very efficient.

Unfortunately, for general cases where inversion symmetry is absent, evaluation of the Z ₂ number becomes difficult. There are mainly three different ways. The first method computes the Z ₂ numbers using the integration of Berry curvature Ω(k) over half of the BZ.^{Reference Fukui and Hatsugai57} In order to do so, one has to set up a fine mesh in k-space and calculate the corresponding quantities for each k point. Since the calculation involves the Berry connection A(k), one has to numerically fix the gauge on the half BZ, which is not easy for the realistic wave functions obtained by first-principles calculations. The second method starts from an artificial system with inversion symmetry, and then smoothly deforms the Hamiltonian toward the realistic one without inversion symmetry. If the energy gap never closes at any point in the BZ during the deformation process, the realistic system must share the same topological nature with the initial reference system, whose Z ₂ number can be easily counted by the parity eigenvalue formula.^{Reference Fu and Kane56} Unfortunately, making sure that the energy gap remains open on the whole BZ is very difficult numerically, especially in 3D. The third method directly calculates the boundary (edge or surface) states. Due to open boundary conditions, the first-principles calculation for the boundary states is numerically heavy. The Wilson loop method that we present here has the following advantages: first, it uses only the periodic bulk system; second, it does not require any gauge-fixing condition—thereby greatly simplifying the calculation; and third, it can be easily applied to general systems with or without inversion symmetry.

As discussed in the previous section, to determine the topological number (either Z or Z ₂), we must determine the key quantity θ(k _y) angle and its 1D evolution. In the single band model, the θ(k _y) angle is simply given as the 1D integration of A _x(k _x,k _y) along the k _x axis. However, for real compounds with multiple bands, its computation needs the Wilson loop method.^{Reference Yu, Qi, Bernevig, Fang and Dai52} Once the θ angle can be computed, all the topological numbers can be determined easily following the previous discussions. For the mathematical details of this method, refer to Reference 52.

The band-dependent phase θ_m(k _y) is exactly what we need to determine the winding number and the corresponding topological numbers, as discussed in the previous section. In practice, for a system with many bands, an arbitrary reference line can be drawn on the cylinder surface, parallel to the k _y-axis, and then a count can be taken of how many times the evolution lines of θ_m(k _y) cross the reference line.^{Reference Yu, Qi, Bernevig, Fang and Dai52} For the calculation of Z ₂ in time-reversal invariant systems, we can further reduce the computational task by considering only the 0 ≤ k _y ≤ π part. As an example, Figure 2 shows the calculated θ_m(k _y) and their evolution for the 3D topological insulator Bi₂Se₃ in the k _z = 0 and π planes.^{Reference Yu, Qi, Bernevig, Fang and Dai52} The θ_m(k _y) lines cross the reference line an odd number of times for the k _z = 0 plane, but an even number of times for the k _z = π plane. Therefore, it is topologically non-trivial with Z ₂ = 1.

Figure 2.

An example of the Wilson loop method for the 3D strong topological insulator Bi₂Se₃ with Z ₂ = 1. The calculated θ_m (k _y) for all occupied bands (m is the band index) are plotted as a function of k _y. The dashed lines are arbitrary reference lines. θ_m (k _y) crosses the reference line an odd number of times in half of the k _y direction (k _y varies from 0 to π or from π to 2π) for the k _z = 0 plane (upper panel), but an even number of times for the k _z = π plane (lower panel). These are re-calculated and re-plotted based on Reference 52.

Finally, in the whole process of the calculations, the advantages of Wannier functions can be used. Particularly, the Maximally localized Wannier functions (MLWFs) have been widely used.^{Reference Marzari and Vanderbilt58,Reference Souza, Marzari and Vanderbilt59} For example, after the self-consistent electronic structure calculations, the MLWF for the energy bands close to the Fermi level can be constructed, which can be reproduced usually to very high accuracy. In this process, however, some of the unoccupied bands need to be included.^{Reference Soluyanov and Vanderbilt51} Once the MLWFs are obtained, calculations can be performed for the eigen wave functions and their inner products for any k-points using the interpolation technique.^{Reference Wang, Yates, Souza and Vanderbilt60,Reference Yates, Wang, Vanderbilt and Souza61} This greatly simplifies the task for the computation of topological numbers. In addition, the MLWFs are also useful for the calculations of boundary states, as discussed later.

Boundary state calculations

After the determination of bulk topological numbers by either of the methods presented previously, the calculation of boundary states of real compounds becomes necessary. This is not just for the purpose of further confirmation of the topological properties, but is also very useful for direct comparison to experiments, such as using angle-resolved photoemission spectroscopy (ARPES) and scanning tunneling microscopy (STM).

In practical terms, two methods are usually used for the calculation of boundary states. First, the boundary states can be directly calculated from first-principles calculations by using a supercell of the object material (called a slab), which is intercalated between two vaccum layers, as is usually done for other surface or edge calculations.^{Reference Yip62} In this way, the atomic details, such as the surface or edge terminations, the atomic relaxations or reconstructions, and possible defects or absorptions, can be carefully treated and studied. This method is straightforward, although it is computationally demanding. Due to the presence of boundary states, we have to use a very large supercell with a thick vacuum region to avoid possible mutual coupling between the boundaries. A convergence test as a function of cell size is usually necessary.

Second, taking advantage of the MLWF as the basis, a Green's function (a tool to solve inhomogenous differential equations with boundary states) can be utilized to calculate the boundary states of a semi-infinite system.^{Reference Harrison63} The MLWF can be used as a basis to construct an effective low-energy Hamiltonian H _MLWF, which can be regarded as a tight-binding model with its hopping parameters determined from first-principle electronic structure calculations. The maximally localized property of MLWF guarantees short range hopping. Using H _MLWF as building blocks, a large Hamiltonian H _semi can be constructed for the semi-infinite system with only one boundary. This H _semi is a block tridiagonal matrix if only the hopping terms are considered between the nearest-neighbor building blocks. Then, the iterative method^64–66 can be used to solve the problem and get the projected Green’s functions on to the boundary. From this boundary Green’s functions, the necessary information can be obtained, such as the charge density of states and the spin density of states at the boundary.^{Reference Zhang, Liu, Qi, Deng, Dai, Zhang and Fang34,Reference Dai, Hughes, Qi, Fang and Zhang67} This procedure is an approximate one, and it cannot treat the atomic details at the boundary precisely. However, it provides the most important information on the gapless boundary states of TIs.

Band inversion mechanism for materials selection

The most intuitive picture and practically useful guideline for searching for TIs is the band inversion mechanism. It involves the energy order switching of low energy electronic bands around certain TRIM points in the BZ. The mechanism was first demonstrated for the HgTe/CdTe quantum well structure^{Reference Bernevig, Hughes and Zhang28} and was later found to be applicable to many other materials, including 2D and 3D TIs^{Reference Ando33,Reference Zhang, Liu, Qi, Dai, Fang and Zhang37,Reference Zhang and Zhang68} and Chern insulators.^{Reference Yu, Zhang, Zhang, Zhang, Dai and Fang23}

As schematically illustrated in Figure 3, we consider two bands (indicated in red and blue, respectively) close to the Fermi level around the Γ point. Far away from the Γ point, the red state is occupied and energetically lower than the blue state. However, close to the Γ point, they have inverted energy order with the blue state lower and occupied at Γ. Assuming that the two bands do not couple when the SOC is absent, they must cross at certain k-points away from Γ. In the presence of the SOC, the two bands will, in general, couple together and open up a gap at the crossing point. As a result, the system becomes an insulator, and most importantly, this gapped state may be topologically non-trivial depending on the coupling term.

Figure 3.

Schematic plot for the band inversion mechanism as discussed in the text. Note. SOC, spin–orbit coupling.

To be more explicit, we consider two parabolic bands with the simplest energy dispersion E _±(k) = ±(−m + k ^{Reference Qi and Zhang2}). The + and – signs are for the blue and red bands, respectively. The band inversion happens when m > 0 is satisfied. Using two bands as a basis, an effective 2×2 Hamiltonian is constructed:

(3)$$H_{eff} \left( k \right) = \left[ {\matrix{ {E_ + \left( k \right)} & {M\left( k \right)^ * } \cr {M\left( k \right)} & {E_ - \left( k \right)} \cr } } \right],$$

where M(k) is the coupling between the two bands and determines the topological properties of the resulting electronic structure. If m < 0 (i.e., without band inversion), in general, nothing is expected but a normal insulator. However, if m > 0 for the band inverted system, we may find many interesting topological phases depending on the choice of M(k). For example, if M(k) = (k _x ± ik _y)ⁿ (n ≥ 1), we will get a Chern insulator with Chern number Z = n (for a 2D system). For time-reversal invariant systems, the two bands should be replaced by four bands with each color (blue and red) representing a Kramers pair of bands, and the model Hamiltonian becomes 4×4. It is block diagonal in the simplest form, such as the Hamiltonian proposed for HgTe/CdTe quantum well.²⁸ If there is a form of M(k) = (k _x ± ik _y)^{Reference Qi and Zhang2n−1}, a TI with Z ₂ = 1 is expected. In the presence of SOC, these forms of the coupling term M(k) are in general allowed and depend on the material details. The band inversion may occur between the s–p, p–p, p–d, or d–f states, which leads to different material classes of TIs.^{Reference Zhang and Zhang68}

The band inversion mechanism discussed previously is useful from a practical view-point, yet it is a crude picture for real compounds. First, this picture comes from the low-energy effective Hamiltonian of a continuous model. If the model is put onto a lattice, we have to be careful about the number of band inversions. Particularly for Z ₂ TIs, if the band inversion happens an even number of times or it happens around an even number of TRIM points, we will return to trivial insulators. Second, real compounds may have many low-energy bands, which must be carefully analyzed. Third, the possible M(k) terms can also be complicated and need careful treatment. In any case, the topology of the whole electronic structure must be determined by first-principles calculations, as presented in previous sections.

Topological electronic materials

First-principles calculations based on DFT have played an important role in studies on topological electronic materials. Many compounds have been theoretically predicted and experimentally confirmed. Although those computational techniques are still not perfect, they have demonstrated their power in materials selection and design, thus greatly motivating the field. It is a hard task to review and summarize the whole field and all known topological materials completely. Instead, here we discuss some recent interesting progress, in which first-principles calculations have demonstrated their predictive power.

Correlated topological insulators

Strong correlation effects can enrich the topological states in condensed matter and lead to “correlated topological states” (CTS), which can be divided into two categories. The first type of CTS can be adiabatically connected to a non-interacting topological state, with their wave functions being single Slater determinants, namely, a matrix form of many electron wave function proposed by J.C. Slater satisfying anti-symmetry requirements and consequently the Pauli principle by changing sign upon exchange of two electrons (or other fermions). This type of CTS is not new, it is a topological state with a renormalized quasi-particle band structure. The other type of CTS cannot be adiabatically connected to any non-interacting topological state. The typical examples in this category are fractional Chern insulators,^{Reference Neupert, Santos, Chamon and Mudry80–Reference Sheng, Gu, Sun and Sheng83} fractional TIs,^{Reference Levin and Stern84,Reference Maciejko, Qi, Karch and Zhang85} quantum spin liquids,^{Reference Meng, Lang, Wessel, Assaad and Muramatsu86} and topological Mott Insulators.^{Reference Pesin and Balents87,Reference Raghu, Qi, Honerkamp and Zhang88} Until now, none of the CTS in the second category have been confirmed in actual materials.

Recently, a special group of mixed valence compounds containing rare-earth elements, including SmB₆, YbB₆, and YbB₁₂, has drawn great attention.^{Reference Dzero, Sun, Galitski and Coleman89–Reference Weng, Zhao, Wang, Fang and Dai91} From the band theory point of view, the “mixed valence” suggests partial occupation of localized 4f bands and implies band inversion between them and dispersive 5d bands, which may lead to non-trivial topological features. Due to the strong Coulomb interaction among the 4f electrons, DFT, which has been successfully applied before, cannot be applied directly to these materials. Advanced techniques, which can treat the strong correlation effects more rigorously, are needed. In addition, it is important to understand how to determine the topological nature of a correlated material, which involves many-body wave functions.

The LDA+Gutzwiller method,^{Reference Deng, Wang, Dai and Fang92} which combines the LDA with the Gutzwiller method^{Reference Gutzwiller93–Reference Gutzwiller95} (a simplified method for projected many-body wave function) for correlated states, is a powerful tool to treat the CTS. This method starts from the LDA Hamiltonian, and supplements it with the interaction term for the localized 4f states, treated using the Gutzwiller-type trial wave functions. The scheme preserves the nice aspect of being variational, however, it is beyond LDA (and also LDA+ U; another modification that takes orbital dependence of the Coulomb and exchange interactions into account^{Reference Anisimov, Zaanen and Andersen96}) because the total energy now relies on the balance between the renormalized kinetic energy of quasi-particle motion and the local interaction energy, which is configuration-dependent. For details of this method, refer to References 92, 97, and 98. We note that including the complete form of local interactions is essential to obtain the correct electronic structure.

The Z ₂ classification of band insulators can also be achieved by looking at the TIs response to interactions with an external electromagnetic field, namely the topological magnetoelectric effect (TME).^{Reference Qi, Li, Zang and Zhang99} A correlated insulator is a TI if the θ-angle defined in TME is π. Wang et al. have developed and simplified the linear response theory for the coefficient of TME.^{Reference Wang and Zhang100} For interacting systems, when the self-energy contains no singularity along the imaginary axis, the formula for the TME coefficient only requires the single particle Green’s function at zero frequency, $\hat g\left( {k,0} \right)$. Because the singular point for the self-energy along the imaginary axis only appears for a Mott insulator with completely localized f-orbitals, which is not the case for Kondo or mixed valence insulators, we can safely apply the previous method to them. Therefore, the way to determine the Z ₂ invariance is by diagonalizing the Hermite matrix $ - \hat g\left( {k,0} \right)^{ - 1} &#x003D; \hat H_0 - \mu _f &#x002B; \hat \Sigma \left( 0 \right)$ and treating the eigenstates with the negative eigenvalues as the “occupied states.” In the Gutzwiller approximation, the low-energy single particle Green’s function of an interacting system can be expressed by the quasi-particle effective Hamiltonian H _eff and the quasi-particle weight $\hat z$ (which in general is a matrix) as,

(4)$$\hat g\left( {k,i\omega } \right) &#x003D; \hat zi\omega - \hat H_{eff} &#x002B; \mu _f &#x002B; \hat g_{ic} \left( {k,i\omega } \right),$$

where the second term describes the incoherent part of the Green’s function, which can be ignored for low frequencies. Therefore, the Z ₂ invariance in the LDA+Gutzwiller method is determined by the occupied eigenstates of the Gutzwiller effective Hamiltonian H _eff, which can be interpreted as the band structure of the “quasi-particles.” Then the Z ₂ invariance can be determined by using the same method developed for the non-interacting insulators.^{Reference Fu, Kane and Mele30,Reference Yu, Qi, Bernevig, Fang and Dai52,Reference Fu and Kane56}

Using the LDA+Gutzwiller method, we have identified that SmB₆ is a topological insulator, with the quasi-particle structures strongly renormalized by the Coulomb interaction.^{Reference Lu, Zhao, Weng, Fang and Dai90} The band inversion happens between the itinerant 5d and the 4f j = 5/2 bands at three X points in its BZ, which leads to three Dirac cones on a generic surface. The topological properties of SmB₆ have been verified recently by various experiments, including transport,^{Reference Wolgast, Kurdak, Sun, Allen, Kim and Fisk101,Reference Li, Xiang, Yu, Asaba, Lawson, Cai, Tinsman, Berkley, Wolgast, Eo, Kim, Kurdak, Allen, Sun, Chen, Wang, Fisk and Li102} ARPES,^{Reference Xu, Shi, Biswas, Matt, Dhaka, Huang, Plumb, Radović, Dil, Pomjakushina, Conder, Amato, Salman, Paul, Mesot, Ding and Shi103–Reference Jiang, Li, Zhang, Sun, Chen, Ye, Xu, Ge, Tan and Niu105} and STM^{Reference Yee, He, Soumyanarayanan, Kim, Fisk and Hoffman106} measurements. Another similar material in this category is YbB₆,^{Reference Weng, Zhao, Wang, Fang and Dai91} for which the band inversion occurs between 5d and the 4f j = 7/2 bands. Since the 4f shell in YbB₆ is nearly fully occupied, the correlation effect as well as the renormalization of the quasi-particle states are much weaker compared to SmB₆. Besides the Z ₂ TIs, the topological crystalline insulator^{Reference Fu54,Reference Hsieh, Lin, Liu, Duan, Bansi and Fu55} can also be found in this class of correlated materials. Applying the LDA+Gutzwiller method again, we found that YbB₁₂^{Reference Weng, Zhao, Wang, Fang and Dai91} is such a correlated topological crystalline insulator characterized by a mirror Chern number of 2, which leads to an even number of Dirac cones on some specific surfaces.

Perspectives and conclusion

There are many other success stories on the selection and design of topological electronic materials, which are not included in this article. As we have mentioned in the introduction, these successes are, on one hand, due to the rapid development in the field of first-principles calculations, and, on the other hand, due to the robustness of topological properties (i.e., its tolerance to error bars). Nevertheless, we are still facing many challenges that require further efforts.

One serious issue is actually an old problem: the LDA or the GGA used in most of the calculations tends to underestimate the bandgap of semiconductors. From the band inversion mechanism, we learn that topological properties crucially depend on the sign of the bandgap, which must be accurately calculated in order to make correct predictions, particularly for narrow gap semiconductors. For example, HgS has been predicted to have an inverted band structure and be a TI based on LDA.^{Reference Virot, Hayn, Richter and van den Brink107} This conclusion was challenged recently by accurate GW calculations,^{Reference Svane, Christensen, Cardona, Chantis, van Schilfgaarde and Kotani108} which suggested that HgS is a normal semiconductor with a positive bandgap. Therefore, from the calculation point of view, when we encounter a compound with small band inversion, we usually need more advanced and accurate methods, like the GW approximation,^{Reference Aryasetiawan and Gunnarsson109} the Hybrid functional,^{Reference Kim and Jordan110,Reference Heyd, Scuseria and Ernzerhof111} or the modified Becke-Johnson approximation proposed recently.^{Reference Tran and Blaha112}

Another problem is even more serious: LDA or GGA, which is basically part of mean-field theory, are qualitatively wrong for strongly correlated systems and should not be used in general for such systems. This is of course a challenging problem and is one of the frontier subjects of modern condensed matter physics. However, as we have discussed, correlated TIs are not only fundamentally important, they are also promising from a practical material point of view; advanced techniques are highly demanded. The computationally cheap and efficient treatment based on the LDA+U scheme may sometimes improve the situation, but in most cases, it is still not good enough. We are happy to see progress in recent developments toward this direction, and the newly developed LDA+Gutzwiller and LDA+DMFT (dynamic mean field theory) methods are being increasingly used to treat correlated TIs. More efforts in this direction are certainly required.

The great successes of studies on TIs are definitely not the end of this field, but rather this should be regarded as an indication of the opening of a much wider field, namely “topological states” or “topological quantum matter.” The discovery of TIs has stimulated the search for newer quantum states and quantum matter, which are topologically non-trivial. For example, if we consider crystal symmetries of solids, we can get other topological states,^{Reference Slager, Mesaros, Juričić and Zaanen113,Reference Hughes, Prodan and Bernevig114} which may be related to the quantization of polarization.^{Reference Jadaun, Xiao, Niu and Banerjee115,Reference Turner, Zhang, Mong and Vishwanath116} The topological semimetal is another interesting state, which is qualitatively different from TIs, and has attracted great recent attention.^{Reference Wan, Turner, Vishwanath and Savrasov4–Reference Wang, Weng, Wu, Dai and Fang11,Reference Liu, Jiang, Zhou, Wang, Zhang, Weng, Prabhakaran, Mo, Peng, Dudin, Kim, Hoesch, Fang, Dai, Shen, Feng, Hussain and Chen42–Reference Borisenko, Gibson, Evtushinsky, Zabolotnyy, Buechner and Cava44} The concept of topology can also be combined with superconductivity, which leads to “tolopogical superconductors”—a promising and challenging field awaiting extensive exploration.^{Reference Qi and Zhang2,Reference Fu and Kane12–Reference Fu and Berg15} For this, however, the predictive power of first-principles calculations is rather limited. For all of these exotic topological quantum states, further explorations for better materials are certainly necessary. We are still in the early stages of possible applications of topological quantum states, and materials are the key to advancing the field.