# Topology of chalcogen chains

Adam Kłosiński<EMAIL_ADDRESS>Institute of Theoretical Physics,
Faculty of Physics, University of Warsaw, Pasteura 5, PL-02093 Warsaw, Poland
Wojciech Brzezicki  International Research Centre MagTop, Institute of Physics
PAS, Aleja Lotników 32/46, PL-02668 Warsaw, Poland Alexander Lau
International Research Centre MagTop, Institute of Physics PAS, Aleja Lotników
32/46, PL-02668 Warsaw, Poland Cliò E. Agrapidis Institute of Theoretical
Physics, Faculty of Physics, University of Warsaw, Pasteura 5, PL-02093
Warsaw, Poland Andrzej M. Oleś  Institute of Theoretical Physics, Jagiellonian
University, Prof. Stanisława Łojasiewicza 11, PL-30348 Kraków, Poland Max
Planck Institute for Solid State Research, Heisenbergstrasse 1, D-70569
Stuttgart, Germany  Jasper van Wezel  Institute for Theoretical Physics
Amsterdam and Delta Institute for Theoretical Physics,
University of Amsterdam, Science Park 904, NL-1098 XH Amsterdam, The
Netherlands  Krzysztof Wohlfeld  Institute of Theoretical Physics, Faculty of
Physics, University of Warsaw, Pasteura 5, PL-02093 Warsaw, Poland

(December 7, 2022)

###### Abstract

We investigate the topological properties of the helical atomic chains
occurring in elemental selenium and tellurium. We postulate a realistic model
that includes spin-orbit interaction and show this to be topologically non-
trivial, with a topological invariant protected by a crystalline symmetry. We
describe the end-states, which are orbitally polarized, with an orbital
density modulation strongly peaked at the edge. Furthermore, we propose a
simplified model that decomposes into three orbital chains, allowing us to
define a topological invariant protected by a crystalline symmetry. We
contrast this result with recent observations made for the orbital Su-
Schrieffer-Heeger model containing a $p$-orbital zigzag chain.

##  I Introduction

Topology has become one of the pillars of modern condensed matter research
over the past decades [1, 2, 3, 4, 5]. Despite the many possible classes of
non-trivial topology identified theoretically [6, 7, 8, 9, 10], and the many
realisations of non-trivial quantum materials identified experimentally [11,
12, 13, 14, 15, 16, 17], suggested direct implementations of the most basic,
prototypical models for topological order are rare. Here, we propose that
elemental chalcogens selenium and tellurium naturally implement an extension
of one of the most fundamental and most well-known topological models in
condensed matter, the Su-Schrieffer-Heeger (SSH) model [18].

The particular model realised in the spiral chains of chalcogen crystals
consists of three independent copies of a period-three SSH model (SSH-3) [19].
This model differs from the standard SSH model and the recently introduced
orbital SSH model [20, 21], by having a three-site rather than two-site
periodicity and lacking exact particle-hole or chiral symmetry [22, 23, 24].
It does, however, have a two-fold rotational symmetry that combines with time-
reversal to guarantee robust end states in the gap between conduction and
valence band. In the chalcogens, the particular charge and orbitally ordered
ground state furthermore causes the topological end states to carry a non-zero
orbital polarization.

The elemental chalcogens have been studied for almost a century [25, 26, 27]
and have been long known to organise into a chiral crystal structure
consisting of weakly coupled spiral chains [25, 26, 24]. Only recently has the
origin of these spiral chains been understood in terms of a combined charge
and orbital instability of a simple-cubic parent structure [22, 28, 29]. In
this instability, a charge density wave develops among quasi-one-dimensional
chains of each of the three possible $p$-orbital orientations. Weak Coulomb
interaction between the chains leads to a relative alignment of the charge
density waves into an emergent three-dimensional structure consisting of
spiral chains. Owing to the orbital origin of the charge density distortions,
the spiral structure is also orbitally ordered [28, 29, 30]. Although the form
of the orbital order depends on the strength of the Coulomb interaction, the
emerging spiral lattice distortions and their origin in combined charge and
orbital order was recently shown to be robust [30].

The topology of the electronic state in the spiral configuration of elemental
chalcogens has not been extensively studied, although there are reports from
electronic structure calculations that it undergoes a topological phase
transition from an insulator to a Weyl semimetal under strain [31, 32].
Additionally, there are recent suggestions that surface states may be
identified on crystals of elemental selenium and tellurium [33], but their
topological origin has not yet been conclusively established. The particular
charge order in the spiral chains, however, is a direct three-site
generalisation of the famous two-site SSH model, while the orbital order is
strongly reminiscent of the so-called orbital SSH model [20, 21]. Here, we
show that the elemental chalcogens are indeed topologically ordered,
exhibiting a direct realisation of some of the most fundamental known
topological models, and giving rise to orbitally polarised end states. The
topological invariant characterising the order and end states coincides with
the one-dimensional line invariant appropriate for time reversal symmetric
systems with a two-fold crystal symmetry [34].

This paper is organized as follows: In Sec. II we introduce the model for a
helical chalcogen chain. We begin by describing the chain geometry (Sec.
II.1). Then we define the model in two different bases, making use of the
helical symmetry (Sec. II.2 and II.3), as well as a simplified version of the
model, for which the bond angle $\alpha=90^{\circ}$ and which we call the
cubic model (Sec. II.4). In the final subsection of Sec. II we discuss the
symmetries of the chalcogen chain (Sec. II.5). In Sec. III we show the
spectrum of the open chalcogen chain obtained using Exact Diagonalisation (ED)
(Sec. III.1), where we observe in-gap states. We also give a description of
the bulk topological invariant (Sec. III.2) and discuss the bulk-boundary
correspondence (Sec. III.3). In the final subsection of this section we
discuss the end states of the chalcogen chain (Sec. III.4). In Sec. IV we
provide a detailed discussion of the topology of the chalcogen chain, going
beyond the invariant defined in Sec. III.2. To that end, in Sec. IV.1 we
analyze the cubic model – a simplified version of the chalcogen chain model,
with the bond angle $\alpha=90^{\circ}$. We show that this simpler model
factorizes into three SSH-3 chains, one of which supports end-states. These
end states are continuously linked to the end-states of the full chalcogen
model, allowing us to relate the topology of the chalcogen chain to that of
the SSH-3 model. In the final subsection of Sec. IV we offer a comparison
between the chalcogen chain model and the orbital-SSH model described in [20,
21]. Finally, in the Appendices we provide a detailed derivation of the end
states in the continuous limit (Appendix A) and a more detailed analysis of
the evolution of the energy spectrum of the chalcogen chain with varying spin-
orbit coupling strength (Appendix B).

##  II Model

###  II.1 Chain geometry and helical symmetry

We previously introduced an electronic model for the helical chains in
elemental chalcogens [30], which contains three $p$-orbitals
$\\{p_{x},p_{y},p_{z}\\}$ – but no spin. Here, we expand this model to include
spin and spin-orbit coupling. Following the conclusions of Ref. [30], we
assume the effects of the weak Coulomb interaction to be negligible. The non-
interacting, electronic model for a single helix then contains three distinct
types of bonds, but it can be mapped onto a model with only a single bond type
using the helical symmetry.

Trigonal selenium and tellurium crystals consist of parallel, weakly-coupled
helical chains arranged on a hexagonal lattice [26, 27, 35]. When viewed from
the top, atoms in each helix collapse onto the vertices of an equilateral
triangle, which demonstrates that the helices are period-three, as shown in
Fig. 1. The period-three helix is characterised by a single bond angle
$\alpha$ (shown in Fig. 2), which can be treated as a free parameter. For
selenium $\alpha_{\text{Se}}\in[102.5^{\circ},105.5^{\circ}]$, while for
tellurium $\alpha_{\text{Te}}\in[102.4^{\circ},103^{\circ}]$. Here we will use
$\alpha=103^{\circ}$, which we call the chalcogen model, though we will also
consider the model for $\alpha=90^{\circ}$, which we call the cubic model.

Figure 1: The geometry of the helical chain: (a) top view, (b) side view. Also
marked on the figure: the helical axis ($\vec{a}$) and the axis normal to the
helical axis ($\vec{r}$).

(a) (b)

Figure 2: The helical chain in three dimensions. The global (a) and the local
(b) orbital bases are shown, with $p_{x}$, $p_{y}$, and $p_{z}$ orbitals
colored red, green and blue, respectively. A shaded gray cylinder is also
marked, whose axis coincides with the helical axis $\vec{a}$ of the chain,
shown in Fig. 1. In both cases the bond angle is $\alpha=103^{\circ}$.

