Revisiting (higher order and crystalline) topology in old models of lead telluride

TL;DR

This paper reexamines a classic PbTe model, revealing it as a higher-order topological crystalline insulator with symmetry-protected boundary modes, and discusses its relevance to real PbTe materials.

Contribution

It demonstrates that the old PbTe model is a higher-order topological crystalline insulator diagnosable by symmetry indicators, linking it to recent predictions for real PbTe.

Findings

The low energy theory corresponds to a higher-order topological crystalline insulator.

Gapless fermions on domain walls are topological boundary modes.

Any symmetric completion of the model shares the same topological phase.

Abstract

In this work, we revisit the model of PbTe presented in [E. Fradkin, E. Dagotto, and D. Boyanovsky, Phys. Rev. Lett. 57, 2967 (1986)]. We show that the low energy theory of this model corresponds to a (higher-order) topological crystalline insulator in space group $Fm\bar{3}m1'$, diagnosable by symmetry indicators. We show that the gapless fermions found on antiphase domain walls are the topological boundary modes of the system, due to a nonvanishing mirror Chern number. Furthermore, we show that any symmetric completion of the model must be in this same topological phase. Finally, we comment on the relationship of this model to realistic PbTe, which has recently been predicted to have a phase which realizes same bulk symmetry indicators.

Tables

Table 1. Table 1: Little group irreps subduced by each of the elementary band representations in our model. The first column gives the name of the elementary band representation (EBR). The subsequent columns give the little group irreps at each of the high symmetry points.

EBR	$\Gamma$	$X$	$L$	$W$
${\left({\overline{E}}_{1g}\right)}_{4a}\uparrow G$	${\overline{\Gamma }}_6$	${\overline{X}}_6$	${\overline{L}}_9$	${\overline{W}}_6$
${\left({\overline{E}}_{1g}\right)}_{4b}\uparrow G$	${\overline{\Gamma }}_6$	${\overline{X}}_6$	${\overline{L}}_8$	${\overline{W}}_7$

Table 2. Table 2: Occupied band irreps for SnTe and PbTe at the high symmetry points. Irreps are listed in order of increasing energy, i.e. those states closest to the fermi level appear at the end of the list. Note that the irreps at Γ Γ \Gamma and X 𝑋 X are identical for the two materials.

$𝐤$	SnTe	PbTe
$\Gamma$	${\overline{\Gamma }}_6,{\overline{\Gamma }}_6,{\overline{\Gamma }}_8,{\overline{\Gamma }}_{11}$	${\overline{\Gamma }}_6,{\overline{\Gamma }}_6,{\overline{\Gamma }}_8,{\overline{\Gamma }}_{11}$
$X$	${\overline{X}}_6,{\overline{X}}_6,{\overline{X}}_8,{\overline{X}}_8,{\overline{X}}_9$	${\overline{X}}_6,{\overline{X}}_6,{\overline{X}}_8,{\overline{X}}_8,{\overline{X}}_9$
$L$	${\overline{L}}_9,{\overline{L}}_8,{\overline{L}}_8,{\overline{L}}_4{\overline{L}}_5,{\overline{L}}_9$	${\overline{L}}_8,{\overline{L}}_9,{\overline{L}}_9,{\overline{L}}_6{\overline{L}}_7,{\overline{L}}_8$
$W$	${\overline{W}}_6,{\overline{W}}_7,{\overline{W}}_7,{\overline{W}}_6,{\overline{W}}_7$	${\overline{W}}_7,{\overline{W}}_6,{\overline{W}}_7,{\overline{W}}_6,{\overline{W}}_6$

Equations

$${\mathbf{t}}_1=\frac{1}{2}(\hat{\mathbf{y}}+\hat{\mathbf{z}}),{\mathbf{t}}_2=\frac{1}{2}(\hat{\mathbf{x}}+\hat{\mathbf{z}}),{\mathbf{t}}_3=\frac{1}{2}(\hat{\mathbf{x}}+\hat{\mathbf{y}}).$$

