# Anomalous dielectric response in insulators with the $\pi$ Zak phase

Yusuke Aihara  These authors contributed equally to this work. Department of
Physics, Tokyo Institute of Technology,2-12-1 Ookayama, Meguro-ku, Tokyo
152-8551, Japan  Motoaki Hirayama  These authors contributed equally to this
work. RIKEN Center for Emergent Matter Science, Wako, Saitama 351-0198, Japan
Shuichi Murakami  Department of Physics, Tokyo Institute of Technology,2-12-1
Ookayama, Meguro-ku, Tokyo 152-8551, Japan  TIES, Tokyo Institute of
Technology, 2-12-1 Ookayama, Meguro-ku, Tokyo 152-8551, Japan

###### Abstract

In various topological phases, nontrivial states appear at the boundaries of
the system. In this paper, we investigate anomalous dielectric response caused
by such states caused by the $\pi$ Zak phase. First, by using the one-
dimensional Su-Schrieffer-Heeger model, we show that, when the system is
insulating and the Zak phase is $\pi$, the polarization suddenly rises to a
large value close to $e/2$, by application of an external electric field. The
$\pi$ Zak phase indicates existence of half-filled edge states, and we
attribute this phenomenon to charge transfer between the edge states at the
two ends of the system. We extend this idea to two- and three-dimensional
insulators with the $\pi$ Zak phase over the Brillouin zone, and find similar
anomalous dielectric response. We also show that diamond and silicon slabs
with (111) surfaces have the $\pi$ Zak phase by ab intio calculations, and
show that this anomalous response survives even surface reconstruction
involving an odd number of original surface unit cells. Another material
example with an anomalous dielectric response is polytetrafluoroethylene
(PTFE), showing plateaus of polarization at $\pm e$ by ab initio calculation,
in agreement with our theory.

###### pacs:

##  I introduction

Topological materials such as topological insulators are currently under
intense investigations due to their unique properties Kane and Mele (2005); Fu
_et al._ (2007); Bernevig  _et al._ (2006). They are called topological
because the $k$-space structure of electronic bands has nontrivial topology.
In the topological phase, nontrivial localized states generally appear at the
boundaries of the system, and this important correspondence is called bulk-
edge correspondence Fu and Kane (2006); Fu  _et al._ (2007). These
edge/surface states often cause distinctive phenomena such as quantum spin
Hall effect.

In this paper we focus on topological phases characterized by the Zak phase.
The Zak phase Zak (1989) is defined as a special case of the Berry phase Berry
(1984). It is useful to determine presence or absence of topological edge
states in systems with time-reversal (TR) and space inversion (SI) symmetries.
In spinless systems, the Zak phase is quantized to be $0$ or $\pi$ under such
symmetries and the systems have topological edge states in the latter case Zak
(1989); Ryu and Hatsugai (2002); Delplace  _et al._ (2011); Kariyado and
Hatsugai (2013); Rhim  _et al._ (2017); van Miert and Ortix (2017); Pletyukhov
_et al._ (2020). The Zak phase is also known to represent the bulk
polarization of the system Resta (1992); King-Smith and Vanderbilt (1993);
Vanderbilt and King-Smith (1993); Resta (1994), and at the same time the
amount of surface charge on a crystalline solid King-Smith and Vanderbilt
(1993); Vanderbilt and King-Smith (1993). In systems with TR and SI
symmetries, the surface charge per unit cell is quantized to be $0$ or $e/2$
modulo $e$ due to the quantization of the Zak phase where $-e$ is the electron
charge. For one-dimensional systems, no edge charge accumulates at both ends
when the Zak is 0 and the charge $e/2$ (mod $e$) accumulates at each edge when
the Zak phase is $\pi$. For two- or three-dimensional insulating systems with
TR and SI symmetries, in which the Zak phase is constant over the Brillouin
zone, no edge charge per unit cell accumulates at both surfaces when the Zak
is 0 and the charge $e/2$ (mod $e$) per unit cell accumulates at both surfaces
when the Zak phase is $\pi$. These huge surface charges do not contradict the
inversion symmetry because an equal amount of charges is accumulated at both
ends. In half-filled one-dimensional insulating systems with the $\pi$ Zak
phase, the edge states at both ends of the system are half-filled under SI
symmetry. These edge states are expected to respond sensitively to an external
field which breaks the symmetry.

In this paper, we apply a static electric field to such insulators with the
$\pi$ Zak phase and calculate their dielectric response. As a result, we find
a sudden uprise of the polarization even for a very weak electric field. First
we show this phenomena in a one-dimensional Su-Schrieffer-Heeger (SSH) model
with a Zak phase equal to $\pi$, which in turn has topological in-gap edge
states of the chain. We also confirm that this anomalous phenomenon is
attributed to the topological edge states. In two- and three-dimensional
insulating systems with $\pi$ Zak phase, we find a similar anomalous behavior
of the polarization when the edge or surface bands are sufficiently flat.
Insulators with the $\pi$ Zak phases have actually been found, such as silicon
and diamond Vanderbilt and King-Smith (1993), and some topological electrides
Hirayama  _et al._ (2017, 2018). For example, diamond and silicon slabs with
the (111) surfaces have the $\pi$ Zak phase and relatively flat midgap states.
We also show such a dielectric response of polytetrafluoroethylene (PTFE) with
plateaus at $P\sim\pm e$ by ab initio calculation, in agreement with our
theory. Our results are expected to be applicable to such materials.

The organization of the paper is as follows. In Sec. II, we first confirm that
the anomalous dielectric response occurs in the one-dimensional SSH model and
also confirm that the phenomenon is caused by the topological edge states. In
Sec. III, we also show similar phenomena in a three-dimensional model. We
discuss several points on realizing anomalous dielectric responses in Sec. IV.
In Sec. V, we show that diamond and silicon slabs with the (111) surfaces have
the $\pi$ Zak phase leading to midgap surface states. We also show an
anomalous dielectric response of polytetrafluoroethylene (PTFE). We then
discuss the relationship between our results and the bulk polarization
described in the modern theory of polarization in Sec. VI. We summarize our
result in Sec. VII. Throughout the paper, we consider systems with weak spin-
orbit coupling, and ignore it. Thus the polarizations and surface charges in
the subsequent results should be doubled, in order to include the spin degree
of freedom.

##  II One-dimensional model (Model I)

In this chapter, we show anomalous dielectric response in the SSH model Su
_et al._ (1979), which we call Model I. This model becomes a one-dimensional
topological insulator by adjusting its parameters. This topological insulator
phase is characterized by the value of the Zak phase equal to $\pi$ and this
model has topological edge states.

###  II.1 Previous research

Figure 1: SSH model. (a) Infinite SSH chain. The shaded square represents a
unit cell. Hopping amplitudes alternate between the intracell hopping $t_{1}$
(red lines) and the intercell hopping $t_{2}$ (blue lines). (b) Finite open
SSH chain. There are $2N$ sites in the chain. (c) Band structure for
$t_{1}=0.5,t_{2}=1$ (red) and $t_{1}=2,t_{2}=1$ (blue). (d) Dielectric
polarization of the SSH model. Here we set $N=20$ and the number of sites is
$2N=40$. The hopping amplitudes are $t_{1}=0.5,t_{2}=1$ for the red line and
$t_{1}=2,t_{2}=1$ for the blue line. In the former case, $P(\varepsilon)$
rapidly changes toward $\pm e/2$ in the vicinity of $\varepsilon=0$.

We first review the definition of the Zak phase and its connection with
electric polarization. For one-dimensional crystalline insulators, the Zak
phase is defined as

$\displaystyle\theta=-i\sum_{n}^{\text{occ.}}\int_{0}^{2\pi/a}dk\bra{u_{n}(k)}\dfrac{\partial}{\partial
k}\ket{u_{n}(k)},$ (1)

where $a$ is the lattice constant, $k$ is the wavenumber, $u_{n}(k)$ is the
periodic part of the Bloch wave function of the $n$th eigenstate, and the
summation is over the occupied states. This Zak phase is defined modulo
$2\pi$. It is quantized as 0 or $\pi$ under the SI symmetry Kariyado and
Hatsugai (2013). When $\theta=\pi$, there exist two degenerated edge states in
the gap when the system preserves chiral symmetry Ryu and Hatsugai (2002).
According to the modern theory of polarization King-Smith and Vanderbilt
(1993); Vanderbilt and King-Smith (1993), the electric polarization $P_{0}$
can be obtained from the Zak phase in the following equation;

$\displaystyle P_{0}=-e\dfrac{\theta}{2\pi}\ ({\rm mod}\ e).$ (2)

This equation indicates that the polarization is $e/2$ (mod $e$) in the system
with the $\pi$ Zak phase. We note that in one-dimensional systems, the
polarization, which is an electric dipole moment per unit length, has the unit
of charge.

In this section, we consider the SSH model. It is a tight-binding model on a
one-dimensional lattice with two sites $A_{i}$ and $B_{i}$ in the $i$th unit
cell and with different intracell ($-t_{1}$) and intercell ($-t_{2}$) hopping
amplitudes between $A$ and $B$ sites as shown in Fig. 1(a). We set $t_{1}$ and
$t_{2}$ to be positive. The distance between the adjacent A site and B site is
set to be $a/2$. The Bloch Hamiltonian is given by

$\displaystyle\mathcal{H}^{\text{bulk}}(k)=\left(\begin{array}[]{cc}0&-t_{1}-t_{2}e^{-ika}\\\
-t_{1}-t_{2}e^{ika}&0\end{array}\right)\ .$ (5)

This model has chiral symmetry, and therefore the spectrum is symmetric with
respect to $E=0$. The energy eigenvalues are given by
$E=\pm\sqrt{t_{1}^{2}+t_{2}^{2}+2t_{1}t_{2}\cos ka}$ and the band gap exists
when $t_{1}\neq t_{2}$ (see Fig. 1(c)). When $t_{1}$ and $t_{2}$ are equal,
the gap closes at $k=\pi/a$.

Here we focus on the case of $t_{1}\neq t_{2}$, where the system is
insulating. We set the Fermi energy to be $E_{f}=0$, i.e. the system is half-
filled. The Zak phase of the SSH model is calculated by substituting the
eigenstates of (5) into (1), and we get

$\displaystyle\theta$ $\displaystyle=\begin{cases}0\ (t_{1}>t_{2})\\\ \pi\
(t_{2}>t_{1})\end{cases}.$ (6)

This value of the Zak phase can be calculated from the parity eigenvalues at
time-reversal invariant momenta (TRIM), $k=0$ and $k=\pi$:

$\displaystyle e^{i\theta}=\prod_{n}^{\text{occ.}}\xi_{n}(k=0)\xi_{n}(k=\pi),$
(7)

where $\xi_{n}(k)$ is the parity eigenvalue of the $n$th band at $k$, and the
product is taken over the occupied bands. In the present case of the SSH
model, from the inversion operator at the center of the unit cell
$\hat{P}=\left(\begin{array}[]{cc}0&1\\\ 1&0\end{array}\right)$, we get