These models do not have translational invariance $T$ by one bond distance
along the helix. However, the helical chain is invariant with respect to a
special translation operator $T\otimes O\otimes U$. This symmetry is a
combination of: the translation $T$ by one bond distance along the helix times
a 120-degree orbital rotation $O$ times a 480-degree spin rotation $U$ around
the axis of the helix $\vec{a}$ shown in Fig. 1 111After a 360-degree rotation
– which, for a 120-degree rotation, corresponds to the special translation by
a full unit cell – the spin basis acquires an extra minus sign. Therefore we
choose $U$ to be a 480-degree spin rotation. With this choice, we have, on the
one hand, $120=480\;\text{mod}\;360$, and on the other $U^{3}=\mathbb{1}$. The
first property is essential in order for the spin-orbit interaction to be
site-independent.. The 120-degree orbital rotation symmetry $O$ is best
visible in Figs. 2(a) and 2(b), where two orbital bases are shown – the basis
before the orbital rotation [see panel (a)] and the basis after the orbital
rotation [see panel (b)].

###  II.2 Period-three model

We choose a coordinate system with the $\hat{z}$ axis along the helical axis
$\vec{a}$ and the $\hat{x}$ axis in the direction of the shortest line
connecting the first site of the helix with its axis. The first site is chosen
arbitrarily. We call this basis the global orbital (GO) basis. The GO basis is
shown in Fig. 2(a). Since the helix is period-three, the unit cell consists of
atoms at positions $\left(r,0,0\right)$,
$\left(-\frac{r}{2},\frac{\sqrt{3}r}{2},a\right)$,
$\left(-\frac{r}{2},-\frac{\sqrt{3}r}{2},2a\right)$, and the translation
vector is $\left(0,0,3a\right)$. On each chain site we have three $p$ orbitals
and spin-orbit coupling, therefore the tight-binding model is given as
follows:

$H=\sum_{i=1}^{3}T_{i}\otimes
h\left(\vec{n}_{i}\right)\otimes\mathbbm{1}_{2}-\lambda\sum_{\alpha=x,y,z}\mathbbm{1}_{N}\otimes
L_{\alpha}\otimes\sigma_{\alpha},$ (1)

where $\sigma_{\alpha}$ are the Pauli matrices describing spin and
$\left(L_{\alpha}\right)_{\beta\gamma}=-i\varepsilon_{\alpha\beta\gamma}$ are
the angular momentum $L=1$ matrices. $\lambda$ is the spin-orbit coupling
parameter. Note that the spin-orbit term is negative, because we consider 2/3
filling. The $T_{i}$ operators are $N\times N$ matrices which describe
nearest-neighbor hopping between sites of the helix of length $N$ (we assume
that $N$ is divisible by three):

$T_{i}=\sum_{j=1}^{N/3-1}\left(\,\left|i+3j\right\rangle\left\langle
i+1+3j\right|+H.c.\,\right),$ (2)

with periodic boundary conditions
$\left|N+1\right\rangle\equiv\left|1\right\rangle$. Finally, $h(\vec{n}_{i})$
are the $3\times 3$ matrices describing hopping between $p_{x}$, $p_{y}$ and
$p_{z}$ orbitals along the bonds. They are given by Slater-Koster rules:

$\displaystyle h\left(\vec{n}\right)$ $\displaystyle=$
$\displaystyle\begin{pmatrix}n_{x}^{2}\delta t+t_{\pi}&-n_{x}n_{y}\delta
t&-n_{x}n_{z}\delta t\\\ -n_{y}n_{x}\delta t&n_{y}^{2}\delta
t+t_{\pi}&-n_{y}n_{z}\delta t\\\ -n_{z}n_{x}\delta t&-n_{z}n_{y}\delta
t&n_{z}^{2}\delta t+t_{\pi}\end{pmatrix}$ (3) $\displaystyle=$
$\displaystyle(t_{\pi}-t_{\sigma})(\vec{n}\cdot\vec{L})^{2}+t_{\sigma}\mathbbm{1}_{3},$

where $\delta t=t_{\sigma}-t_{\pi}$ and $t_{\pi}$($t_{\sigma}$) are the
bonding amplitudes of the $\pi$ ($\sigma$) bonds [37]. Since the helix is
period-three, we have three non-equivalent bond directions denoted as
normalized vectors $\vec{n}_{i}$. They can be derived from the atomic
positions as functions of the bond angle $\alpha$. They are given by:

$\displaystyle\vec{n}_{1}$ $\displaystyle=$
$\displaystyle\frac{1}{g\left(\alpha\right)}\left(-\frac{3}{2},\frac{\sqrt{3}}{2},\sqrt{3}f\left(\alpha\right)\right),$
(4) $\displaystyle\vec{n}_{2}$ $\displaystyle=$
$\displaystyle\frac{1}{g\left(\alpha\right)}\left(0,-\sqrt{3},\sqrt{3}f\left(\alpha\right)\right),$
(5) $\displaystyle\vec{n}_{3}$ $\displaystyle=$
$\displaystyle\frac{1}{g\left(\alpha\right)}\left(\frac{3}{2},\frac{\sqrt{3}}{2},\sqrt{3}f\left(\alpha\right)\right),$
(6)

with
$f\left(\alpha\right)=\sqrt{\frac{3}{4}\left(\cos\frac{\alpha}{2}\right)^{-2}-1}$,
and $g\left(\alpha\right)=\frac{3}{2}\left(\cos\frac{\alpha}{2}\right)^{-1}$.
Alternatively, they can be expressed by the ratio of the helix pitch to the
distance of the chain atoms from the helical axis $r$, namely $q=a/r$. The
value of $q$ is related to the bond angle $\alpha$ by

$\cos\alpha=\frac{3-2q^{2}}{6+2q^{2}}\,.$ (7)

Note that $\frac{\pi}{3}<\alpha<\pi$ while $0<q<\infty$.

(a) (b)

Figure 3: Band structure of the chalcogen model described by Eq. (1) (with
$\alpha=103^{\circ}$, $t_{\pi}/t_{\sigma}$ = -1/3) with periodic boundary
conditions: (a) without spin-orbit coupling, $\lambda/t_{\sigma}=0$, and (b)
with spin-orbit coupling $\lambda/t_{\sigma}=0.2$. Occupied bands are colored
blue, unoccupied bands are colored gray.

The band structure for the model of Eq. (1) is shown in Fig. 3. It exhibits
band gaps at $1/3$ and $2/3$ filling, which are substantial for realistic
values of the spin-orbit coupling parameter $\lambda$, as discussed in Sec.
III.1 below.

###  II.3  Period-one model

Note that the model of Eq. (1) has a symmetry whenever $N$ is divisible by
three described by the operator $Q$,

$Q=T\otimes O\otimes
U=T\otimes\exp\left[i\frac{2\pi}{3}L_{z}\right]\otimes\exp\left[i\frac{8\pi}{3}\frac{\sigma_{z}}{2}\right].$
(8)

This represents a cyclic shift by one site along the helix combined with a
120-degree rotation of orbitals and a 480-degree rotation of spins. $O$ and
$U$ are the orbital and spin rotations discussed in Sec. II.1 and $T$ is:

$T=\sum_{r=1}^{N-1}\left(\,\left|r\right\rangle\left\langle
r+1\right|+\left|N\right\rangle\left\langle 1\right|\,\right).$ (9)

Here, $r$ labels the sites along the chain.

From $\left[H,Q\right]=0$ and the choice of $N\in 3\mathbb{N}$ we see that in
the eigenbasis of $Q$ the Hamiltonian $H$ is decomposed into $N$ $6\times$6
diagonal blocks. The easiest way to perform such decomposition, however, is by
defining a unitary transformation,

$V=\sum_{r=0}^{N-1}\left|r+1\right\rangle\left\langle
r+1\right|\otimes\exp\\!\left[i\frac{2\pi
r}{3}L_{z}\right]\otimes\exp\\!\left[i\frac{8\pi
r}{3}\frac{\sigma_{z}}{2}\right],$ (10)

that acts locally in each unit cell. Transforming the Hamiltonian to the new
basis, which we will call the local orbital (LO) basis, we get

$\displaystyle H^{\prime}$ $\displaystyle\\!=\\!$ $\displaystyle\left(T\otimes
h\left(\vec{n}_{1}\right)\exp\\!\left[-i\frac{2\pi}{3}L_{z}\right]\\!\otimes\\!\exp\\!\left[-i\frac{8\pi}{3}\frac{\sigma_{z}}{2}\right]\\!+\\!H.c.\right)$
(11) $\displaystyle-$
$\displaystyle\lambda\sum_{\alpha=x,y,z}\mathbbm{1}_{N}\otimes
L_{\alpha}\otimes\sigma_{\alpha},$

where $H^{\prime}=VHV^{\dagger}$, with $H$ defined in Eq. (1). In the LO
basis, the Hamiltonian has a single-site unit cell. Now we can easily go to
$k-$space by inserting the eigenvalues of $T$ to get:

$\displaystyle H^{\prime}(k)$ $\displaystyle\\!=\\!$
$\displaystyle\left(\\!e^{ik}h\\!\left(\vec{n}_{1}\right)\exp\\!\left[-i\frac{2\pi}{3}L_{z}\right]\\!\otimes\exp\\!\left[-i\frac{8\pi}{3}\frac{\sigma_{z}}{2}\right]\\!+H.c.\\!\right)$
(12) $\displaystyle-$
$\displaystyle\lambda\sum_{\alpha=x,y,z}L_{\alpha}\otimes\sigma_{\alpha}.$

Note that this Hamiltonian can also be put in a form that does not depend on
the choice of the coordinate system by substituting

$\displaystyle L_{z}$ $\displaystyle\to$ $\displaystyle\vec{m}\cdot\vec{L},$
$\displaystyle\sigma_{z}$ $\displaystyle\to$
$\displaystyle\vec{m}\cdot\vec{\sigma}$ (13)

in $H^{\prime}$. In the above, $\vec{m}$ is a normalized vector along the axis
of the helix and $h\\!\left(\vec{n}_{1}\right)$ is already given in
coordinate-independent form in Eq. (3).

###  II.4 The cubic model ($\alpha=90^{\circ}$)

Figure 4: The helical chain for $\alpha=90^{\circ}$. At each site of the chain
the GCO basis is shown (see Sec. II.4), with $p_{x}$, $p_{y}$, and $p_{z}$
orbitals colored red, green and blue respectively. A shaded gray cylinder is
also marked, whose axis coincides with the helical axis $\vec{a}$ of the
chain, shown in Fig. 1. The cubic environment is made apparent.

The model becomes especially simple for $\alpha=90^{\circ}$ (the cubic model
case), when a different basis choice is convenient. In the period-1 case we
use the coordinate independent form of $H^{\prime}(k)$ and set the direction
of the helix as $\vec{m}=(1,1,1)/\sqrt{3}$ and of the first bond as
$\vec{n}_{1}\\!=\\!(0,0,1)$. This is a valid basis choice for any $\alpha$,
and we call this basis choice the local cubic orbital (LCO) basis, as it makes
the cubic model Hamiltonian especially simple. In this basis, the orbital
hopping part of the Hamiltonian becomes:

$\displaystyle
h\left(\vec{n}_{1}\right)\exp\left[-i\frac{2\pi}{3}\vec{m}\cdot\vec{L}\right]=\begin{pmatrix}0&t_{\pi}&0\\\
0&0&t_{\pi}\\\ t_{\sigma}&0&0\end{pmatrix}.$ (14)

Transforming back to the period-3 model, for $\alpha=90^{\circ}$ we get Eqs.
(1)-(3) with

$\displaystyle\vec{n}_{1}$ $\displaystyle=$ $\displaystyle(0,0,1),$ (15)
$\displaystyle\vec{n}_{2}$ $\displaystyle=$ $\displaystyle(1,0,0),$ (16)
$\displaystyle\vec{n}_{3}$ $\displaystyle=$ $\displaystyle(0,1,0).$ (17)

We call the LCO basis transformed back to the period-three model the global
cubic orbital (GCO) basis. It is shown in Fig. 4 for the $\alpha=90^{\circ}$
case.

###  II.5 Symmetries of the model

Besides the $T\otimes O\otimes U$ helical symmetry the Hamiltonian has two
other important symmetries:

(i) Time-reversal symmetry $\Theta=i{\cal K}\sigma_{y}$ ($\cal K$ is complex
conjugation), which leads to Kramers degeneracies at time-reversal invariant
momenta since $\Theta^{2}=-1$, and

$\Theta H^{\prime}(k)\Theta^{-1}=H^{\prime}(-k)^{\star}.$ (18)

(ii) $\mathcal{R}$ – a 180-degree rotation symmetry around the axis $\vec{r}$
normal to the helical axis and going through one of the atoms (see Fig. 1). In
the LO basis the form of $\mathcal{R}$ is

$\mathcal{R}=\exp\left[i\pi
L_{x}\right]\otimes\exp\left[i\pi\frac{\sigma_{x}}{2}\right],$ (19)

and it acts on the Hamiltonian as

$\mathcal{R}H^{\prime}(k)\mathcal{R}^{-1}=H^{\prime}(-k).$ (20)

As a further consequence of these two symmetries, there exist an anti-unitary
operator $\mathcal{V}=\Theta\mathcal{R}$ such that

$\mathcal{V}H^{\prime}(k)\mathcal{V}^{-1}=H^{\prime}(k)^{\star},$ (21)

and $\mathcal{V}^{2}=\mathbb{1}$. Consequently, there exists a basis in which
the Hamiltonian is purely real.

##  III Results

###  III.1 Exact diagonalization of an open chain

We begin with the exact diagonalisation (ED) of the chalcogen model
Hamiltonian in Eq. (1) with open boundary conditions. We use the chain length
$N=501$. The spectrum is presented in Fig. 5. We find that in a finite
chalcogen chain with weak spin-orbit interaction four states separated from
the continuum appear in each of the two gaps. Their presence suggests a
nontrivial topology of the chain – this will be investigated in the following
subsections. The results of Fig. 5 also show that as the spin-orbit coupling
$\lambda$ increases, the upper gap – where the Fermi level is located – closes
for $\lambda_{\text{crit.}}/t_{\sigma}\approx 0.57$ (the value of
$\lambda/t_{\sigma}$ for selenium and tellurium is $0.15-0.33$). The in-gap
states disappear after the gap closing, indicating a topological phase
transition.

The in-gap states in the lower gap disappear into the bulk at around
$\lambda/t_{\sigma}=0.8$, while the gap remains open even for very large
$\lambda/t_{\sigma}$. We will disregard the lower gap in our analysis as it
lies deep below the Fermi energy. A more detailed discussion of both gaps is
presented in Appendix B.

Figure 5:  The spectrum of the chalcogen model in Eq. (1)
($\alpha=103^{\circ}$, $t_{\pi}/t_{\sigma}=-1/3$) on a finite open chain of
length $N=501$, for different values of $\lambda/t_{\sigma}$. The end states
are colored red, the occupied bulk states are colored blue and the unoccupied
bulk states are colored gray.

###  III.2 Bulk topological properties

In this subsection we examine the bulk properties and the topological
invariant. The system has a spatial symmetry $\cal R$ that can be regarded as
a spinful inversion. Therefore, inspired by Ref. [38], we search for the
$\mathbb{Z}^{\geq}$ invariant by looking at the spectrum of the Wilson loop
operator ${\cal W}$. Its matrix elements can be calculated using a gauge-
invariant formula,

${\cal W}_{nm}=\left\langle E_{n}(k_{1})\right|{\cal
P}_{n_{f}}(k_{2})\dots\cal P_{n_{f}}(k_{N})\left|E_{m}(k_{1})\right\rangle,$
(22)

where $k_{j}=2\pi j/N$, while $\left|E_{m}(k_{j})\right\rangle$ are the
eigenstates of $H^{\prime}(k)$ in Eq. (12), and $N$ should be taken large
enough to assure convergence. Notice that only the filled bands are considered
in the band labels $m$ and $n$. The operators ${\cal P}_{n_{f}}(k)$ are
projectors on the occupied subspace for a given $k$:

${\cal P}_{n_{f}}(k)=\sum_{n=1}^{n_{f}}\left|E_{n}(k)\right\rangle\left\langle
E_{n}(k)\right|.$ (23)

Here, $n_{f}$ is the index of the occupied state corresponding to the largest
eigenvalue (energy).

Due to the symmetry $\cal R$, the spectrum of the $\cal W$ operator consists
of the values: $+1$, $-1$ or pairs of complex conjugate numbers of unit
modulus [38]. The number of $-1$ appearances we denote as $N_{(-1)}$. Note
that if $\cal R$ was spinless, as in Ref. [38], then $N_{(-1)}$ would be equal
to the difference in number of occupied states for which ${\cal R}=-1$ between
$k=0$ and $k=\pi$ and it would be a topological invariant. In the present
case, however, such difference gives zero because in each Kramers doublet we
always find two opposite ${\cal R}$ eigenvalues.

From the general classification involving mirror symmetries [8], and more
specific work of Ref. [34], it follows that we can have only one non-trivial
phase. Therefore we define a $\mathbb{Z}_{2}$ topological invariant as

$\nu=\frac{1}{2}N_{(-1)}\mod 2.$ (24)