$$H_{\mathrm{FDB}}=m{\tau }_z+t{\tau }_y\sum _{\mu =x,y,z}{\sigma }_{\mu }\cos \left(\frac{k_{\mu }}{2}\right),$$

$$\Delta (C_{3,111})$$

$$\Delta (I)$$

$$\Delta (T)$$

$$\Delta (C_{4x})={\tau }_z{\sigma }_y\exp (-i\pi /4{\sigma }_x).$$

$$\Delta (C_{4x})=\exp (-i\frac{\pi }{4}{\sigma }_x).$$

$$H_{\mathrm{SOC}}=t{\tau }_y\sum _{\mu ,\nu ,\lambda =x,y,z}\left({ϵ}^{\mu \nu \lambda }{\sigma }_{\mu }\sin \frac{k_{\nu }}{2}\sin k_{\lambda }\right).$$

$$H_{\mathrm{hop}}=\sum _{\mu =x,y,z}\left[{\delta }_1{\tau }_z\left(1+\cos k_{\mu }\right)+{\delta }_2{\tau }_x\cos \frac{k_{\mu }}{2}\right].$$

$$H_{\mathrm{iFDB}}=m{\tau }_z+H_{\mathrm{SOC}}+H_{\mathrm{hop}}.$$

$$\nu ={\kappa }_1-2{\kappa }_4\mathrm{mod}8,$$