$\displaystyle\xi(k=0)=1,\
\xi(k=\pi)=\left\\{\begin{array}[]{ll}1&t_{1}>t_{2}\\\
-1&t_{2}>t_{1}\end{array}\right.$ (10)

and we reproduce the result (6). This also means that the bands are inverted
at $k=\pi$ when $t_{1}$ is changed across the gap closing at $t_{1}=t_{2}$.

From Eq. (2), Eq. (6) indicates that an electric charge $e/2\ ({\rm mod}\ e)$
accumulates at both ends when $t_{2}>t_{1}$ meaning that each of the
topological edge states are half filled King-Smith and Vanderbilt (1993);
Vanderbilt and King-Smith (1993). In order to have such topological edge
states due to the $\pi$ Zak phase, the chain should be fully covered by the
chosen unit cells. In this way, in order to see bulk-edge correspondence for
the Zak phase, the termination of the crystal is related with the choice of
the unit cell with which the Zak phase is calculated. We note that appearance
of the polarization $P_{0}$ does not contradict the inversion symmetry because
$P_{0}\equiv e/2\equiv-e/2\ $(mod $e$).

###  II.2 Dielectric response

We next consider the dielectric polarization in a finite SSH chain (see Fig. 1
(b)). The chain consists of $2N$ sites, and each site is represented as $\sl
A_{i}$ and $\sl B_{i}(i=1,2,\cdots N)$ as shown in the figure. The Hamiltonian
of the finite SSH chain is given by

$\displaystyle H_{0}$
$\displaystyle=-t_{1}\sum_{n=1}^{N}a^{\dagger}_{n}b_{n}-t_{2}\sum_{n=1}^{N-1}a^{\dagger}_{n+1}b_{n}+{\sl
h.c.},$ (11)

where $a^{\dagger}_{i}$ and $b^{\dagger}_{i}$ are creation operators at the
$\sl A_{i}$ and $\sl B_{i}$ sites respectively. The Hamiltonian with static
electric field $\varepsilon$ is given by

$\displaystyle H$ $\displaystyle=H_{0}+H_{\varepsilon},$ (12) $\displaystyle
H_{\varepsilon}$ $\displaystyle=\varepsilon
ea\sum_{n=1}^{N}\left[\left(-\dfrac{3}{4}-\dfrac{N}{2}+n\right)a^{\dagger}_{n}a_{n}\right.$
$\displaystyle\left.+\left(-\dfrac{1}{4}-\dfrac{N}{2}+n\right)b^{\dagger}_{n}b_{n}\right].$
(13)

The dielectric polarization is given by ;

$\displaystyle
P(\varepsilon)=-\dfrac{1}{a(N-1/2)}\sum_{n=1}^{N}\dfrac{\partial
E_{n}}{\partial\varepsilon},$ (14)

where $E_{n}$ with $n=1,2,\cdots,N$ are the energy eigenvalues of (12) below
$E_{f}=0$, and $a(N-1/2)$ is the length of the chain. Note that unlike the
electric polarization $P_{0}$ determined from the bulk wave function, the
dielectric polarization $P(\varepsilon)$ is defined for a chain with a finite
length $N$, and is dependent on $N$. The results for two cases of $t_{1}=0.5$,
$t_{2}=1$ and $t_{1}=2$, $t_{2}=1$ are shown in Fig. 1(d), with the band
structures given in Fig. 1(c).

The result shows that the dielectric polarization takes an anomalously large
value of $P(\varepsilon)=\pm e/2$ even for a very weak electric field when
$t_{1}<t_{2}$. This behavior is dramatically different from that in the other
case $t_{1}>t_{2}$, where $P(\varepsilon)$ linearly increases with a small
slope. We note that in the previous work Combes  _et al._ (2016) the
polarization of the Rice-Mele model, which is the SSH model with a staggered
potential added, is calculated with the same formalism, and its relation to
the Zak phase is discussed. Nonetheless, an abrupt change of the polarization
as a function of the electric field in the SSH model, which we find in this
paper, is not studied in Ref. Combes  _et al._ , 2016. Differences between the
present work and Ref. Combes  _et al._ , 2016 will be discussed in Sec. VI.

The $\pi$ Zak phase indicates presence of topological edge states, and we
attribute this anomalous dielectric response to these edge states. Because of
the chiral symmetry, the edge state at each end of the chain is at the zero
energy, when hybridization between the edge states at the two ends is
neglected. These edge states are half filled. Thus, at an either edge, there
appears an edge charge equal to $e/2$ as shown in Fig. 2(a). We attribute the
anomalous dielectric response for the insulating phase with the $\pi$ Zak
phase to these charges at the edge states. When a small electric field
$\varepsilon$ is applied, the two edge states will have different energies as
shown in Fig. 2(b), giving rise to an immediate transfer of charges from one
end to the other. This concept is verified using an effective two-site model
(Fig. 2(c)). This model consists of two sites, each of which corresponds to
one of the two edge sites of the SSH model. The parameter $t^{\ast}$
represents an effective hopping amplitude between the two sites, L and R,
representing the left and right ends of the chain, respectively. Let
$\ket{R}=(1,0)^{T}\left[\ket{L}=(0,1)^{T}]\right)$ denote the state where the
electron is at the R [L] site. Then, the effective Hamiltonian with the
electric field takes the form

$\displaystyle H^{\ast}=\left(\begin{array}[]{cc}\varepsilon
ea(2N-1)/4&-t^{\ast}\\\ -t^{\ast}&-\varepsilon ea(2N-1)/4\end{array}\right)\
.$ (17)

We put one electron into the model. The effective dielectric polarization is
given by

$\displaystyle P^{\ast}(\varepsilon)$
$\displaystyle=\frac{e}{2}\left[\braket{R}{-}\braket{-}{R}-\braket{L}{-}\braket{-}{L}\right]$
$\displaystyle=-\text{sgn}(\varepsilon)\dfrac{e}{2\sqrt{(\frac{2t^{\ast}}{(N-1/2)ae\varepsilon})^{2}+1}}\
,$ (18)

where $\ket{-}$ is the eigenstate of Eq. (17) with the negative eigenvalue.
The hopping parameter $t^{\ast}$ is determined by fitting Eq. (18) for
$P^{\ast}(\varepsilon)$ to our numerical result $P(\varepsilon)$ within the
region of a weak electric field. The comparisons between (14) and (17) are
shown in Fig. 2(d). We can see the two graphs fit well even for a strong
electric field and evolve in the same way as $t_{2}/t_{1}$ changes. This
result confirms that the anomalous dielectric response is caused by the
topological edge states. We also note that this response around
$\varepsilon\sim 0$ becomes sharper for a longer chain, because in Eq. (18)
$t^{*}$ becomes smaller and $N-1/2$ becomes larger for a longer chain.

Figure 2: Dielectric response of the SSH model (Model I). (a)(b) Schematic
picture of the anomalous dielectric polarization for $t_{2}\gg t_{1}$. In such
a case, the levels in each dimer formed by the strong bonds are split by the
hopping $t_{2}$, and there appears an isolated state at each end of the chain.
The orange circles represent electrons and red semicircle represents a half of
an electron. In (a), the electric field $\varepsilon$ is zero, and the states
at each end is half-filled. In (b), the electric field is applied, and there
appears a difference between the levels at the two ends, and a half of an
electron charge is carried from one end to the other. (c) Schematic picture of
the effective two-site model. (d-1)-(d-4) Dielectric polarization for the SSH
model (red lines) and effective two-site model (blue lines). The hopping
amplitudes are taken as $t_{1}=1$, (d-1) $t_{2}=1.2$, (d-2) $t_{2}=1.3$, (d-3)
$t_{2}=1.4$, and (d-4) $t_{2}=1.5$.

##  III three-dimensional model (Model II)

The discussion in the previous section on a one-dimensional model can be
extended to higher-dimesions. In this section, we show that a jump of the
dielectric response also appears in a simple three-dimensional tight-binding
model on the diamond lattice with anisotropic hopping amplitudes, which we
call Model II. In addition, we show other examples showing the similar effect,
two two-dimensional models and one three-dimensional model in Appendices A and
C, respectively. In three-dimensional systems, the polarization, which is an
electric dipole moment per unit volume, has the unit of charge divided by an
area.

The diamond lattice and its slab with (111) surfaces are shown in Fig. 3 (a)
and (b). The lattice consists of two sublattices ${\rm A}$ and ${\rm B}$. Let
$a$ denote the distance between two neighboring ${\rm A}$ sites. The nearest
bond vectors $\bm{d}_{i=1,2,3,4}$ from an A site to the adjacent B sites are
given by

$\displaystyle\bm{d}_{1}$
$\displaystyle=\frac{\sqrt{6}a}{4}(0,0,1),\bm{d}_{2}=\frac{\sqrt{6}a}{4}(-\frac{2\sqrt{2}}{3},0,-\frac{1}{3}),$
$\displaystyle\bm{d}_{3}$
$\displaystyle=\frac{\sqrt{6}a}{4}(\frac{\sqrt{2}}{3},\frac{{\sqrt{6}}}{3},-\frac{1}{3}),\
\bm{d}_{4}=\frac{\sqrt{6}a}{4}(\frac{\sqrt{2}}{3},-\frac{{\sqrt{6}}}{3},-\frac{1}{3}),$

as shown in Fig. 3.

Figure 3: Model II on the diamond lattice. (a) Diamond lattice. The lattice
constant is $a$ and the nearest bond vectors from site A to B are
$\bm{d}_{i=1,2,3,4}$. (b) Slab of the model on the diamond lattice with a
(111) surface. The red dotted line shows the unit cell. (c) Schematic figure
of the model with anisotropic hopping amplitude (Model II). The hopping
amplitude $-t_{2}$ on the bond along (111) direction is different from that
along other directions, $-t_{1}$. (e) Dielectric polarization of a slab of the
Model IV on a diamond lattice with (111) surfaces. The polarization
$P(\varepsilon)$ rapidly changes toward $\pm e/(2S)$ in the vicinity of
$\varepsilon=0$. $S$ is the area of the surface unit cell. Here we set $N=10$
and the number of sites in the unit cell is $20$. The hopping amplitudes are
$(t_{1},t_{2})=(1,4)$.

###  III.1 Zak phase and dielectric polarization for three-dimensional systems

First, we review the Zak phase and its connection with electric polarization
in three-dimensional systems. When we take the integration path along the
reciprocal lattice vector $\bm{G}_{\perp}$, the Zak phase is defined as

$\displaystyle\theta(\bm{k}_{\parallel})=-i\sum_{n}^{\text{occ.}}\int_{0}^{|\bm{G}_{\perp}|}dk_{\perp}\bra{u_{n}(\bm{k})}\dfrac{\partial}{\partial
k_{\perp}}\ket{u_{n}(\bm{k})},$ (20)

where $\bm{k}=k_{\perp}\bm{n}+\bm{k}_{\parallel}$,
$\bm{n}=\bm{G}_{\perp}/|\bm{G}_{\perp}|$ and
${\bm{k}}_{\parallel}\perp{\bm{G}}_{\perp}$. The Zak phase is related with the
polarization for the systems with surface perpendicular to $\bm{G}_{\perp}$.
This phase is quantized as 0 or $\pi$ modulo $2\pi$ under both the TR and SI
symmetries Kariyado and Hatsugai (2013). In the region of the surface
Brillouin zone where the Zak phase is $\pi$, there exist degenerated surface
states at zero energy when the chiral symmetry is present Ryu and Hatsugai
(2002). The electric polarization $P_{0}$ along ${\bm{n}}$ can be obtained
from the Zak phase as shown in the following equation;

$\displaystyle
P_{0}=-e\int_{\bm{k}_{\parallel}\in\text{2DBZ}}\dfrac{d^{2}\bm{k}_{\parallel}}{(2\pi)^{3}}\theta(\bm{k}_{\parallel})\
({\rm mod}\ e/S),$ (21)

where $S$ is the area of surface unit cell, and the integral is over the two-
dimensional Brillouin zone for $\bm{k}_{\parallel}$.

###  III.2 Model II: Anisotropic tight-binding model on the diamond lattice

We consider a model as shown in Fig. 3 (c). The system has anisotropic hopping
amplitudes; $-t_{2}$ is the hopping amplitude along the bonds along the (111)
direction, and $-t_{1}$ is that along the other directions. We set $t_{1}$ and
$t_{2}$ positive. The Hamiltonian is given by

$\displaystyle H_{0}=-\sum_{\langle ij\rangle}t_{ij}c^{\dagger}_{i}c_{j}\ ,$
(22)

where $t_{ij}$ is the hopping amplitude from site $i$ to $j$. First, we
consider a bulk system to calculate the Zak phase. The bulk Bloch Hamiltonian
is given by

$\displaystyle H^{\text{bulk}}(\bm{k})$
$\displaystyle=\left(\begin{array}[]{cc}0&R(\bm{k})\\\
R^{\ast}(\bm{k})&0\end{array}\right),$ (25) $\displaystyle R(\bm{k})$
$\displaystyle=-t_{1}e^{-i{\bm{k}}\cdot{\bm{d}_{1}}}-t_{2}(e^{-{\bm{k}}\cdot{\bm{d}_{2}}}+e^{-i{\bm{k}}\cdot{\bm{d}_{3}}}+e^{-{\bm{k}}\cdot{\bm{d}_{4}}}),$
(26)