This expression is equivalent to the invariant of Ref. [34], and is valid as
long as $N_{(-1)}$ is even. This is always satisfied because of Kramers
degeneracy, which implies that end states in our model come in pairs, and
therefore $N_{(-1)}$ has to be even. For the chalcogen model we find that the
invariant takes values $\nu=1$ and $\nu=0$ at $2/3$-filling for
$\lambda<\lambda_{\text{crit.}}$ and $\lambda>\lambda_{\text{crit.}}$
respectively. The critical value $\lambda_{\text{crit.}}$ coincides with the
gap closing in Fig. 5.

Finally, we would like to point out that in one dimension, the bulk-boundary
correspondence appears robustly only when there is a symmetry that binds the
end-states to zero energy, so either chirality or particle-hole symmetry. In
any other case, end-states that appear in the gap can be shifted up or down in
energy by a local potential and forced to join the bulk continuum. There are,
however, results indicating that a non-trivial Zak phase of the bulk bands
leads to charge accumulation at the ends of the system which sometimes
manifests itself as the presence of end states [39, 40]. In our case, the
appearance of the end-states in Fig. 5 is fully consistent with the invariant
$\nu$ in Eq. (24). However, in order to better understand what conditions lead
to their presence, in the next section we will discuss the bulk-boundary
correspondence in our model.

###  III.3 Bulk-boundary correspondence

In order to establish the robustness of the end states, we performed an
expansion of the model defined by Eq. (12) around the point of the upper gap
closing. There are two values of pseudo-momentum $k$ for which the gap closes
(see Fig. 3), symmetric with respect to the point $k=0$. For now we will focus
on the gap closing at pseudo-momentum $k_{0}>0$.

In order to expand the Hamiltonian around the gap closing point, we first
apply a basis change $u_{k_{0},\lambda_{\text{crit.}}}$ to the eigenbasis of
the chalcogen model Hamiltonian in Eq. (12) at $k=k_{0}$ and
$\lambda=\lambda_{\text{crit.}}$. Next, we apply a projection $P$ of the
resulting Hamiltonian onto the two eigenstates forming the Dirac cone at
$k_{0}$ to obtain an effective $2\times 2$ Hamiltonian

$h(k,\lambda)=Pu_{k_{0},\lambda_{\text{crit.}}}H^{\prime}(k)\,u_{k_{0},\lambda_{\text{crit.}}}^{\dagger}P.$
(25)

By doing this we find that the eigenbasis of the model at the gap-closing
point $(k_{0},\lambda_{\text{crit.}})$ is exactly the basis in which the
entire Hamiltonian
$u_{k_{0},\lambda_{\text{crit.}}}H^{\prime}(k)u_{k_{0},\lambda_{\text{crit.}}}^{\dagger}$
is real and whose existence is guaranteed by symmetry (see Sec. II.5).
Consequently, $h(k,\lambda)$ is also real and we can write it as

$h(k,\lambda)=s(k-k_{0},\lambda)\sigma_{0}+f(k-k_{0},\lambda)\sigma_{x}+g(k-k_{0},\lambda)\sigma_{z},$
(26)

where $s(0,\lambda_{\text{crit.}})=E_{0}$ (the gap closing energy) and
$f(0,\lambda_{\text{crit.}})=g(0,\lambda_{\text{crit.}})=0$. After subtracting
the scalar term, the effective $2\times 2$ Hamiltonian

$\tilde{h}(k,\lambda)=h(k,\lambda)-s(k-k_{0},\lambda)\sigma_{0},$ (27)

becomes chiral.

The standard bulk-boundary correspondence guarantees that the chiral
Hamiltonian $\tilde{h}(k,\lambda)$ leads to zero-energy end states upon
interfacing with a trivial phase. For $h(k,\lambda)$, however, this is not
necessarily the case. This is due to the fact that the scalar term, while
irrelevant for the topology of the bulk, does not commute with
$\tilde{h}(k,\lambda)$ on an open chain and can therefore affect the end
states. Nevertheless, a detailed calculation shows that in this case the end
states survive in the gap (see Appendix A for details).

###  III.4 Orbital polarization of the end states

Figure 6: The end-state orbital density for a chain cut across two $y$ bonds:
(a) for the chalcogen model of Eq. (1) ($\alpha=103^{\circ}$,
$t_{\sigma}/t_{\pi}=-1/3$) with
$\lambda/t_{\sigma}=0.2<\lambda_{\text{crit.}}/t_{\sigma}$, (b) for the cubic
model ($\alpha=90^{\circ}$, $t_{\sigma}/t_{\pi}=-1/3$) with
$\lambda/t_{\sigma}=0$, see Eq. (28) in Sec. IV.1. In both panels we use the
GCO basis defined in Sec. II.4 (see Fig. 4). Also shown in panel (b): a
schematic picture of the end-state charge density along the helix for the
cubic model. The calculation was performed for a chain of length $N=213$
sites.

We have analyzed the character of the end states by looking at their
eigenvectors obtained using ED on an open chain (see Fig. 6) of size $N=213$.
The most transparent representation is obtained in the GCO basis. Looking at
Fig. 6(a), we see that the charge density of the end states decays
exponentially away from the edge and that the end states have a strong orbital
polarisation. The latter is dictated by the bulk bond pattern and carries a
predominantly $p_{y}$ orbital character. Since in our calculation the break in
the chain occurs where a $p_{y}$-type dimer would occur in an infinite chain.

One can identify this result as one closely related to the result presented in
Fig. 6(b). It shows the end states of the cubic model ($\alpha=90^{\circ}$) in
the absence of spin-orbit coupling ($\lambda=0$). In this case the end states
have a purely $p_{y}$ orbital character and the exponentially decaying end
state is modulated by a period-three, 1-1-0 charge density wave (CDW). The
physical understanding of this profound result is studied in detail in Sec.
IV.1.

We conclude that the spin-orbit interaction and the increase in bond angle
change the nature of the end states only slightly, introducing a small
admixture of the other two $p$-orbitals and weakly breaking the 1-1-0 CDW
modulation in the $p_{y}$ orbital density. However, the density of the $p_{y}$
orbital in the end state is still almost zero on every third site. The non-
zero $p_{x}$ and $p_{z}$ orbital density is due to the fact that both relevant
terms – the finite spin-orbit coupling $\lambda$ and the deviation from the
90-degree bond angle $\alpha$ – introduce coupling between the $p$ orbitals,
which is absent in the model without these terms.

##  IV Discussion

###  IV.1 Cubic model ($\alpha=90^{\circ}$) with $\lambda=0$

In what follows we discuss the physical origin of the topological invariant
and the orbital polarization of the end states. To this end we turn our
attention to the cubic model, which is considerably simpler than the chalcogen
model. We uncover a mapping of the cubic model for $\lambda=0$ onto a SSH-3
model with a strong-weak-weak band pattern. In the following subsections we
discuss the implications of this finding by showing: (i) its relevance for the
chalcogen model, and (ii) its relation to the (orbital) SSH model.

We consider a simplified cubic ($\alpha=90^{\circ}$) model:

$H_{0}\equiv H\Big{|}_{\lambda=0,\;\alpha=90^{\circ}},$ (28)

where $H$ is defined in Eq. (1). Here it is useful to use the GCO basis, in
which the Hamiltonian takes a particularly simple form (see Sec. II.4 for
details). In this basis the $x$, $y$ and $z$ axes point along the bonds of a
simple cubic lattice and the helix bonds correspond to consecutive bonds on a
simple cubic lattice (see Fig. 4).

The Hamiltonian of Eq. (28) separates into two independent spin channels,

$H_{0}=\bigoplus_{\sigma}H_{0}^{\sigma}.$ (29)

Furthermore, each spin block separates into three independent orbital
channels,

$H_{0}^{\sigma}=\bigoplus_{\gamma=1}^{3}H_{0}^{\sigma,\gamma},$ (30)

each described by a SSH-3 chain. The chains corresponding to different
orbitals are shifted by one lattice distance from one another. Their $k$-space
Hamiltonians (in the GCO basis) are given by,

$H_{0}^{\sigma,1}(k)=\begin{pmatrix}0&t_{\pi}&t_{\sigma}e^{-ik}\\\
t_{\pi}&0&t_{\pi}\\\ t_{\sigma}e^{ik}&t_{\pi}&0\end{pmatrix},$ (31)
$H_{0}^{\sigma,2}(k)=\begin{pmatrix}0&t_{\sigma}&t_{\pi}e^{-ik}\\\
t_{\sigma}&0&t_{\pi}\\\ t_{\pi}e^{ik}&t_{\pi}&0\end{pmatrix},$ (32)
$H_{0}^{\sigma,3}(k)=\begin{pmatrix}0&t_{\pi}&t_{\pi}e^{-ik}\\\
t_{\pi}&0&t_{\sigma}\\\ t_{\pi}e^{ik}&t_{\sigma}&0\end{pmatrix}.$ (33)