$${\kappa }_1^{\mathrm{FDB}}=\{\begin{array}{ll}0,& m>0\\ 4,& m<0\end{array}$$

$$H_{\mathrm{FDB}}(L+\mathbf{k})\approx m{\tau }_z+\sum _{\mu =x,y,z}\frac{t}{2}{\tau }_y{\sigma }_{\mu }{k}_{\mu }.$$

$$H_{\mathrm{iFDB}}(L+\mathbf{k})$$

$$+\sqrt{3}{\delta }_2{k}_c{\tau }_x+2T\sqrt{3}\left(k_a{\tau }_y{s}_a+k_b{\tau }_y{s}_b\right)$$

Full text

Introduction.

One of the most striking features of topological insulators is the presence of protected gapless modes at surfaces, interfaces, and defectsSu et al. (1979); Teo and Kane (2010). The best-known examples of this are the helical modes at the boundary of a two-dimensional “strong” topological insulatorKane and Mele (2005); Bernevig et al. (2006); König et al. (2007), and the single Dirac fermion at the two-dimensional boundary of a three-dimensional “strong” topological insulatorFu et al. (2007); Xia et al. (2009). Both of these surface modes are protected by time-reversal symmetry alone. It has recently been appreciated that, with additional crystalline symmetries, more exotic topologically protected boundary features may emerge, such as the multiple Dirac fermions at symmetric boundaries of a mirror Chern insulatorTeo et al. (2008); Hsieh et al. (2012); Tanaka et al. (2012). Even more surprising, “higher-order TIs” with (roto)-inversion symmetries may feature topologically protected “hinge” modes which propagate on boundaries of two (or more) fewer dimensions than the bulkFu (2011); Wieder et al. (2018a); Song et al. (2017); Benalcazar et al. (2017a, b); Schindler et al. (2018a); Po et al. (2017); Kruthoff et al. (2017); Bradlyn et al. (2017); Khalaf et al. (2018); Song et al. (2018).

In light of these recent developments, we may be tempted to take a fresh look at the old observation that a simplified tight-binding model for PbTe was found to host four Dirac fermions on a two-dimensional antiphase boundaryFradkin et al. (1986); Fradkin (1987); Boyanovsky et al. (1987). It was realized early on that this effective model did not quite respect the symmetries of the crystalWilczek (1987); Haldane (1988), and the even number of Dirac cones disqualified this model from being a strong TI regardlessTchernyshyov (2000); Fu and Kane (2007). Furthermore, realistic PbTe is known to have zero mirror Chern numberHsieh et al. (2012), and nevertheless the effective model does not possess the requisite mirror symmetry. It is natural to ask, then: are the domain wall fermions of Refs. Fradkin et al., 1986; Fradkin, 1987; Boyanovsky et al., 1987 a signifier of any topological crystalline insulating (TCI) phase?

In this work, we answer this question in the affirmative. First, we review the model of Ref. Fradkin et al., 1986 (hereafter referred to as the FDB model), and show how to modify it in order to respect the symmetries of the cubic space group $Fm\bar{3}m1^{\prime}$ (225)111the final $1^{\prime}$ denotes the presence of time-reversal symmetryAroyo (2016). In particular, we will show that the low-energy effective models obtained from either the FDB model or our symmetric model are unitarily equivalent. Next, we will show that our improved tight-binding model captures the transition between a trivial insulator and a mirror Chern and higher-order topological insulator (HOTI), with symmetry-indicated topological index $\nu=4\in\mathbb{Z}_{8}$. Because the dynamics of this transition are captured by the low energy effective model, we will show that the four Dirac cones in the FDB model at an antiphase domain wall are the boundary modes of this mirror and higher-order topological insulator. We will show that any symmetric completion of the FDB model has $\nu=4$, with mirror Chern number $\nu_{m_{1\bar{1}0}}=2$ and higher-order “$S_{4}$” invariant $\delta_{S_{4}}=1$Khalaf et al. (2018); Song et al. (2018); Po et al. (2017); Vergniory et al. (2019) (Here $m_{1\bar{1}0}$ denotes the mirror about $\mathbf{\hat{x}}-\mathbf{\hat{y}}$, and $S_{4}$ signifies a fourfold roto-inversion). Finally, we will connect these results to the recent predictionVergniory et al. (2019); Zhang et al. (2019); Tang et al. (2019) that PbTe may be, in some cases, a TCI/HOTI. In doing so, we will see the importance of careful structural determination for finding small-gap topological materials.

Effective tight-binding model

Let us start by reviewing the FDB model, exploring its shortcomings, and constructing an improved model with the same phenomenology. Let us start as did the authors of Ref. Fradkin et al., 1986 by noting that PbTe has a rocksalt structure, with the symmetries of space group $Fm\bar{3}m1^{\prime}$ (225), the centrosymmetric, symmorphic space group with a face-centered cubic Bravais lattice and octahedral point group. We take as a basis for the Bravais lattice (we set the lattice constant $a=1$ for simplicity)

[TABLE]

The rocksalt structure has Te atoms located at the $4a$ Wyckoff position, with reduced coordinates $\mathbf{q}_{a}=(0,0,0)$, and Pb atoms at the $4b$ Wyckoff position with reduced coordinates $\mathbf{q}_{b}=(1/2,1/2,1/2)$. The point group of PbTe is generated by a threefold rotation $C_{3,111}$ about the body diagonal of the unit cell, a fourfold rotation $C_{4x}$ about the $\mathbf{\hat{x}}=\mathbf{t}_{2}+\mathbf{t_{3}}-\mathbf{t}_{1}$ axis, spatial inversion $I$, and time-reversal symmetry $T$. The original model of FDB consisted of spin-1/2 s orbitals on the $4a$ and $4b$ Wyckoff positions as a proxy for the Te and Pb atoms (whether we start from s or p orbitals, the sign of the inversion matrix will change, but all topological properties remain invariant). They added a staggered on-site potential taking opposite values $\pm m$ on the Pb and Te sites, and a nearest-neighbor spin-dependent hopping. Letting $\vec{\tau}$ be a set of Pauli matrices acting in the orbital (Te,Pb) basis, and letting $\vec{\sigma}$ be the set of Pauli matrices acting on spin, we can write the Bloch Hamiltonian for the FDB model as

[TABLE]

where we have taken the liberty of restoring a neglected factor of $i$ in the hopping term to restore time-reversal symmetryHaldane (1988)222While this makes our Eq. (2) formally different from the original FDB model, we view the lack of time-reversal symmetry in Ref. Fradkin et al. (1986) to be a simple oversight. Eq. (2) is manifestly time-reversal, inversion, and $C_{3,111}$ symmetric; each of these symmetries acts trivially in the orbital $\vec{\tau}$ space, and as a rotation (the identity rotation for inversion) on the spin degrees of freedom. Precisely, we have for these symmetries $\{g\}$ that $\Delta(g)^{-1}H_{\mathrm{FDB}(\mathbf{k})}\Delta(g)=H_{\mathrm{FDB}}(g\mathbf{k})$, with the matrix representation $\Delta$ specified by

[TABLE]

The spectrum of this Hamiltonian consists of two sets of doubly-degenerate bands (due to $IT$ symmetry) separated by a spectral gap which is smallest at the $L$ point $(1/2,1/2,1/2)$ and given by $\delta E_{L}=2m$. Depending on the sign of $m$, there is a band inversion at the $L$ point: for $m>0$ the valence bands carry the representation $\bar{L}_{9}$ of the little group of $L$ (with inversion eigenvalues $(-1,-1)$), while for $m<0$ the valence bands carry the representation $\bar{L}_{8}$ of the little group of $L$ (with inversion eigenvalues $(+1,+1)$). Because there are four $L$ points in the FCC Brillouin zone, this is not a $\mathbb{Z}_{2}$ nontrivial TIFu and Kane (2007), and so we must look for nontrivial TCI invariants.

Here, however, we run into a problem: The Hamiltonian $H_{\mathrm{FDB}}$ is $C_{4x}$-symmetric, but the matrix representative of $C_{4x}$ can be seen to be

[TABLE]

This has the unfortunate property that $\Delta(C_{4x})^{4}=+1$, rather than $-1$ as needed for a double-valued representation. Repairing this by multiplying by a factor of $\sqrt{i}$ is futile, as then $\Delta(C_{4x})$ and $\Delta(T)$ no longer commute. Thus, the FDB Hamiltonian does not have the symmetries of the space group $Fm\bar{3}m1^{\prime}$.

To repair the symmetries, we seek hopping terms which vanish at the $L$ point, are nonvanishing everywhere else, and transform in the representation given by Eqs. (3–5) along with

[TABLE]

In this way, we will replicate the band inversion at $L$ in our symmetric model. Let us first fix the spin-orbit coupling term. Noting that the matrices $\tau_{x,y}\sigma_{\mu}$ transform in a pseudovector representation under roto-inversions, we need them to be multiplied by functions of $\mathbf{k}$ which are also pseudovectors. Combining this observation with time-reversal invariance and the boundary conditions on the Bloch Hamiltonian, we find that the simplest choice of SOC term which vanishes at $L$ is

[TABLE]

For now, let us overlook the long (5-th nearest neighbor) range of this coupling in light of its mathematical simplicity. We may be tempted to take $m\tau_{z}+H_{\mathrm{SOC}}$ as our improved Hamiltonian, however $H_{\mathrm{SOC}}$ vanishes along the whole $\Gamma-L$ line, rather than just at the $L$ point. To remedy this, we can add two additional spin-independent hopping terms

[TABLE]

We take for our full improved Hamiltonian

[TABLE]

when $m=0$ with $t$ and $\delta_{2}$ nonzero, this model is gapless only at the $L$ point. for $m\neq 0$, a nonzero $\delta_{1}$ ensures that spectrum of this Hamiltonian is gapped. We see that there is thus an insulator-to-insulator transition driven by an inversion of bands at the $L$ point, just as in the original $FDB$ model. We show the bulk spectrum for positive and negative values of $m$ in Fig. 1a and 1b.

Topological properties

Let us now examine the topological nature of this band-inversion transition. We will start by analyzing the band representations in the model following Refs. Bradlyn et al. (2017); Cano et al. (2018); Vergniory et al. (2019). Following the notation of the Bilbao Crystallographic ServerAroyo et al. (2011, 2006a, 2006b); Elcoro et al. (2017), our model is induced from orbitals transforming under the $\bar{E}_{1g}$ representation of the point group $O_{h}$ on the $4a$ site, and a second set of orbitals transforming in the same representation on the $4b$ site. The four bands in our model thus transform under the composite band representation $(\bar{E}_{1g})_{4a}\uparrow G\oplus(\bar{E}_{1g})_{4b}\uparrow G$. These two elementary band representations subduce representations of the little group $G_{\mathbf{k}}$ at each of the high symmetry momenta $\Gamma,X,W,$ and $L$ in the Brillouin zone; we summarize these representations in Table 1 below. In the atomic limit of $H_{iFDB}$, we have $m\gg t,\delta_{1},\delta_{2}$, and so the valence bands of our model transform in the $(\bar{E}_{4b}\uparrow G)$ elementary band representation, with occupied little group representations $\bar{\Gamma}_{6},\bar{X}_{6},\bar{W}_{6},\bar{L}_{9}$; this can be checked explicitly using the representation $\Delta$ defined in Eqs. (3–5) and (7) and taking into account the boundary conditions $|u_{n\mathbf{k}+n_{i}\mathbf{g}_{i}}\rangle=(\tau_{z})^{\sum_{i}n_{i}}|u_{n}\mathbf{k}\rangle$, where $\{\mathbf{g}_{i}\}$ is a basis for the reciprocal lattice, and $\{n_{i}\}$ are integers; this expresses the “tight-binding gauge” boundary conditionsVanderbilt (2018). Upon inverting bands by taking $t\gtrsim\delta_{2}\gg\delta_{1}>0>m,|m|\ll t$, we see that the occupied band irreps at $\Gamma,X$ and $W$ have not changed. At $L$ however, the wavefunctions now transform under the $\bar{L}_{8}$ irrep. Examining the full table of elementary band representations for $Fm\bar{3}m1^{\prime}$Bilbao Crystallogr. Server (2017), we find that the collection $(\bar{\Gamma}_{6},\bar{X}_{6},\bar{W}_{6},\bar{L}_{8})$ of occupied little group representations cannot be subduced by an integer sum or difference of elementary band representations; we thus conclude that this phase of our model is a symmetry-indicated, stable topological crystalline insulator. Going further, we can attempt to express the irrep multiplicities in this model as a rational sum of those in EBRs. Reading off the denominator of the rational coefficients, we find that symmetry-indicated TCIs in this space group are classified by an index $\nu\in\mathbb{Z}_{8}$. The index $\nu$ can be expressed as Po et al. (2017); Khalaf et al. (2018); Song et al. (2018)

[TABLE]

where $4\kappa_{1}$ is the sum of occupied band inversion eigenvalues, and $2\sqrt{2}\kappa_{4}$ is the sum of occupied band $IC_{4z}$ eigenvalues. Note that both $\kappa_{1}$ and $\kappa_{4}$ are integersKhalaf et al. (2018). For our model, we have $\nu=4$ in the topological phase, and $\nu=0$ in the trivial phase. As shown in Refs. Khalaf et al., 2018; Song et al., 2018, a minimal model for a TCI with $\nu=4$ arises through a “stacking” (in Hilbert space) of four $\mathbb{Z}_{2}$ topological insulators. We thus expect to find protected gapless states on symmetric surfaces of this model, as we will discuss further below.

Furthermore, note that while the original FDB model is not $C_{4z}$ symmetric, it is inversion symmetric. Thus, the sum of occupied band inversion eigenvalues $\kappa^{\mathrm{FDB}}_{1}$ can be computed for that model, and is given by

[TABLE]

Next, note that $IC_{4z}$ is not in the little group of any of the $L$ points. Therefore, in any symmetric completion of the FDB model the occupied band $IC_{4z}$ eigenvalues do not change as a function of $m$, and so neither does $\kappa_{4}$. Thus, in any symmetric completion of the FDB model, the index $\nu$ must change by $4$ as the sign of $m$ changes. Assuming additionally that the $m>0$ phase is connected to the (unobstructedBradlyn et al. (2017)) atomic limit, we deduce that the band inversion in the FDB model becomes, when cubic symmetry is enforced, the transition between phases with $\nu=0$ and $\nu=4$.

As discussed in detail in Refs. Song et al., 2018; Khalaf et al., 2018, the value of $\nu$ does not uniquely determine the topological phase of a system in space group $Fm\bar{3}m1^{\prime}$. In particular, with $\nu=4$, there are two possible types of topological phase: the first has a mirror Chern number $\nu_{m_{z}}=4\mod 8$ associated with the $m_{z}$ mirror symmetry, while the second has both a mirror Chern number $\nu_{m_{1\bar{1}0}}=2\mod 8$ associated with the diagonal $\hat{\mathbf{x}}-\hat{\mathbf{y}}$ mirror symmetry, as well as a non-vanishing higher-order topological index. It is this latter phase which describes our current model.

We can make these statements more precise by examining the low-energy $\mathbf{k}\cdot\mathbf{p}$ theory for the topological transition in both $H_{FDB}$ and $H_{iFDB}$. Starting with the original FDB model, we find by expanding Eq. (2) that

[TABLE]

On the other hand, performing the same expansion of $H_{\mathrm{iFDB}}$, setting $\delta_{2}=t$, and defining the rotated coordinates $(k_{a},k_{b},k_{c})$ and spin matrices $(s_{a},s_{b},s_{c})$ 333Explicitly, these rotations are given by \left(\begin{array}[]{c}k_{1}\\ k_{2}\\ k_{3}\\ \end{array}\right)=\left(\begin{array}[]{ccc}-\frac{1}{2\sqrt{6}}&-\frac{1}{2\sqrt{2}}&\frac{1}{\sqrt{3}}\\ -\frac{1}{2\sqrt{6}}&\frac{1}{2\sqrt{2}}&\frac{1}{\sqrt{3}}\\ \frac{1}{\sqrt{6}}&0&\frac{1}{\sqrt{3}}\\ \end{array}\right)\left(\begin{array}[]{c}k_{a}\\ k_{b}\\ k_{c}\\ \end{array}\right) and \left(\begin{array}[]{c}s_{x}\\ s_{y}\\ s_{z}\\ \end{array}\right)=\left(\begin{array}[]{ccc}\frac{1}{\sqrt{2}}&\frac{\sqrt{3}+\sqrt{6}}{6}&\frac{\sqrt{2}-1}{2\sqrt{3}}\\ -\frac{1}{\sqrt{2}}&\frac{\sqrt{3}+\sqrt{6}}{6}&\frac{\sqrt{2}-1}{2\sqrt{3}}\\ 0&\frac{1}{\sqrt{6}}-\frac{1}{\sqrt{3}}&\frac{1}{\sqrt{6}}+\frac{1}{\sqrt{3}}\\ \end{array}\right)\left(\begin{array}[]{c}s_{a}\\ s_{b}\\ s_{c}\\ \end{array}\right)

yields to quadratic order

[TABLE]

Up to a choice of basis for the Dirac matrices, this is the BHZ model Hamiltonian for a topological insulatorBernevig et al. (2006); Bernevig and Hughes (2013) – note that because of our choice of boundary conditions and our expansion about the $L$ point, inversion symmetry is represented by $\Delta_{L}(I)=\tau_{z}$ in the $\mathbf{k}\cdot\mathbf{p}$ expansion. Eq. (14) is also equivalent to Eq. (13) if we take $\delta_{2}=2T=t/2$. Since there are $4$ L points, we see that this is a TCI rather than a TI transition. Furthermore, note that the plane $k_{x}=k_{y}$ corresponds to the plane $k_{b}=0$, and is invariant under $m_{1\bar{1}0}$; this symmetry is represented at the $L$ point by the matrix $\Delta_{L}(m_{1\bar{1}0})=\exp(i\pi/(2\sqrt{2})(s_{a}))$. Restricting $H_{\mathrm{iFDB}}$ to the mirror plane, we find that the Hamiltonian is block diagonal in the basis of $m_{1\bar{1}0}$ eigenstates, and describes a Chern insulator transition in each mirror subspace. Since there are two $L$ points in this mirror plane, we thus deduce, following Ref. Hsieh et al., 2012, that this model corresponds to a mirror Chern insulator with $\nu_{m_{1\bar{1}0}}=2$. To verify this, we extract the mirror Chern number from the flow of hybrid Wannier charge centersTaherinejad et al. (2014); Gresch et al. (2017), i.e. from the Wilson loopYu et al. (2011); Alexandradinata et al. (2014). We show in Fig. 1c the $k_{2}$ directed Wilson loop for the occupied bands as a function of $k_{3}$. We see that in the $k_{1}=0$ plane the Wilson loop phases exhibit a nontrivial winding with winding number $2$; since in the FCC Brillouin zone this is the only mirror invariant planeHsieh et al. (2012), this signifies the mirror Chern number $\nu_{m_{1\bar{1}0}}=2$Alexandradinata et al. (2016).

As pointed out in Refs. Song et al., 2018; Schindler et al., 2018a; Khalaf et al., 2018, models in space group $Fm\bar{3}m1^{\prime}$ with $\nu=4$ and $\nu_{m_{1\bar{1}0}}=2$ also have a nontrivial higher order index $\delta_{S_{4}}=1$ protected by fourfold rotoinversion symmetry. This phase has four Dirac cones on a $\mathbf{\hat{z}}$-normal surface, thus explaining the four Dirac fermions on an antiphase boundary found in Refs.Fradkin et al. (1986); Fradkin (1987); Boyanovsky et al. (1987): At the level of the tight-binding model, an antiphase domain wall is simply a boundary between the trivial ($m>0$) and topological ($m<0$) phase of the model. These domain wall fermions are indeed of topological origin, and are symmetry protected in any symmetric extension of the FDB model. We show the spectrum of the topological phase of Eq. (10) on a $\mathbf{\hat{z}}$-normal slab in Fig. 1d

Finally, note that the original low-energy FDB Hamiltonian Eq. (13) at linear order has an accidental mirror symmetry $\tilde{\Delta}_{L}(m_{1\bar{1}0})=\tau_{z}\exp(i\pi/(2\sqrt{2})(s_{x}-s_{y}))$, which leaves the plane $k_{x}=k_{y}$ invariant. While this symmetry is inconsistent with the inversion symmetry of the full model (and so will be broken by crystal-symmetry preserving perturbations, including higher-order terms in the $\mathbf{k}\cdot\mathbf{p}$ expansion), it explains why the authors of Ref. Fradkin et al., 1986; Fradkin, 1987; Boyanovsky et al., 1987 were able to find domain wall fermions in their model.

Ab-initio results

In the previous sections, we have seen how symmetrized completions of the FDB model of PbTe yield higher order topological and mirror Chern insulators. Under certain conditions, the band structure for realistic PbTe as computed with ab-initio methods realizes this same $\nu=4$ symmetry-indicated phase. This can be seen by analyzing the material catalogues of Refs. Vergniory et al., 2019; Zhang et al., 2019, which both report the value of $\nu=4$ for PbTe. For confirmation, we have recomputed the ab-initio band structure of PbTe using Density Functional Theory (DFT)Hohenberg and Kohn (1964); Kohn and Sham (1965) as implemented in the Vienna Ab initio Simulation Package (VASP)Kresse and Furthmueller (1996); Kresse and Hafner (1993). We use the structural parameters as reported in Ref. Noda et al., . The interaction between ion cores and valence electrons was treated by the projector augmented-wave methodKresse and Joubert (1999), the generalized gradient approximation (GGA) for the exchange-correlation potential with the Perdew-Burke-Ernkzerhof for solids parameterization Perdew et al. (1996) and spin-orbit coupling was taken into account by the second variation methodHobbs et al. (2000). A Monkhorst-Pack centered at $\Gamma$ k-point grid of (11$\times$11$\times$11) for reciprocal space integration and 500 eV energy cutoff of the plane-wave expansion have been used. We show the band structure in Fig. 2, with an inset highlighting the rather small gap at $L$. Employing the VASPtoTrace toolBilbao Crystallogr. Server (2018); Vergniory et al. (2019), we compute the little group irreps of the occupied bands at the high symmetry points, shown in Table 2; we give the irreps of SnTe as well for comparison. By using Eq. (11), we see that $\nu=4$ for both SnTe and PbTe. Note, in fact, that the irrep labels for SnTe and PbTe differ only in a shift of the origin of the system by $(\frac{1}{2}\frac{1}{2}\frac{1}{2})$ Furthermore, the topological transition to $\nu=4$ in the realistic material is driven by a band inversion of the irreps at $L$ relative to $W$, just as in the FDB model. To fully determine the topological phase, we evaluate the mirror Chern number of the occupied bands using Z2PackGresch et al. (2017). We find a mirror Chern number $\nu_{m_{110}}=\nu_{m_{1\bar{1}0}}=2$, just as in the iFDB model444These two mirror chern numbers are equal due to the cubic symmetry of the system. (See Fig. 3).

However, it is well-accepted that the mirror Chern number $\nu_{m_{1\bar{1}0}}$ in PbTe is zero under ambient experimental conditions, while it is $2$ for SnTeHsieh et al. (2012). To reconcile this with the nontrivial $\nu=4$ topological index, we note that in addition to the structure used for the ab initio calculations here and in Refs. Vergniory et al., 2019; Zhang et al., 2019, PbTe has 41 other entries in the ICSDFachinformationszentrum Karlsruhe (2016) in the space group $Fm\bar{3}m1^{\prime}$Vergniory et al. (2019). A DFT analysis of other structures (for instance, the structure reported in Ref. Sternberg et al., 1982) yields a trivial index $\nu=0$ due to a band de-inversion at $L$, in agreement with the experimental findings. This highlights the fact that for small band gap insulators, one must be cautious in extracting the band topology from ab initio calculations; for PbTe in particular, the failure of semilocal DFT to correctly produce the (correct sign of the) experimental band gap in certain cases has been noted previouslyHummer et al. (2007).

To investigate this systematically, we have computed the band structures and topological index $\nu$ for all $42$ entries of PbTe in the ICSD, using the PBE functional. The input parameters for these compounds differ only in the reported lattice constant $a_{0}$, which range between $6.157$Å and $6.543$Å. For the six reported structures with $a_{0}\leq 6.44$Å, PBE predicts $\nu=4$; for the remaining with larger lattice constants we find $\nu=0$. In the Supplementary Material we give a table summarizing our DFT calculations555See Supplemental Material at https://journals.aps.org/prmaterials/supplemental/10.1103/PhysRevMaterials.3.041202/suppmat.pdf for the topological index for all reported entries of PbTe.. Taken optimistically, This shows that PbTe is very close to a topological phase transition, which may be tunable as a function of external parameters such as hydrostatic pressure.

Conclusion

We have revisited the effective model of PbTe as presented in Refs. Fradkin et al., 1986; Fradkin, 1987; Boyanovsky et al., 1987. We have shown that the domain wall fermions in the FDB model, long derided as nontopological, are signatures of the topological surface states present in any symmetric completion of the model, protected by mirror and fourfold rotoinversion symmetries. Furthermore, we show that ab initio calculations reveal that some of the reported structures of realistic PbTe are in this same symmetry-indicated class of materials, at least within the GGA. This shows that PbTe is an ideal platform for exploring structurally-tunable topological behavior. Finally, while within the context of our effective model there is no difference between an antiphase domain wall and a domain wall with the vacuum, this is not necessarily true in a more realistic system. Given the recent focus on defect response of higher-order topological insulatorsSlager et al. (2015); Queiroz et al. (2018); Wieder et al. (2018b); Benalcazar et al. (2018), it would be interesting to examine this more carefully for both SnTe and PbTe structural variations in future work.Schindler et al. (2018b); Robredo et al. (2019)

Acknowledgements.

Acknowledgements:

BB and MGV thank E. Fradkin for helpful discussions, and for pointing out Refs. Fradkin (1987); Boyanovsky et al. (1987). BB additionally acknowledges discussions with B. J. Wieder. This work was initiated at the Aspen Center for Physics, which is supported by National Science Foundation grant PHY-1607611. MGV acknowledges the IS2016-75862-P national project of the Spanish MINECO.