where $\bm{k}=(k_{x},k_{y},k_{z})$. The band gap is open when $t_{2}>3t_{1}$
Takahashi and Murakami (2013), and it is around $E=0$. We set the Fermi energy
to be $E_{f}=0$. Now we calculate the Zak phase $\theta(\bm{k}_{\parallel})$
along ${\bm{G}}_{\perp}=(0,0,\sqrt{6}\pi/a)$ for
$\bm{k}_{\parallel}=(k_{x},k_{y})$. Because of the time-reversal and inversion
symmetries, the Zak phase is quantized as $0$ or $\pi$ modulo $2\pi$. As an
example, we calculate $\theta_{z}(\bm{k}_{\parallel}=0)$ from the products of
the parity eigenstates at the $\Gamma$ ($\bm{k}=0$) and the $L$
($\bm{k}=(0,0,\sqrt{6}\pi/(2a)$) points in the following way Hughes  _et al._
(2011)

$\displaystyle
e^{i\theta_{z}(\bm{k}_{\parallel}=0)}=\prod_{n}^{\text{occ.}}\xi_{n}(\Gamma)\xi_{n}(L)=-1,$
(27)

because these two points are time-reversal invariant momenta on the
$\bm{k}_{\parallel}=0$ line. In this calculation, the unit cell is taken to be
a pair of the A and B sites displaced by ${\bm{d}_{2}}$, as determined from
the crystal termination in this case, and the inversion center is in the
middle of the unit cell. Therefore this equation indicates
$\theta_{z}(\bm{k}_{\parallel})\equiv\pi$ for any $\bm{k}_{\parallel}$ because
the Zak phase is constant in the two-dimensional Brillouin zone for
${\bm{k}}_{\parallel}$ when the band gap is open.

The Hamiltonian of the slab with the (111) surfaces (Figs. 3 (b) and (c)) is
given by

$\displaystyle H_{0}(\bm{k}_{\parallel})$
$\displaystyle=-\left[t_{2}\sum_{n=1}^{N-1}e^{i\tilde{\bm{d}_{1}}\cdot{\bm{k}}_{\parallel}}a^{\dagger}_{\bm{k}_{\parallel},n+1}b_{\bm{k}_{\parallel},n}\right.$
$\displaystyle\ \ \
+\left.t_{1}\sum_{n=1}^{N}\left(e^{i\tilde{\bm{d}_{2}}\cdot{\bm{k}_{\parallel}}}+e^{i\tilde{\bm{d}_{3}}\cdot{\bm{k}_{\parallel}}}+e^{i\tilde{\bm{d}_{4}}\cdot{\bm{k}_{\parallel}}}\right)a^{\dagger}_{\bm{k}_{\parallel},n}b_{\bm{k}_{\parallel},n}\right]$
$\displaystyle\ \ \ +h.c.,$ (28)

where $a_{\bm{k}_{\parallel},n}^{\dagger}(b_{\bm{k}_{\parallel},n}^{\dagger})$
is the creation operator of an electron with Bloch wave vector
$\bm{k}_{\parallel}$ at the $n$th A(B) site in the unit cell and
$\tilde{\bm{d}_{i}}=(d_{ix},d_{iy})$. When we apply an electric field of
strength $\varepsilon$ along the $z$ axis, the Hamiltonian is given by

$\displaystyle H(\bm{k}_{\parallel})$
$\displaystyle=H_{0}(\bm{k}_{\parallel})+H_{\varepsilon}(\bm{k}_{\parallel}),$
(29) $\displaystyle H_{\varepsilon}(\bm{k}_{\parallel})$
$\displaystyle=\dfrac{\sqrt{6}}{4}ae\varepsilon\sum_{n=1}^{N}\left[\left(\dfrac{4}{3}n-\dfrac{2}{3}N-\dfrac{5}{6}\right)a^{\dagger}_{\bm{k}_{\parallel},n}a_{\bm{k}_{\parallel},n}\right.$
$\displaystyle\ \ \ \ \ \ \
\left.+\left(\dfrac{4}{3}n-\dfrac{2}{3}N-\dfrac{1}{2}\right)b^{\dagger}_{\bm{k}_{\parallel},n}b_{\bm{k}_{\parallel},n}\right],$
(30)

The dielectric polarization is given by

$\displaystyle
P(\varepsilon)=-\dfrac{1}{a(4N/3-1)}\sum_{n=1}^{N}\int_{\bm{k}_{\parallel}\in\text{2DBZ}}\dfrac{d^{2}\bm{k}_{\parallel}}{(2\pi)^{2}}\dfrac{\partial
E_{n}(\bm{k}_{\parallel})}{\partial\varepsilon},$ (31)

where the integral is taken over the two-dimensional Brillouin for
${\bm{k}}_{\parallel}$, $E_{n}(\bm{k}_{\parallel})$ is the $n$th lowest
eigenvalue of (29), and $a(4N/3-1)$ is the thickness of the system. The
calculation results are shown in Fig. 3(e). The dielectric polarization takes
a value $P(\varepsilon)\sim\pm e/2S$ even for a very weak electric field,
where $S$ is an area of the surface unit cell. This indicates that a half of
an electron per surface unit cell is accumulated on each surface whenever
$t_{2}>3t_{1}$, and the response of this surface charge gives rise to the
anomalous behavior of the polarization. Thus, the anomalous dielectric
response occurs in this system when it is insulating.

##  IV Remarks on various effects on the anomalous dielectric responses

We have discussed several models in one, and three dimensions with the $\pi$
Zak phase, and demonstrated that there is an abrupt jump in the dielectric
response. Because they are idealized models, we make several remarks in order
to see how the anomalous dielectric responses in general systems including in
real materials.

###  IV.1 Effect of surface-state dispersion

Here, we discuss flatness of surface-state dispersions and its relation to
anomalous dielectric response. In the models in Secs. II and III and in
Appendices A and C, the edge/surface states are exactly at the zero energy and
are dispersionless, which is due to chiral symmetry and particle-hole symmetry
in these models. Under the chiral symmetry the energy band becomes symmetric
with respect to $E\leftrightarrow-E$, which means that the energy eigenvalues
satisfy $E_{n,\mathbf{k}}=-E_{m,\mathbf{k}}$. Under the particle-hole symmetry
the energy eigenvalues satisfy $E_{n,\mathbf{k}}=-E_{m,\mathbf{-k}}$. Because
of the inversion symmetry of the sysem, both the chiral symmetry and the
particle-hole symmetry are preserved here, which makes the topological
boundary states to be exactly at the zero energy. Therefore, as we discussed
so far, their response under a small electric field becomes anomalous, since
even for a small electric field, the energies of all the states on one surface
become positive, while those on the other surface become negative, and all the
surface charges are transferred at the same time. In contrast, general systems
do not have chiral or particle-hole symmetry, and the edge/surface states have
dispersions. This edge-/surface-state dispersion will make the dielectric
response less abrupt, because upon the application of the electric field, not
all the charges in these edge/surface states are carried from one end to the
other at the same time. Thus depending on the dispersion of the edge/surface
states, the response becomes less anomalous. We show several examples for
dispersions of surface states in real materials with the $\pi$ Zak phase in
the next section.

Thus far we have seen that the polarization jumps to $e/(2S)$ by a weak
electric field when the surface-state dispersion is flat. In this discussion
we do not take into account electron correlation. If we take electron-electron
interaction into account, the dielectric response may change. When a weak
electric field is applied, the electrons at the surface may segregate and form
a superstructure, so that the correlation energy may overcome the potential
difference of the surface states at the both ends of the slab. Such kinds of
formation of charge segregation will make the dielectric response less abrupt.
We do not argue this phenomena any further here, because it should be
sensitive to various details of the system such as the form and the strength
of the electron-electron interaction, dispersion of the surface states, and so
on.

###  IV.2 Effect of the depolarization field

So far we have ignored an influence of the depolarization field. When an
electric field is applied and a polarization is induced, the induced
polarization generates a depolarization field. This depolarization field
depends on the geometry of the system. Particularly when the shape of the
system is an ellipsoid or a large thin slab, the influence can be easily
calculated because the depolarization field is uniform, represented by the
depolarization factor.

The electric field $\varepsilon$ used in the discussion so far represents the
total electric field felt by the electrons. It is a sum of the external
electric field $\varepsilon_{\text{ext}}$ and the depolarization field
$\varepsilon_{\text{dep}}$. The depolarization field
$\varepsilon_{\text{dep}}$ is expressed as
$\varepsilon_{\text{dep}}=-LP/\epsilon_{0}$, where $P$ is the dielectric
polarization, $L$ ($0\leq L\leq 1$) is the depolarization factor, and
$\epsilon_{0}$ is the permittivity of vacuum. Therefore, the total electric
field $\varepsilon$ takes the form;

$\displaystyle\varepsilon=\varepsilon_{\text{ext}}+\varepsilon_{\text{dep}}=\varepsilon_{\text{ext}}-LP(\varepsilon)/\epsilon_{0},$
(32)

Here, the dependence of $P(\varepsilon)$ on $\varepsilon$ has been discussed
so far in this present paper. Now, we focus on the case with the Zak phase
equal to $\pi$ for any value of $\bm{k}_{\|}$, and for simplicity we assume
that the topological surface states have flat dispersion. Then, the response
of the polarization $P$ to the electric field $\varepsilon$ is given by

$\displaystyle P=\frac{e}{2S}\text{sgn}(\varepsilon),$ (33)

where $S$ is an area of the surface unit cell, as shown in Fig. 4(a). From
this equation, and Eq. (32), the dielectric response $P$ versus the external
electric field $\varepsilon_{\text{ext}}$ is shown in Fig. 4(b). Thus, the
slope around $\varepsilon_{\text{ext}}\sim 0$ is $\epsilon_{0}/L$.

The present system corresponds to a dielectric with the dielectric constant
$\epsilon=\infty$, because a nonzero polarization arises in response to an
infinitesimal electric field. As a comparison, let us consider a uniform
dielectric with a dielectric constant $\epsilon$. Then by including the effect
of the depolarization field, with the depolarization factor $L$, its
dielectric response is given by

$\displaystyle P=(\epsilon-\epsilon_{0})\varepsilon,\
\varepsilon=\varepsilon_{\text{ext}}-\frac{LP}{\epsilon_{0}}$ (34)

Thus the ratio of the polarization to the external electric field
$\varepsilon_{\text{ext}}$, is given by
$P/\varepsilon_{\text{ext}}=\frac{\epsilon_{0}(\epsilon-\epsilon_{0})}{L(\epsilon-\epsilon_{0})+\epsilon_{0}}$.
By comparison with our result in Fig. 4(b), the present system indeed
corresponds to the dielectric with $\epsilon=\infty$.

Figure 4: Effect of depolarization field for an ideal case, where the Zak
phase equal to $\pi$ and the dispersion of the topological surface states is
flat. (a) Polarization $P$ versus the electric field $\varepsilon$. (b)
Polarization $P$ versus the external electric field $\varepsilon_{\rm ext}$.

###  IV.3 Dynamical Dielectric Response in an insulator with the $\pi$ Zak
phase

So far we have discussed how anomalous dielectric response appears in a static
case. We here discuss dynamical response; we suppose there is no electric
field when the time $T$ is in $T<0$, and then we apply a constant electric
field when $T\geq 0$. To study dynamical response we study the following
model, representing the two states at the two ends of the one-dimensional
system similar to Eq. (17):

$H=\left(\begin{array}[]{cc}V(T)&-t^{*}\\\ -t^{*}&-V(T)\end{array}\right).$
(35)

Here, we set $\ket{R}=(1,0)^{T}\left[\ket{L}=(0,1)^{T}]\right)$ to denote the
state at the right (left) end, similarly to Eq. (17), as shown in Fig. 5(a).
In addition, $-t^{*}$ ($t^{*}>0$) is the hopping between the two states at the
two ends, and $\pm V(T)$ are the energies of these two states, due to the
external electric field. Here $t^{*}$ is an exponentially decaying function of
the chain length. When $T<0$, the external electric field is zero and
$V(T)=0$, and when $T>0$, the external electric field is applied and
$V(T)=\bar{\varepsilon}$, where $\varepsilon$ is a constant proportional to
the electric field (Fig. 5(b)). We set the Fermi energy $E_{F}=0$ so that only
one state is occupied. At the time $T=0$, the occupied state is expressed as

$\psi(T=0)=\frac{1}{\sqrt{2}}\left(\begin{array}[]{c}1\\\
1\end{array}\right).$ (36)

The evolution of this state is described by the time-dependent Schrödinger
equation

$i\hbar\frac{\partial}{\partial T}\psi(T)=H\psi(T),$ (37)

Through a straightforward calculation, we obtain

$\psi(T)=\frac{1}{2\sqrt{2}A}\left(\begin{array}[]{c}e^{-iAT/\hbar}(A-t^{*}+\bar{\varepsilon})+e^{iAT/\hbar}(A+t^{*}-\bar{\varepsilon})\\\
e^{-iAT/\hbar}(A-t^{*}-\bar{\varepsilon})+e^{iAT/\hbar}(A+t^{*}+\bar{\varepsilon})\end{array}\right),$
(38)

for $T>0$, where $A=\sqrt{\bar{\varepsilon}^{2}+t^{*2}}$. Thus the probablity
for the electron to reside at the left (right) end of the chain is calculated
as

$n_{R(L)}(T)=\frac{1}{2}\left(1\mp\frac{2t^{*}\bar{\varepsilon}}{A^{2}}\sin^{2}\frac{2AT}{\hbar}\right)$
(39)

Therefore, the dynamical polarization $P^{*}_{dyn.}$ is given by

$P^{*}_{dyn.}(T)=-\frac{e}{2}(n_{R}-n_{L})=\frac{e}{2}\frac{2t^{*}\bar{\varepsilon}}{A^{2}}\sin^{2}\frac{2AT}{\hbar}$
(40)

On the other hand, the static polarization is also calculated in response to
the static electric field, similarly to Sec. IIB. By putting
$V\equiv\bar{\varepsilon}$, the occupied state is given by

$\Psi=\frac{1}{\sqrt{2A(A+\bar{\epsilon})}}\left(\begin{array}[]{c}t^{*}\\\
A+\bar{\epsilon}\end{array}\right),$ (41)

Thus the probablity for the electron to reside at the left (right) end of the
chain is calculated as

$n_{R(L)}=\frac{1}{2}\left(1\mp\frac{\bar{\varepsilon}}{A}\right)$ (42)

Therefore, the static response of the polarization $P^{*}_{sta.}$ is given by

$P^{*}_{sta.}=-\frac{e}{2}(n_{R}-n_{L})=\frac{e}{2}\frac{\bar{\varepsilon}}{A}$
(43)

The dynamical and static polarizations are shown in Fig. 5, for
$t^{*}=10^{-2}$ and $\bar{\varepsilon}=10^{-3},10^{-2},10^{-1}$.

The dynamical polarization $P^{*}_{dyn.}$ undergoes Rabi oscillations. When
$\bar{\varepsilon}=10^{-3}\ll t^{*}$, the polarization oscillates with its
time average almost equal to the static value $\langle
P_{dyn.}^{*}(T)\rangle\sim
P^{*}_{sta.}\sim\frac{e}{2}\frac{\bar{\varepsilon}}{t^{*}}$ On the other hand,
when $\bar{\varepsilon}=10^{-1}\gg t^{*}$, the dynamical polarization with the
Rabi oscillation is much smaller than the static value. This Rabi oscillation
appears due to the interference between two eigenstates. Therefore, when we
include dissipation to other degrees of freedom such as phonons, the state
showing the Rabi oscillation will gradually relax to the lowest-energy state
of the static Hamiltonian, and this oscillating polarization $P^{*}_{dyn.}$
will gradually approach the value of the static one $P^{*}_{sta.}$.
Nonetheless, in the present context, the overlap between the two states at the
two ends of the chain is small, and the dissipation from the higher-energy
state to the ground state is very slow. Thus it may take a long time to
approach the static value of the polarization, which shows anomalous behaviors
proposed in the previous Sections.