Thus, each of the three SSH-3 chains has a strong bond followed by two weak
bonds, i.e., a strong-weak-weak bond pattern. This is also visible in the
dimer pattern presented in Fig. 7. Note that the properties of the SSH-3 model
are quite different than those of the standard SSH model, as the above
Hamiltonians are not chiral-symmetric and are gapless at half-filling.
Therefore we cannot get an insulating phase protected by a winding number in
contrast to the standard SSH model (see Sec. IV.3).

From the point of view of the crystalline symmetry all three models of Eqs.
(31), (32), and (33) are inversion-symmetric, being a manifestation of the
$\mathcal{R}$ symmetry present in the helix, but only in the first SSH-3 model
the unit cell is compatible with the inversion symmetry taking a form of

${\cal I}=\begin{pmatrix}0&0&1\\\ 0&1&0\\\ 1&0&0\end{pmatrix},$ (34)

and yielding

${\cal I}H_{0}^{\sigma,1}(k){\cal I}^{-1}=H_{0}^{\sigma,1}(-k).$ (35)

Now we can use the inversion-invariant of Ref. [38] to characterize the
insulating phases of $H_{0}^{\sigma,1}(k)$ at $n=1/3$ and $n=2/3$ filling. The
invariant is given by the difference of number of occupied bands with
inversion eigenvalue $-1$ between $k=0$ and $k=\pi$. If there are no
degeneracy at high-symmetry points, this comes down to:

$\nu^{\prime}=\sum_{n=1}^{n_{f}}|\langle E_{n}(0)\,{\cal
I}\,|E_{n}(0)\rangle-\langle E_{n}(\pi)|\,{\cal I}\,|E_{n}(\pi)\rangle|,$ (36)

where $\nu^{\prime}$ is a positive-integer topological invariant and $n_{f}$
is the index of the topmost occupied band. For our choice of $t_{\sigma}$ and
$t_{\pi}$, we find that $\nu^{\prime}=1$ both for $n_{f}=1$ ($n=1/3$ filling)
and $n_{f}=2$ ($n=2/3$ filling). Thus, according to the arguments given in
Ref. [39], we can expect that the open-chain version supports end states,
which we discuss below.

Since the models $H_{0}^{\sigma,\gamma}(k)$ are related by a cyclic shift by
one lattice site, no matter how we cut the helical chain of the cubic model we
should observe a pair of end states (in each spin channel). These end states
are nothing other than the end states of the inversion-symmetric SSH-3 chain
$H_{0}^{\sigma,1}(k)$.

Figure 7: Dimers forming in the cubic model with $\lambda=0$ in Eq. (28). The
helical chain is defined in a cubic environment. Each dimer type, indicated by
a different colour, appears on every third bond. For each dimer type a
distinct orbital tunneling amplitude is amplified, so that each orbital
experiences a $t_{\sigma}$-$t_{\pi}$-$t_{\pi}$ strong-weak-weak tunneling
amplitude when moving along the chain. This gives an SSH-3 model for each
orbital flavor.

By analogy to the ordinary SSH model, it is possible to derive the exact
analytic form of the end-states of the model $H_{0}^{\sigma,1}$. It turns out
that their energy is equal to $\epsilon_{\pm}=\pm t_{\pi}$, so they satisfy
the real-space Schrödinger equation of the form

$(H_{0}^{\sigma,1}-\epsilon_{\pm})|\psi_{\pm}\rangle=0,$ (37)

with the eigenvectors of the left end states given by

$\langle
i|\psi_{-}^{L}\rangle=\begin{pmatrix}\eta&\eta&0&\eta^{2}&\eta^{2}&0&\eta^{3}&\eta^{3}&0&\cdots\end{pmatrix},$
(38)

and

$\langle
i|\psi_{+}^{L}\rangle=\begin{pmatrix}-\eta&\eta&0&\eta^{2}&-\eta^{2}&0&-\eta^{3}&\eta^{3}&0&\cdots\end{pmatrix},$
(39)

where $\eta$ stands for the ratio of the hopping amplitudes,
$\eta=|t_{\pi}/t_{\sigma}|$. This shows that these states can be normalized
only for $\eta<1$ (this marks the boundary of the topological phase). The
right end-states share the same energies and can be obtained from the left
ones by applying inversion symmetry $\cal I$, so effectively reversing the
order of coefficients of the above vectors. These states are the exponentially
decaying 1-1-0 CDW end states of Fig. 6(b) 222For a finite chain the vectors
$\\{\ket{\psi^{L}_{\pm}}\\}$ satisfy
$\left(H_{0}^{\sigma,1}-\epsilon_{\pm}\right)\ket{\psi^{L}_{\pm}}=\left|\xi_{\pm}\right>,$
(40) where $\left|\xi\right>$ is a vector with all coefficients but the last
equal to 0 and the last coefficient
$|\bra{N}\ket{\xi}|=|\left(|t_{\pi}|^{N+1}/|t_{\sigma}|^{N}\right)|$ –
exponentially suppressed with increasing system size.. We have verified that
$\ket{\psi^{L}_{\pm}}$ summed over both spin channels agree with the orbital
densities shown in Fig. 6(b), as obtained from exact diagonalization of
$H_{0}$ in Eq. (28).

###  IV.2 Relation between cubic and chalcogen models

Figure 8: The spectrum of the cubic model ($\alpha=90^{\circ}$,
$t_{\pi}/t_{\sigma}=-1/3$) on a finite open chain of Eq. (41) for chain length
$N=501$, for different values of $\lambda/t_{\sigma}$. The end states are
colored red, the occupied bulk states are colored blue and the unoccupied bulk
states are colored gray.

We now turn our attention to a more complex model than that in the previous
subsection. While we still keep $\alpha=90^{\circ}$, we consider a finite
value of the spin-orbit coupling $\lambda$. The Hamiltonian of this cubic
model is

$H_{\text{cubic}}(\lambda)\equiv H\Big{|}_{\alpha=90^{\circ}},$ (41)

where $H$ is given by Eq. (1). In Fig. 8 we plot the spectrum of the cubic
model (41) on a finite open chain with $\alpha=90^{\circ}$ as a function of
$\lambda/t_{\sigma}$. The result is very similar to the one obtained for the
chalcogen model (Fig. 5). For any value of $\lambda$ below the gap closing at
$\lambda_{\text{crit.}}/t_{\sigma}\approx 0.68$ we observe end states.
Moreover, the end states with and without spin-orbit coupling are
adiabatically connected in both gaps.

The similarity of Fig. 8 to Fig. 5, as well as that of Fig. 6(a) to Fig. 6(b),
shows that the difference between the realistic chalcogen model of Eq. (1) and
the simplified cubic model of Eq. (41) is purely quantitative. Consequently,
we conclude that the bulk-edge correspondence seen in the cubic model of Eq.
(28) with $\lambda=0$ (see also Fig. 7) using the factorization into three
SSH-3 models, holds also for small but finite deviations from
$\alpha=90^{\circ}$ bond angle and for $\lambda_{\text{crit.}}>\lambda>0$.
This means that the simple picture of the three SSH-3 models on a helix
captures the topological properties of the more realistic chalcogen model.
Even though the SSH-3 invariant in Eq. (36) present in the simpler cubic model
without spin-orbit coupling is no longer well-defined in the chalcogen model,
there exists a more general invariant given by Eq. (24) that characterizes the
topological phase of both models. The end states in both models are
adiabatically linked and the end states of the chalcogen model are most easily
explained using the simple physics of the SSH-3 model.

###  IV.3 Comparison with the orbital SSH model

Figure 9: Dimers forming in the orbital SSH model. The zig-zag chain is
defined on a simple square lattice. Each dimer type, indicated by a different
color, appears on every second bond. For each dimer type a distinct orbital
tunneling amplitude is amplified, so that each orbital experiences a strong-
weak tunneling amplitude when moving along the chain. This gives an SSH model
for each orbital flavor.

In this section we study electrons in a right-angle zigzag chain on a plane
with two active orbitals – the in-plane $p$ orbitals $p_{x}$ and $p_{y}$ [20,
21], which are parallel to the two zigzag bonds. This Hamiltonian has two
lattice sites in the unit cell, $a$ and $b$. We consider the filling of 1/2.
This model is the orbital SSH model discussed in Refs. [20, 21], and is given
by the following Hamiltonian

$H_{\text{SSH}}=\sum_{i=1}^{2}\tilde{T}_{i}\otimes\tilde{h}_{i}\otimes\mathbbm{1}_{2},$
(42)

where $\tilde{T}_{i}$ operators are $N\times N$ matrices which describe
nearest-neighbor hopping between sites of a zigzag chain of length $N$
(assuming that $N$ is divisible by $2$):

$\tilde{T}_{i}=\sum_{j=1}^{N/2-1}\left(\,\left|i+2j\right\rangle\left\langle
i+1+2j\right|+H.c.\,\right),$ (43)

with periodic boundary conditions
$\left|N+1\right\rangle\equiv\left|1\right\rangle$. The $2\times 2$ matrices
$\tilde{h}_{i}$ describe hopping between $p_{x}$ and $p_{y}$ orbitals along
the two distinct bonds:

$\tilde{h}_{1}=\left(\begin{array}[]{cc}t_{\sigma}&0\\\
0&t_{\pi}\end{array}\right),$ (44)
$\tilde{h}_{2}=\left(\begin{array}[]{cc}t_{\pi}&0\\\
0&t_{\sigma}\end{array}\right).$ (45)

The orbital SSH model of Eq. (42) is topologically non-trivial with four in-
gap states. The topological phase is easy to understand by analogy to the SSH
model [18]. First, note that the spin channels separate and one needs to
consider only one spin. Second, note that the $\propto t_{\sigma}$ tunneling
processes in Eqs. (44)-(45) link $p_{x}$ orbitals along $x$ bonds and $p_{y}$
orbitals along $y$ bonds, while the $t_{\pi}$ tunneling processes link the
remaining orbitals on each bond. Since there are no processes mixing orbital
flavors, this gives rise to two copies of a half-filled SSH model, one in each
orbital channel. The tunneling amplitudes in each channel follow a weak-strong
pattern ($\left|t_{\pi}/t_{\sigma}\right|<1$), with the patterns shifted by
one lattice distance from each other. The emerging dimer pattern is
illustrated in Fig. 9.

In the SSH model we find end states in a finite chain obtained by cutting an
infinite chain across two strong bonds (dimers), that is bonds with a larger
absolute value of tunneling amplitude – in our case $t_{\sigma}$. In other
words, the presence of end states is dependent on the cut, i.e., if one cuts
across two weak bonds, no end states are found. In Fig. 9 we see that no
matter how we cut the chain in the orbital SSH model we will always cut two
dimers – it is just the orbital character of these dimers that depends on the
cut. We should therefore always find end states in the orbital SSH chain.

The analogy with two decoupled and shifted SSH chains (per spin channel)
allows us to understand the four end states in the orbital SSH model. They
stem from the bulk invariant – the winding number of the SSH model. This is in
contrast to the cubic model at $\lambda=0$ of Eq. (28) (and, therefore, also
the chalcogen model). For the SSH model the invariant comes from the
fundamental chiral symmetry and is therefore robust, unaffected by lattice
deformations. In the cubic model, on the other hand, we find an underlying
SSH-3 model which is gapless at half-filling and cannot exhibit end-states
protected by chiral symmetry. Instead the end states are protected by the
crystalline symmetry (see above).

##  V Conclusion

In this work, we have studied the topological properties of the single-element
chalcogen crystals selenium and tellurium. To this end, we have concentrated
our efforts on investigating an electronic model with active $p$ orbital
degrees of freedom in helical chains – since such chains form the backbone of
the chalcogen crystal structure. While we leave the coupling between the
chalcogen chains for a later study, this is not expected to alter the main
findings of this work on the existence of topological end states. Rather, the
coupling is expected to allow the end states to bundle into coherent surface
states.

Our work has two main results: First, we have shown that the occupied
electronic bands of a realistic model for the chalcogen chains carry a
$\mathbb{Z}_{2}$ topological invariant. The invariant is non-trivial for
realistic values of spin-orbit coupling and becomes trivial only for
unrealistically large values. The invariant is related to the number of
particular eigenvalues of the Wilson loop operator $\nu$ [38], and is
protected by the 180∘ rotation symmetry $\mathcal{R}$ around an axis $\vec{r}$
normal to the helical chain. The end states arising in the topologically non-
trivial phase of the open chain are orbitally polarised, i.e., depending on
the particular cut of the helical chain, they have a predominantly $p_{x}$,
$p_{y}$ or $p_{z}$ orbital character.

Second, we have demonstrated that the onset of the nontrivial topology can be
understood in terms of a simplified cubic model. This model also describes a
helical chain with active $p$ orbital degrees of freedom, but, unlike for the
chalcogens, with the bond angle within the helix being 90∘ and with vanishing
spin-orbit coupling. While the topological properties of the cubic model
without spin-orbit coupling are equivalent to the ones found in the realistic
chalcogen chain, its topology can also be simply understood as stemming from
the topological properties of an SSH-3 model. The SSH-3 model shows a strong-
weak-weak bond pattern in one of the three $p$ orbital flavors ($p_{x}$ or
$p_{y}$ or $p_{z}$), its topological invariant $\nu^{\prime}$ is protected by
the inversion symmetry $\mathcal{I}$, and the end states are entirely
orbitally-polarised.

Note that the SSH-3 model is not only different from the standard SSH model
[18], but also from its orbital variant [21]. The latter contains a set of two
standard SSH chains in each spin channel, one per each $p$-orbital channel.
Furthermore, both the orbital SSH model and the standard SSH model have a
winding-number topological invariant protected by chiral symmetry.

Note added. – After this work was completed, we became aware of a recent study
of the electronic structure of Se and Te [42]. There, it is suggested that the
one-dimensional Se/Te chains are chiral SSH chains. In contrast, we show that
the topology of the chalcogen chain stems from the exact Lau-Brink-Ortix
invariant [34]. It can also be understood by analogy to the SSH-3 model, which
is topologically distinct from the chiral SSH model.

## Acknowledgements

We kindly acknowledge support by the National Science Centre (NCN, Poland)
under Projects No. 2016/22/E/ST3/00560, 2019/34/E/ST3/00404, and
2021/43/B/ST3/02166. The work is supported by the Foundation for Polish
Science through the IRA Programme co-financed by EU within SG OP Programme.
A.L. acknowledges support from a Marie Skłodowska-Curie Individual Fellowship
under grant MagTopCSL (ID 101029345). A.M.O. is grateful for support via the
Alexander von Humboldt Foundation Fellowship [43] (Humboldt-Forschungspreis).

For the purpose of Open Access, the authors have applied a CC-BY public
copyright licence to any Author Accepted Manuscript (AAM) version arising from
this submission.

W.B. and A.L. have contributed equally to this paper.

##  Appendix A End states in the continuous limit

In this section we investigate whether the chirality-breaking scalar term in
the chalcogen chain model linearized around the upper gap closing [Eq. (26)]
affects the in-gap states. To this end we use the standard adiabatic argument
saying that interfacing a non-trivial phase of $h(k,\lambda)$ with the vaccum
is equivalent to going from the non-trivial to trivial phase by changing the
parameters of the model. This involves passing through the gapless point. The
linear expansion of $h(k,\lambda)$ around the gap closing point is given by

$\displaystyle\begin{split}h_{\text{lin.}}&=h(k_{0}+\delta
k,\lambda_{\text{crit.}}+\delta\lambda)\\\
&=E_{0}+\delta\lambda\left(a_{1}+a_{2}\sigma_{x}+a_{3}\sigma_{z}\right)+\delta
k\left(b_{1}+b_{2}\sigma_{z}\right)+\\\ &+{\cal O}\left(\delta
k^{2},\delta\lambda^{2}\right),\end{split}$ (46)

where $a_{i}$ and $b_{i}$ are real numbers. To study the real-space interface
between trivial and non-trivial phases we work in the continuum setting:

$\displaystyle\delta k$ $\displaystyle=$ $\displaystyle-i\partial_{x},$ (47)
$\displaystyle\delta\lambda$ $\displaystyle=$ $\displaystyle\lambda(x),$

and we search for solutions of the Schrodinger equation of the form:

$h_{{\rm lin.}}\left(\begin{array}[]{c}\psi_{a}(x)\\\
\psi_{b}(x)\end{array}\right)=E\left(\begin{array}[]{c}\psi_{a}(x)\\\
\psi_{b}(x)\end{array}\right),$ (48)

or equivalently:

$\displaystyle\begin{split}&-i\;M\partial_{x}\left(\begin{array}[]{c}\psi_{a}(x)\\\
\psi_{b}(x)\end{array}\right)+R\lambda(x)\left(\begin{array}[]{c}\psi_{a}(x)\\\
\psi_{b}(x)\end{array}\right)\\\
&=(E-E_{0})\left(\begin{array}[]{c}\psi_{a}(x)\\\
\psi_{b}(x)\end{array}\right),\end{split}$ (49)

with

$M=\left(\begin{array}[]{cc}b_{1}+b_{2}&0\\\ 0&b_{1}-b_{2}\end{array}\right),$
(50)

and

$R=\left(\begin{array}[]{cc}a_{1}+a_{3}&a_{2}\\\
a_{2}&a_{1}-a_{3}\end{array}\right).$ (51)