Figure 5: Dynamics of polarization under the electric field. (a) Schematic
figure of the setup. (b) External electric field used in the calculation,
which is turned on at $T=0$. (c)-(e) Time dependence of the dynamical
polarization wih $t^{*}=10^{-2}$ and $\bar{\varepsilon}=$ (c) $10^{-3}$, (d)
$10^{-2}$, and (e) $10^{-1}$. The solid lines show the dynamic response
$P^{*}_{dyn.}(T)$ and the broken lines show the static response
$P^{*}_{sta.}$.

##  V Proposals for real materials

In this section, we propose some materials whose Zak phase has a quantized
value of $\pi$ along some direction. In such cases, there exist topological
in-gap states, similar to the in-gap states at the two ends of the SSH model
within the topological regime. To search for candidate materials, we point out
that the $\pi$ Zak phase typically means existence of dangling bonds. One can
see it from the simple example of the SSH model. The $\pi$ Zak phase and
resulting in-gap states appear when the bond within each dimer is cut at the
end of the system. This is nothing but a dangling bond. It is also similar in
two and three dimensions. The $\pi$ value of the Zak phase is realized when
the Wannier orbitals are “cut into halves” at the surfaces. In covalent
crystals, where the Wannier orbitals are in the middle of the bonds, this
typically leads to dangling bonds. On the other hand, if the Wannier orbitals
are at the atomic sites, one cannot cut the crystal at the center of the
Wannier orbitals. Thus, realization of a surface supporting the in-gap surface
states with $\pi$ Zak phase is not straightforward from the viewpoint of
stability of materials.

Here we propose three classes of materials, which circumvent this obstacle
against realizing the in-gap surface states due to the $\pi$ Zak phase. One is
a cubic semiconductor with the (111) surface, such as diamond and silicon. In
this case, the surface has in-gap surface states due to the $\pi$ Zak phase,
if the surface reconstruction is absent. In reality, dangling bonds exist at
the surface, leading to surface reconstructions. Nevertheless, as we show
later, when the surface reconstruction leads to the formation of surface
superstructure with enlarged surface unit cell by $M$ times, where $M$ is an
integer, the Zak phase becomes $M\pi$ modulo $2\pi$. Therefore, when $M$ is
odd, the Zak phase remains $\pi$ and the topological in-gap states survive.

Another example is a carbon nanotube. In this case, the covalent bond of the
$sp_{2}$ hybrid orbital gives rise to the $\pi$ Zak phase. In nanotubes,
existence of a gapless state and symmetry of the system change depending on
the chirality $(m,n)$. Here, we focus on a zigzag nanotube with chirality
$(m,0)$ with $m$ being an odd number other than a multiple of 3, so that it is
insulating and the Zak phase is $\pi$.

The other material is Sc2C, one of the topological electrides Hirayama  _et
al._ (2018). It is a layered material, and the Wannier orbitals reside between
the atomic layers; therefore, one can create the surface by cutting the system
between the atomic layers, leading to the $\pi$ Zak phase.

###  V.1 Diamond and silicon

In this section, we show that diamond and silicon have the $\pi$ Zak phase
along the [111] direction over the entire two-dimensional Brillouin zone along
the (111) surface, similarly to Sec. III. In fact, this was first pointed out
in Ref. Vanderbilt and King-Smith (1993). In Ref. Vanderbilt and King-Smith
(1993), various cubic semiconductors are studied, and it was found that their
(111) and ($\bar{1}\bar{1}\bar{1}$) surfaces have one half of an electron per
1$\times$ 1 surface unit cell (excluding the spin degree of freedom),
corresponding to the $\pi$ Zak phase.

In the following, we show how the topological surface states in diamond and
silicon appear, based on the generalized gradient approximation (GGA) of the
density functional theory (DFT) and discuss physical origins for the
topological nature of the bands. The computational condition is shown in
Appendix D. Their band structures are shown in Fig. 6. Although diamond has a
larger gap than silicon, both have qualitatively the same band structure. The
valence and conduction bands near the Fermi level originate from $sp_{3}$
hybrid orbitals. The energy gap originates from covalent bonds of the $sp_{3}$
hybrid orbitals.

Table 1: Zak phase of diamond and silicon. The Zak phase is $\pi$ along [111] and along [11$\bar{1}$] $\bm{n}$ | Zak phase |  number of dangling bond
---|---|---
[111] | $\pi$ | $1$
[11$\bar{1}$] | $\pi$ | $3$
[001] | $0$ | $2$
[110] | $0$ | $2$
Figure 6: Band structures of (a) diamond and (b) silicon. Figure 7: Band
structures for a slab of diamond along (111) plane. The unit cell of the slab
contains 60 atoms. (a) Slab with periodic boundary condition. The system has
no surface, and therefore surface states do not appear. (b) Slab with (111)
surfaces. There is one dangling bond per surface unit cell, and there are two
degenerated surface states in the gap. The dispersion of the surface states is
almost flat. (c) Slab with (11$\overline{\rm 1}$) surfaces. There are three
dangling bonds per surface unit cell, and there are six surface states in the
gap. They are doubly degenerate. Dispersion of edge states in this direction
is large compared with that in the (111) direction in (b). Figure 8: Changes
in the band structure of diamond with increasing the lattice constant. We set
the lattice constant to be the original one times (a) $1.0$, (b) $1.2$, (c)
$1.3$ and (d) 1.4. The gap closing and band inversion occurs at $\Gamma$ in
(c).

As we explain here, the Zak phase is related with dangling bonds at the
surface. The Zak phases of these materials along [111], [11$\overline{\rm
1}$], [001] and [110] directions are calculated by ab intio calculation. Table
1 shows the calculation results of the Zak phase for diamond and silicon. It
shows that the Zak phases along [111] and along [11$\overline{\rm 1}$]
directions are $\pi$ modulo $2\pi$, while those along [001] and along [110]
directions are $0$ modulo $2\pi$. In these covalent crystals, the value of the
Zak phase $\theta=\pi$ ($\theta=0$) coincides with an odd number (an even
number) of dangling bonds per surface unit cell on each surface. This is
because one dangling bond per surface unit cell supports one surface mode,
which adds $\pi$ to the Zak phase. In the following, we confirm this for the
diamond lattice. The energy bands of diamond slabs with (111) and
(11$\overline{\rm 1}$) surfaces are shown in Fig. 7. Figure 7 (a) shows the
band structure when a periodic boundary condition is imposed. Here, dangling
bonds do not exist and there is no surface state. On the other hand, Fig. 7
(b) shows the band structure for a slab with (111) surfaces. There is a
dangling bond per unit cell, supporting one surface mode on each surface, Thus
the two surfaces of the slab support two surface modes, as Fig. 7 (b) shows.
On the other hand, for the (11$\bar{1}$) surfaces (see Fig. 7 (c)), there are
three dangling bonds per unit cell on each surface, and as a result, there
exist six surface modes in the gap.

Here, we discuss flatness of surface-state dispersions and its relation to
anomalous dielectric response. As we discussed in Sec. VA, the models
discussed in Secs. II-IV have dispersionless topological surface states due to
chiral symmetry, and an anomalous dielectric response is expected. In
contrast, real materials do not have chiral symmetry. As a result, in silicon
and diamond, the surface states protected by the $\pi$ Zak phase have
dispersions. In such cases, the two surface modes on the two surfaces are half
filled and are metallic. Under a small electric field, the surface remains
metallic; in other words, there is still an overlap in energy between the
surface states on one surface and those on the other surface. Therefore, an
abrupt charge transfer does not take place, and the polarization does not
jump. The strength of the electric field has a threshold fro a charge transfer
from one surface to the other, leading to the uprise of the polarization. In
silicon, we see from Figs. 7 (b)(c) that the band structure of the edge states
are almost flat in the slab with (111) surfaces, while far from flat with
(11$\bar{1}$) surfaces. Therefore, the anomalous dielectric response is
expected in the slab with (111) surfaces.

Since any system in the atomic limit is always in a topologically trivial
phase, topological phase transition should occur when the interatomic distance
in diamond is extended from the topological insulator phase with a $\pi$ Zak
phase. We show the change of band structure in diamond with various lattice
constant in Fig. 8. We set the lattice constant to be the original one times
$1.0,1.2,1.3$ and $1.4$. These figures show that the band inversion occurs
when the lattice constant is extended by $1.3$ times. If we further enlarge
the lattice constant, it becomes a filling-enforced semimetal, because
threefold degenerate states at the $\Gamma$ point are at the Fermi energy, and
they bridge between the valence and conduction bands.

###  V.2 Surface reconstruction

Thus far, we have ignored the surface reconstruction. Since the electron
energy in a dangling bond is higher than that in the valence band, the atoms
on the surface are generally displaced so as to reduce the number of dangling
bonds. It is called surface reconstruction. The effect of surface
reconstruction has been discussed in Sec. IIID in Ref. Vanderbilt and King-
Smith (1993); when the surface unit cell is multiplied by $M$ times by surface
reconstruction, where $M$ is an integer, the polarization is defined modulo
$e/(MS)$ instead of $e/S$. Here we discuss the Zak phase and resulting
topological surface states after surface reconstruction. It is well known that
the silicon (111) surface forms $7\times 7$ structure Schlier and Farnsworth
(1959). In this structure, $49$ original surface unit cells form one new unit
cell. The Zak phase is multiplied by 49 times, because the surface Brillouin
zone is folded down by a factor of $1/49$. We note that the integrand in the
formula of the Zak phase remains the same because the Zak phase is a bulk
quantity independent of the details of the surface. Thus as a result the Zak
phase changes from $\theta=\pi$ to $\theta=49\pi\equiv\pi$ modulo $2\pi$, and
the bulk polarization for a weak electric field is expected to be $e/(2S\times
49)=e/(98S)$ mod $e/(49S)$ per original surface unit cell. Here $49S$ is the
area of the surface superstructure. In general, when $M$ unit cells form a
superstructure, the Zak phase $\theta$ is given by

$\displaystyle\theta$ $\displaystyle\equiv
M\pi\equiv\begin{cases}\pi~{}(\text{mod}\ 2\pi)&(\text{$M=$ odd}),\\\
0~{}(\text{mod}\ 2\pi)&(\text{$M=$ even}).\end{cases}$ (44)

Here, $\theta$ does not depend on the wave vector $\bm{k}_{\parallel}$ because
we consider only insulators. This result agrees with a counting argument of
dangling bonds. When $M=$ even, all dangling bonds form pairs with covalent
bonds, and there are no surface states in the gap, and the Zak phase becomes
zero. Thus, the dielectric response is normal. Meanwhile when $M$ is odd, the
electric polarization converges to $P=e/(2MS)~{}\text{mod}\ e/(MS)$, where
$MS$ is a surface area of the superstructure.

The previous studies show the correspondence between the survival of surface
states and the parity of the number of unit cells constituting the
superstructure in silicon Takagi and Okada (2008); Smeu  _et al._ (2012). In
particular in Ref. Smeu  _et al._ , 2012, the density of states of the silicon
(111) surfaces with surface reconstuctions $2\times 2$,
$\sqrt{3}\times\sqrt{3}$, $5\times 5$ and $7\times 7$ is shown. Among these
four cases, only for the $2\times 2$ surface, where $M=2\times 2=4$ is even,
the surface states have a gap, whereas in the other three cases with odd $M$,
the surface is metallic, suggesting existence of half-filled midgap surface
states. These results are in perfect agreement with our results in this paper.

Figure 9: Topological edge states in a zigzag carbon nanotube. (a) Crystal
structure of the zigzag carbon nanotube with the (7,0) chirality. (b)
Electronic band structure of the (7,0) nanotube. The Fermi energy is at $E=0$.
(c) Eigen energies of the (7,0) nanotube with a finite length, containing 10
unit cells, i.e. 280 carbon atoms. (d) Distribution of the electron density
for the eigenstate at the Fermi level of the (7,0) nanotube with finite
length. It is localized at the edge of the nanotube.

###  V.3 Carbon nanotube

Next, we show the result of a carbon nanotube. While a carbon nanotube of the
armchair type is always semimetallic, that of the zigzag type has a gap when
the chirality is $(3i+1,0)$ or $(3i+2,0)$ $(i=1,2,...)$. As we limit ourselves
on insulating systems, we focus on the $(m,0)$ zigzag nanotube, with $m$ not
an integer multiple of three. As the number of dangling bonds in this nanotube
is $m$, we restrict ourselves to the case with odd $m$, so that the Zak phase
is $\pi$. Figure 9(a) shows a structure of zigzag nanotube with a chirality of
(7,0). The $p_{z}$ orbitals forming the $\pi$ bonds have a small gap around
the Fermi level (Fig. 9(b)). Because the blue region in Fig. 9(a) has about
twice as much transfer as the red region, the Wannier center of the $\pi$
covalent bond are located at the blue region. Therefore, the $p_{z}$ orbital
forming the $\pi$ bonds has no polarization in the unit cell corresponding to
the zigzag edge. On the other hand, since $\sigma$ covalent bonds are cut at
the $(m,0)$ zigzag edge, the Zak phase is $\pi$ when $m$ is an odd integer.
Figure 9(c) shows the eigen energies of a finite nanotube with 10 unit cells
in length. Half-filled edge states exist at the Fermi level. These are
topological states originated from the $\pi$ Zak phase and are localized at
the edges (Fig. 9(d)).

###  V.4 Topological electrides

In our previous work Hirayama  _et al._ (2018), we proposed that the
electrides can be good candidates for topological materials. Among such
topological electrides, we proposed that Sc2C is an insulator with $\pi$ Zak
phase along the [111] direction. It is a layered material, and is a two-
dimensional eletride, where the Wannier orbitals reside between the atomic
layers along the (111) plane. It can be cleaved along the (111) plane, which
cuts the Wannier orbitals into half. This is similar to the SSH model on the
topological regime, and it leads to in-gap surface states, which is a remnant
of the electronic states between the atomic layers. Since existence of
electronic states at interstitial regions is a hallmark of the electrides,
this characteristic of the electrides is the key to realize the surface with
the $\pi$ Zak phase.

The $\pi$ Zak phase in this topological electride shows that the notion of the
Zak phase has wider applicability than counting of dangling bonds at the
surface. In Sec VIA, we have discussed that in covalent crystals, the number
of dangling bonds being even and odd corresponds to the Zak phase 0 or $\pi$,
respectively. On the other hand, in Sc2C, we can make the (111) surface
without appearance of dangling bonds. Even in such crystals we can adopt the
notion of the Zak phase, and indeed we can see appearance of midgap
topological surface states on the (111) surface of Sc2C.

###  V.5 Ab initio calculation for dielectric response

Figure 10: Dielectric response of the PTFE. (a) Structure of the PTFE. (b)
Electronic band structure of the PTFE. The Fermi energy is at $0$. (c)
Covalent bonds of the PTFE. (d) Enegy eigenvalues of [C2F4]8. (e) Dielectric
response of the finite PTFE. $P$ is the polarization. (f) LUMO of [C2F4]8
under the electric field $\varepsilon=9$ (GV/m).

Here we show the ab initio calculation for the dielectric response. First, we
calculate the dielectric response of PTFE [C2F4]n. Figure 10(a) shows the
crystal structure of the PTFE. The one-dimensional band structure of the PTFE
has a large band gap (Fig. 10(b)). Such a large gap is suitable for the
anomalous dielectric response (see Figs. 2(a)(b)). The C atoms are connected
by the $\sigma$-bonds originating from the $sp^{3}$ orbitals (Fig. 10(c)).
Because the two orbitals per each site form covalent bondings, this system is
regarded as a superposition of two SSH models with dimers at different
positions, as shown in Fig 10(c). Therefore, the Zak phase is always $\pi$
regardless of the choice of the unit cell.

Next, we calculate the energy eigenvalues of a finite system of the PTFE (Fig.
10(d)). The topological edge states originating from the $\pi$ Zak phase exist
at the Fermi level ($E=0$ in Fig. 10(d)) . Next we apply the electric field
along the C chain to finite systems and calculate the polarization (Fig
10(e)). We here discuss a condition for the system length to exhibit the
anomalous dielectric response. If the length of the PTFE is too short, the
hybridization between edges become large, which makes the dielectric response
less sharp (Fig. 2(c)). On the other hand, if the length of the PTFE is too
long, the energy loss of the Coulomb interaction causes a gradual charge
transfer in the system to locally satisfy charge neutrality. Therefore, in a
realistic system with a finite bulk gap, the length of the system should be
moderate in order to exhibit the anomalous response. From the above reasons,
we calculate the dielectric response of [C2F4]6, [C2F4]8, and [C2F4]12 in Fig
10(e), where the magnitude of ${\bm{P}}$ is shown. The electric response
becomes flat around $P\sim\pm e$ in all the three cases. It perfectly agrees
with Fig. 2(d), considering the spin degeneracy ($=2$). Moreover, the slope of
the dielectric response becomes sharp for a longer chain, which also agrees
with our theory. We show the lowest unoccupied molecular orbital (LUMO) of
[C2F4]8 under the electric field $\varepsilon=9$ (GV/m) in Fig 10(f), which
corresponds to one of the topological edge states originating from the $\pi$
Zak phase. When the electric field becomes even larger, the electrons in the
bulk region participate in screening the polarization, and the dielectric
response becomes linear again (see Fig. 2(d)).

We also calculate the dielectric response of a (111) diamond slab with
$1\times 1$ surface structure discussed in Sec. V.1. We find that the
polarization does not show an anomalous response but has a linear dependence
on the electric field. We attribute the absence of the anomalous behavior of
the dielectric response to the weak ${\bm{k}}$ dependence of the topological
surface states in diamond (111) surface with $1\times 1$ structure (Fig.
7(b)), as discussed in Sec. IV.1. On the other hand, in the silicon with
realistic $7\times 7$ surface reconstruction, the anomalous dielectri response
is more likely because of the two reasons. One is because the topological
surface bands are almost flat, and the other is because the larger unit cell
leads to the smaller energy increase of Coulomb repulsion per two-dimensional
unit cell associated with the charge transfer, making the bulk region less
likely to participate in the anomalous response. Therefore, nano-sized silicon
with the surface reconstruction would be promising to realize the anomalous
response.

##  VI Bulk physics versus surface physics

In this section, we discuss how the anomalous dielectric response is related
with the bulk polarization, described in the modern theory of polarization. In
the modern theory of polarization Resta (1992); King-Smith and Vanderbilt
(1993); Vanderbilt and King-Smith (1993); Resta (1994), the notion of the bulk
polarization is introduced, and it is formulated in term of the Zak phase, in
the absence of an electric field. A dielectric response to an electric field
is also discussed along the same line Nunes and Vanderbilt (1994); Nunes and
Gonze (2001); Combes  _et al._ (2016). In these works, the dielectric response
is studied as a bulk effect, which is considered to be valid in the
thermodynamic limit.

We need some care in comparing our results with those in the previous works on
bulk dielectric responses Nunes and Vanderbilt (1994); Nunes and Gonze (2001);
Souza  _et al._ (2002); Umari and Pasquarello (2002); Combes  _et al._ (2016).
In particular, the dielectric response is formulated by constructing the
Wannier-Stark ladder perturbatively for a weak electric field in Ref. Combes
_et al._ , 2016, and a similar result on dielectric response is obtained in
Ref. Swiecicki and Sipe, 2014 as well. In Ref. Combes  _et al._ , 2016, this
formalism is applied to the Rice-Mele model; because the SSH model is a
special case of the Rice-Mele model, we can directly compare those results
with ours in this paper. As we explain later in detail, in the results in Ref.
Combes  _et al._ , 2016 for the SSH model with the $\pi$ Zak phase, there is
no divergence in the dielectric response $\chi\equiv\frac{dP}{d\varepsilon}$,
which seems to contradict our result.

This discrepancy is related with an interplay between the bulk and surface
physics. To explain this, we focus on degrees of freedom in the choices of the
crystal terminations and those of the unit cell Vanderbilt and King-Smith
(1993); Chen and Lee (2011). Let us consider 1D systems with inversion
symmetry, such as the SSH model in Sec. II. In this case, as mentioned
earlier, the unit cell should be invariant under inversion symmetry. There are
two choices of the inversion center, displaced by a half of the primitive
translation vector, and correspondingly, we have two choices of unit cells. It
is shown in Figs. 11 (a-1) and (b-1). This correspondence between the
inversion center and the unit cell is required for the relationship between
the Zak phase and inversion parities in Eq. (7). By changing the choice of the
unit cell, the parity at $k=\pi$ changes sign, and the value of the Zak phase
is switched between $0$ and $\pi$. In Ref. Vanderbilt and King-Smith, 1993, it
is discussed that the crystal termination is related to the choice of the unit
cell so that the unit cell, with displacement by crystal translation vectors,
should fully cover the whole system with surface termination. Thus, there are
two ways of terminations for the SSH model, shown in Figs. 11 (a-2) and (b-2),
and they correspond to the choices of the unit cells shown in Figs. 11 (a-1)
and (b-1), respectively. By employing this correspondence, the relationship
between the $0$ ($\pi$) Zak phase and absence (presence) of the in-gap
boundary states is guaranteed.

We now turn to the Rice-Mele model studied in Refs. Combes  _et al._ , 2016;
Vanderbilt and King-Smith, 1993. In Ref. Combes  _et al._ , 2016, this model
is defined on a lattice in Fig. 1(a), and is described by a $2\times 2$ Bloch
Hamiltonian

$H(k)=\left(\begin{array}[]{cc}\Delta&2t(\cos\frac{k}{2}-i\delta\sin\frac{k}{2})\\\
2t(\cos\frac{k}{2}+i\delta\sin\frac{k}{2})&-\Delta\end{array}\right)$ (45)

in the basis of A and B sublattices, where $t(1\pm\delta)$ represents the
alternating hopping, and $\Delta$ is a staggered on-site potential. If we
start from our SSH model in Eq. (5), this Rice-Mele model is realized by
putting $t_{1}=-t(1+\delta)$, $t_{2}=-t(1-\delta)$, and by adding a staggered
on-site potential $\pm\Delta$ at the A and B sublattices. In Ref. Combes  _et
al._ , 2016, the spontaneous magnetization and susceptibility is studied for
the values of parameters $t=1$, $\Delta=0.6\cos\theta$ and
$\delta=0.6\sin\theta$ with a parameter $\theta$. We note that $\theta=3\pi/2$
corresponds to the SSH model with the $\pi$ Zak phase. It is found that the
dielectric response $\chi\equiv\frac{dP}{d\varepsilon}$ has no singularity
even at $\theta=3\pi/2$, and that the bulk spontaneous polarization, defined
modulo $e$, has no anomaly at $\theta=3\pi/2$. It means that there is no
anomalous dielectric response in the bulk polarization even at
$\theta=3\pi/2$. On the other hand, the polarization for a finite chain is
also calculated in Refs. Combes  _et al._ , 2016; Vanderbilt and King-Smith,
1993, and it has a jump between $e/2$ and $-e/2$ at $\theta=3\pi/2$. This
implies that the system is at the verge between the states with polarizations
$+e/2$ and $-e/2$, and that a tiny electric field can easily drive the system
to obtain polarization $\pm e/2$. Thus, the bulk polarization is defined only
modulo $e$, while the polarization for a finite chain is not modulo $e$. The
latter can be determined without ambiguity once the finite-length chain,
including its terminations, is fixed. The bulk polarization and the
polarization for a finite chain are identical only modulo $e$ (if we ignore
small finite-size effect, which converges to zero for a long-chain limit).
Thus the jump of the polarization between $-e/2$ and $e/2$ never appears in
the bulk polarization, because of its “mod $e$” nature. In other words, in the
SSH model, the anomalous dielectric response discussed in the present paper
depends on the terminations, and in the case $|t_{2}|>|t_{1}|$, it appears
only in the termination in Fig. 11(a-2) but not in the termination in Fig.
11(b-2). Thus it is natural that the bulk polarization does not have a
singularity for the SSH model with the $\pi$ Zak phase (i.e. $\theta=3\pi/2$),
where we have shown the anomalous dielectric response.

Figure 11: Choice of the unit cell and corresponding crystal termination in
the SSH model. (a-1) and (b-1) show two choices of unit cell shown as a grey
box. (a-2) and (b-2) are the SSH chain with terminations corresponding to
(a-1) and (b-1), respectively. In each case, the star represents the
corresponding inversion center.

##  VII Conclusion

To summarize, we have shown anomalous dielectric response in insulators with
the $\pi$ Zak phase. We demonstrated it in various tight-binding models in
one, two and three dimensions. In each model, the polarization suddenly rises
to a large value close to $e/2$ per edge or surface unit cell by application
of an external electric field when the Zak phase is $\pi$ over the entire
Brillouin zone. By using effective models, we have confirmed that the
topological surface states protected by the $\pi$ Zak phase give rise to this
phenomenon. We also show, by using ab-intio calculation, that diamond and
silicon slabs with (111) surfaces have the $\pi$ Zak phase over the entire
Brillouin zone and relatively flat midgap states. Our results are expected to
be applicable to these materials even when surface reconstruction occurs. We
also show such a dielectric response of the PTFE with plateaus at $P\sim\pm e$
by ab initio calculation, in agreement with our theory.

We have also discussed that this anomalous dielectric response do not appear
for the bulk polarization described by the modern theory of polarization,
because this effect depends on the termination of the system, and therefore is
attributed to the surface effect.

###### Acknowledgements.

The authors gratefully acknowledge useful discussions with D. Vanderbilt. This
work was supported by JSPS KAKENHI Grant No. JP18H03678, by MEXT KAKENHI Grant
No. JP20H04633, by CREST, JST (No. JP-MJCR14F1), and by the MEXT Elements
Strategy Initiative to Form Core Research Center (TIES), Grant Number
JPMXP0112101001.

##  Appendix A  Two-dimensional models (Models III and IV)

In Sec. II, we have seen that the dielectric polarization in the SSH model
abruptly increases up to $\pm e/2$ in response to a weak electric field when
the Zak phase is $\pi$. We expect similar phenomena in two-dimensional
insulators. To confirm that, we consider two tight-binding models on the
honeycomb lattice: the Model III with anisotropic hopping amplitude and the
Model IV with three orbitals per lattice site. In both models, the lattice
consists of two sublattices A and B. Let $a$ denote the distance between the
two nearest A sites, and we set the nearest bond vectors from the A site to
the adjacent B sites $\bm{d}_{i=1,2,3}$ to be
$\bm{d}_{1}=\frac{a}{\sqrt{3}}(0,1),\bm{d}_{2}=\frac{a}{\sqrt{3}}(-\frac{\sqrt{3}}{2},-\frac{1}{2}),\bm{d}_{3}=\frac{a}{\sqrt{3}}(\frac{\sqrt{3}}{2},-\frac{1}{2})$
as shown in Fig. 12(a). We note that in two-dimensional systems, the
polarization, which is an electric dipole moment per unit area, has the unit
of charge divided by a length.