We look for solutions of the homogenous equation, that is we look for a state
whose energy is equal to the gap closing energy $E_{0}$. We use the ansatz

$\displaystyle\psi_{a}(x)=C_{a}\exp\left(\gamma\int\lambda(x)dx\right)$ (52)
$\displaystyle\psi_{b}(x)=C_{b}\exp\left(\gamma\int\lambda(x)dx\right).$ (53)

The homogenous equation now reads:

$\lambda(x)\exp\left(\gamma\int\lambda(x)dx\right)A\left(\begin{array}[]{c}C_{a}\\\
C_{b}\end{array}\right)=0,$ (54)

with

$A=\left(\begin{array}[]{cc}-ib_{+}\gamma+a_{+}&a_{2}\\\
a_{2}&-ib_{-}\gamma+a_{-})\end{array}\right).$ (55)

If $\lambda(x)\neq 0$, a solution exists for such $\gamma$ that
$\text{Det}(A)=0$. In the case of the chalcogen chain with realistic
parameters, that is $\alpha=103^{\circ}$, $t_{\sigma}=1$, and
$t_{\pi}/t_{\sigma}=-1/3$, we found two solutions: $\gamma_{1}=-0.582+0.056i$,
and $\gamma_{2}=0.582+0.056i$, i.e., $\gamma_{2}=-\gamma_{1}^{*}$.

Depending on the shape of the domain wall, only one of these solutions is
physical and we can set the other to zero so that $\Psi_{i}(\pm\infty)=0$. For
instance, choosing

$\displaystyle\lambda(x)$ $\displaystyle=$
$\displaystyle\Lambda\;\tanh\left(\frac{x}{w}\right),$ (56)
$\displaystyle\int\lambda(x)dx$ $\displaystyle=$ $\displaystyle\Lambda
w\log\cosh\left(\frac{x}{w}\right)>0\,,$

where for $x\to-\infty$ we have a topological phase with
$\lambda=\lambda_{\text{crit.}}-\Lambda$, and for $x\to+\infty$ a trivial
phase with $\lambda=\lambda_{\text{crit.}}+\Lambda$, we set
$C_{a,2}=C_{b,2}=0$. Then the solution is

$\Psi_{1}(x)=\left(\begin{array}[]{c}C_{a,1}\\\
C_{b,1}\end{array}\right)\exp\left(\gamma_{1}w\Lambda\log\cosh\frac{x}{w}\right).$
(57)

We obtain the localization length of the domain-wall state by taking the limit
of large $x$,

$\Psi_{1}(x\to\pm\infty)\propto e^{\left|x\right|\Lambda\gamma_{1}},$ (58)

hence the localisation length is

$\xi=\frac{1}{\Lambda\left|\text{Re}(\gamma_{1})\right|}\,.$ (59)

This characterizes one end state. A similar analysis of the other gap closing
point, for $k_{0}^{\prime}<0$, yields another.

Since the homogenous solution describes a state with energy $E_{0}$ – that is
in the middle of the gap for all values of $\lambda(x)<\lambda_{\text{crit.}}$
– and whose localization length is finite, the conclusion from the above
calculation is that in the vicinity of the gap closing the localized end state
will survive in the gap even in the presence of the scalar term which breaks
the chiral symmetry. Let us reiterate that this conclusion applies only in the
region around the gap closing, where the linearized model is valid. In
particular, the pinning of end states to the gap for all
$\lambda(x)<\lambda_{\text{crit.}}$ is a feature which is in contrast with the
dispersion of the end states in the upper gap observed in the full model (see
Fig. 5). The reason for this is that the linearized model is a 2-band
$p$-orbital model and consequently the spin-orbit interaction in this model is
much simpler than that in the full model, namely,

$H_{\text{soc},\text{lin.}}\propto-\lambda L_{z}\otimes\sigma_{z}.$ (60)

As a result, the linearized model does not capture the evolution of the end-
state energy. It does, however, show that the end states survive near the gap
closing, and taken together with the results of exact diagonalization (see
Fig. 5) it allows us to conclude that the end states do indeed survive for
$\lambda(x)<\lambda_{\text{crit.}}$.

##  Appendix B DOS as a function of $\lambda$

Figure 10: The spectrum of the cubic model of Eq. (41) ($\alpha=90^{\circ}$,
$t_{\sigma}/t_{\pi}=-1/3$) for both positive and negative values of
$\lambda/t_{\sigma}$. The end states are colored red, the occupied bulk states
are colored blue and the unoccupied bulk states are colored gray. The
$\lambda/t_{\sigma}<0$ part of the spectrum corresponds to the filling 1/3,
while the $\lambda/t_{\sigma}\geq 0$ part of the spectrum to the filling 2/3 –
as indicated by the coloring of bands. The inset shows the spin-orbit coupling
eigenvalues in the $p$ shell. A symmetry is visible in the DOS upon changing
the sign of the spin-orbit interaction, the filling from 1/3 to 2/3 or vice
versa, and performing an electron-hole transformation.

In this Appendix, we provide a simple understanding of the evolution of DOS of
the cubic model of Eq. (41) ($\alpha=90^{\circ}$) with varying spin-orbit
interaction strength $\lambda$. Even though for simplicity we only consider
the cubic case, the results of this section are also valid for the chalcogen
model of Eq. (1) ($\alpha=103^{\circ}$).

An observation can be made about Fig. 8 – the spectrum of the cubic model on
an open chain for $\lambda=0$, given by Eq. (28), is symmetric with respect to
$E/t_{\sigma}=0$. However, this model is not chiral – this is a feature of the
open chain alone. Because of this, for $\lambda=0$, the spectrum around 1/3
filling looks the same as the spectrum around 2/3 filling. For $\lambda\neq
0$, however, the physics in the two cases becomes different. While in the
upper gap we observe a gap closing, in the lower gap we do not. This asymmetry
is the consequence of the fact that spin-orbit interaction breaks particle-
hole symmetry. The large $\lambda/t_{\sigma}$ limit, as well as the gap
closing, can be better understood by considering the orbital character of the
bands.

The key observation is that, generally, and in particular at the gap closing
momentum $k=\pi/3$, the top and bottom bands have a predominantly $l_{z}=\pm
1$ character, while the middle band has a predominantly $l_{z}=0$ character,
where $l_{z}$ is the projection of the electron angular momentum onto the $z$
axis of the orbital basis [30]. Consequently, the $j=1/2$ pseudospin states
are concentrated mostly in the middle band, while the $j=3/2$ pseudospin
states dominate the top and bottom bands. Because the spin-orbit energy is
positive for the $j=1/2$ pseudospin states and negative for the $j=3/2$
pseudospin states, we observe an upwards trend in the energy evolution of the
middle band, which is not observed in the top and bottom bands. What follows
is that there occurs a crossing of the two top bands but no crossing of the
bottom two bands, a features which is visible in Fig. 5 and in Fig. 8.

Finally, note that if the filling were 1/3 instead of 2/3, the spin-orbit
interaction would change sign and the DOS in Figs. 5 and 8 would flip around
the $E/t_{\sigma}=0$ axis in addition to the Fermi level moving to the lower
gap—see Fig. 10. Consequently, the topological properties for the 1/3 filling
are the same as for the 2/3 filling.

## References

  * Moore [2010] J. E. Moore, The birth of topological insulators, Nature 464, 194 (2010).
  * Hasan and Kane [2010] M. Z. Hasan and C. L. Kane, Colloquium: Topological insulators, Rev. Mod. Phys. 82, 3045 (2010).
  * Qi and Zhang [2011] X.-L. Qi and S.-C. Zhang, Topological insulators and superconductors, Rev. Mod. Phys. 83, 1057 (2011).
  * Hasan and Moore [2011] M. Z. Hasan and J. E. Moore, Three-dimensional topological insulators, Ann. Rev. Cond. Matt. Phys. 2, 55 (2011).
  * Tokura  _et al._ [2019] Y. Tokura, K. Yasuda, and A. Tsukazaki, Magnetic topological insulators, Nat. Rev. Phys. 1, 126 (2019).
  * Altland and Zirnbauer [1997] A. Altland and M. R. Zirnbauer, Nonstandard symmetry classes in mesoscopic normal-superconducting hybrid structures, Phys. Rev. B 55, 1142 (1997).
  * Chiu  _et al._ [2013] C.-K. Chiu, H. Yao, and S. Ryu, Classification of topological insulators and superconductors in the presence of reflection symmetry, Phys. Rev. B 88, 075142 (2013).
  * Chiu and Schnyder [2014] C.-K. Chiu and A. P. Schnyder, Classification of reflection-symmetry-protected topological semimetals and nodal superconductors, Phys. Rev. B 90, 205136 (2014).
  * Shiozaki and Sato [2014] K. Shiozaki and M. Sato, Topology of crystalline insulators and superconductors, Phys. Rev. B 90, 165114 (2014).
  * Shiozaki  _et al._ [2016] K. Shiozaki, M. Sato, and K. Gomi, Topology of nonsymmorphic crystalline insulators and superconductors, Phys. Rev. B 93, 195413 (2016).
  * König  _et al._ [2007] M. König, S. Wiedmann, C. Brüne, A. Roth, H. Buhmann, L. W. Molenkamp, X.-L. Qi, and S.-C. Zhang, Quantum spin Hall insulator state in HgTe quantum wells, Science 318, 766 (2007).
  * Dziawa  _et al._ [2012] P. Dziawa, B. J. Kowalski, K. Dybko, R. Buczko, A. Szczerbakow, M. Szot, E. Łusakowska, T. Balasubramanian, B. M. Wojek, M. H. Berntsen, O. Tjernberg, and T. Story, Topological crystalline insulator states in Pb1-xSnxSe, Nat. Mater. 11, 1023 (2012).
  * Tanaka  _et al._ [2012] Y. Tanaka, Z. Ren, T. Sato, K. Nakayama, S. Souma, T. Takahashi, K. Segawa, and Y. Ando, Experimental realization of a topological crystalline insulator in SnTe, Nat. Phys. 8, 800 (2012).
  * Mourik  _et al._ [2012] V. Mourik, K. Zuo, S. M. Frolov, S. R. Plissard, E. P. A. M. Bakkers, and L. P. Kouwenhoven, Signatures of Majorana fermions in hybrid superconductor-semiconductor nanowire devices, Science 336, 1003 (2012).
  * Wang  _et al._ [2013] Z. Wang, H. Weng, Q. Wu, X. Dai, and Z. Fang, Three-dimensional dirac semimetal and quantum transport in Cd3As2, Phys. Rev. B 88, 125427 (2013).
  * Lv  _et al._ [2015] B. Q. Lv, H. M. Weng, B. B. Fu, X. P. Wang, H. Miao, J. Ma, P. Richard, X. C. Huang, L. X. Zhao, G. F. Chen, Z. Fang, X. Dai, T. Qian, and H. Ding, Experimental discovery of Weyl semimetal TaAs, Phys. Rev. X 5, 031013 (2015).
  * Xiao  _et al._ [2018] D. Xiao, J. Jiang, J.-H. Shin, W. Wang, F. Wang, Y.-F. Zhao, C. Liu, W. Wu, M. H. W. Chan, N. Samarth, and C.-Z. Chang, Realization of the axion insulator state in quantum anomalous Hall sandwich heterostructures, Phys. Rev. Lett. 120, 056801 (2018).
  * Heeger  _et al._ [1988] A. J. Heeger, S. Kivelson, J. R. Schrieffer, and W. P. Su, Solitons in conducting polymers, Rev. Mod. Phys. 60, 781 (1988).
  * He and Chien [2020] Y. He and C.-C. Chien, Non-Hermitian generalizations of extended Su–Schrieffer–Heeger models, J. Phys. Cond. Matter 33, 085501 (2020).
  * St-Jean  _et al._ [2017] P. St-Jean, V. Goblot, E. Galopin, A. Lemaître, T. Ozawa, L. Le Gratiet, I. Sagnes, J. Bloch, and A. Amo, Lasing in topological edge states of a one-dimensional lattice, Nat. Photon. 11, 651 (2017).
  * Sun  _et al._ [2020] G. Sun, W.-L. You, and T. Zhou, Topological phases of spinless $p$-orbital fermions in zigzag optical lattices, Phys. Rev. A 102, 063301 (2020).
  * Fukutome [1984] H. Fukutome, Charge density wave-soliton model for Se and Te, Prog. Theor. Phys. 71, 1 (1984).
  * Shimoi and Fukutome [1992] Y. Shimoi and H. Fukutome, Vector charge density wave model of metallic and trigonal Te, Prog. Theor. Phys. 87, 307 (1992).
  * Nakamura and Ikawa [2002] K. Nakamura and A. Ikawa, Vibrational properties of a regular helical Se chain, Phys. Rev. B 66, 024306 (2002).
  * Von Hippel [1948] A. Von Hippel, Structure and conductivity in the VI b group of the periodic system, J. Chem. Phys. 16, 372 (1948).
  * Reitz [1957] J. R. Reitz, Electronic band structure of selenium and tellurium, Phys. Rev. 105, 1233 (1957).
  * Olechna and Knox [1965] D. J. Olechna and R. S. Knox, Energy-band structure of selenium chains, Phys. Rev. 140, A986 (1965).
  * Silva  _et al._ [2018] A. Silva, J. Henke, and J. van Wezel, Elemental chalcogens as a minimal model for combined charge and orbital order, Phys. Rev. B 97, 045151 (2018).
  * Silva and van Wezel [2018] A. Silva and J. van Wezel, The simple-cubic structure of elemental Polonium and its relation to combined charge and orbital order in other elemental chalcogens, SciPost Phys. 4, 028 (2018).
  * Kłosiński  _et al._ [2021] A. Kłosiński, A. M. Oleś, C. E. Agrapidis, J. van Wezel, and K. Wohlfeld, Chalcogenic orbital density waves in the weak-and strong-coupling limit, Phys. Rev. B 103, 235123 (2021).
  * Agapito  _et al._ [2013] L. A. Agapito, N. Kioussis, W. A. Goddard, and N. P. Ong, Novel family of chiral-based topological insulators: Elemental tellurium under strain, Phys. Rev. Lett. 110, 176401 (2013).
  * Ideue  _et al._ [2019] T. Ideue, M. Hirayama, H. Taiko, T. Takahashi, M. Murase, T. Miyake, S. Murakami, T. Sasagawa, and Y. Iwasa, Pressure-induced topological phase transition in noncentrosymmetric elemental tellurium, PNAS 116, 25530 (2019).
  * Pal  _et al._ [2013] A. Pal, S. N. Shirodkar, S. Gohil, S. Ghosh, U. V. Waghmare, and P. Ayyub, Multiferroic behavior in elemental selenium below 40 K: effect of electronic topology, Sci. Rep. 3, 1 (2013).
  * Lau  _et al._ [2016] A. Lau, J. van den Brink, and C. Ortix, Topological mirror insulators in one dimension, Phys. Rev. B 94, 165164 (2016).
  * Matsui [2014] M. Matsui, Role of interchain interaction in determining the band gap of trigonal selenium: A density functional theory study with a linear combination of Bloch orbitals, J. Phys. Chem. C 118, 19294 (2014).
  * Note [1] After a 360-degree rotation – which, for a 120-degree rotation, corresponds to the special translation by a full unit cell – the spin basis acquires an extra minus sign. Therefore we choose $U$ to be a 480-degree spin rotation. With this choice, we have, on the one hand, $120=480\tmspace+{.2777em}\text{mod}\tmspace+{.2777em}360$, and on the other $U^{3}=\mathbb{1}$. The first property is essential in order for the spin-orbit interaction to be site-independent.
  * Slater and Koster [1954] J. C. Slater and G. F. Koster, Simplified LCAO method for the periodic potential problem, Phys. Rev. 94, 1498 (1954).
  * Alexandradinata  _et al._ [2014] A. Alexandradinata, X. Dai, and B. A. Bernevig, Wilson-loop characterization of inversion-symmetric topological insulators, Phys. Rev. B 89, 155114 (2014).
  * Rhim  _et al._ [2017] J.-W. Rhim, J. Behrends, and J. H. Bardarson, Bulk-boundary correspondence from the intercellular Zak phase, Phys. Rev. B 95, 035421 (2017).
  * Silva and van Wezel [2021] A. Silva and J. van Wezel, Topological states between inversion symmetric atomic insulators, SciPost Phys. 10, 137 (2021).
  * Note [2] For a finite chain the vectors $\\{\ket{\psi^{L}_{\pm}}\\}$ satisfy  $\left(H_{0}^{\sigma,1}-\epsilon_{\pm}\right)\ket{\psi^{L}_{\pm}}=\left|\xi_{\pm}\right>,$ (61) where $\left|\xi\right>$ is a vector with all coefficients but the last equal to 0 and the last coefficient $|\bra{N}\ket{\xi}|=|\left(|t_{\pi}|^{N+1}/|t_{\sigma}|^{N}\right)|$ – exponentially suppressed with increasing system size.
  * Zhang  _et al._ [2022] R. Zhang, J. Deng, Y. Sun, Z. Fang, Z. Guo, and Z. Wang, Large shift current, $\pi$ Zak phase and unconventional nature in Se and Te, arXiv:2211.04116 (2022).
  * [43] https://www.humboldt-foundation.de/web/humboldt-preis.html.