###  A.1 Zak phase and dielectric polarization for two-dimensional systems

First, we review the Zak phase and its connection with electric polarization
in two-dimensional insulators. The Zak phase with an integration path along a
reciprocal lattice vector $\bm{G}_{\perp}$ is defined as

$\displaystyle\theta(k_{\parallel})=-i\sum_{n}^{\text{occ.}}\int_{0}^{|\bm{G}_{\perp}|}dk_{\perp}\bra{u_{n}(\bm{k})}\dfrac{\partial}{\partial
k_{\perp}}\ket{u_{n}(\bm{k})},$ (46)

where $\bm{k}=k_{\perp}\bm{n}+k_{\parallel}\bm{l}$ and
$\bm{n}={\bm{G}}_{\perp}/|\bm{G}_{\perp}|,\bm{n}\perp\bm{l},|{\bm{l}}|=1$.
This Zak phase is related with a polarization in a ribbon geometry with its
edges along $\bm{l}$, defined at $k_{\parallel}$ in the surface Brillouin
zone. The Zak phase $\theta(k_{\parallel})$ is quantized as 0 or $\pi$ modulo
$2\pi$ under both the TR and SI symmetries Kariyado and Hatsugai (2013).
Therefore, in insulators, this quantized value of $\theta(k_{\parallel})$ is
constant for any $k_{\parallel}$. In the region of the surface Brillouin zone
where $\theta(k_{\parallel})=\pi$, there exist degenerated edge states at zero
energy when the system has chiral symmetry Ryu and Hatsugai (2002). The
electric polarization $P_{0}$ across the ribbon can be obtained from

$\displaystyle
P_{0}=-e\int_{0}^{2\pi/a}\dfrac{dk_{\parallel}}{(2\pi)^{2}}\theta(k_{\parallel})\
({\rm mod}\ e/a),$ (47)

where $a$ is the lattice constant along the edge parallel to ${\bm{l}}$. This
equation indicates that the polarization takes a large value of $e/(2a)$ in
the system with the $\pi$ Zak phase for any $k_{\parallel}$. This nonzero
value of the polarization $P_{0}$ does not contradict the inversion symmetry
of the system, because $P\equiv e/(2a)\equiv-e/(2a)\ $(mod $e/a$).

###  A.2 Model III: Anisotropic tight-binding model on the honeycomb lattice

In this section, we consider an anisotropic tight-binding model on the
honeycomb lattice. The hopping amplitudes between the nearest-neighbor sites
are anisotropic. The hopping amplitude along $\bm{d}_{1}$ is $-t_{2}$ and
those along other directions are $-t_{1}$, as shown in Fig. 12 (c). $t_{1}$
and $t_{2}$ are set to be positive. The Hamiltonian is given by

$\displaystyle H_{0}=-\sum_{\langle ij\rangle}t_{ij}c^{\dagger}_{i}c_{j}\ ,$
(48)

where $t_{ij}$ is the hopping amplitude from site $i$ to $j$. First, the bulk
Hamiltonian is given by

$\displaystyle H^{\text{bulk}}(\bm{k})$
$\displaystyle=\left(\begin{array}[]{cc}0&R(\bm{k})\\\
R^{\ast}(\bm{k})&0\end{array}\right)$ (51) $\displaystyle\ R(\bm{k})$
$\displaystyle=-t_{2}e^{-i{\bm{k}}\cdot{\bm{d}_{1}}}-t_{1}(e^{-{\bm{k}}\cdot{\bm{d}_{2}}}+e^{-i{\bm{k}}\cdot{\bm{d}_{3}}}).$
(52)

The energy eigenvalues are given by

$\displaystyle E_{\bm{k}}^{2}$
$\displaystyle=\left(t_{2}+2t_{1}\cos\frac{k_{x}a}{2}\cos\frac{\sqrt{3}k_{y}a}{2}\right)^{2}$
$\displaystyle+4t_{1}^{2}\cos^{2}\frac{k_{x}a}{2}\sin^{2}\frac{\sqrt{3}k_{y}a}{2}.$
(53)

A band gap exists when $t_{2}>2t_{1}$ (see Fig. 12(e)) Dietl  _et al._ (2008);
Takahashi and Murakami (2013). Thus, we can define the dielectric polarization
only when $t_{2}>2t_{1}$. By taking $\bm{n}=(0,1)$ and $\bm{l}=(1,0)$, the Zak
phase along the $y$-axis is given by

$\displaystyle\theta_{y}{(k_{x})}=\begin{cases}\pi\ \left(\cos
k_{x}a<\dfrac{t_{2}}{2t_{1}}\right)\\\ 0\ (\text{otherwise})\end{cases}.$ (54)

When the band gap is open, i.e. $t_{2}>2t_{1}$, the Zak phase equals to $\pi$
in the whole Brillouin zone. The $\pi$ Zak phase represents the polarization
$e/(2a)({\rm mod}\ e/a)$ along the $y$-axis Delplace  _et al._ (2011). Here,
thanks to the SI symmetry, the Zak phase $\theta_{y}(k_{x})$ with $k_{x}$
being one of the one-dimensional TRIM, i.e. $k_{x}=0$ and $k_{x}=\pi/a$, can
be obtained from parity eigenvalues of the occupied states at TRIM. Namely,

$\displaystyle e^{i\theta_{y}(k_{x}=0)}$
$\displaystyle=\prod_{n}^{\text{occ.}}\xi_{n}(0,0)\xi_{n}\left(0,\frac{2\pi}{\sqrt{3}a}\right),$
(55) $\displaystyle e^{i\theta_{y}\left(k_{x}=\frac{\pi}{a}\right)}$
$\displaystyle=\prod_{n}^{\text{occ.}}\xi_{n}\left(\frac{\pi}{a},\frac{\pi}{\sqrt{3}a}\right)\xi_{n}\left(\frac{\pi}{a},-\frac{\pi}{\sqrt{3}a}\right),$
(56)

where $\xi_{n}({\bm{k}})$ is the parity eigenvalue of the $n$th eigenstate. It
reproduces the result in Eq. (54). Here, the unit cell for the bulk used in
the calculation of $\xi_{n}$ is taken as a pair of a A-site and a B-site,
displaced by ${\bm{d}_{2}}$, so that the the unit cell covers the whole ribbon
by translating the unit cell by the translation vectors, and the inversion
center is located at the midpoint of the unit cell.

Figure 12: Two-dimensional models III and IV. (a) Honeycomb lattice. The
lattice constant is $a$ and the nearest-neighbor bond vectors from the A site
to the B site are $\bm{d}_{i=1,2,3}$. (b) Zigzag-edge nano ribbon. The red
dotted line shows the unit cell in the ribbon. (c) Schematic figure of the
model with anisotropic hopping amplitudes (Model III). (d) Schematic figure of
the model with three orbitals per lattice site (Model IV). (e) Band structure
for the Model II with $t_{1}=1$, $t_{2}=3$. The bands are shown along high-
symmetry lines shown in red in the two-dimensional Brillouin zone in the right
panel. (f) Dielectric polarization $P(\varepsilon)$ of the model III in a
nano-ribbon with zigzag edges. $P(\varepsilon)$ rapidly changes toward $\pm
e/2a$ in the vicinity of $\varepsilon=0$. Here we set $N=10$ and the number of
sites in the unit cell is $20$. The hopping amplitudes $(t_{1},t_{2})$ are
$(1,3)$ for the red line and (1,10) for the blue line. As the anisotropy
increases, the change in the polarization becomes sharper.

We next study the dielectric polarization for a ribbon with zigzag edges as
shown in Figs. 12 (b) and (c). There are $2N$ sites in the unit cell of the
zigzag-edge nano ribbon. Let $\sl A_{i}$ and $\sl B_{i}$ $(i=1,2,\cdots,N)$
denote the $2N$ lattice sites in the unit cell as shown in Fig. 12 (b). The
Hamiltonian of this system is given by

$\displaystyle H_{0}(k_{x})$
$\displaystyle=-t_{1}\sum_{n=2}^{N}e^{-i\tilde{d}_{1}k_{x}}b^{\dagger}_{k_{x},n-1}a_{k_{x},n}$
$\displaystyle\ \
-t_{2}\sum_{n=1}^{N}\left(e^{-i\tilde{d}_{2}k_{x}}+e^{-i\tilde{d}_{3}k_{x}}\right)b^{\dagger}_{k_{x},n}a_{k_{x},n}+{\sl
h.c.}$ (57)

where $a_{k_{x},n}^{\dagger}(b_{k_{x},n}^{\dagger})$ and
$a_{k_{x},n}(b_{k_{x},n})$ are the creation and annihilation operators of an
electron with Bloch wave vector $k_{x}$ at the $A_{n}$ ($B_{n}$) site in the
unit cell, and $\tilde{d}_{i}~{}(i=1,2,3)$ is the $x$ component of
$\bm{d}_{i}$. When an electric field of strength $\varepsilon$ along the $y$
axis is applied, the Hamiltonian is given by

$\displaystyle H(k_{x})$ $\displaystyle=H_{0}(k_{x})+H_{\varepsilon}(k_{x}),$
(58) $\displaystyle H_{\varepsilon}(k_{x})$
$\displaystyle=\dfrac{ae\varepsilon}{\sqrt{3}}\sum_{n=1}^{N}\left[\left(\dfrac{3}{2}n-\dfrac{3}{4}N-1\right)a^{\dagger}_{k_{x},n}a_{k_{x},n}\right.$
$\displaystyle\ \ \ \ \ \ \
\left.+\left(\dfrac{3}{2}n-\dfrac{3}{4}N-\dfrac{1}{2}\right)b^{\dagger}_{k_{x},n}b_{k_{x},n}\right],$
(59)

Then the dielectric polarization is given by

$\displaystyle
P(\varepsilon)=-\dfrac{2}{\sqrt{3}a(N-2/3)}\sum_{n=1}^{N}\int\dfrac{dk_{x}}{2\pi}\dfrac{\partial
E_{n}(k_{x})}{\partial\varepsilon},$ (60)

where $E_{n}(k_{x})$ is the $n$th lowest eigenvalue of $H(k_{x})$, and
$\sqrt{3}a(N-2/3)$ is the width of the ribbon. We set the Fermi energy to be
$E_{f}=0$. The calculation results are shown in Fig. 12(f). It shows that when
$t_{2}>2t_{1}$, the dielectric polarization takes a value
$P(\varepsilon)\sim\pm e/(2a)$ even for a very weak electric field. This
indicates that a half of an electron per unit cell of an edge is accumulated
at the edge whenever $t_{2}>2t_{1}$. Thus, the anomalous dielectric response
occur in this system when it is insulating.

Figure 13: Dielectric polarization of the Model III (red lines) and that of
the effective two-site model (blue lines). We set $N=10$ for the Model II. The
hopping amplitudes are taken as $t_{1}=1$, (a) $t_{2}=2.5$, (b) $t_{2}=3$, (c)
$t_{2}=3.5$ and (d) $t_{2}=4$.

We attribute this phenomenon to topological edge states. To confirm this, we
use an effective model which describes the present model at a fixed value of
$k_{x}$, similar to Sec. II. It consists of two sites, each of which
corresponds to one of the edge sites of the ribbon. We call the two sites $L$
and $R$, and the state with the electron occupying the R[L] site is
represented as $\ket{R}=(1,0)^{T}\left[\ket{L}=(0,1)^{T}]\right)$ (see Fig.
12(b)). Then, the effective Hamiltonian is given by

$\displaystyle
H^{\ast}(k_{x})=\left(\begin{array}[]{cc}{\sqrt{3}(N-2/3)a\varepsilon}/{4}&-t^{\ast}(k_{x})\\\
-t^{\ast}(k_{x})&-{\sqrt{3}(N-2/3)a\varepsilon}/{4}\end{array}\right)\ ,$ (63)

where $t^{\ast}(k_{x})$ is the effective hopping amplitude. Since
$t^{\ast}(k_{x})$ depends on the wave vector $k_{x}$, it cannot be determined
in the same way as in the previous section. Instead, we analytically calculate
the edge states of the model, by applying the method in Ref. Wakabayashi  _et
al._ (2010) to our anisotropic model. The edge states when $t_{2}>2t_{1}$ is
given by

$\displaystyle\begin{split}\ket{\Psi^{1}_{\rm
edge}(k_{x})}&=C_{N}\sum_{n=1}^{N}(-1)^{n}\sinh[\eta_{k}(N+1-n)a^{\dagger}_{k_{x},n}\ket{0},\\\
\ket{\Psi^{2}_{\rm
edge}(k_{x})}&=C_{N}\sum_{n=1}^{N}(-1)^{N+1-n}\sinh[\eta_{k}n]a^{\dagger}_{k_{x},n}\ket{0},\end{split}$
(64)

where

$\displaystyle|C_{N}|^{2}=2[\sinh(\eta_{k}N)\cosh(\eta_{k}(N+1))/\sinh(\eta_{k})-N]$
(65)

is a normalization constant. Here, $\eta_{k}$ is a non-zero solution of

$\displaystyle\begin{cases}-t_{1}\sinh\eta_{k_{x}}N-g_{k}\sinh\eta_{k_{x}}(N+1)=0,\
\ \text{for}\ 0<|ak_{x}|<\pi,\\\
-t_{1}\sinh\eta_{k_{x}}N+g_{k}\sinh\eta_{k_{x}}(N+1)=0,\ \ \text{for}\
\pi<|ak_{x}|<2\pi,\end{cases}$ (66)

where $g_{k}=-2t_{1}\cos\frac{ka}{2}$, and we assume $t_{2}>2t_{1}$ and $N\gg
1$. The derivation of this equation is shown in Appendex B. One can evaluate
the effective hopping amplitude $t^{\ast}(k_{x})$ as

$\displaystyle t^{\ast}(k_{x})$
$\displaystyle=-\bra{L(k_{x})}H_{0}(k_{x})\ket{R(k)}$
$\displaystyle=(-1)^{N+1}2N|C_{N}|^{2}\left[t_{1}\left(\cosh(\eta_{k}N)-\dfrac{\sinh(\eta_{k}N)}{\sinh(\eta_{k})}\right)\right.$
$\displaystyle\left.+t_{2}\cos\dfrac{k_{x}}{2}\left(\cosh(\eta_{k}N)-\dfrac{\sinh(\eta_{k}(N+1))}{\sinh(\eta_{k})}+1\right)\right].$
(67)

The effective dielectric polarization is given by

$\displaystyle P^{\ast}(\varepsilon)$
$\displaystyle=\int_{0}^{2\pi/a}\dfrac{dk_{x}}{2\pi}p^{\ast}(\varepsilon,k_{x})$
(68) $\displaystyle p^{\ast}(\varepsilon,k_{x})$
$\displaystyle=\frac{e}{2}\left[\braket{L}{-}\braket{-}{L}-\braket{R}{-}\braket{-}{R}\right]$
$\displaystyle=\text{sgn}(\varepsilon)\dfrac{e}{2\sqrt{(2t^{\ast}(k_{x})/(3N/2-1)ae\varepsilon)^{2}+1}},$
(69)

where $\ket{-}$ is the eigenstate of (63) at $k_{x}$ with a negative
eigenvalue. We plot the polarization of this effective model in Fig. 13. It
agrees well with the result from the Model II, meaning that our scenario of
the anomalous polarization is valid.

###  A.3 Model IV: Model with three orbitals per lattice site on the honeycomb
lattice

Next, we consider a tight-binding model with three orbitals, $s,$ $p_{x}$ and
$p_{y}$, per lattice site as shown in Fig. 12 (d). The Hamiltonian is given by

$\displaystyle H_{0}=-\sum_{\langle
ij\rangle}\sum_{\alpha,\beta=s,x,y}t^{ij}_{\alpha\beta}c^{\dagger}_{i\alpha}c_{j\beta}\
,$ (70)

where $s$, $x$ and $y$ represent $s,p_{x}$ and $p_{y}$ orbitals respectively.
$c^{\dagger}_{i\alpha}$ and $c_{i\alpha}$ are creation and annihilation
operators for an electron in the $\alpha$ orbital at the $i$ site. The hopping
integral $t^{ij}_{\alpha\beta}$ is determined from the Slater-Koster
parameters Slater and Koster (1954) $V_{ss\sigma},V_{sp\sigma},V_{pp\sigma}$
and $V_{pp\pi}$ as shown in Table 2. In the same way as in the previous
section, we can obtain the bulk Hamiltonian $H^{bulk}(\bm{k})$ in a $6\times
6$ matrix form. Using its eigenstates, we numerically confirmed that the gap
is open around $E=0$ and the Zak phase along the $y$-axis becomes $\pi$ in the
whole Brillouin zone when
$(V_{ss\sigma},V_{sp_{\sigma}},V_{pp\sigma},V_{pp\pi})=(4.43,-3.79,-5.66,1.83)$.
These values are taken from those for carbon atoms in diamond obtained by
fitting from the values of the first-principle calculation Laref  _et al._
(1998) in the unit of eV. We set the Fermi energy to be $E_{f}=0$. Therefore,
an anomalous dielectric response is expected with these Slater-Koster
parameters.

Table 2: Expression of $t^{ij}_{\alpha\beta}$ for the model IV and the model V
Hattori  _et al._ (2017). The table shows the hopping amplitude
$t^{ij}_{\alpha\beta}$ from the $\alpha$ orbital at the site $i$ to the
$\beta$ orbital at the site $j~{}(\alpha,\beta=s,p_{x},p_{y},p_{z}).$ Here
$l_{ij}$, $m_{ij}$ and $n_{ij}$ are the $x,y$ and $z$ components of the
direction cosines measured from site $i$ to $j$. Other components such as
$t^{sy}_{ij},t^{yy}_{ij}$ and $t^{yz}_{ij}$ are obtained from the components
in the table by a replacement with the correspondence, $x\leftrightarrow
l_{ij}$, $y\leftrightarrow m_{ij}$, and $z\leftrightarrow n_{ij}$.
$t^{ij}_{\alpha\beta}$
---
$t_{ij}^{ss}$ |  $l_{ij}V_{ss}$
$t_{ij}^{sx}=-t_{ij}^{xs}$ |  $l_{ij}V_{sp}$
$t_{ij}^{xx}$ |  $l^{2}_{ij}V_{pp\sigma}+(1-l^{2}_{ij})V_{pp\pi}$
$t_{ij}^{xy}$ |  $l_{ij}m_{ij}(V_{pp\sigma}-V_{pp\pi})$
Figure 14: Dielectric polarization of the Model IV in a zigzag-edge nano
ribbon. The model is on the honeycomb lattice with anisotropic hopping
amplitudes. $P(\varepsilon)$ rapidly changes toward $\pm e/(2a)$ in the
vicinity of $\varepsilon=0$. Here we set $N=10$ and the number of sites in the
unit cell of the ribbon is $20$. The Slater-Koster parameters
$(V_{ss},V_{sp},V_{pp\sigma},V_{pp\pi})$ = $(4.43,-3.79,-5.66,1.83)$.

The Hamiltonian of this system is given by

$\displaystyle H(k)$ $\displaystyle=H_{0}(k)+H_{\varepsilon}(k),$ (71)
$\displaystyle H_{\varepsilon}$
$\displaystyle=\dfrac{ae\varepsilon}{\sqrt{3}}\sum_{\alpha}\sum_{n=1}^{N}\left[\left(\dfrac{3}{2}n-\dfrac{3}{4}N-1\right)a^{\dagger}_{k,\alpha,n}a_{k,\alpha,n}\right.$
$\displaystyle\ \ \ \ \ \ \
\left.+\left(\dfrac{3}{2}n-\dfrac{3}{4}N-\dfrac{1}{2}\right)b^{\dagger}_{k,\alpha,n}b_{k,\alpha,n}\right],$
(72)

where $a_{k,\alpha,n}^{\dagger}(b_{k,\alpha,n}^{\dagger})$ is the creation
operator of an electron with the Bloch wave vector $k$ in the $\alpha$ orbital
at the $n$th A(B) site. $H_{0}(k)$ is the Bloch form of the Hamiltonian
$H_{0}$ with the Bloch wavevector $k$ along the $x$ axis. The dielectric
polarization is given by

$\displaystyle
P(\varepsilon)=-\dfrac{2}{\sqrt{3}a(N-2/3)}\sum_{n=1}^{3N}\int_{0}^{2\pi/a}\dfrac{dk}{2\pi}\dfrac{\partial
E_{n}(k)}{\partial\varepsilon}$ (73)

where $E_{n}(k)$ is the $n$th lowest eigenvalue of (71). The calculation
results are shown in Fig. 14. The dielectric polarization takes a value
$P(\varepsilon)\sim\pm e/(2a)$ for a very weak electric field. Thus, an
anomalous dielectric response occurs even in this case, and it is also
attributed to midgap edge states as is similar to the models I and II.

##  Appendix B Analytical derivation of the edge states for Model III

In this section, we show the derivation of the eigenstates of (57). In Ref.
Wakabayashi  _et al._ (2010), wavefunctions for the edge states of Model II in
an isotropic case is analytically derived. Here we extend this method to an
anisotropic case. According to Ref. Wakabayashi  _et al._ (2010), eigenstates
of the model on the anisotropic honeycomb lattice in a nano-ribbon are
expressed as

$\displaystyle\ket{\Psi(k_{x})}=\sum_{n=1}^{N}\left[a_{n}(k_{x})\ket{A_{n}(k_{x})}+b_{n}(k_{x})\ket{B_{n}(k_{x})}\right],$
(74)

where $k_{x}$ is the wave number along the edge, and
$a_{n}(k_{x})[b_{n}(k_{x})]$ denotes the coefficient for the wave function at
$A_{n}[B_{n}]$ site, $\ket{A_{n}(k_{x})}[\ket{B_{n}(k_{x})}]$. By using Eqs.
(57) and (74), the matrix elements of the Hamiltonian are given by

$\displaystyle\begin{split}\bra{B_{n}(k_{x})}H_{0}(k_{x})\ket{A_{n}(k_{x})}&=-2t_{1}\cos\dfrac{ka}{2}\equiv
g_{k}\\\
\bra{B_{n-1}(k_{x})}H_{0}(k_{x})\ket{A_{n}(k_{x})}&=-t_{2}.\end{split}$ (75)

Here, we impose a boundary condition

$\displaystyle a_{N+1}(k_{x})=b_{0}(k_{x})=0.$ (76)

The Schr$\ddot{\text{o}}$dinger equation
$H_{0}(k_{x})\ket{\Psi(k_{x})}=E_{k_{x}}\ket{\Psi(k_{x})}$ gives the following
equations;

$\displaystyle\begin{split}&g_{k}a_{n}(k_{x})-t_{2}a_{n+1}(k_{x})=E_{k_{x}}b_{n}(k_{x}),\\\
&g_{k}b_{n}(k_{x})-t_{2}b_{n-1}(k_{x})=E_{k_{x}}a_{n}(k_{x}).\\\ \end{split}$
(77)

Because it can be expressed as an eigenvalue problem of a $2N\times 2N$
matrix, it has $2N$ solutions. We assume the following forms for
$a_{n}(k_{x})$ and $b_{n}(k_{x})$:

$\displaystyle\begin{split}a_{n}(k_{x})&=Ae^{ipan}+Be^{-ipan},\\\
b_{n}(k_{x})&=Ce^{ipan}+De^{-ipan}.\end{split}$ (78)

From the boundary condition, the coefficients $B$ and $D$ are determined as
$B=-Az^{2}$ and $D=-C$ with $z=e^{ip(N+1)}$. Substituting these values into
Eq. (77), we get

$\displaystyle M\begin{pmatrix}A\\\ C\end{pmatrix}=0,\
M=\begin{pmatrix}m_{11}&m_{12}\\\ m_{21}&m_{22}\end{pmatrix},$ (79)

where the matrix elements are

$\displaystyle\begin{split}m_{11}&=-E_{k_{x}}(e^{ipan}-z^{2}e^{-ipan}),\\\
m_{12}&=e^{ipan}(g_{k}-t_{1}e^{-ipa})-e^{-ipan}(g_{k}-t_{1}e^{ipa}),\\\
m_{21}&=e^{ipan}(g_{k}-t_{1}e^{ipa})-e^{-ipan}z^{2}(g_{k}-t_{1}e^{-ipa}),\\\
m_{22}&=-E_{k_{x}}(e^{ipan}-e^{-ipan}).\end{split}$ (80)

Equation (79) has nontrivial solutions only when $\det M=0$. Then, we get the
following two equations;

$\displaystyle E_{k_{x}}^{2}$
$\displaystyle=(g_{k}-t_{2}e^{ipa})(g_{k}-t_{2}e^{-ipa}),$ (81) $\displaystyle
E_{k_{x}}^{2}(1+z^{2})$
$\displaystyle=(g_{k}-t_{2}e^{-ipa})^{2}z^{2}+(g_{k}-t_{2}e^{ipa})^{2}.$ (82)

By eliminating $E_{k_{x}}$ from these equations, we obtain an equation for
$p$:

$\displaystyle F(p,N)\equiv-t_{2}\sin paN+g_{k}\sin pa(N+1)=0.$ (83)

This equation holds when $pa=0,\pm\pi$. However, this corresponds to the
trivial solutions $a_{n}=b_{n}=0$ for all $n$, and should be discarded. As
shown in Ref. Wakabayashi  _et al._ (2010), there are only $2N-2$ solutions
which correspond to the bulk states for $p$ when

$\displaystyle|g_{k}|\leq\dfrac{N}{N+1}t_{2}.$ (84)

The remaining two solutions are expected to correspond to the edge states.
These edge states can be obtained by analytical continuation as

$\displaystyle pa=\begin{cases}&\pi\pm i\eta_{k_{x}}=p_{\pi},\
k^{L}_{c}a<|k_{x}a|<\pi,\\\ &0\pm i\eta_{k_{x}}=p_{0},\
\pi<|k_{x}a|<k^{R}_{c}a,\end{cases}$ (85)

where the wavevectors $k^{L}_{c}$ and $k^{R}_{c}$ are defined as

$\displaystyle\begin{cases}\cos\dfrac{k^{L}_{c}a}{2}&=\dfrac{t_{2}}{2t_{1}}\dfrac{N}{N+1},\\\
\cos\dfrac{k^{R}_{c}a}{2}&=-\dfrac{t_{2}}{2t_{1}}\dfrac{N}{N+1}.\end{cases}$
(86)

Substituting (85) into (83), we get

$\displaystyle\begin{cases}-t_{1}\sinh\eta_{k_{x}}N-g_{k}\sinh\eta_{k_{x}}(N+1)=0,\
\ \text{for}\ p_{\pi},\\\
-t_{1}\sinh\eta_{k_{x}}N+g_{k}\sinh\eta_{k_{x}}(N+1)=0,\ \ \text{for}\
p_{0}.\end{cases}$ (87)

There exists one solution for $\eta_{k_{x}}$ each equation in Eq. (87) when
Eq. (84) is satisfied. When $t_{2}>2t_{1}$ and in the limit of large $N$, Eq.
(87) is always satisfied and the two edge states are given by Eq. (85).

Since we focus on the edge states, we consider only the solution of (87).
Then, the coefficients $a_{n}$ and $b_{n}$ are given by

$\displaystyle a_{n}(k_{x})$
$\displaystyle=\pm(-1)^{N+1-n}C_{N}\sinh\eta_{k_{x}}(N+1-n),$ (88)
$\displaystyle b_{n}(k_{x})$ $\displaystyle=(-1)^{n}C_{N}\sinh\eta_{k_{x}}n,$
(89)

where $C_{N}$ is the normalization factor expressed in (65). The two edge
states $\ket{\Psi_{\rm edge}^{+}(k_{x})}$ and $\ket{\Psi_{\rm
edge}^{-}(k_{x})}$ are obtained by substituting these into (74). These edge
states $\ket{\Psi_{\rm edge}^{\pm}(k_{x})}$ are not localized at one edge, but
are symmetric and antisymmetric linear combinations of edge states at an
either end. Therefore, to obtain approximate edge states at one edge, we take
linear combinations as follows:

$\displaystyle\ket{\Psi^{1}_{\rm edge}(k_{x})}$
$\displaystyle=\left(\ket{\Psi_{edge}^{+}(k_{x})}+\ket{\Psi_{edge}^{-}(k_{x})}\right)/\sqrt{2},$
(90) $\displaystyle\ket{\Psi^{2}_{\rm edge}(k_{x})}$
$\displaystyle=\left(\ket{\Psi_{edge}^{+}(k_{x})}-\ket{\Psi_{edge}^{-}(k_{x})}\right)/\sqrt{2}.$
(91)

These are given in Eq. (64) and represent states localized at one end.

##  Appendix C Three-dimensional Model V: Model with four orbitals per lattice
site on the diamond lattice

In this section, we show the anomalous dielectric response in a three-
dimensional tight-binding model (Model V) on the diamond lattice, with four
orbitals per unit cell. As is similar to Model II, we use the diamond lattice
and its slab with (111) surfaces are shown in Fig. 3 (a) and (b). Then, we
consider a tight-binding model with four orbitals, $s,p_{x},p_{y}$ and
$p_{z}$, per site on the diamond lattice as shown in Fig. 15 (d). The
Hamiltonian is

$\displaystyle H_{0}=-\sum_{\langle
ij\rangle}\sum_{\alpha,\beta=s,x,y,z}t^{ij}_{\alpha\beta}c^{\dagger}_{i\alpha}c_{j\beta}\
,$ (92)

where $s$, $x$, $y$ and $z$ represent $s,p_{x}$,$p_{y}$ and $p_{z}$ orbitals,
respectively. Here, $c^{\dagger}_{i\alpha}$ and $c_{i\alpha}$ are creation and
annihilation operators for an electron with the $\alpha$ orbital at $i$ site.
The hopping integrals $t^{ij}_{\alpha\beta}$ are taken as the Slater-Koster
form shown in Table 2. We can obtain the bulk Hamiltonian $H^{bulk}(\bm{k})$
as an $8\times 8$ matrix form. Using its eigenstates, we numerically confirmed
that the gap is open and that the Zak phase along the $z$-axis becomes $\pi$
in the whole Brillouin zone when
$(V_{ss\sigma},V_{sp_{\sigma}},V_{pp\sigma},V_{pp\pi})=(4.43,-3.79,-5.66,1.83)$,
which are obtained by fitting from the values of the first-principle
calculation with the unit of eV Laref  _et al._ (1998). We take the Fermi
energy to be $E_{f}=0$. Therefore, the anomalous dielectric response is
expected when the Slater-Koster parameters take these values.

The Hamiltonian of the slab with the (111) surfaces is given by

$\displaystyle H(\bm{k}_{\parallel})$
$\displaystyle=H_{0}(\bm{k}_{\parallel})+H_{\varepsilon}(\bm{k}_{\parallel}),$
(93) $\displaystyle H_{\varepsilon}(\bm{k}_{\parallel})$ $\displaystyle=$
$\displaystyle\dfrac{\sqrt{6}}{4}a$ $\displaystyle
e\varepsilon\sum_{n=1}^{N}\sum_{\alpha=s,x,y,z}\left[\left(\dfrac{4}{3}n-\dfrac{2}{3}N-\dfrac{5}{6}\right)a^{\dagger}_{\bm{k}_{\parallel},\alpha,n}a_{\bm{k}_{\parallel},\alpha,n}\right.$
$\displaystyle\ \ \ \ \ \ \
\left.+\left(\dfrac{4}{3}n-\dfrac{2}{3}N-\dfrac{1}{2}\right)b^{\dagger}_{\bm{k}_{\parallel},\alpha,n}b_{\bm{k}_{\parallel},\alpha,n}\right],$
(94)

where
$a_{\bm{k}_{\parallel},\alpha,n}^{\dagger}(b_{\bm{k}_{\parallel},\alpha,n}^{\dagger})$
is a creation operator of an electron with a Bloch wave vector
$\bm{k}_{\parallel}$ in the $\alpha$ orbital of the $n$th A(B) site in the
unit cell. The dielectric polarization is given by

$\displaystyle
P(\varepsilon)=-\dfrac{1}{a(4N/3-1)}\sum_{n=1}^{4N}\int_{\bm{k}_{\parallel}\in\text{2DBZ}}\dfrac{d^{2}\bm{k}_{\parallel}}{(2\pi)^{2}}\dfrac{\partial
E_{n}(\bm{k}_{\parallel})}{\partial\varepsilon},$ (95)

where $E_{n}(\bm{k}_{\parallel})$ is the $n$th lowest eigenvalue of (93). The
calculation results are shown in Fig. 15 (f). The dielectric polarization
takes a value $P(\varepsilon)=\pm e/(2S)$ even for a very weak electric field.
This indicates that a half of an electrons are accumulated at the surface per
unit cell. Thus, the anomalous dielectric response occurs in this system when
it is insulating.

Figure 15: Model V on the diamond lattice. (a) Schematic figure of the model
with four orbitals $s,p_{x},p_{y}$ and $p_{z}$ orbitals per lattice site
(Model V). (b) Dielectric polarization of a slab of the Model V on a diamond
lattice with the (111) surfaces. $P(\varepsilon)$ rapidly changes toward $\pm
e/(2S)$ in the vicinity of $\varepsilon=0$, where $S$ is an area of the
surface unit cell. Here we set $N=10$ and the number of sites in the unit cell
is $20$. The Slater-Koster parameters are
$(V_{ss},V_{sp},V_{pp\sigma},V_{pp\pi})$ = $(4.43,-3.79,-5.66,1.83)$.

##  Appendix D Computational conditions

The electronic structures of diamond, silicon, carbon nanotube, and PTFE are
obtained from the GGA of the DFT. We use the ab initio code OpenMX based on
localized basis functions and norm-conserving pseudopotentials. We use the
Perdew-Burke-Ernzerhof (PBE) functional in the GGA Perdew  _et al._ (1996). We
employ the $8\times 8\times 8$ and $8\times 8\times 1$ $\bm{k}$-point mesh for
the bulk and the slab calculation for diamond and silicon, respectively, and
the $1\times 1\times 16$ $\bm{k}$-point mesh for carbon nanotube and PTFE. The
valence orbital set is $s2p2d1$ for C, $s2p3d1$ for Si, and $s2p2d1$ for F.
The energy cutoff for the numerical integrations is 150 Ry.

## References

  * Kane and Mele (2005) C. L. Kane and E. J. Mele, Phys. Rev. Lett. 95, 146802 (2005).
  * Fu  _et al._ (2007) L. Fu, C. L. Kane, and E. J. Mele, Phys. Rev. Lett. 98, 106803 (2007).
  * Bernevig  _et al._ (2006) B. A. Bernevig, T. L. Hughes, and S.-C. Zhang, Science 314, 1757 (2006).
  * Fu and Kane (2006) L. Fu and C. L. Kane, Physical Review B 74, 195312 (2006).
  * Zak (1989) J. Zak, Phy. Rev. Let 62, 2747 (1989).
  * Berry (1984) M. V. Berry, Proc. R. Soc. Lond. A 392, 45 (1984).
  * Ryu and Hatsugai (2002) S. Ryu and Y. Hatsugai, Phy. Rev. Lett. 89, 077002 (2002).
  * Delplace  _et al._ (2011) P. Delplace, D. Ullmo, and G. Montambaux, Phys. Rev. B 84, 195452 (2011).
  * Kariyado and Hatsugai (2013) T. Kariyado and Y. Hatsugai, Phys. Rev. B 88, 245126 (2013).
  * Rhim  _et al._ (2017) J.-W. Rhim, J. Behrends, and J. H. Bardarson, Phys. Rev. B 95, 035421 (2017).
  * van Miert and Ortix (2017) G. van Miert and C. Ortix, Phys. Rev. B 96, 235130 (2017).
  * Pletyukhov  _et al._ (2020) M. Pletyukhov, D. M. Kennes, J. Klinovaja, D. Loss, and H. Schoeller, Phys. Rev. B 101, 161106 (2020).
  * Resta (1992) R. Resta, Ferroelectrics 136, 51 (1992).
  * King-Smith and Vanderbilt (1993) R. D. King-Smith and D. Vanderbilt, Phys. Rev. B 47, 1651 (1993).
  * Vanderbilt and King-Smith (1993) D. Vanderbilt and R. D. King-Smith, Phys. Rev. B 48, 4442 (1993).
  * Resta (1994) R. Resta, Rev. Mod. Phys. 66, 899 (1994).
  * Hirayama  _et al._ (2017) M. Hirayama, R. Okugawa, T. Miyake, and S. Murakami, Nat. Commun 8, 14022 (2017).
  * Hirayama  _et al._ (2018) M. Hirayama, S. Matsuishi, H. Hosono, and S. Murakami, Phys. Rev. X 8, 031067 (2018).
  * Su  _et al._ (1979) W. Su, J. Schrieffer, and A. J. Heeger, Phys. Rev. Lett. 42, 1698 (1979).
  * Combes  _et al._ (2016) F. Combes, M. Trescher, F. Piéchon, and J.-N. Fuchs, Phys. Rev. B 94, 155109 (2016).
  * Takahashi and Murakami (2013) R. Takahashi and S. Murakami, Phys. Rev. B 88, 235303 (2013).
  * Hughes  _et al._ (2011) T. L. Hughes, E. Prodan, and B. A. Bernevig, Phys. Rev. B 83, 245132 (2011).
  * Schlier and Farnsworth (1959) R. Schlier and H. Farnsworth, J. Chem. Phys 30, 917 (1959).
  * Takagi and Okada (2008) Y. Takagi and S. Okada, Surf. Sci. 602, 2876 (2008).
  * Smeu  _et al._ (2012) M. Smeu, H. Guo, W. Ji, and R. A. Wolkow, Phys. Rev. B 85, 195315 (2012).
  * Nunes and Vanderbilt (1994) R. W. Nunes and D. Vanderbilt, Phys. Rev. Lett. 73, 712 (1994).
  * Nunes and Gonze (2001) R. W. Nunes and X. Gonze, Phys. Rev. B 63, 155107 (2001).
  * Souza  _et al._ (2002) I. Souza, J. Íñiguez, and D. Vanderbilt, Phys. Rev. Lett. 89, 117602 (2002).
  * Umari and Pasquarello (2002) P. Umari and A. Pasquarello, Phys. Rev. Lett. 89, 157602 (2002).
  * Swiecicki and Sipe (2014) S. D. Swiecicki and J. E. Sipe, Phys. Rev. B 90, 125115 (2014).
  * Chen and Lee (2011) K.-T. Chen and P. A. Lee, Phys. Rev. B 84, 113111 (2011).
  * Dietl  _et al._ (2008) P. Dietl, F. Piéchon, and G. Montambaux, Phys. Rev. Lett. 100, 236405 (2008).
  * Wakabayashi  _et al._ (2010) K. Wakabayashi, K.-i. Sasaki, T. Nakanishi, and T. Enoki, Sci. Technol. Adv. Mater. 11, 054504 (2010).
  * Slater and Koster (1954) J. C. Slater and G. F. Koster, Phys. Rev. 94, 1498 (1954).
  * Laref  _et al._ (1998) A. Laref, B. Bouhafs, H. Aourag, and N. Bouarissa, J. Phys: Condens. Matter 10, 3195 (1998).
  * Hattori  _et al._ (2017) A. Hattori, S. Tanaya, K. Yada, M. Araidai, M. Sato, Y. Hatsugai, K. Shiraishi, and Y. Tanaka, J. Phys: Condens. Matter 29, 115302 (2017).
  * Perdew  _et al._ (1996) J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996).