# Emergent non-Hermitian boundary contributions to charge pumping and electric
polarization

K. Kyriakou  K. Moulopoulos  University of Cyprus, Department of Physics,
1678, Nicosia, Cyprus

###### Abstract

The phenomenon of charge pumping and the modern theory of electric
polarization are reconsidered by analytically taking into account emergent
non-Hermitian contributions. These are accounted for through the use of an
extended definition of the velocity operator and are determined by means of a
dynamic Hellmann–Feynman theorem (DHFT) that we derive here for the first
time. The DHFT introduces generalized Berry curvatures and it is valid for
calculating observables non-perturbatively, hence with results valid to all
orders of the external fields. By using the extended velocity operator we
rigorously show how the charge pumping is linked up with the boundaries of the
material (with the non-Hermiticity being essential for this connection), and
by means of the DHFT we show that the well-known topological quantization of
the pumped charge breaks down due to a non-integrable Aharonov-Anandan phase
in driven non-equilibrium processes whenever the periodic gauge cannot be
applied to the Floquet-Bloch states. Likewise, we show that the electronic
polarization change has an additional non-Hermitian contribution, which is
overlooked in the modern theory of electric polarization. The non-Hermitian
contribution is by definition a bulk quantity that may equally be evaluated as
a boundary quantity due to a symmetric structure that allows the bulk
integration to be transformed into a boundary one. This non-Hermitian
contribution is very sensitive to the realistic boundary conditions imposed on
the wavefunctions and it is therefore expected to be significant in biased
insulators where charge accumulation over their boundaries is present during
the process that causes the polarization change. Finally, we show how a well-
defined surface-charge theorem can be formulated in terms of the boundary non-
Hermitian contribution.

###### pacs:

##  I Introduction

Boundary effects, combined with the involvement of dynamics and topological
concepts have become the cornerstones of the quantum mechanical description of
topological condensed matter systems Roy and Harper (2017); Caio  _et al._
(2015); Kitagawa  _et al._ (2010). The interplay between bulk and boundary
contributions within a quantum system Qi  _et al._ (2006); Wang  _et al._
(2017); Hatsugai and Fukui (2016); Wang  _et al._ (2016, 2018), seems to be
far less clear than the counterpart in a classical system. A simple quantum
theoretical framework that, at the same time takes into account the dynamics,
the bulk and boundary contributions of an observable, together with the use of
topological concepts, seems to be absent. In the theoretical studies of charge
pumping and electric polarization, dynamics and bulk contributions, as well as
topological concepts, are always explicitly involved. It would therefore be a
desirable theoretical advance, if a simple theoretical method could be given
that takes into account all involved aspects, together with explicit boundary
contributions, in a unified formulation. The present article proposes
precisely such an advance in an analytically careful and rigorous manner.

The central quantity to be calculated in charge pumping studies is the
electrons’ collective center of mass displacement
${\displaystyle\Delta\\!\left\langle\mathbf{r}\right\rangle_{\text{coll}}/N}$
at zero voltage-bias (absence of mean force) after a time period $T$ that the
Hamiltonian has gone a periodic change, where $N$ is the number of electrons
that have been displaced. The adiabatic (equilibrium) charge pumping was first
studied by Thouless in the seminal work of Ref. Thouless, 1983 and serves as a
simple yet fundamental example of topology in quantum systems. It shares the
same topological origin as the integrally quantized Hall conductivity Thouless
_et al._ (1982); Kohmoto (1985) and may thus be regarded as a dynamical
version of the Integer Quantum Hall Effect (IQHE). This will actually show up
in our derivation of the dynamic Hellmann–Feynman theorem (DHFT), where the
IQHE and Thouless charge pumping will appear as parts of a single result; this
connection appeared analytically in earlier unpublished work Kyriakou and
Moulopoulos (2015) but without the emergent non-Hermitian contribution that is
the central focus of the present work.

The hallmark of this effect of (equilibrium) charge pumping is the transport
of a precisely quantized amount of charge during an adiabatic cycle in the
parameter space. This remarkable phenomenon, has been demonstrated by
experiments in a one-dimensional chain of ultra-cold atoms trapped in an
optical lattice Lu  _et al._ (2016); Nakajima  _et al._ (2016); Lohse  _et
al._ (2016). On the other hand, driven non-equilibrium (non-adiabatic)
particle processes using the Floquet-Bloch bands $\varepsilon_{a}(\mathbf{k})$
and the counterpart Floquet states, were first studied in Ref. Shih and Niu,
1994 where it was shown that the quantization of the particle transport breaks
down due to emergence of band gaps at the quasienergy bands. Similar findings
were presented in a recent study Privitera  _et al._ (2018), where, by a
careful Floquet analysis of a closed, clean, and non-interacting driven Rice-
Mele model in the thermodynamic limit, it was found that the pumped charge
deviates from the topologically quantized value by an exponentially small
amount, provided that the initial state coincides with the lowest-energy
Floquet state. Specifically, in Refs Shih and Niu, 1994; Privitera  _et al._ ,
2018 a Hellmann-Feynman theorem for time-periodic systemsSambe (1973) has been
used (derived from the eigenvalue equation of the so-called Floquet
Hamiltonian) that, however, does not take into account non-Hermitian, Berry
curvature and Aharonov-Anandan phase contributions; in this respect, charge
pumping has been studied solely in terms of the spectral properties of Floquet
quasienergy alone. Along these lines but in a different context (by studying
the Floquet operator spectra) Ref. Kitagawa  _et al._ , 2010 has provided a
topological characterization of periodically driven quantum systems in terms
of the quasienergy winding alone.

In the development of the modern theory of electric polarization Resta (1992);
King-Smith and Vanderbilt (1993); Ortiz and Martin (1994); Resta (1994); Resta
and Vanderbilt (2007); Resta (2010); Vanderbilt (2018) for periodic and
extended systems, the charge pumping process has played a central role. One of
the differences between charge pumping and polarization change processes, is
the voltage-bias (non-zero mean force) that is present in the latter process
while is absent in the former. Polarization changes ${\Delta\mathbf{P}}$ have
been theoretically realized to be more fundamental quantities than the
“absolute” bulk polarization ${\mathbf{P}}$ itself, and they can be quantified
by using the Berry phase of the electronic wavefunctions. The modern theory
avoids addressing the “absolute” polarization of a given equilibrium state,
quite in agreement with the experiments, which invariably measure polarization
differences. The polarization difference is defined in Ref. Resta and
Vanderbilt, 2007 as the time-integrated transient macroscopic current that
flows through the insulating sample during the switching process

$\Delta\mathbf{P}=\int_{0}^{T}\mathbf{J}(t)dt,$ (1)

provided that the time-dependent Hamiltonian remains insulating at all times.
Although Eq. (1) is the fundamental equation of modern theory, its
differential form

$\frac{d\mathbf{P}(t)}{dt}=\mathbf{J}(t),$ (2)

has the same structure as the one that gives the classical polarization
current with respect to bound charges
${\displaystyle\frac{d\mathbf{P_{b}}(t)}{dt}=\mathbf{J}_{b}(t)}$. On the other
hand, by having in mind that in a periodic band insulator the distinction
between “free” (de-localized) and “bound” (localized) charges is ambiguous,
one can realize that the polarization difference ${\Delta\mathbf{P}}$ given by
Eq. (1) must generally have an extra contribution due to de-localized (bound)
Bloch states that may be regarded as the “free” charge states contribution. We
show that such a contribution indeed exists and is captured by an emerging
non-Hermitian term that accounts for such delocalized states. Therefore, we
propose that the collective induced electric polarization difference
$\Delta\mathbf{P}$ must be defined (assuming non-interacting particles and by
deploying the many-body Slater determinant wavefunction) as

$\Delta\mathbf{P}=\frac{1}{V}\sum\limits_{i}q_{i}\
\Delta\left\langle\psi_{i}(t)|\,\mathbf{r}\,|\psi_{i}(t)\right\rangle,$ (3)

where $q_{i}$ and
${\displaystyle\Delta\left\langle\psi_{i}(t)|\,\mathbf{r}\,|\psi_{i}(t)\right\rangle}$
are the charge and the displacement of each particle respectively, and $V$ is
the volume of the material. In Ref. Kyriakou and Moulopoulos, 2018 it was
shown that, for a periodic and extended system in the thermodynamic limit,
even though each electron’s position expectation value
${\displaystyle\left\langle\psi_{i}(t)|\,\mathbf{r}\,|\psi_{i}(t)\right\rangle}$
is an undefined quantity, its displacement
${\displaystyle\Delta\left\langle\psi_{i}(t)|\,\mathbf{r}\,|\psi_{i}(t)\right\rangle}$
is a well-defined quantity provided that one uses the extended velocity
operator definition

$\Delta\left\langle\psi_{i}(t)|\,\mathbf{r}\,|\psi_{i}(t)\right\rangle=\int_{0}^{T}\\!\\!\left\langle\psi_{i}(t)|\,\mathbf{v}_{ext}\,|\psi_{i}(t)\right\rangle
dt,$ (4)

where the expectation values are assumed to be evaluated over the entire
volume $V$ of the system. The extended velocity operator expectation value
${\left\langle\psi_{i}(t)|\,\mathbf{v}_{ext}\,|\psi_{i}(t)\right\rangle}$ has
two separate contributions (originating from the additive structure of the
extended operator itself, see below Eq. (7)); one is due to the standard
velocity operator
${\displaystyle\mathbf{v}=\frac{i}{\hbar}\left[H(\mathbf{r},t),\mathbf{r}\right]}$,
and another due to an emerging non-Hermitian contribution that is attributed
to the boundary velocity operator defined as
${\mathbf{v}_{b}=\frac{i}{\hbar}\\!\left({H(\mathbf{r},t)}^{+}-H(\mathbf{r},t)\right)\\!\mathbf{r}}$.
The boundary velocity expectation value has the physical meaning of a position
weighted probability flux through the face (boundary) of the material. For
bulk localized states when the wavefunctions are zero over the material
boundaries, the boundary velocity expectation value turns to zero; on the
other hand, bulk delocalized Bloch states give non-zero boundary velocity
expectation value. In this framework, Eq. (3) together with Eq. (4) may be
regarded as an extended theoretical method within the modern theory of
polarization, that simultaneously takes into account both the “bound” charge
contribution (due to localized states) and the “free” charge contribution (due
to delocalized states) described by this emergent non-Hermitian effect, and as
such it has not yet been considered in the modern theory. (We shall see later
in Eq. (II) that such a contribution is necessary for a complete study of
polarization changes and the associated induced currents.)

In order to calculate exactly the boundary and the extended velocity operator
expectation values (without using time-dependent perturbation theory),
quantities that, as we shall see, are involved in the charge pumping and
polarization change respectively, we have derived in the present work a
dynamical Hellmann-Feynman theorem (DHFT). This theorem can be applied to
arbitrary time-dependent states with arbitrary time-dependent parameters, and
it is valid for calculating observables up to all orders of the external
fields. The DHFT explicitly involves generalized Berry curvatures as well as
boundary contributions due to an emerging non-Hermitian effect term. By using
this DHFT, one is enabled to study simultaneously in the simplest manner
possible the dynamics of the system, the emerging non-Hermitian contributions
and the involvement of topology, in a unified formalism. The DHFT can also be
employed to non-equilibrium driven quantum processes, i.e. to Floquet states,
where the quest for topological invariants is still in an ongoing research
state Kitagawa  _et al._ (2010); Dehghani  _et al._ (2015); Lindner  _et al._
(2011); Carpentier  _et al._ (2015); Rudner  _et al._ (2013); Budich  _et al._
(2017); Oka and Aoki (2009); Privitera  _et al._ (2018); Wang  _et al._
(2015); Kolodrubetz  _et al._ (2018).

With respect to the new DHFT with the inclusion of the non-Hermitian boundary
effect (that is here incorporated for the first time into the general
foundations of the theory of pumping and polarization), it is useful to
provide some further details on how this has been accomplished: we first re-
examine the Thouless adiabatic particle transport (for an insulator with a
fully occupied valence band) by using the extended velocity operator
definition and the adiabatic limit of the DHFT without using (first-order)
time-dependent perturbation theory. Instead of assuming a sliding potential
${V(\mathbf{r}-\upsilon\,t\,\bm{e}_{x})}$ like in the original problem
Thouless (1983), we assume a microscopic electric field ${\displaystyle
E_{x}^{(micro)}(x,t)=E(t)cos(2\pi x/a)}$ along $\bm{e}_{x}$ direction. The
microscopic electric field is periodic in space (with period equal to the
lattice constant $a$) as well as over time (with period equal to $T$) that (i)
originates from a scalar potential that crucially breaks the bulk mirror
symmetry of the Hamiltonian along $\bm{e}_{x}$ direction and (ii) gives zero
macroscopic (spatially averaged) electric field, therefore no voltage-bias is
present during the pumping. By using the emerging non-Hermitian effect
boundary velocity operator ${\mathbf{v}_{b}}$, we show that, for delocalized
Bloch states, each electron’s displacement over the material boundaries can be
calculated by means of
${\displaystyle\Delta\left\langle{\mathbf{r}}\right\rangle\\!|_{bound}=\int_{0}^{T}\\!\\!\left\langle\mathbf{v}_{b}\right\rangle
dt}$. Given now that each electron’s ground eigenstate evolves in time
adiabatically, and by using the DHFT, we find that the collective non-
Hermitian electrons’ displacement in the $x$ direction and over the boundaries
${\sum_{i}^{N}\Delta\left\langle x_{i}\right\rangle\\!|_{bound}\ /\ {N}}$ is
quantized in units of lattice constant due to topology. Hence, we demonstrate
this way that, in the adiabatic limit, the non-Hermitian boundary term is
needed in order to obtain the standard topological invariant (in fact, an
omission of this would lead to a paradox in topological arguments - see
discussion after Eq. (58)). We then relax adiabaticity and study the non-
equilibrium driven charge pumping. We assume delocalized Floquet-Bloch states
Gómez-León and Platero (2013) for the electrons’ motion, and use the extended
velocity operator definition together with the non-adiabatic form of DHFT, in
order to evaluate the non-Hermitian collective displacement of the centre of
mass of the electrons ${\sum_{i}^{N}\Delta\left\langle
x_{i}\right\rangle\\!|_{bound}\ /\ {N}}$. We find that, for a fully occupied
Floquet-Bloch band $\varepsilon_{a}$, the quantization of the pumped charge
(per cycle) breaks down due to a non-trivial (non-integrable) Aharonov-Anandan
phase Aharonov and Anandan (1987). Particularly, in non-equilibrium driven
processes the particle transport deviates from the topologically quantized
values as the difference of the Aharonov-Anandan phases at the edges of the
Brillouin zone. This driven particle transport breaking, emerges whenever the
Floquet-Bloch states ${\Phi_{\varepsilon_{a}}(\mathbf{r},t,\mathbf{k})}$ of
the electrons cannot satisfy the requirements of a periodic gauge along the
direction of the externally applied electric field, that is
${\Phi_{\varepsilon_{a}}(\mathbf{r},t,\mathbf{k})\neq\Phi_{\varepsilon_{a}}(\mathbf{r},t,\mathbf{k}+G_{x}\bm{e}_{x})}$,
where the electric field is assumed to be in $\bm{e}_{x}$ direction and
$G_{x}$ is any reciprocal lattice vector in $\bm{e}_{x}$.

Subsequently, we extend the standard outcomes of the modern theory of
polarization by evaluating the electronic polarization change by means of Eqs.
(3) – (4) and by using the adiabatic form of DHFT. In this framework, we
rigorously take into account both contributions, the ones due to localized
states as well as the ones due to delocalized states. By using the DHFT we
calculate the collective electrons’ displacement
${\sum_{i}^{N}\Delta\left\langle\mathbf{r}_{i}\right\rangle}$ as a continuous
$\mathbf{k}$-space formula, which in turn provides the induced electronic
polarization change $\Delta\mathbf{P}$ according to Eq. (3). Surprisingly, we
find that the adiabatically induced electronic polarization change, has two
extra, overlooked, contributions which are absent from the modern theory King-
Smith and Vanderbilt (1993); Ortiz and Martin (1994); Resta (1994); Resta and
Vanderbilt (2007); Resta (2010); Vanderbilt (2018). The one is attributed to a
Berry curvature and is not zero only within ferromagnetic insulators where
time-reversal symmetry is broken, and the other, is a non-Hermitian boundary
contribution that captures the delocalization of the wavefunctions. The non-
Hermitian contribution is by definition a bulk quantity that, however, may
equivalently be evaluated as a boundary quantity due to a symmetry that allows
the bulk integration to be transformed into a boundary integration. When one
takes into account only bulk localized states, that is, the wavefunctions are
zero over the material boundaries, the overlooked non-Hermitian contribution
always goes to zero. This non-Hermitian contribution is very sensitive to the
realistic boundary conditions of the wavefunctions during the time-periodic
process that causes the polarization change. We therefore expect it to be
significant in biased insulators with electron charge accumulation over their
boundaries during the process that causes the polarization change (i.e. during
the spontaneous polarization change experiments in a ferroelectric material,
or during the induced polarization change experiments in a piezoelectric
material).

By using the non-Hermitian contribution to the polarization change, we show
how one can extend the surface-charge theorem of Ref. Vanderbilt and King-
Smith, 1993. Specifically, we argue that a macroscopic surface charge density
change theorem can be formulated by means of

$\Delta\sigma_{\text{surf}}=\Delta\mathbf{P_{b}}\cdot\hat{\mathbf{e}}_{a},$
(5)

where, ${\Delta\mathbf{P_{b}}}$ is the polarization change that is attributed
only to boundary properties of the electrons’ wavefunctions and
$\hat{\mathbf{e}}_{a}$ is the unit vector normal to the surface of the
material. When only bulk localized states are taken into account, the surface
charge change $\Delta\sigma_{\text{surf}}$ Eq. (5) turns to zero by
definition, as expected.

We have organized the article as follows. In Sec. II we review the extended
and the boundary velocity operator definitions. In Sec. III we present the
dynamic Hellmann-Feynman theorem (DHFT) that we analytically derive. By then
using the boundary velocity operator definition and the DHFT, we study the
adiabatic charge pumping in Sec. IV and the counterpart driven and non-
equilibrium charge pumping in Sec. V. We then reconsider and re-study the
adiabatic electronic polarization change in Sec. VI and finally, we summarize
and conclude in Sec. VII. Some very useful limiting forms of the dynamic
Hellmann-Feynman theorem are derived in Appendix B, while some helpful
discrete symmetries of the curvatures are given in Appendix C. It is important
to also note that, for completeness, the recently discussed boundary integral
theory of emergent non-Hermiticities, widely overlooked throughout the history
of quantum mechanics (and first noted in ref. 40), is developed in an
alternative way (and within the formulation developed in the present article)
in Appendix A.

##  II Extended velocity operator

The extended velocity operator ${\mathbf{v}_{ext}}$ was first introduced in
Ref. Kyriakou and Moulopoulos, 2018 in order to take into account an anomaly
of the position operator $\mathbf{r}$ that creates a boundary, non-Hermitian
effect term, in the Ehrenfest theorem. Namely, by defining the velocity
operator as the operator that its expectation value is always equal with the
rate of change of the position operator expectation value, it turns out that
the velocity must be defined in an extended theoretical framework

$\left\langle\mathbf{v}_{ext}\right\rangle=\frac{d}{dt}\\!\left\langle{\mathbf{r}}\right\rangle,$
(6)

where the extended velocity operator has two terms

$\mathbf{v}_{ext}=\mathbf{v}+\mathbf{v}_{b},$ (7)

and the expectation values in Eq. (6) are evaluated over the entire volume $V$
of the system.

The first term

$\mathbf{v}=\frac{i}{\hbar}\left[H(\mathbf{r},t),\mathbf{r}\right]$ (8)

is the standard velocity operator, while the second term

$\mathbf{v}_{b}=\frac{i}{\hbar}\\!\left({H(\mathbf{r},t)}^{+}-H(\mathbf{r},t)\right)\\!\mathbf{r}$
(9)

is a new object that we have termed the boundary velocity operator.

The introduction of this new operator $\mathbf{v}_{b}$ is rather naturally
motivated by Refs. Esteve, 1986, 2002 and its expectation value is not zero
only whenever the position operator becomes anomalous due to the non-Hermitian
effect, in which case there are paradoxes first noted in Ref. Hill, 1973.
Specifically, the position operator becomes an anomalous operator whenever the
state $\mathbf{r}\Psi(\mathbf{r},t)$ does not lie within the given Hilbert
space of the Hermitian Hamiltonian ${H(\mathbf{r},t)}$. As consequence,
whenever periodic boundary conditions (PBCs) are adopted for the wavefunction
$\Psi(\mathbf{r},t)$ the boundary velocity expectation value
${\left\langle\mathbf{v}_{b}\right\rangle}$ is not zero. Although the
expectation value of the boundary velocity operator Eq. (9) given by

$\left\langle\mathbf{v}_{b}\right\rangle=\frac{i}{\hbar}(\left\langle
H(\mathbf{r},t)\Psi(t)|\,\mathbf{r}\Psi(t)\right\rangle-\left\langle\Psi(t)|\,H(\mathbf{r},t)\mathbf{r}\Psi(t)\right\rangle)$
(10)

is by definition a bulk quantity, through a space-volume integration in
position representation (assuming for the moment a 3D system), it can always
and equivalently be evaluated as a boundary quantity due to the symmetric
structure of the integrands that allows an integration by parts. (Everything
however carries out in an obvious manner to other dimensionalities as well,
with the results that follow still being valid but with appropriate
adjustments, with boundaries being lines (edges) in a 2D and points (ends) in
a 1D system).

In this respect, by working in position representation, for real scalar and
vector potentials and after a straightforward integration by parts, the
expectation value of Eq. (9) can be given in the form

$\left\langle\mathbf{v}_{b}\right\rangle=-\oiint_{S}\mathbf{r}\,(\mathbf{J}_{pr}(\mathbf{r},t)\\!\cdot\\!d\mathbf{S})+\frac{i\hbar}{2m}\oiint_{S}\left|\Psi(\mathbf{r},t)\right|^{2}\\!d\mathbf{S}$
(11)

with $S$ being the terminated boundary surface of the system where the
boundary conditions are imposed, and
${\mathbf{J}_{pr}(\mathbf{r},t)\\!=\\!\text{Re}[\Psi(\mathbf{r},t)^{\dagger}\mathbf{v}\,\Psi(\mathbf{r},t)]}$
is the standard local probability current density. It is to be noted that,
apart from the imaginary part (which will always be cancelled out in the total
velocity expectation value as we shall see below), the boundary velocity
expectation value ${\left\langle\mathbf{v}_{b}\right\rangle}$ has the meaning
of a position weighted probability flux through the face of the material.

The general form of Eq. (11) can be further reduced for periodic systems.
Specifically, by assuming a Bloch wavefunction
${\displaystyle\Psi_{\mathbf{k}}(\mathbf{r},t)=\frac{1}{\sqrt{N}}\,e^{\displaystyle
i\mathbf{k.r}}u_{\mathbf{k}}(\mathbf{r},t)}$ and cell normalization convention
${\displaystyle\left\langle\Psi_{\mathbf{k}}(t)|\Psi_{\mathbf{k}}(t)\right\rangle=\left\langle
u_{\mathbf{k}}(t)|u_{\mathbf{k}}(t)\right\rangle_{cell}=1}$ (where $N$ is the
total number of the unit cells enclosed within the volume $V$ of the system),
then, Eq. (10) truncates into the form

$\left\langle\mathbf{v}_{b}\right\rangle=\frac{i}{\hbar}\left(\left\langle
H_{k}\,u_{\mathbf{k}}(t)|\,\mathbf{r}\,u_{\mathbf{k}}(t)\right\rangle_{cell}-\left\langle
u_{\mathbf{k}}(t)|\,H_{k}\,\mathbf{r}\,u_{\mathbf{k}}(t)\right\rangle_{cell}\right)$
(12)

where ${H_{k}}$ is given by
${H_{k}=e^{\displaystyle-i\mathbf{k.r}}H(\mathbf{r},t)e^{\displaystyle
i\mathbf{k.r}}}$. In deriving Eq. (12) we have used the fact that
${u_{\mathbf{k}}(\mathbf{r},t)}$ lies within the domain of definition of the
Hamiltonian ${H_{k}}$, that is ${\left\langle
H_{k}\,u_{\mathbf{k}}(t)|\,u_{\mathbf{k}}(t)\right\rangle-\left\langle
u_{\mathbf{k}}(t)|\,H_{k}\,u_{\mathbf{k}}(t)\right\rangle=0}$. By then
exploiting the symmetry of the integrands and performing integration by parts,
Eq. (12) takes the simplified form

$\left\langle\mathbf{v}_{b}\right\rangle=-\oiint_{cell}\\!\\!\\!\mathbf{r}\,(\,\mathbf{J}_{pr}(\mathbf{r},t,\mathbf{k})\\!\cdot\\!d\mathbf{S}\,)$
(13)

that is valid for periodic systems.

The first term of Eq. (11) can be seen as a position-weighted probability flux
through the boundaries of the system, while the second and purely imaginary
part, cancels a possible imaginary remnant part of the standard velocity
operator expectation which is given by

$\left\langle\mathbf{v}\right\rangle=\iiint_{V}\mathbf{J}_{pr}(\mathbf{r},t)dV-\frac{i\hbar}{2m}\oiint_{S}\left|\Psi(\mathbf{r},t)\right|^{2}\\!d\mathbf{S}$
(14)

(by adding Eq. (11) and Eq. (14), that is taking
${\left\langle\mathbf{v}_{ext}\right\rangle=\left\langle\mathbf{v}\right\rangle+\left\langle\mathbf{v}_{b}\right\rangle}$,
we see that ${\left\langle\mathbf{v}_{ext}\right\rangle}$ is always a real
quantity, as should be expected from Eq. (6)). Incidentally, these elementary
facts (never been pointed out in the literature to the best of our knowledge)
exemplify the general inadequacy (for certain boundary conditions) of the
standard velocity operator in basic quantum mechanics; most importantly for
the present work, they display the necessity of defining and using the
augmented velocity operator $\mathbf{v}_{ext}$ as a first example (with more
demonstrations of this necessity to be discussed in the more focused solid
state applications that follow). A point that we wish to make is that this
need for the use of $\mathbf{v}_{ext}$ is required for serious investigations,
especially if one wants to have an analytical machinery, that can address
issues with possible subtleties (such as the ones that we will confront with
in the present work).

The extended velocity operator guarantees that, each Hamiltonian’s stationary
eigenstate indexed by $n$, will always produce zero displacement
(irrespectively of the static potentials) for the electron

$\frac{d}{dt}\\!\left\langle{\mathbf{r}}\right\rangle_{n}=\left\langle\mathbf{v}_{ext}\right\rangle_{n}=\left\langle\mathbf{v}\right\rangle_{n}+\left\langle\mathbf{v}_{b}\right\rangle_{n}=0,$
(15)

as expected from the trivial fact that the position operator expectation value
$\left\langle{\mathbf{r}}\right\rangle_{n}$ is a static quantity with respect
to any stationary state (since the phase factors containing time multiply to
1). However, for a stationary and extended plane wave state of a free electron
of mass $m$ with well-defined momentum $\hbar\mathbf{k}$ in a finite volume
$V$, Eq. (15) is violated (hence a paradox) whenever the boundary velocity
expectation value ${\left\langle\mathbf{v}_{b}\right\rangle_{n}}$ is not taken
into account (as the standard velocity
${\left\langle\mathbf{v}_{n}\right\rangle_{n}}=\hbar\mathbf{k}/m$ is generally
nonzero). This apparent paradox is bypassed within the extended velocity
operator definition, as it turns out that the boundary velocity contributes a
quantity of equal magnitude and opposite sign to the bulk electrons’ velocity
$\left\langle\mathbf{v}\right\rangle_{n}$ resulting in zero displacement
${\displaystyle\Delta{\left\langle\mathbf{r}\right\rangle}_{n}=0}$ at every
instant $t$ for the assumed stationary state. In this framework, Eq. (15) can
be used to form a bulk-boundary correspondence in a general sense for every
stationary state, namely
${\left\langle\mathbf{v}\right\rangle_{n}=-\left\langle\mathbf{v}_{b}\right\rangle_{n}}$
; this is an example, therefore, of a bulk formulation that properly takes
into account boundary currents that are rigorously related to the bulk band
structure. We will see in what follows that similar line of reasoning applied
to (and extending) the Hellmann-Feynman theorem will lead to the same
conclusion about various systems of modern interests (see discussions after
Eq. (37)).

By taking into account Eqs. 6 \- 7 and Eq. 11, we can formally define a non-
Hermitian displacement of the electron through the sample boundaries
${\Delta\left\langle{\mathbf{r}}\right\rangle\\!|_{bound}}$ by using the
identification

$\left\langle\mathbf{v}_{b}\right\rangle=\frac{d}{dt}\\!\left\langle{\mathbf{r}}\right\rangle\\!|_{bound},$
(16)

which results into

$\displaystyle\Delta\left\langle{\mathbf{r}}\right\rangle\\!|_{bound}=\int_{0}^{t}\left\langle\mathbf{v}_{b}\right\rangle
dt^{\prime}.$ (17)

With the aid of Eq. (7) – (9), the extended velocity operator can be recast in
the form

$\mathbf{v}_{ext}=\frac{i}{\hbar}({H(\mathbf{r},t)}^{+}\mathbf{r}-\mathbf{r}H(\mathbf{r},t))$
(18)

which, by direct application of Eq. (7) with respect to the Hamiltonian’s
instantaneous eigenstates labelled by $n$ and $m$, gives the off-diagonal
position matrix elements

$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\displaystyle\frac{i}{\hbar}(E_{m}(t)-E_{n}(t))\left\langle
m(t)|\mathbf{r}|n(t)\right\rangle$ $\displaystyle=$ $\displaystyle\left\langle
m(t)|\mathbf{v}|n(t)\right\rangle$ $\displaystyle+$ $\displaystyle\left\langle
m(t)|\mathbf{v}_{b}|n(t)\right\rangle$ (19)

that are proportional to the electron’s transition dipole moment. We see,
therefore, that the emission and absorption of photons can be rigorously
related with boundary properties owing to the off-diagonal boundary velocity
matrix elements. It is thus essential to point out that, Eq. (II), without the
boundary velocity matrix elements taken into account, is generally inadequate
to be employed as a starting point formula - as often done in the modern
theory (i.e. in Ref. Vanderbilt, 2018) - in order to theoretically study the
electronic polarization change and the associated induced current, especially
if subtle issues need to be addressed for all types of realistic boundary
conditions.

##  III Derivation of the Dynamic Hellmann-Feynman theorem (DHFT)

We consider a real vector parameter ${\mathbf{R}}$ that has an arbitrary time-
dependence (without any adiabatic approximation involved), namely,
${\mathbf{R}=\bm{\mathcal{R}}(t,\mathbf{R_{o}})}$ where ${\mathbf{R_{o}}}$ is
the initial value of the parameter satisfying
${\mathbf{R_{o}}=\bm{\mathcal{R}}(0,\mathbf{R_{o}})}$. Therefore, the
parameter satisfies the general equation of motion
${\displaystyle\mathbf{R}=\mathbf{R_{o}}+\int_{0}^{\displaystyle
t}\frac{\partial\mathbf{R}}{\partial t^{\prime}}dt^{\prime}}$, and its time
derivative is given by ${\displaystyle\frac{\partial\mathbf{R}}{\partial
t}=\frac{\partial\bm{\mathcal{R}}(t,\mathbf{R_{o}})}{\partial t}}$. The
theorem that we are about to prove is for a continuous vector parameter
${\mathbf{R}}$, the initial value ${\mathbf{R_{o}}}$ been therefore assumed to
have continuous values as well. The Hamiltonian of the system
${H(t,\mathbf{R})}$, apart from the implicit time-dependence (via the
parameters) may also have an arbitrary explicit time-dependence (something
that is often not discussed). The derivation that is given below owes its
existence to the Hamiltonian being the generator of time evolution of quantum
states. We provide the derivation for a single-particle state while the
generalization to a many-particle system is straightforward. A particle’s
motion is generally encoded in its normalized time-dependent state
${\left|\Psi(t,\mathbf{R})\right\rangle}$ which evolves either by the time-
dependent Schrödinger equation for non-relativistic and spinless particle, or
by the time-dependent Dirac equation for spinful particle. We assume for
simplicity a one particle quantum system. The motion of the particle is
described by a general state, not necessarily an eigenstate of the Hamiltonian
nor a localized state (such as a narrow wave packet). The system is assumed to
be closed
${\left\langle\Psi(t,\mathbf{R})|\Psi(t,\mathbf{R})\right\rangle=1}$, and the
quantum state time evolution is determined by the time-dependent equation

$i\hbar\frac{d}{dt}\left|\Psi(t,\mathbf{R})\right\rangle=H(t,\mathbf{R})\left|\Psi(t,\mathbf{R})\right\rangle,$
(20)

where the Hamiltonian is either of Schrödinger or Dirac type. The time
derivative in Eq. (20) is the total time derivative given by

$\frac{d}{dt}=\frac{\partial}{\partial t}+\frac{\partial\mathbf{R}}{\partial
t}\\!\cdot\\!\bm{\nabla_{\mathbf{R}}}$ (21)

where
${\displaystyle\bm{\nabla_{\mathbf{R}}}=\sum_{i=1}^{3}\mathbf{e}_{i}\frac{\partial}{\partial
R_{i}}}$. The initial value of the parameter ${\mathbf{R_{o}}}$ that
implicitly enters Eq. (20) can be used to label the quantum states
${\left|\Psi(t,\mathbf{R})\right\rangle}$. The expectation value of the
Hamiltonian

$\left\langle\Psi(t,\mathbf{R})\right|H(t,\mathbf{R})\left.\\!\Psi(t,\mathbf{R})\right\rangle=E(t,\mathbf{R})$
(22)

can be viewed as the instantaneous time-dependent “energy” of the particle
${E(t,\mathbf{R})}$. Differentiation with respect to the parameter
${\mathbf{R}}$ of both sides of Eq. (22) gives

$\left\langle\Psi\,\right|\left(\bm{\nabla_{\mathbf{R}}}H\right)\left.\Psi\right\rangle=\bm{\nabla_{\mathbf{R}}}E-\left\langle\bm{\nabla_{\mathbf{R}}}\Psi\,|\,H\Psi\right\rangle-\left\langle\Psi|\,H\,\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle.$
(23)

Taking now into account that the parameter gradient operator
${\bm{\nabla_{\mathbf{R}}}}$ can generally be an anomalous operator, that is,
the states ${\left|\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle}$ may not fall
within the domain of the Hermitian Hamiltonian, something that can be
expressed by the non-trivial inequality ${\left\langle
H\Psi\,|\,\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle=\left\langle\Psi\,|\,H^{+}\,\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle\neq\left\langle\Psi\,|\,H\,\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle}$,
we recast Eq. (23) into the form

$\displaystyle\left\langle\Psi\,\right|\left(\bm{\nabla_{\mathbf{R}}}H\right)\left.\Psi\right\rangle$
$\displaystyle=$
$\displaystyle\bm{\nabla_{\mathbf{R}}}E-\left\langle\bm{\nabla_{\mathbf{R}}}\Psi\,|\,H\Psi\right\rangle-\left\langle
H\Psi|\,\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle$ (24)
$\displaystyle+\left\langle\Psi\,|\,\left(H^{+}-H\right)\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle.$

By then using Eq. (20) in Eq. (24) we find

$\displaystyle\left\langle\Psi\,\right|\left(\bm{\nabla_{\mathbf{R}}}H\right)\left.\Psi\right\rangle$
$\displaystyle=$
$\displaystyle\bm{\nabla_{\mathbf{R}}}E\,+\,\left\langle\Psi\,|\left(H^{+}-H\right)\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle$
$\displaystyle-i\hbar\left(\left\langle\bm{\nabla_{\mathbf{R}}}\Psi\,|\,\frac{d\Psi}{dt}\right\rangle-\left\langle\frac{d\Psi}{dt}|\,\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle\right).$

Applying then the total time derivative Eq. (21) on Eq. (III) we obtain

$\displaystyle\left\langle\Psi\,\right|\left(\bm{\nabla_{\mathbf{R}}}H\right)\left.\Psi\right\rangle$
$\displaystyle=$
$\displaystyle\bm{\nabla_{\mathbf{R}}}E\,+\,\left\langle\Psi\,|\left(H^{+}-H\right)\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle$
$\displaystyle-i\hbar\left(\left\langle\bm{\nabla_{\mathbf{R}}}\Psi\,|\,\frac{\partial\Psi}{\partial
t}\right\rangle-\left\langle\frac{\partial\Psi}{\partial
t}|\,\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle\right)$
$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!-i\hbar\left(\left\langle\bm{\nabla_{\mathbf{R}}}\Psi\,|\,\frac{\partial\mathbf{R}}{\partial
t}\\!\cdot\\!\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle-\left\langle\frac{\partial\mathbf{R}}{\partial
t}\\!\cdot\\!\bm{\nabla_{\mathbf{R}}}\Psi|\,\bm{\nabla_{\mathbf{R}}}\Psi\right\rangle\right).$

Using then the vector identity
${\displaystyle\frac{\partial\mathbf{R}}{\partial
t}\times(\mathbf{A}\times\mathbf{B})=\mathbf{A}\left(\frac{\partial\mathbf{R}}{\partial
t}\\!\cdot\\!\mathbf{B}\right)-\left(\frac{\partial\mathbf{R}}{\partial
t}\\!\cdot\\!\mathbf{A}\right)\mathbf{B}}$ in the last term on the right hand
side of Eq. (III), we finally obtain our new dynamic Hellmann-Feynman theorem
(DHFT), namely

$\displaystyle\left\langle\bm{O}(t,\mathbf{R})\right\rangle$ $\displaystyle=$
$\displaystyle\bm{\nabla_{\mathbf{R}}}E(t,\mathbf{R})\,+\,\bm{\mathcal{S}}(t,\mathbf{R})$
(27)
$\displaystyle-\,\hbar\,\bm{\mathcal{E}}(t,\mathbf{R})-\hbar\,\frac{\partial\mathbf{R}}{\partial
t}\\!\times\bm{\mathcal{B}}(t,\mathbf{R}),$

where

$\left\langle\bm{O}(t,\mathbf{R})\right\rangle=\left\langle\Psi(t,\mathbf{R})\,\right|\left(\bm{\nabla_{\mathbf{R}}}H(t,\mathbf{R})\right)\left.\Psi(t,\mathbf{R})\right\rangle$
(28)

is the observable in quest, whereas

$\bm{\mathcal{B}}(t,\mathbf{R})=i\left\langle\bm{\nabla_{\mathbf{R}}}\Psi(t,\mathbf{R})\right|\times\left|\bm{\nabla_{\mathbf{R}}}\Psi(t,\mathbf{R})\right\rangle$
(29)

and

$\displaystyle\bm{\mathcal{E}}(t,\mathbf{R})$ $\displaystyle=$ $\displaystyle
i\left\langle\bm{\nabla_{\mathbf{R}}}\Psi(t,\mathbf{R})\,|\,\frac{\partial\Psi(t,\mathbf{R})}{\partial
t}\right\rangle$ (30)
$\displaystyle-i\left\langle\frac{\partial\Psi(t,\mathbf{R})}{\partial
t}|\,\bm{\nabla_{\mathbf{R}}}\Psi(t,\mathbf{R})\right\rangle$

are generalized Berry curvatures (of both ‘magnetic’ and ‘electric’ type) in
the ${t\\!\times\\!\mathbf{R}}$ space, each of them being a real quantity as
manifest from their definition.

The ${\bm{\mathcal{S}}(t,\mathbf{R})}$ term is an emergent non-Hermitian
contribution to the above theorem which is by definition a bulk quantity.
However, as shown in Appendix A, by working in position representation and
assuming real scalar and vector potentials, the
${\bm{\mathcal{S}}(t,\mathbf{R})}$ term is always transformed (due to
symmetric form of the integrands) and evaluated by means of a boundary
integral over the system’s boundaries and is given (assuming a 3D system,
although this can trivially be generalized) by

$\displaystyle\bm{\mathcal{S}}(t,\mathbf{R})$ $\displaystyle=$
$\displaystyle\left\langle\Psi(t,\mathbf{R})\,|\,\left(H(t,\mathbf{R})^{+}-H(t,\mathbf{R})\right)\bm{\nabla_{\mathbf{R}}}\Psi(t,\mathbf{R})\right\rangle$
$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!=\frac{i\hbar}{2}\oiint_{S}\mathbf{n}\\!\cdot\\!\left(\,(\mathbf{v}\,\Psi(t,\mathbf{R}))^{+}+\Psi(t,\mathbf{R})^{+}\,\mathbf{v}\,\right)\\!\,\bm{\nabla_{\mathbf{R}}}\Psi(t,\mathbf{R})\,dS,$

where
${\displaystyle\mathbf{v}=\frac{i}{\hbar}\left[H(t,\mathbf{R}),\mathbf{r}\right]}$
is the standard velocity operator and ${\mathbf{n}}$ is the unit vector that
is locally normal to the surface that encloses the system. As shown in
Appendix A, the boundary integral form of the non-Hermitian contribution
${\bm{\mathcal{S}}(t,\mathbf{R})}$ given by Eq. (III), is valid for every
energy scale of the particle, that is, either the time evolution is governed
by a (i) spinless non-relativistic Hamiltonian, or (ii) a spinful Dirac
Hamiltonian, or a (iii) a spinful Hamiltonian in the non-relativistic limit.
In all above cases, the non-Hermitian contribution
${\bm{\mathcal{S}}(t,\mathbf{R})}$ is extremely sensitive to the realistic
boundary conditions on the wavefunctions. Moreover,
${\bm{\mathcal{S}}(t,\mathbf{R})}$ is a real quantity whenever the operator
${\bm{\nabla_{\mathbf{R}}}H}$ behaves as a Hermitian operator. This can be
shown by starting from a closed system assumption
${\left\langle\Psi(t,\mathbf{R})\right|H\left.\\!\Psi(t,\mathbf{R})\right\rangle-\left\langle
H\Psi(t,\mathbf{R})\right|\left.\\!\Psi(t,\mathbf{R})\right\rangle=0}$, and
then differentiate both sides with respect to the parameter ${\mathbf{R}}$
giving

$\bm{\mathcal{S}}(t,\mathbf{R})=\bm{\mathcal{S}}(t,\mathbf{R})^{\displaystyle*}+\left(\,\left\langle\Psi\,\right|\left(\bm{\nabla_{\mathbf{R}}}H\right)\left.\Psi\right\rangle-\left\langle\left(\bm{\nabla_{\mathbf{R}}}H\right)\,\Psi\,\right|\left.\\!\Psi\right\rangle\,\right).$
(32)

Interestingly, whenever the assumed Hamiltonian does not have explicit
parameter dependence, but the considered states do depend on the parameters
(as i.e. it happens for the original $H$ in a periodic system (before passing
to the Bloch Hamiltonian), which does not depend on the Bloch wavevector k,
although its eigenstates, the standard Bloch states, do), then, the non-
Hermitian contribution is always a real quantity; in the example given above,
the non-Hermitian contribution evaluated with respect to a Bloch state by
identifying ${\mathbf{R}=\mathbf{k}}$ is always a real quantity. Eq. (32) will
be crucial in properties such as the surface-charge theorem discussed towards
the end of the paper (see Eq. (100)).

Due to the structure of the DHFT, the last term on the right side of Eq. (27)
gives a contribution that is transverse to the direction of the parameter
variation ${\partial\mathbf{R}}$, namely,
${\displaystyle\left\langle\bm{O}(t,\mathbf{R})\right\rangle_{tran}=-\hbar\,\frac{\partial\mathbf{R}}{\partial
t}\\!\times\bm{\mathcal{B}}(t,\mathbf{R})}$ that satisfies
${\left\langle\bm{O}(t,\mathbf{R})\right\rangle_{tran}\cdot\partial\mathbf{R}=0}$.
This occurs in IQHE, but this term is also expected to contribute
significantly to non-linear quantum processes whenever the assumed parameter
varies with time, i.e within non-linear materials where longitudinal as well
as transverse polarization changes occurs simultaneously.

###  III.1 Application of the DHFT to a state with ansatz form
$\displaystyle\left|\Psi(t,\mathbf{R})\right\rangle=e^{\displaystyle
i\Theta(t,\mathbf{R})}\left|\Phi(t,\mathbf{R})\right\rangle$.

When the considered quantum state, assumed to be a solution of Eq.(20), has
the ansatz form

$\left|\Psi(t,\mathbf{R})\right\rangle=e^{\displaystyle
i\Theta(t,\mathbf{R})}\left|\Phi(t,\mathbf{R})\right\rangle,$ (33)

then, as shown in Appendix B, the DHFT Eq. (27) retains its structure form but
all involved quantities are now evaluated with respect to the quantum state
$\left|\Phi(t,\mathbf{R})\right\rangle$, namely

$\displaystyle\left\langle\bm{O}(t,\mathbf{R})\right\rangle$ $\displaystyle=$
$\displaystyle\bm{\nabla_{\mathbf{R}}}E_{\Phi}(t,\mathbf{R})\,+\,\bm{\mathcal{S}}_{\Phi}(t,\mathbf{R})$
(34)
$\displaystyle\\!\\!\\!\\!\\!\\!-\hbar\,\frac{\partial\mathbf{R}}{\partial
t}\\!\times\bm{\mathcal{B}}_{\Phi}(t,\mathbf{R})\,-\,\hbar\,\bm{\mathcal{E}}_{\Phi}(t,\mathbf{R}).$

In other words, the observables in quest
${\left\langle\bm{O}(t,\mathbf{R})\right\rangle}$ given by Eq. (28) are
insensitive to phases ${e^{\displaystyle i\Theta(t,\mathbf{R})}}$ that do not
depend explicitly on the position coordinate ${\mathbf{r}}$. Such simple
applications of the DHFT are described below:

For an adiabatically evolved eigenstate of a time-dependent Hamiltonian, we
replace in Eq.(33),
${\displaystyle\Theta(t,\mathbf{R})\mapsto\Theta_{n}(t,\mathbf{R})}$ and
${\left|\Phi(t,\mathbf{R})\right\rangle\mapsto\left|n(t,\mathbf{R})\right\rangle}$,
where ${\Theta_{n}(t,\mathbf{R})}$ is the total (sum of dynamic and geometric)
phase and ${\left|n(t,\mathbf{R})\right\rangle}$ is the instantaneous
eigenstate of the Hamiltonian satisfying the eigenvalue equation
${H(t,\mathbf{R})\left|n(t,\mathbf{R})\right\rangle=E_{n}(t,\mathbf{R})\left|n(t,\mathbf{R})\right\rangle.}$

For a stationary eigenstate of a static Hamiltonian where the assumed
parameter is static ${\displaystyle\frac{\partial\mathbf{R}(t)}{\partial
t}=0}$, we replace in Eq.(33),
${\displaystyle\Theta(t,\mathbf{R})\mapsto\Theta_{n}(t,\mathbf{R})}$ and
${\left|\Phi(t,\mathbf{R})\right\rangle\mapsto\left|n(\mathbf{R})\right\rangle}$,
where ${\Theta_{n}(t,\mathbf{R})}$ is the dynamic phase and
${\left|n(\mathbf{R})\right\rangle}$ is the instantaneous eigenstate of the
Hamiltonian satisfying the eigenvalue equation
${H(\mathbf{R})\left|n(\mathbf{R})\right\rangle=E_{n}(\mathbf{R})\left|n(\mathbf{R})\right\rangle.}$
The Berry curvature ${\bm{\mathcal{E}}_{n}(\mathbf{R})}$ is zero by definition
according to Eq. (137) owing to ${\left|n(\mathbf{R})\right\rangle}$ that does
not have explicit time dependence, whereas the Berry curvature
${\bm{\mathcal{B}}_{n}(\mathbf{R})}$ is not zero but it does not contribute
into the DHFT Eq. (34) due to zero velocity of the parameter.

For a non-equilibrium and nonadiabatic quantum state where the assumed
parameter may have arbitrary high magnitude of velocity
${\displaystyle\frac{\partial\mathbf{R}(t)}{\partial t},}$ that, however,
evolves in time by a time-periodic Hamiltonian
${H(t+T,\mathbf{R}(t+T))=H(t,\mathbf{R}(t))}$, where
${\mathbf{R}(t+T)=\mathbf{R}(t)}$ and ${T}$ is the period of driving, we
replace in Eq.(33)
${\displaystyle\Theta(t,\mathbf{R})\mapsto-\frac{\varepsilon_{a}\,t}{\hbar}}$
and
${\left|\Phi(t,\mathbf{R})\right\rangle\mapsto\left|\Phi_{a}(t,\mathbf{R})\right\rangle}$,
where ${\left|\Phi_{a}(t,\mathbf{R})\right\rangle}$ is the Floquet mode which
is periodic in time
${\left|\Phi_{a}(t,\mathbf{R})\right\rangle=\left|\Phi_{a}(t+T,\mathbf{R})\right\rangle}$
and satisfies the eigenvalue equation
${\displaystyle(H(t,\mathbf{R})-i\hbar\frac{d}{dt}\,)\Phi_{a}(\mathbf{r},t,\mathbf{k})=\varepsilon_{a}\,\Phi_{a}(\mathbf{r},t,\mathbf{k}),}$
with ${\varepsilon_{a}}$ being the static quasienergy. In conclusion, we argue
that the DHFT Eq. (34) can be applied to general topological band theory.

###  III.2 The need for the non-Hermitian term in the DHFT

The simplest example where one can ascertain the necessity of presence of the
boundary non-Hermitian term ${\bm{\mathcal{S}}(t,\mathbf{R})}$ is a free
electron motion in 1D. In this respect, we assume the Hamiltonian
${\displaystyle H=-\frac{\hbar^{2}}{2m}\frac{d^{2}}{dx^{2}}}$ which does not
have any parameter dependence. We consider an electron’s motion described by
the eigenstate ${\displaystyle\Phi_{k}(x)=\frac{1}{\sqrt{L}}e^{\displaystyle
ikx}}$ of the Hamiltonian with well defined energy ${\displaystyle
E(k)=\frac{\hbar^{2}k^{2}}{2m}}$, where $L$ is the length of the system. We
identify the static electron’s 1D wave vector $k$ (assumed to take continuous
values) as the parameter ${R\equiv k}$ and examine the necessity of a nonzero
${S(k)}$.

Specifically, by employing the standard static HF theorem without the non-
Hermitian term, that is, the 1D analog of Eq. (148) derived in Appendix B.3
without the boundary term ${S(k)}$, will eventually lead to a paradox since
${\displaystyle\left\langle\frac{dH}{dk}\right\rangle=0}$ while
${\displaystyle\frac{dE(k)}{dk}=\frac{d}{dk}\left(\frac{\hbar^{2}k^{2}}{2m}\right)=\frac{\hbar^{2}k}{m}\neq
0}$. The paradox is resolved whenever one takes into account the boundary non-
Hermitian term ${S(k)}$. Specifically, the 1D analogue of Eq. (148) is given
by

$0=\frac{dE(k)}{dk}+S(k),$ (35)

where ${S(k)}$ is in 1D truncated to a two point formula (cf. Eq. (138)) given
by

$S(k)=\frac{i\hbar}{2}\left[\left(\,(\text{v}\,\Phi_{k}(x))^{\displaystyle*}+\Phi_{k}(x)^{\displaystyle*}\,\text{v}\,\right)\\!\,\frac{d\Phi_{k}(x)}{dk}\right]_{x}^{x+L}.$
(36)

We substitute in Eq. (36) the standard velocity operator
${\displaystyle\text{v}=-\frac{i\hbar}{m}\frac{\partial}{\partial x}}$, as
well as $\text{v}\,\Phi_{k}(x)=\displaystyle\frac{1}{\sqrt{L}}\frac{\hbar
k}{m}e^{\displaystyle ikx}=\frac{\hbar k}{m}\Phi_{k}(x)$ and
${\displaystyle\frac{d\Phi_{k}(x)}{dk}=\frac{1}{\sqrt{L}}i\,x\,e^{\displaystyle
ikx}=i\,x\,\Phi_{k}(x)}$. Eq. (36) then takes the value

$\displaystyle S(k)$ $\displaystyle=$
$\displaystyle\frac{i\hbar}{2}\left[\left(\,\frac{\hbar
k}{m}\Phi_{k}(x)^{\displaystyle*}-\frac{i\hbar}{m}\Phi_{k}(x)^{\displaystyle*}\,\frac{\partial}{\partial
x}\,\right)\\!\,i\,x\,\Phi_{k}(x)\right]_{x}^{x+L}$ $\displaystyle=$
$\displaystyle\frac{i\hbar}{2}\left[\frac{i\,\hbar\,k}{m}|\Phi_{k}(x)|^{2}x+\frac{\hbar}{m}|\Phi_{k}(x)|^{2}\right.$
$\displaystyle\left.+i\,x\,\Phi_{k}^{\displaystyle*}\,(\text{v}\,\Phi_{k}(x))\frac{}{}\right]_{x}^{x+L}$
$\displaystyle=$
$\displaystyle\frac{i\hbar}{2}\left[2\frac{i\,\hbar\,k}{m}|\Phi_{k}(x)|^{2}x+\frac{\hbar}{m}|\Phi_{k}(x)|^{2}\right]_{x}^{x+L}$
$\displaystyle=$
$\displaystyle\frac{i\hbar^{2}}{2mL}\left[2\,i\,k\,x+1\frac{}{}\right]_{x}^{x+L}$
$\displaystyle=$
$\displaystyle\frac{i\hbar^{2}}{2mL}(2\,i\,k\,L)=-\frac{\hbar^{2}k}{m},$

which together with ${\displaystyle\frac{dE(k)}{dk}=\frac{\hbar^{2}k}{m}}$
verifies Eq. (35) and guarantees the validity of the static HF theorem.
Whenever, therefore, the momentum derivative operator
${\displaystyle\frac{d}{dk}}$ becomes anomalous, the boundary non-Hermitian
contribution cannot be excluded from the standard, static HF theorem; a
similar kind of anomaly as in Eq. (35), was noticed in Ref. Esteve  _et al._ ,
2010 but, notably, without indicating any usage into conventional band theory.
The above non-Hermitian contribution satisfies the antisymmetry
${S(k)=-S(-k)}$, which, as shown in Appendix C, is a general property of band
theory (recall that the non-Hermitian contribution is always evaluated with
respect to a Bloch state, Eq. (III)) when the Hamiltonian satisfies mirror
symmetry along $x$.

###  III.3 The non-Hermitian contribution as a probe of the dispersion
relation

If we identify the wave vector as a static parameter entering Eq. (34), that
is ${\mathbf{R}=\mathbf{k}}$, then, whenever the initial static Hamiltonian of
a bulk periodic system does not depend explicitly on the wave vector
${\bm{\nabla_{\mathbf{k}}}H(\mathbf{k})=0}$, since the wave vector is only
involved only in the Bloch states ${\Psi_{n}(\mathbf{r},\mathbf{k})}$, one can
in principle use the non-Hermitian boundary term
${\bm{\mathcal{S}}_{n}(\mathbf{k})}$ (where now the label $n$ refers to
consideration of a particular band) as a probe for detecting information about
the energy dispersion relation (through the boundary behavior of the Bloch
functions), namely

$\bm{\nabla_{\mathbf{k}}}E_{n}(\mathbf{k})=-\bm{\mathcal{S}}_{n}(\mathbf{k}).$
(37)

Consequently, by defining the bulk localized eigenstates of the Hamiltonian as
the ones that satisfy ${\bm{\mathcal{S}}_{n}(\mathbf{k})=0}$ for every crystal
momentum ${\mathbf{k}}$ in the Brillouin zone (which in the simplest scenario
occurs when the Bloch wavefunction and all of its derivatives are zero over
the boundaries of the material), then, it is evident from Eq. (37) that these
kind of states will have flat band dispersion relations
${\bm{\nabla_{\mathbf{k}}}E_{n}(\mathbf{k})=0}$, indicating the existence of a
large degree of degeneracy. Notably, flat bands are known to be formed in
conventional band theory of periodic (perfectly ordered) systems Leykam  _et
al._ (2018); Röntgen  _et al._ (2018); Maimaiti  _et al._ (2017); Bodyfelt
_et al._ (2014); Danieli  _et al._ (2015); Ramachandran  _et al._ (2017) and
are attributed to compact localized states that are created due to destructive
interference, but, their formation still lacks a comprehensive generating
principle Maimaiti  _et al._ (2017). In this framework we argue that, whenever
the non-Hermitian contribution of a given band is zero, it signals the
creation of a flat band.

Moreover, by employing Eq. (37), one can theoretically explain the changing of
the dispersion relation slope in relation to the localization of the involved
states. For example, in the Integer Quantum Hall Effect, the bulk localized
states are described by flat bands (Landau levels) $E_{n}(\mathbf{k})\equiv
E_{n}=\textit{constant}$, while the states near the boundaries have dispersion
relations with non-zero slopes that are approximately linear Halperin (1982)
(this linearity being with respect to guiding center coordinates in real
space, which however for the simple Landau problem in a Landau gauge is
directly proportional to a wave vector value k, hence the linearity being
again in k). This kind of behavior, can qute possibly be theoretically
explained when the non-Hermitian boundary term is taken into account.
Moreover, in Dirac and Weyl semimetals Armitage  _et al._ (2018); Vafek and
Vishwanath (2014); Yan and Felser (2017); Jia  _et al._ (2016); Burkov (2018),
linear dispersion relations in k are found near the vicinity of energy-
degenerate points in the Brillouin zone which represent the bulk electrons’
physics near these degeneracy points in k-space. Although the electrons’
motions are not fully relativistic (apart from correction due to spin-orbit
coupling), these are theoretically explained by employing various model
Hamiltonians ${H(\mathbf{k})}$ which are effectively (formally) relativistic.
We propose that, these emerging relativistic linear dispersion relations that
are accompanied by exotic boundary behaviors such as Fermi arcs, can
alternatively be explained by employing the non-Hermitian boundary term.
Specifically, by taking a Taylor expansion of the non-Hermitian boundary term
${\bm{\mathcal{S}}_{n}(\mathbf{k})}$ near a degeneracy point
${\mathbf{k}_{o}}$ we find

$\displaystyle\bm{\nabla_{\mathbf{k}}}E_{n}(\mathbf{k})$ $\displaystyle=$
$\displaystyle-\bm{\mathcal{S}}_{n}(\mathbf{k}_{o})-\left((\mathbf{k}-\mathbf{k}_{o})\\!\cdot\\!\bm{\nabla_{\mathbf{k}}}\right)\bm{\mathcal{S}}_{n}(\mathbf{k})\arrowvert_{\mathbf{k}=\mathbf{k}_{o}}$
(38) $\displaystyle+O(\mathbf{k}^{2})+\dots.$

With the dispersion relation ${E_{n}(\mathbf{k})}$ profile experimentally
obtained, one can extract information about the boundary term
${\bm{\mathcal{S}}_{n}(\mathbf{k})}$ by keeping only up to first order terms
in the Taylor expansion Eq. (38). By then also employing the Bloch functions
${\Psi_{n}(\mathbf{r},\mathbf{k})}$ in Eq. (138), further information about
the behavior of the wavefunction at the boundary can be extracted.

##  IV Thouless adiabatic charge pump

Let us now re-examine the adiabatic charge pumping of non-interacting spinless
electrons in a slightly different context than the original method of Ref
Thouless, 1983. Rather than assuming a sliding potential
${V((\mathbf{r}-\mathbf{u}\,t)\,\bm{e}_{x})}$, we assume an electric field
along $\bm{e}_{x}$ direction which is periodic in space (over $x$ coordinate
with period equal to the lattice constant) as well as over time $t$ (with
elementary time period equal to $T$). By then using the non-Hermitian boundary
velocity operator and the adiabatic limit of DHFT, we here evaluate the
collective, center of mass, electrons’ displacement in the $\bm{e}_{x}$
direction and over the boundaries ${\sum_{i}^{N}\Delta\left\langle
x_{i}\right\rangle\\!|_{bound}\ /\ {N}}$ without using (first-order) time-
dependent perturbation theory.

By attributing the electric field to a scalar potential only, each spinless
electron’s wavefunction evolves in time by the Hamiltonian

$H(\mathbf{r},t)=\frac{1}{2m}{\mathbf{p}^{2}+e\phi(x,t)}+V_{crys}(\mathbf{r})$
(39)

where $m$ and $-e$ are the mass and charge of the electron respectively,
$\mathbf{p}$ is the canonical momentum operator, and $V_{crys}(\mathbf{r})$ is
the crystal potential. The applied scalar potential is given by

$\phi(x,t)=-\frac{a}{2\pi}E(t)sin(2\pi x/a),$ (40)

where $a$ is the primitive cell length in $x$ direction, and satisfies
${\phi(x+na,t)=\phi(x,t)}$ with $n$ an integer, as well as,
${\phi(x,T)=\phi(x,0)=0}$ and
${\left.\displaystyle\frac{\partial\phi(x,t)}{\partial
t}\right|_{t=0}=\left.\displaystyle\frac{\partial\phi(x,t)}{\partial
t}\right|_{t=T}=0}$. Therefore, the Hamiltonian Eq. (39) is periodic with
respect to the elementary time period $T$ satisfying
${H(\mathbf{r},0)=H(\mathbf{r},T)}$ and ${\left.\displaystyle\frac{\partial
H(\mathbf{r},t)}{\partial t}\right|_{t=0}=\left.\displaystyle\frac{\partial
H(\mathbf{r},t)}{\partial t}\right|_{t=T}=0}$ , and it is also periodic over
position coordinates $H(\mathbf{r}+\mathbf{R},t)=H(\mathbf{r},t)$, where
$\mathbf{R}$ is a Bravais lattice vector.

Provided that the crystal potential $V_{crys}(\mathbf{r})$ is centromeric, the
applied scalar potential crucially breaks the space-inversion symmetry of the
Hamiltonian as well as the mirror symmetry along the $\bm{e}_{x}$ direction
due to ${\phi(x,t)\neq\phi(-x,t)}$. This results to the Berry curvature’s
(defined by Eq. (30) and assumed to be evaluated with respect to an
instantaneous eigenstate of the Hamiltonian) breaking of antisymmetry
${\mathcal{E}_{n,\,x}(t,k_{x},k_{y},k_{z})\neq-\mathcal{E}_{n,\,x}(t,-k_{x},k_{y},k_{z})}$,
which in turn, gives non-zero collective charge pumping. Analogous mirror
symmetry breaking over $x$ direction is involved into the case of the Thouless
sliding scalar potential ${V(x-\upsilon\,t,y,z)\neq V(-x-\upsilon\,t,y,z)}$.
Due to the applied scalar potential, each electron is subject to an applied
microscopic electric field
${\displaystyle\mathbf{E}^{(micro)}(x,t)=-\bm{\nabla}\phi(x,t)=E(t)cos(2\pi
x/a)\,\bm{e}_{x}}$ along $\bm{e}_{x}$ direction that, is periodic over
position coordinate $x$ (with period equal to the lattice constant $a$) and
over time ($E(T)=E(0)=0$). The applied microscopic electric field gives zero
macroscopic (spatially averaged) electric field, therefore no voltage-bias is
present during the pumping. Crucially, as we will show later, if an
intermediate (middle-) time period symmetry is present, that is,
$E(T)=E(0)=E(T/2)$, the charge pumping turns into zero due to time reversal
symmetry.

Each electron’s state evolves in time according to the time-dependent
Schrödinger equation (TDSE)

$i\hbar\frac{d}{dt}\left|\Psi(t)\right\rangle=\left(\frac{1}{2m}{\mathbf{p}^{2}+e\phi(x,t)}+V_{crys}(\mathbf{r})\right)\left|\Psi(t)\right\rangle.$
(41)

We assume that each electron is initially in a non-degenerate ground state,
${\left|\Psi(t_{o})\right\rangle\equiv\left|\psi_{n}(t_{o})\right\rangle}$,
where $n$ labels the band. The ground state energy band is for all moments
${0\leq t\leq T}$ separated by a gap from the upper bands. Each electron’s
ground state wavefunction satisfies periodic boundary conditions over the
material’s boundaries

$\psi_{n}(\mathbf{r}+\mathbf{L},t_{o})=\psi_{n}(\mathbf{r},t_{o}),$

where $\mathbf{L}$ denotes the (vector) length of the material in each
direction. Assuming that the electric field ${\mathbf{E}(x,t)}$ is slowly
changing with time, meaning ${T\rightarrow\infty}$, each ground state evolves
in time adiabatically

$\Psi_{n}(\mathbf{r},t)=e^{\displaystyle
i\Theta_{n}(t)}\psi_{n}(\mathbf{r},t),$

where ${\Theta_{n}(t)}$ is the total phase of the wavefunction (the sum of the
dynamic and the geometric adiabatic phase). Each wavefunction
${\psi_{n}(\mathbf{r},t)}$ satisfies the instantaneous eigenvalue equation

$\left(\frac{1}{2m}{\mathbf{p}^{2}+e\phi(x,t)}+V_{crys}(\mathbf{r})\right)\psi_{n}(\mathbf{r},t)=E_{n}\,\psi_{n}(\mathbf{r},t).$
(42)

By following Ref. 9, we assume that the electron’s wavefunction has the form

$\psi_{n}(\mathbf{r},t,\mathbf{k})=e^{\displaystyle
i\mathbf{k}\\!\cdot\\!\mathbf{r}}u_{n}(\mathbf{r},t,\mathbf{k}),$ (43)

and the ${u_{n}(\mathbf{r},t,\mathbf{k})}$ satisfies periodic boundary
conditions
$u_{n}(\mathbf{r}+\mathbf{L},t,\mathbf{k})=u_{n}(\mathbf{r},t,\mathbf{k})$,
which constrains the allowed wave vector ${\mathbf{k}}$-values by the
condition ${e^{\displaystyle i\mathbf{k}\\!\cdot\\!\mathbf{L}}=1}$. In order
to apply the DHFT given by Eq. (27) we need to define a continuous parameter
${\mathbf{R}}$ and its corresponding velocity
${\displaystyle\frac{\partial\mathbf{R}}{\partial t}}$. The parameter
definition is made by using the assumed form of the quantum state. Namely, by
having in mind that, in the presence of translational invariance of the
Hamiltonian, the wave vector ${\mathbf{k}}$ is a conserved quantity of the
dynamics, we define the parameter ${\mathbf{R}}$ as the static wave vector
${\mathbf{R}\equiv\mathbf{k}}$ which has zero corresponding velocity
${\displaystyle\frac{\partial\mathbf{R}}{\partial
t}\equiv\displaystyle\frac{\partial\mathbf{k}}{\partial t}=0}$. In order for
the application of the theorem to be valid, the parameter ${\mathbf{k}}$ has
to be a continuous one. Therefore, we assume that the length of the system
${\mathbf{L}}$ is infinite, which results to infinitesimal minimum spacing
between the allowed values of the wave vectors ${\Delta\mathbf{k}\rightarrow
0}$. By then noting that the Hamiltonian ${H(\mathbf{r},t)}$ does not depend
on any parameter

$\bm{\nabla_{\mathbf{k}}}H(\mathbf{r},t)=0,$

application of the adiabatic form of DHFT with respect to the state
${\psi_{n}(\mathbf{r},t,\mathbf{k})}$, namely, Eq. (B.2) derived in Appendix
B.2 gives

$0=\bm{\nabla_{\mathbf{k}}}E_{n}(t,\mathbf{k})\,+\,\bm{\mathcal{S}}_{n}(t,\mathbf{k})\,-\,\hbar\,\bm{\mathcal{E}}_{n}(t,\mathbf{k}),$
(44)

where

$\displaystyle\bm{\mathcal{E}}_{n}(t,\mathbf{k})$ $\displaystyle=$
$\displaystyle
i\left\langle\bm{\nabla_{\mathbf{k}}}\psi_{n}(t,\mathbf{k})\,|\,\frac{\partial\psi_{n}(t,\mathbf{k})}{\partial
t}\right\rangle$
$\displaystyle-i\left\langle\frac{\partial\psi_{n}(t,\mathbf{k})}{\partial
t}|\,\bm{\nabla_{\mathbf{k}}}\psi_{n}(t,\mathbf{R})\right\rangle$

is the earlier defined electric field type of Berry curvature, and

$\displaystyle\bm{\mathcal{S}}_{n}(t,\mathbf{k})$ $\displaystyle=$
$\displaystyle\left\langle\psi_{n}(t,\mathbf{k})\,|\,\left(H(\mathbf{r},t)^{+}-H(\mathbf{r},t)\right)\bm{\nabla_{\mathbf{k}}}\psi_{n}(t,\mathbf{k})\right\rangle$
(46) $\displaystyle=$
$\displaystyle\frac{i\hbar}{2}\oiint_{S}\mathbf{n}\\!\cdot\\!\left(\,(\mathbf{v}\,\psi_{n})^{+}+\psi_{n}^{+}\,\mathbf{v}\,\right)\\!\,\bm{\nabla_{\mathbf{k}}}\psi_{n}\,dS$

is the non-Hermitian contribution, where $\mathbf{n}$ is the unit vector that
is locally normal to the surface $S$.

By then using Eq. (43) into Eq. (46), that is, by substituting the identity

$\bm{\nabla}_{\mathbf{k}}\psi_{n}(\mathbf{r},t,\mathbf{k})=i\,\mathbf{r}\,\psi_{n}(\mathbf{r},t,\mathbf{k})\,+\,e^{\displaystyle
i\mathbf{k}\\!\cdot\\!\mathbf{r}}\,\bm{\nabla}_{\mathbf{k}}u_{n}(\mathbf{r},t,\mathbf{k}),$

we find an expression that relates the non-Hermitian contributions

$\bm{\mathcal{S}}_{n}(t,\mathbf{k})=\hbar\left\langle\,\mathbf{v}_{b}(t,\mathbf{k})\,\right\rangle_{n}\,+\,\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k}),$
(47)

where ${\left\langle\,\mathbf{v}_{b}(t,\mathbf{k})\,\right\rangle_{n}}$ is the
non-Hermitian boundary velocity expectation value (that is non-zero when the
position operator $\mathbf{r}$ becomes anomalous) given by Eq. (10).
${\bm{\mathcal{S}}_{n}(t,\mathbf{k})}$ and
${\left\langle\,\mathbf{v}_{b}(t,\mathbf{k})\,\right\rangle_{n}}$ are
evaluated with respect to the Bloch state
${\psi_{n}(\mathbf{r},t,\mathbf{k})}$, while
${\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k})}$ is a non-Hermitian boundary
contribution that is evaluated with respect to the instantaneous cell periodic
states ${u_{n}(\mathbf{r},t,\mathbf{k})}$ (whose value is not zero when the
momentum gradient operator $\bm{\nabla_{\mathbf{k}}}$ becomes anomalous with
respect to the ${u_{n}(\mathbf{r},t,\mathbf{k})}$) and is given by

$\displaystyle\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k})$ $\displaystyle=$
$\displaystyle\left\langle
u_{n}(t,\mathbf{k})\,|\,\left(H_{k}(t,\mathbf{k})^{+}-H_{k}(t,\mathbf{k})\right)\bm{\nabla_{\mathbf{k}}}u_{n}(t,\mathbf{k})\right\rangle$
(48) $\displaystyle=$
$\displaystyle\frac{i\hbar}{2}\oiint_{S}\mathbf{n}\\!\cdot\\!\left(\,(\mathbf{v}_{k}\,u_{n})^{+}+u_{n}^{+}\,\mathbf{v}_{k}\,\right)\\!\,\bm{\nabla_{\mathbf{k}}}u_{n}\,dS,$

where
${\displaystyle\mathbf{v}_{k}=\frac{i}{\hbar}\left[H_{k}(\mathbf{r},t,\mathbf{k}),\mathbf{r}\right]}$
is the standard velocity operator,
${H_{k}(\mathbf{r},t,\mathbf{k})=e^{\displaystyle
i\mathbf{k}\\!\cdot\\!\mathbf{r}}H(\mathbf{r},t)e^{\displaystyle
i\mathbf{k}\\!\cdot\\!\mathbf{r}}}$, and ${u_{n}\equiv
u_{n}(\mathbf{r},t,\mathbf{k})}$.

By substituting Eq. (47) into Eq. (44) we find the non-Hermitian boundary
velocity of the electron given by

$\left\langle\mathbf{v}_{b}(t,\mathbf{k})\right\rangle_{n}=-\frac{1}{\hbar}\bm{\nabla_{\mathbf{k}}}E_{n}(t,\mathbf{k})\,-\,\frac{1}{\hbar}\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k})\,+\,\bm{\mathcal{E}}_{n}(t,\mathbf{k}).$
(49)

Because we have assumed periodic boundary conditions,
${\psi_{n}(\mathbf{r}+\mathbf{L},t,\mathbf{k})=\psi_{n}(\mathbf{r},t,\mathbf{k})}$
as well as
${u_{n}(\mathbf{r}+\mathbf{L},t,\mathbf{k})=u_{n}(\mathbf{r},t,\mathbf{k})}$,
the quantum state ${\bm{\nabla}_{\mathbf{k}}u_{n}(\mathbf{r},t,\mathbf{k})}$
is periodic with respect to space coordinates
${\bm{\nabla}_{\mathbf{k}}u_{n}(\mathbf{r}+\mathbf{L},t,\mathbf{k})=\bm{\nabla}_{\mathbf{k}}u_{n}(\mathbf{r},t,\mathbf{k})}$,
it therefore belongs within the domain of definition of the Hamiltonian
${H_{k}(\mathbf{r},t,\mathbf{k})}$ which results into

$\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k})=0.$ (50)

This can be shown by expanding the cell periodic functions in a Fourier series
over all reciprocal lattice vectors ${\mathbf{G}}$, that is,
${u_{n}(\mathbf{r},t,\mathbf{k})=\sum_{\mathbf{G}}C_{n}(t,\mathbf{k},\mathbf{G})\
e^{\displaystyle-i\mathbf{G}\\!\cdot\\!\mathbf{r}}}$, which gives
${\bm{\nabla}_{\mathbf{k}}u_{n}(\mathbf{r},t,\mathbf{k})=\sum_{\mathbf{G}}\bm{\nabla}_{\mathbf{k}}C_{n}(t,\mathbf{k},\mathbf{G})\
e^{\displaystyle-i\mathbf{G}\\!\cdot\\!\mathbf{r}}}$. Substitution of
${u_{n}(\mathbf{r},t,\mathbf{k})}$ and
${\bm{\nabla}_{\mathbf{k}}u_{n}(\mathbf{r},t,\mathbf{k})}$ in the definition
of the non-Hermitian term, namely
${\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k})=\left\langle
H_{k}\,u_{n}\,|\,\bm{\nabla_{\mathbf{k}}}u_{n}\right\rangle-\left\langle
u_{n}\,|\,H_{k}\,\bm{\nabla_{\mathbf{k}}}u_{n}\right\rangle}$, results in
${\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k})=0}$ (for spinful electrons,
although, when spin-Hall effect is present, this non-Hermitian term may not be
zero).

On the other hand, the quantum state
${\mathbf{r}\,\psi_{n}(\mathbf{r},t,\mathbf{k})}$ lies outside of the domain
of definition of periodic wavefunctions and leaves a residue in the DHFT
resulting in

$\left\langle\mathbf{v}_{b}(t,\mathbf{k})\right\rangle_{n}\neq 0.$ (51)

Eq. (47) together with Eq. (50) indicates that the non-Hermitian term
${\bm{\mathcal{S}}_{n}(t,\mathbf{k})}$, evaluated with respect to the Bloch
state ${\psi_{n}(\mathbf{r},t,\mathbf{k})}$, is never zero when periodic
boundary conditions are adopted for the wavefunctions. To the best of our
knowledge, it is the first time that the particle transport is explicitly
related with a non-Hermitian boundary velocity
${\left\langle\mathbf{v}_{b}(t,\mathbf{k})\right\rangle_{n}}$, together with
the electric field type of Berry curvature
${\bm{\mathcal{E}}_{n}(t,\mathbf{k})}$ that originates from the DHFT.

The electron’s non-Hermitian boundary velocity in $x$ direction is therefore
given by

$\left\langle\text{v}_{b}(t,\mathbf{k})\right\rangle_{n,\,x}=-\frac{1}{\hbar}\frac{\partial
E_{n}(t,\mathbf{k})}{\partial k_{x}}\,+\,\mathcal{E}_{n,\,x}(t,\mathbf{k})$
(52)

where the electric curvature is given from

$\displaystyle\mathcal{E}_{n,\,x}(t,\mathbf{k})$ $\displaystyle=$
$\displaystyle i\left\langle\frac{\partial\psi_{n}(t,\mathbf{k})}{\partial
k_{x}}\,|\,\frac{\partial\psi_{n}(t,\mathbf{k})}{\partial t}\right\rangle$
(53)
$\displaystyle-i\left\langle\frac{\partial\psi_{n}(t,\mathbf{k})}{\partial
t}|\,\frac{\partial\psi_{n}(t,\mathbf{k})}{\partial k_{x}}\right\rangle$

By defining the collective electrons’ displacement (for a fully occupied band)
through the sample boundaries in the $x$ direction, as
${\sum_{i}^{N}\Delta\left\langle
x_{i}\right\rangle_{n}\\!|_{bound}=\sum\limits_{\mathbf{k}}\int\limits_{0}^{T}\\!\left\langle\text{v}_{b}(t,\mathbf{k})\right\rangle_{n,\,x}dt}$,
then, for a cubic crystal with lattice constant $a$ and in the thermodynamic
limit, this displacement is given by

$\displaystyle\sum_{i}^{N}\Delta\left\langle x_{i}\right\rangle\\!|_{bound}$
$\displaystyle=$
$\displaystyle-\frac{V}{(2\pi)^{3}}\int\limits_{0}^{T}dt\iiint\limits_{BZ}\upsilon_{n,x}(t,\mathbf{k})d^{3}k$
$\displaystyle+\frac{V}{(2\pi)^{3}}\int\limits_{0}^{T}dt\iiint\limits_{BZ}\mathcal{E}_{n,\,x}(t,\mathbf{k})d^{3}k.$

where
${\displaystyle\upsilon_{n,x}(t,\mathbf{k})=\frac{1}{\hbar}\frac{\partial
E_{n}(t,\mathbf{k})}{\partial k_{x}}}$ is the $x$ component of the group
velocity.

The group velocity contribution is always truncated into a boundary
contribution over the Brillouin edges

$\int\limits_{0}^{T}dt\iiint\limits_{BZ}\bm{\upsilon}_{n}(t,\mathbf{k})d^{3}k=\frac{1}{\hbar}\int\limits_{0}^{T}dt\oiint\limits_{BZ}E_{n}(t,\mathbf{k})\,\mathbf{n}_{k}\,d^{2}k,$

where ${\mathbf{n}_{k}}$ is the unit vector that is locally normal to the
Brillouin zone boundaries, without any gauge phase ambiguities involved. On
the contrary, the electric field type of Berry curvature term of Eq. (IV) is
not generally amenable to further simplifications (unless some symmetries are
operating, which is thoroughly discussed next).

###  IV.1 Symmetry Considerations

If the time-dependent Hamiltonian is space-inversion symmetric, then,
${E_{n}(t,\mathbf{k})=E_{n}(t,-\mathbf{k})}$ and the group velocity
${\displaystyle\bm{\upsilon}_{n}(t,\mathbf{k})=\frac{1}{\hbar}\bm{\nabla_{\mathbf{k}}}E_{n}(t,\mathbf{k})}$
will satisfy
${\displaystyle\bm{\upsilon}_{n}(t,\mathbf{k})=-\frac{1}{\hbar}\bm{\nabla_{(-\mathbf{k})}}E_{n}(t,-\mathbf{k})=-\bm{\upsilon}_{n}(t,-\mathbf{k})}$,
resulting in
${\displaystyle\iiint_{BZ}\\!\\!\bm{\upsilon}_{n}(t,\mathbf{k})d^{3}k=0}$,
which means that the first term on the right side of Eq. (IV) integrates to
zero. Similarly, if the Hamiltonian is mirror symmetric along $x$, then,
${E_{n}(t,-k_{x},k_{y},k_{z})=E_{n}(t,k_{x},k_{y},k_{z})}$ resulting in
${\displaystyle\int\limits_{-\pi/a}^{+\pi/a}\frac{\partial
E_{n}(t,\mathbf{k})}{\partial k_{x}}\,dk_{x}=0}$, which also indicates that
the first term on the second side of Eq. (IV) integrates to zero. On the other
hand, due to the assumed scalar potential Eq. (40), space-inversion symmetry
and mirror symmetry along $x$ are broken, therefore, the first term on the
second side of Eq. (IV) cannot be zero due to space inversion or mirror
symmetry. If the scalar potential satisfies
${\phi(x,-t)^{\displaystyle*}=\phi(x,t)}$, then, the Hamiltonian is time-
reversal symmetric and the energy satisfies
${E_{n}(t,\mathbf{k})=E_{n}(-t,-\mathbf{k})}$, resulting in
${\displaystyle\bm{\upsilon}_{n}(t,\mathbf{k})=-\bm{\upsilon}_{n}(-t,-\mathbf{k})}$
for the group velocity. When the eigenvalue of the “energy” $E_{n}$ in the
instantaneous eigenvalue equation Eq. (42) does not have explicit time-
dependence ${\displaystyle\left\langle\psi_{n})\right|\frac{\partial
H(\mathbf{r},t)}{\partial
t}\left|\psi_{n}(t,\mathbf{k})\right\rangle=\frac{\partial
E_{n}(t,\mathbf{k})}{\partial t}=0}$, then, the system does not absorb energy
from the environment. In this framework the group velocity does not have
explicit time-dependence, thus, due to the assumed time-reversal symmetry it
satisfies
${\displaystyle\bm{\upsilon}_{n}(\mathbf{k})=-\bm{\upsilon}_{n}(-\mathbf{k})}$,
and as a result, the group velocity on the second side of Eq. (IV) integrates
to zero. On the other hand, when the energy has explicit time-dependence the
group velocity on the second side of Eq. (IV) does not necessarily integrate
to zero unless an extra, intermediate (middle) time period symmetry of the
Hamiltonian is present, namely,
${H(\mathbf{r},t)=H(\mathbf{r},t+T)=H(\mathbf{r},t+T/2)}$. The latter symmetry
indicates that the wavefunctions differs by a phase
${\psi_{n}(\mathbf{r},0)=e^{\displaystyle
i\Theta_{1k}}\psi_{n}(\mathbf{r},T)=e^{\displaystyle
i\Theta_{2k}}\psi_{n}(\mathbf{r},T/2)}$, which guarantees that the integration
${\int\limits_{0}^{T}dt\iiint\limits_{BZ}\bm{\upsilon}_{n}(t,\mathbf{k})d^{3}k=\int\limits_{-T/2}^{T/2}dt\iiint\limits_{BZ}\bm{\upsilon}_{n}(t+T/2,\mathbf{k})d^{3}k}$
turns into zero due to (i)
${\bm{\upsilon}_{n}(t+T/2,\mathbf{k})=\bm{\upsilon}_{n}(t,\mathbf{k})}$ and
(ii)
${\displaystyle\bm{\upsilon}_{n}(t,\mathbf{k})=-\bm{\upsilon}_{n}(-t,-\mathbf{k})}$.
Counter to a time-reversal symmetric Hamiltonian, when the assumed (real)
scalar potential breaks time-reversal symmetry (i.e. owing to ${E(-t)\neq
E(t)}$), then, the group velocity of the time-reversal breaking state in Eq.
(IV) does not necessarily integrate to zero (unless we are assuming compact
localized states with a flat bands).

On the other hand, when the translation invariant Hamiltonian has translation-
inversion symmetry (see Appendix C), namely, it remains invariant under the
inversions ${(i,\mathbf{k})\mapsto(-i,-\mathbf{k})}$, this symmetry dictates
that ${E_{n}(t,\mathbf{k})=E_{n}(t,-\mathbf{k})}$ and the group velocity
integrates to zero, hence, for the assumed scalar potential Eq. (40) the first
term on the right side of Eq. (IV) is zero.

Interestingly, if an additional confinement potential ${V_{c}(y,z)}$ that
localize the electrons near to the $x$ axis is assumed in the Hamiltonian Eq.
(39), then, the translation symmetry is broken along $y$ and $z$ coordinates
and the energy band of the 3D system becomes ${E_{n}\equiv E_{n}(t,k_{x})}$,
therefore, ${E_{n}(t,-k_{x})=E_{n}(t,k_{x})}$ due to time-reversal symmetry.
This is precisely the assumption that was made in Ref Thouless, 1983 where the
adiabatic particle transport of a 3D system was found quantized.

As ofr the Berry curvature ${\mathcal{E}_{n,\,x}(t,\mathbf{k})}$, as shown in
Appendix C, it does not satisfy the antisymmetry relation
${\mathcal{E}_{n,\,x}(t,k_{x},k_{y},k_{z})=-\mathcal{E}_{n,\,x}(t,-k_{x},k_{y},k_{z})}$
due to the mirror symmetry that is broken along $x$, therefore, the
integration
${\displaystyle\int\limits_{-\pi/a}^{+\pi/a}\mathcal{E}_{n,\,x}(t,\mathbf{k})\,dk_{x}}$
does not turn into zero. On the other hand, if time-reversal symmetry is
present
${\mathcal{E}_{n,\,x}(t,\mathbf{k})=-\mathcal{E}_{n,\,x}(-t,-\mathbf{k})}$
together with an assumed intermediate time period symmetry of the Hamiltonian
${H(\mathbf{r},t)=H(\mathbf{r},t+T)=H(\mathbf{r},t+T/2)}$, then,
${\mathcal{E}_{n,\,x}(t,\mathbf{k})=\mathcal{E}_{n,\,x}(t+T/2,\mathbf{k})}$
and the Berry curvature integrates to zero according to
$\int\limits_{0}^{T}dt\iiint\limits_{BZ}\mathcal{E}_{n,\,x}(t,\mathbf{k})d^{3}k=\int\limits_{-T/2}^{T/2}dt\iiint\limits_{BZ}\mathcal{E}_{n,\,x}(t+T/2,\mathbf{k})d^{3}k=\int\limits_{-T/2}^{T/2}dt\iiint\limits_{BZ}\mathcal{E}_{n,\,x}(t,\mathbf{k})d^{3}k=0$.
(For consequences of all this see below.)

###  IV.2 Topological arguments

The integration of the Berry curvature ${\mathcal{E}_{n,\,x}(t,\mathbf{k})}$
over the manifold spanned by ${dk_{x}\,dt}$ is quantized due to topology and
it is equal to

$\left(\,\int\limits_{0}^{T}\int\limits_{-\pi/a}^{+\pi/a}\mathcal{E}_{n,\,x}(t,\mathbf{k})\,dt\,dk_{x}\right)=2\pi\,C_{1}^{(n)},\
\ \ \ \ \ C_{1}^{(n)}\in\mathbb{Z}$ (55)

where ${C_{1}^{(n)}}$ is the first Chern number Nakahara (2003). The first
Chern number defines the mapping from the parameter space ${(t,\mathbf{k})}$
to the complex projective space of normalized Bloch states
${\Psi_{n}(\mathbf{r},t,\mathbf{k})}$. Non-zero Chern number indicates that
the mapping is topologically non-trivial. Specifically, due to the assumed
periodic boundary conditions at every instant,
${\psi_{n}(\mathbf{r}+\mathbf{L},t,\mathbf{k})=\psi_{n}(\mathbf{r},t,\mathbf{k})}$
as well as
${u_{n}(\mathbf{r}+\mathbf{L},t,\mathbf{k})=u_{n}(\mathbf{r},t,\mathbf{k})}$,
it implies that the Bloch states are generally periodic over the parameter
${\mathbf{k}}$ up to a phase
$\displaystyle\psi_{n}(\mathbf{r},t,\mathbf{k}+\mathbf{G})=e^{\displaystyle
i\Theta_{\\!{}_{G}}}\,\psi_{n}(\mathbf{r},t,\mathbf{k})\thinspace$ where
${\mathbf{G}}$ is a reciprocal lattice vector. Similarly, because we have
assumed time-periodic Hamiltonian, the Bloch states are periodic over ${t}$ up
to a phase ${\displaystyle\psi_{n}(\mathbf{r},t+T,\mathbf{k})=e^{\displaystyle
i\Theta_{{}_{T}}}\,\psi_{n}(\mathbf{r},t,\mathbf{k}).}$ The last two equations
imply that the Bloch states ${\Psi_{n}(\mathbf{r},t,\mathbf{k})}$ are not in
general single-valued quantities in the parameter space ${t\times\mathbf{k}}$.
Non-zero first Chern number ${C_{1}^{(n)}}$ indicates that the Bloch states
${\Psi_{n}(\mathbf{r},t,\mathbf{k})}$ cannot be everywhere single-valued in
the parameter space domain ${0\leq t\leq T}$ and
${\displaystyle-\frac{\pi}{a}\leq k_{x}\leq\frac{\pi}{a}}$.

###  IV.3 Quantization of the non-Hermitian displacement

Assuming, either a time-reversal symmetric Hamiltonian and a system that does
not absorb energy from the environment, or, bulk compact localized states with
flat bands, then, by employing Eq. (55) in Eq. (IV), it turns out that the
collective electrons’ displacement through the system’s boundaries in the $x$
direction is given by

$\displaystyle\sum_{i}^{N}\Delta\left\langle
x_{i}\right\rangle\\!|_{Boundary}$ $\displaystyle=$
$\displaystyle\frac{V}{(2\pi)^{2}}\,C_{1}^{(n)}\int\limits_{-\pi/a}^{+\pi/a}\\!\\!\\!dk_{y}\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\\!\\!\\!dk_{z}$
(56) $\displaystyle=$ $\displaystyle\,C_{1}^{(n)}N\alpha,$

where $N$ is the total number of unit cells of the material. Eq. (56) implies
that, after one period $T$ of adiabatic driving, the displacement of the
center of mass of the electrons in the $x$ direction and over the boundaries,
is quantized in units of lattice constant $a$ and is given by

$\displaystyle\frac{1}{N}\sum_{i}^{N}\Delta\left\langle
x_{i}\right\rangle\\!|_{bound}=\,C_{1}^{(n)}\,\alpha,$ (57)

meaning that the center of mass performs a quantized rigid displacement over
the boundaries along the $x$ direction, or equivalently, the electrons behaves
as an incompressible fluid in $\mathbf{k}$-space. Accordingly, the collective
electrons’ displacement over $y$ and $z$ directions are zero

$\sum_{i}^{N}\Delta\left\langle
y_{i}\right\rangle\\!|_{bound}=\sum_{i}^{N}\Delta\left\langle
z_{i}\right\rangle\\!|_{bound}=0,$ (58)

due to the mirror inversion symmetries (assumed to be present) that implies
${\mathcal{E}_{n,\,y}(t,k_{x},k_{y},k_{z})=-\mathcal{E}_{n,\,y}(t,k_{x},-k_{y},k_{z})}$
as well as
${\mathcal{E}_{n,\,z}(t,k_{x},k_{y},k_{z})=-\mathcal{E}_{n,\,z}(t,k_{x},k_{y},-k_{z})}$
which in turn gives zero displacements.

We conclude that, if the non-Hermitian term
${\bm{\mathcal{S}}_{n}(t,\mathbf{k})}$ is not taken into account in the
adiabatic form of DHFT Eq. (44), then, this will eventually lead to an
apparent paradox since the curvature flux must always give zero value,
${0=\int\limits_{0}^{T}\int\limits_{-\pi/a}^{+\pi/a}\mathcal{E}_{n,\,x}(t,\mathbf{k})\,dt\,dk_{x}}$
due to the above assumed symmetries. On the other hand, by keeping track of
the emerging non-Hermitian contribution
${\bm{\mathcal{S}}_{n}(t,\mathbf{k})}$, the apparent paradox is resolved, and,
owing to ${\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k})=0}$ (due to the
cell-periodic wavefunctions that satisfy periodic boundary conditions over the
material boundaries at every instant), the flux of the non-Hermitian term
${\bm{\mathcal{S}}_{n}(t,\mathbf{k})=\hbar\left\langle\,\mathbf{v}_{b}\,\right\rangle_{n}}$
(due to Eq. (47)) turns out a topological invariant.

##  V Non-adiabatic Charge Pump

We now employ the Floquet states and counterpart Floquet-Bloch bands
$\varepsilon_{a}(\mathbf{k})$ in order to theoretically study the non-
equilibrium and driven charge pumping with the DHFT. Each spinless electron’s
wavefunction evolves in time by the Hamiltonian, ${\displaystyle
H(\mathbf{r},t)=\frac{1}{2m}{\mathbf{p}^{2}+e\phi(x,t)}+V_{crys}(\mathbf{r}),}$
where $m$ and $e$ are the mass and charge of the electron respectively,
$\mathbf{p}$ is the canonical momentum operator, and $V_{crys}(\mathbf{r})$ is
the crystal potential. The applied scalar potential is given by

$\phi(x,t)=\frac{a}{2\pi}E(t)\,cos(2\pi x/a+\theta),$ (59)

where $a$ is the primitive cell length in $x$ direction, ${E(t)}$ is the
amplitude of the microscopic electric field, and $\theta$ is a dimensionless
constant. The scalar potential satisfies ${\phi(x+na,t)=\phi(x,t)}$ where $n$
is an integer, as well as ${\phi(x,T)=\phi(x,0)=0}$. Whenever $\theta$ is zero
(or more generally $\theta=2\pi n$), the Hamiltonian has mirror symmetry along
$\bm{e}_{x}$ direction, otherwise the mirror symmetry is broken. The
macroscopic (spatially averaged) electric field is zero, thus, no voltage-bias
is present during the pumping. The time-dependent amplitude of the electric
field $E(t)$ is periodic ${E(t+T)=E(t)}$, but, has arbitrary large velocity
$\displaystyle\frac{dE(t)}{dt}$, therefore, no adiabatic time-evolution can be
implied. Instead, we use the non-equilibrium Floquet states and apply the DHFT
of Eq. (34).

Each electron’s wavefunction evolves in time according to the TDSE

$i\hbar\frac{d}{dt}\Psi(\mathbf{r},t,\mathbf{k})=H(\mathbf{r},t)\,\Psi(\mathbf{r},t,\mathbf{k}),$
(60)

and the solution of Eq. (60) is considered to be a non-equilibrium Floquet
state
${\displaystyle\Psi_{a}(\mathbf{r},t,\mathbf{k})=e^{-\displaystyle\frac{i}{\hbar}\varepsilon_{a}(\mathbf{k})\,t}\Phi_{a}(\mathbf{r},t,\mathbf{k})}$.
The Floquet modes ${\Phi_{a}(\mathbf{r},t,\mathbf{k})}$ are periodic in time
${\Phi_{a}(\mathbf{r},t,\mathbf{k})=\Phi_{a}(\mathbf{r},t+T,\mathbf{k})}$ and
satisfy the eigenvalue equation

$\displaystyle\left(H(\mathbf{r},t)-i\hbar\frac{d}{dt}\right)\Phi_{a}(\mathbf{r},t,\mathbf{k})$
$\displaystyle=$ $\displaystyle
H_{F}(\mathbf{r},t)\Phi_{a}(\mathbf{r},t,\mathbf{k})$ (61) $\displaystyle=$
$\displaystyle\varepsilon_{a}(\mathbf{k})\Phi_{a}(\mathbf{r},t,\mathbf{k}),$

where ${H_{F}(\mathbf{r},t)}$ is the so-called Floquet Hamiltonian, and
${\varepsilon_{a}(\mathbf{k})}$ is the (static) quasienergy that is restricted
into the interval
${\displaystyle-\pi\hbar/T\leq\varepsilon_{a}(\mathbf{k})\leq\pi\hbar/T}$
called first Floquet-zone which contains all the physically non-equivalent
quantum states. Furthermore, we assume that each Floquet mode is expressed
with respect to the Floquet-Bloch state Gómez-León and Platero (2013)

$\Phi_{a}(\mathbf{r},t,\mathbf{k})=e^{\displaystyle
i\mathbf{k}\\!\cdot\\!\mathbf{r}}u_{a}(\mathbf{r},t,\mathbf{k}),$ (62)

where ${u_{a}(\mathbf{r},t,\mathbf{k})}$ is periodic both in ${\mathbf{r}}$
and $t$. Whenever the Floquet Hamiltonian ${H_{F}(\mathbf{r},t)}$ is time-
reversal invariant the quasienergy spectrum is symmetric with respect to
crystal momentum inversion
${\varepsilon_{a}(\mathbf{k})=\varepsilon_{a}(-\mathbf{k})}$. Moreover, when
the Floquet Hamiltonian satisfies mirror symmetry, i.e. along $x$, then, for
the specific direction it satisfies
${\varepsilon_{a}(k_{x},k_{y},k_{z})=\varepsilon_{a}(-k_{x},k_{y},k_{z})}$.

The parameter definition is made by using the assumed form of the quantum
state. Due to the presence of translational invariance of the Hamiltonian, the
wave vector ${\mathbf{k}}$ is a conserved quantity of the dynamics. We
therefore define the parameter ${\mathbf{R}}$ as the static wave vector
${\mathbf{R}\equiv\mathbf{k}}$ which has zero corresponding velocity
${\displaystyle\frac{\partial\mathbf{R}}{\partial
t}\equiv\displaystyle\frac{\partial\mathbf{k}}{\partial t}=0}$.

The Hamiltonian ${H(\mathbf{r},t)}$ does not depend on any parameter
${\bm{\nabla_{\mathbf{k}}}H(\mathbf{r},t)=0}$, therefore, application of the
DHFT Eq. (34) with respect to the Floquet state,
${\displaystyle\Psi_{a}(\mathbf{r},t,\mathbf{k})=e^{-\displaystyle\frac{i}{\hbar}\varepsilon_{a}(\mathbf{k})\,t}\Phi_{a}(\mathbf{r},t,\mathbf{k})}$,
gives

$0=\bm{\nabla_{\mathbf{k}}}E_{a}(t,\mathbf{k})\,+\,\bm{\mathcal{S}}_{a}(t,\mathbf{k})\,-\,\hbar\,\bm{\mathcal{E}}_{a}(t,\mathbf{k}),$
(63)

where the “energy” is given by

$\displaystyle E_{a}(t,\mathbf{k})$ $\displaystyle=$
$\displaystyle\left\langle\Phi_{a}(t,\mathbf{k})\right|H(t)\left|\Phi_{a}(t,\mathbf{k})\right\rangle$
(64) $\displaystyle=$
$\displaystyle\varepsilon_{a}(\mathbf{k})+i\hbar\left\langle\Phi_{a}(t,\mathbf{k})\,|\,\frac{d\Phi_{a}(t,\mathbf{k})}{dt}\right\rangle,$

the non-Hermitian term from

$\bm{\mathcal{S}}_{a}(t,\mathbf{k})=\left\langle\Psi_{a}(t,\mathbf{k})|\left(H(\mathbf{r},t)^{+}-H(\mathbf{r},t)\right)\bm{\nabla_{\mathbf{k}}}\Psi_{a}(t,\mathbf{k})\right\rangle,$
(65)

and the electric field type of curvature by

$\displaystyle\bm{\mathcal{E}}_{a}(t,\mathbf{k})$ $\displaystyle=$
$\displaystyle
i\left\langle\bm{\nabla_{\mathbf{k}}}\Phi_{a}(t,\mathbf{k})\,|\,\frac{\partial\Phi_{a}(t,\mathbf{k})}{\partial
t}\right\rangle$
$\displaystyle-i\left\langle\frac{\partial\Phi_{a}(t,\mathbf{k})}{\partial
t}|\,\bm{\nabla_{\mathbf{k}}}\Phi_{a}(t,\mathbf{R})\right\rangle.$

By then using the identity

$\displaystyle\bm{\nabla}_{\mathbf{k}}\Psi_{a}(\mathbf{r},t,\mathbf{k})$
$\displaystyle=$
$\displaystyle-\frac{i}{\hbar}\bm{\nabla_{\mathbf{k}}}\varepsilon_{a}(\mathbf{k})\
t\ \Psi_{a}(\mathbf{r},t,\mathbf{k})\ +\
i\,\mathbf{r}\,\Psi_{a}(\mathbf{r},t,\mathbf{k})$
$\displaystyle+\,e^{-\displaystyle\frac{i}{\hbar}\varepsilon_{a}(\mathbf{k})\,t}e^{\displaystyle
i\mathbf{k}\\!\cdot\\!\mathbf{r}}\,\bm{\nabla}_{\mathbf{k}}u_{a}(\mathbf{r},t,\mathbf{k}),$

we find the equation that relates the non-Hermitian terms

$\bm{\mathcal{S}}_{a}(t,\mathbf{k})=\hbar\left\langle\,\mathbf{v}_{b}(t,\mathbf{k})\,\right\rangle_{a}\,+\,\widetilde{\bm{\mathcal{S}}_{a}}(t,\mathbf{k}),$
(66)

where ${\left\langle\,\mathbf{v}_{b}(t,\mathbf{k})\,\right\rangle_{a}}$ is the
boundary velocity evaluated with respect to the Floquet mode
${\Phi_{a}(\mathbf{r},t,\mathbf{k})}$ and
${\widetilde{\bm{\mathcal{S}}_{a}}(t,\mathbf{k})}$ is evaluated with respect
to the Floquet-Bloch state

$\widetilde{\bm{\mathcal{S}}_{a}}(t,\mathbf{k})=\left\langle
u_{a}(t,\mathbf{k})|\left(H_{k}(\mathbf{r},t)^{+}-H_{k}(\mathbf{r},t)\right)\bm{\nabla_{\mathbf{k}}}u_{a}(t,\mathbf{k})\right\rangle.$
(67)

By substituting Eq. (66) into Eq. (63) we find the boundary velocity of the
electron

$\left\langle\mathbf{v}_{b}(t,\mathbf{k})\right\rangle_{a}=-\frac{1}{\hbar}\bm{\nabla_{\mathbf{k}}}E_{a}(t,\mathbf{k})\,-\,\frac{1}{\hbar}\widetilde{\bm{\mathcal{S}}_{a}}(t,\mathbf{k})\,+\,\bm{\mathcal{E}}_{a}(t,\mathbf{k}),$
(68)

and, by substituting Eq. (64) into Eq. (68), we finally find the non-Hermitian
boundary velocity in $x$ direction for a single Floquet mode that is given by

$\displaystyle\left\langle\text{v}_{b}(t,\mathbf{k})\right\rangle_{a,\,x}$
$\displaystyle=$
$\displaystyle-\frac{1}{\hbar}\frac{\partial\varepsilon_{a}(\mathbf{k})}{\partial
k_{x}}\ -\ i\frac{\partial}{\partial
k_{x}}\left\langle\\!\Phi_{a}(t,\mathbf{k})\,|\,\frac{d\Phi_{a}(t,\mathbf{k})}{dt}\\!\right\rangle$
(69)
$\displaystyle\,-\,\frac{1}{\hbar}\widetilde{\mathcal{S}_{a,\,x}}(t,\mathbf{k})\,+\,\mathcal{E}_{a,\,x}(t,\mathbf{k}).$

In this respect, each electron’s displacement over the material’s boundaries
in the $x$ direction after the time of period $T$ is given by

$\displaystyle\Delta\left\langle x(\mathbf{k})\right\rangle_{a}\\!|_{bound}$
$\displaystyle=$
$\displaystyle-\frac{1}{\hbar}\frac{\partial\varepsilon_{a}(\mathbf{k})}{\partial
k_{x}}T\ -\ \frac{\partial\beta_{a}(\mathbf{k})}{\partial k_{x}}$
$\displaystyle-\,\frac{1}{\hbar}\int_{0}^{T}\\!\\!\\!\widetilde{\mathcal{S}_{a,\,x}}(t,\mathbf{k})dt\
+\ \int_{0}^{T}\\!\\!\\!\mathcal{E}_{a,\,x}(t,\mathbf{k})dt$

where, interestingly, the quantity ${\beta_{a}(\mathbf{k})}$ showing up is the
non-adiabatic Aharonov-Anandan phase of the Floquet mode

$\beta_{a}(\mathbf{k})=\int_{0}^{T}i\left\langle\Phi_{a}(t,\mathbf{k})\,|\,\frac{d\Phi_{a}(t,\mathbf{k})}{dt}\right\rangle
dt.$ (71)

For a spinless electron and time-reversal symmetric Floquet Hamiltonian, the
Floquet modes satisfy
${\Phi_{a}(\mathbf{r},t,\mathbf{k})=\Phi_{a}(\mathbf{r},-t,-\mathbf{k})^{\displaystyle*}}$
(up to a constant ${U(1)}$ phase ${e^{\displaystyle i\lambda}}$ which we
neglected for simplicity). In this respect, it is easy to show that the
Floquet modes satisfy ${\displaystyle
i\left\langle\Phi_{a}(t,\mathbf{k})\,|\,\frac{d\Phi_{a}(t,\mathbf{k})}{dt}\right\rangle=i\left\langle\Phi_{a}(-t,-\mathbf{k})\,|\,\frac{d\Phi_{a}(-t,-\mathbf{k})}{d(-t)}\right\rangle,}$
which results to an Aharonov-Anandan phase that is symmetric under the
inversion of crystal momentum
${\beta_{a}(\mathbf{k})=\beta_{a}(-\mathbf{k})}$. Analogously, for a space-
inversion symmetric Floquet Hamiltonian
(${\Phi_{a}(\mathbf{r},t,\mathbf{k})=\Phi_{a}(-\mathbf{r},t,-\mathbf{k})}$ up
to a constant ${U(1)}$ phase), the Aharonov-Anandan phase satisfies
${\beta_{a}(\mathbf{k})=\beta_{a}(-\mathbf{k})}$, while, for a mirror-
symmetric Floquet Hamiltonian, i.e. along $x$, it satisfies
${\beta_{a}(k_{x},k_{y},k_{z})=\beta_{a}(-k_{x},k_{y},k_{z})}$.

Interestingly, whenever the Floquet Hamiltonian satisfies time-reversal
symmetry but the mirror symmetry is broken along $x$, then,
$\varepsilon_{a}(k_{x},k_{y},k_{z})-\varepsilon_{a}(-k_{x},-k_{y},-k_{z})=0$
(due to time-reversal symmetry) while
${\varepsilon_{a}(k_{x},k_{y},k_{z})-\varepsilon_{a}(-k_{x},k_{y},k_{z})\neq
0}$ (due to breaking of the mirror symmetry), for every crystal momentum
within the first Floquet-Brillouin zone. This implies that there is at least
one discontinuity of the Floquet quasienergy band (which extends in the
interval
${\displaystyle-\pi\hbar/T\leq\varepsilon_{a}(\mathbf{k})\leq\pi\hbar/T}$)
while traversing the $k_{x}$ momenta of the first Brillouin zone, which
results to an apparent winding (along the quasienergy axe) of the Floquet
band. Provided that there is no gap in the first Floquet-Brillouin zone, each
apparent winding (discontinuity) contributes a quantized amount of
${\displaystyle|\Delta\varepsilon_{a}(\mathbf{k})|=2\pi\hbar/T}$ resulting
into
${\displaystyle\int_{BZ}\frac{\partial\varepsilon_{a}(\mathbf{k})}{\partial
k_{x}}dk_{x}=n\,2\pi\hbar/T}$ where $n$ is the integer number of windings.

Similarly, when the Floquet Hamiltonian satisfies time-reversal symmetry but
the mirror symmetry is broken along $x$, then, the Aharonov-Anandan phase
satisfies $\beta_{a}(k_{x},k_{y},k_{z})-\beta_{a}(-k_{x},-k_{y},-k_{z})=0$
(due to time-reversal symmetry), but,
${\beta_{a}(k_{x},k_{y},k_{z})-\beta_{a}(-k_{x},k_{y},k_{z})\neq 0}$ (due to
breaking of the mirror symmetry) for every crystal momentum within the first
Floquet-Brillouin zone. This implies that there is at least one discontinuity
of the Aharonov-Anandan phase ${\beta_{a}(\mathbf{k})}$ with unspecified jump
while traversing the $k_{x}$ momenta of the first Brillouin zone.
Specifically, although the Aharonov-Anandan phase values (modulo $2\pi$) lie
in the interval ${0\leq\beta_{a}(\mathbf{k})\leq 2\pi}$, the jump is not
necessarily equal to $2\pi$ because the phase is not coming out as a product
of an eigenvalue equation, in contrast to the Floquet energy
${\varepsilon_{a}(\mathbf{k})}$ 111Each Floquet energy gives the total phase
of the quantum state after the time period of $T$, that is
$\displaystyle\frac{1}{\hbar}\varepsilon_{a}(\mathbf{k})T$, therefore, all
different quasienergy eigenvalues define the projective Hilbert space of rays
of the single valued Floquet modes which is the eigenvalue of the Floquet
Hamiltonian. In this respect, for mirror-breaking (along $x$) Floquet
Hamiltonians, we cannot naively assume that
${\displaystyle\int_{BZ}\frac{\partial\beta_{a}(\mathbf{k})}{\partial
k_{x}}dk_{x}=0}$, and the Aharonov-Anandan phase must be taken into account
when calculating the collective non-Hermitian displacement of the electrons.
We conclude that, for mirror-breaking Floquet Hamiltonians, the Aharonov-
Anandan phase behaves as a non-integrable quantity with respect to crystal
momentum, and as a consequence it contributes to charge pumping an amount
proportional to its non quantized winding.

###  V.1 Charge pumping from the lowest-energy Floquet state

We assume that each independent electron is initially in the lowest-energy
non-equilibrium Floquet state
${\displaystyle\Psi_{a}(\mathbf{r},t,\mathbf{k})=e^{-\displaystyle\frac{i}{\hbar}\varepsilon_{a}(\mathbf{k})\,t}\Phi_{a}(\mathbf{r},t,\mathbf{k})}$,
therefore, the period-averaged “energy” $\displaystyle
E_{a}(\mathbf{k})=\frac{1}{T}\int_{0}^{T}\\!\\!\left\langle\Phi_{a}(t,\mathbf{k})\right|H(t)\left|\Phi_{a}(t,\mathbf{k})\right\rangle
dt$ is assumed to have the lowest value. By employing Eq. (69) and for a fully
occupied quasienergy band, the non-Hermiticity-related electrons’ displacement
in the $x$ direction and over the boundaries of the material (in the
thermodynamic limit) is given by

$\displaystyle\frac{(2\pi)^{3}}{V}\sum_{i}^{N}\Delta\left\langle
x_{i}\right\rangle\\!|_{bound}\,$ $\displaystyle=$
$\displaystyle-\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{y}\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{z}\left(\frac{T}{\hbar}\int\limits_{-\pi/a}^{+\pi/a}\frac{\partial\varepsilon_{a}(\mathbf{k})}{\partial
k_{x}}\\!dk_{x}\right)\ -\
\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{y}\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{z}\left(\int\limits_{-\pi/a}^{+\pi/a}\frac{\partial\beta_{a}(\mathbf{k})}{\partial
k_{x}}\\!dk_{x}\right)$
$\displaystyle\\!\\!\\!-\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{y}\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{z}\left(\,\int\limits_{0}^{T}\int\limits_{-\pi/a}^{+\pi/a}\widetilde{\mathcal{S}_{a,\,x}}(t,\mathbf{k})\,dt\,dk_{x}\right)\
+\
\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{y}\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{z}\left(\,\int\limits_{0}^{T}\int\limits_{-\pi/a}^{+\pi/a}\mathcal{E}_{a,\,x}(t,\mathbf{k})\,dt\,dk_{x}\right)\\!\\!.$

We are assuming periodic boundary conditions,
${\psi_{a}(\mathbf{r}+\mathbf{L},t,\mathbf{k})=\psi_{n}(\mathbf{r},t,\mathbf{k})}$
as well as
${u_{a}(\mathbf{r}+\mathbf{L},t,\mathbf{k})=u_{n}(\mathbf{r},t,\mathbf{k})}$,
thus, the non-Hermitian contribution is in this case zero

$\widetilde{\mathcal{S}_{a,x}}(t,\mathbf{k})=0.$ (73)

####  V.1.1 Mirror-symmetric Floquet Hamiltonian

For a mirror-symmetric Floquet Hamiltonian along $x$ we assume ${\theta=0}$ in
Eq. (59). By employing the eigenvalue equation Eq. (61), the quasienergy and
the Aharonov-Anandan phase satisfies,
${\varepsilon_{a}(k_{x},k_{y},k_{z})=\varepsilon_{a}(-k_{x},k_{y},k_{z})}$ and
${\beta_{a}(k_{x},k_{y},k_{z})=\beta_{a}(-k_{x},k_{y},k_{z})}$ respectively;
while the Berry curvature satisfies
${\mathcal{E}_{a,\,x}(t,k_{x},k_{y},k_{z})=-\mathcal{E}_{a,\,x}(t,-k_{x},k_{y},k_{z})}$
as shown in Appendix C. In this respect, the collective non-Hermitian pumped
charge per cycle is trivially zero,

$\displaystyle\frac{1}{N}\sum_{i}^{N}\Delta\left\langle
x_{i}\right\rangle\\!|_{bound}\,=0.$

####  V.1.2 Mirror-breaking Floquet Hamiltonian

For a mirror-symmetric along $x$ Floquet Hamiltonian we assume ${\theta\neq
2\pi n}$ in Eq. (59). Due to mirror symmetry breaking,
$\varepsilon_{a}(k_{x},k_{y},k_{z})\neq\varepsilon_{a}(-k_{x},k_{y},k_{z})$,
as well as $\beta_{a}(k_{x},k_{y},k_{z})\neq\beta_{a}(-k_{x},k_{y},k_{z})$ and
$\mathcal{E}_{a,\,x}(t,k_{x},k_{y},k_{z})\neq-\mathcal{E}_{a,\,x}(t,-k_{x},k_{y},k_{z})$.
The collective, non-adiabatic and formally originating from non-Hermiticity
pumped charge of the center of mass of the electrons is therefore given by

$\displaystyle\frac{1}{N}\sum_{i}^{N}\Delta\left\langle
x_{i}\right\rangle\\!|_{bound}\,=-\frac{a^{3}\,T}{(2\pi)^{3}\,\hbar}\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{y}\\!\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{z}\\!\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\frac{\partial\varepsilon_{a}(\mathbf{k})}{\partial
k_{x}}\\!dk_{x}\ -\
\frac{a^{3}}{(2\pi)^{3}}\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{y}\\!\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{z}\\!\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\frac{\partial\beta_{a}(\mathbf{k})}{\partial
k_{x}}\\!dk_{x}\ +\ C_{1}^{(n)}\,\alpha,$

where we have used that
${\left(\,\int\limits_{0}^{T}\int\limits_{-\pi/a}^{+\pi/a}\mathcal{E}_{a,\,x}(t,\mathbf{k})\,dt\,dk_{x}\right)=2\pi\,C_{1}^{(a)}}$
with ${C_{1}^{(a)}}$ the first Chern number. Non-zero first Chern number
${C_{1}^{(a)}}$ indicates that the Floquet states
${\Phi_{a}(\mathbf{r},t,\mathbf{k})}$ cannot be everywhere single-valued in
the parameter space domain ${0\leq t\leq T}$ and
${\displaystyle-\frac{\pi}{a}\leq k_{x}\leq\frac{\pi}{a}}$. Namely, although
the Floquet states are periodic with respect to time by definition,
${\Phi_{a}(\mathbf{r},t,\mathbf{k})=\Phi_{a}(\mathbf{r},t+T,\mathbf{k})}$,
non-zero Chern number indicates that they don’t satisfy the periodic gauge
along $x$, that is,
${\Phi_{a}(\mathbf{r},t,\mathbf{k})\neq\Phi_{a}(\mathbf{r},t,\mathbf{k}+G_{x}\,\mathbf{e}_{x})}$
where $G_{x}$ is any reciprocal lattice vector in $\mathbf{e}_{x}$ direction.

If we assume that there is no gap present in the first Floquet-Brillouin zone,
then, the collective, non-adiabatic and non-Hermiticity-related pumped charge
of the center of mass of the electrons is given by

$\displaystyle\frac{1}{N}\sum_{i}^{N}\Delta\left\langle
x_{i}\right\rangle\\!|_{bound}\,$ $\displaystyle=$ $\displaystyle-n\,a\ +\
C_{1}^{(n)}\alpha$
$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!-\frac{a^{3}}{(2\pi)^{3}}\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{y}\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{z}\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\frac{\partial\beta_{a}(\mathbf{k})}{\partial
k_{x}}\\!dk_{x},$

where the first term gives a quantized contribution in terms of the integer
number $n$ of windings of the Floquet-Brillouin zone and the last term gives a
generically non-quantized contribution in terms of the non-integrable
Aharonov-Anandan phase. Interestingly, a naive application of the periodic
gauge on the Floquet modes along $\mathbf{e}_{x}$ direction, that is,
${\Phi_{a}(\mathbf{r},t,\mathbf{k})=\Phi_{a}(\mathbf{r},t,\mathbf{k}+G_{x}\mathbf{e}_{x})}$,
results to the cancellation of the Aharonov-Anandan contribution to the pumped
charge

$\displaystyle\int\limits_{-\pi/a}^{+\pi/a}\frac{\partial\beta_{a}(\mathbf{k})}{\partial
k_{x}}\\!dk_{x}$ $\displaystyle=$
$\displaystyle\beta_{a}(\mathbf{k}+G_{x}\mathbf{e}_{x})-\beta_{a}(\mathbf{k})$
$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!=i\int_{0}^{T}\left\langle\Phi_{a}(t,\mathbf{k}+G_{x}\mathbf{e}_{x})\,|\,\frac{d\Phi_{a}(t,\mathbf{k}+G_{x}\mathbf{e}_{x})}{dt}\right\rangle
dt$
$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!-i\int_{0}^{T}\left\langle\Phi_{a}(t,\mathbf{k})\,|\,\frac{d\Phi_{a}(t,\mathbf{k})}{dt}\right\rangle
dt$
$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!=0.$

In a slightly different context, presumably a more intuitive one, by employing
Eq. (68), the pumped charge can be expressed in terms of the dynamical phase
as

$\displaystyle\frac{1}{N}\sum_{i}^{N}\Delta\left\langle
x_{i}\right\rangle\\!|_{bound}\,$ $\displaystyle=$ $\displaystyle
C_{1}^{(n)}\,\alpha$
$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!+\frac{a^{3}}{(2\pi)^{3}}\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{y}\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{z}\\!\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\frac{\partial
D_{a}(\mathbf{k})}{\partial k_{x}}\\!dk_{x},$

where the dynamical phase is given by $\displaystyle
D_{a}(\mathbf{k})=-\frac{1}{\hbar}\int_{0}^{T}\\!\\!\left\langle\Phi_{a}(t,\mathbf{k})\right|H(t)\left|\Phi_{a}(t,\mathbf{k})\right\rangle
dt$. If the Hamiltonian is mirror-inversion symmetric along $x$, then,
${C_{1}^{(n)}=0}$, and, due to the dynamical phase that satisfies
$\displaystyle D_{a}(k_{x},k_{y},k_{z})=D_{a}(-k_{x},k_{y},k_{z})$, the pumped
charge turns into zero. On the other hand, for the present mirror-reversal
breaking case $\displaystyle D_{a}(k_{x},k_{y},k_{z})\neq
D_{a}(-k_{x},k_{y},k_{z})$, we don’t a priori have any topological argument
that the difference of the dynamical phase over the edges of the first
Brillouin zone should be quantized. In conclusion, for a mirror-breaking
Hamiltonian, the quantization of the pumped charge breaks down due to a non-
integrable Aharonov-Anandan phase even if no gaps exist in the Floquet-
Brillouin zone, provided that the insulator is a topologically nontrivial one,
where the periodic gauge cannot be taken (the well-known obstruction in
satisfying it).

###  V.2 Trivialization of the non-adiabatic displacement

The boundary velocity Eq. (69) can be simplified by assuming that,
${\displaystyle\frac{\partial}{\partial k_{x}}\frac{\partial}{\partial
t}\left|\Phi_{a}(t,\mathbf{k})\right\rangle-\frac{\partial}{\partial
t}\frac{\partial}{\partial
k_{x}}\left|\Phi_{a}(t,\mathbf{k})\right\rangle=0}$, which is the
integrability condition that guarantees that the Floquet modes
${\left|\Phi_{a}(t,\mathbf{k})\right\rangle}$ are analytic (single valued)
functions of $t$ and $k_{x}$ in the ${t\\!\times\\!k_{x}}$ space. Therefore,
either the material is a conventional band insulator if we take into account
the entire ${t\\!\times\\!k_{x}}$ space, or we are assuming a fiber bundle
theory with patches. In this framework, the Berry curvature takes the form

$\displaystyle\mathcal{E}_{a,x}(t,\mathbf{k})$ $\displaystyle=$ $\displaystyle
i\frac{\partial}{\partial
k_{x}}\left\langle\Phi_{a}(t,\mathbf{k})\,|\,\frac{\partial\Phi_{a}(t,\mathbf{k})}{\partial
t}\right\rangle$ $\displaystyle-i\frac{\partial}{\partial
t}\left\langle\Phi_{a}(t,\mathbf{k})|\,\frac{\partial\Phi_{a}(t,\mathbf{k})}{\partial
k_{x}}\right\rangle$

and the boundary velocity Eq. (69) truncates into

$\displaystyle\left\langle\text{v}_{b}(t,\mathbf{k})\right\rangle_{a,\,x}$
$\displaystyle=$
$\displaystyle-\frac{1}{\hbar}\frac{\partial\varepsilon_{a}(\mathbf{k})}{\partial
k_{x}}\ -\ i\frac{\partial}{\partial
t}\left\langle\\!\Phi_{a}(t,\mathbf{k})\,|\,\frac{\partial\Phi_{a}(t,\mathbf{k})}{\partial
k_{x}}\\!\right\rangle$ (76)
$\displaystyle\,-\,\frac{1}{\hbar}\widetilde{\mathcal{S}_{a,\,x}}(t,\mathbf{k}).$

Therefore, by assuming that ${\widetilde{\mathcal{S}_{a,x}}(t,\mathbf{k})=0}$,
the center of mass, non-adiabatic displacement over the boundaries Eq. (V.1)
takes the form

$\displaystyle\frac{1}{N}\sum_{i}^{N}\Delta\left\langle
x_{i}\right\rangle\\!|_{bound}\,=-\frac{a^{3}\,T}{(2\pi)^{3}\,\hbar}\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{y}\\!\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{z}\\!\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\frac{\partial\varepsilon_{a}(\mathbf{k})}{\partial
k_{x}}\\!dk_{x}\ -\
\frac{a^{3}}{(2\pi)^{3}}\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{y}\\!\\!\\!\int\limits_{-\pi/a}^{+\pi/a}\\!dk_{z}\\!\int\limits_{0}^{T}\frac{\gamma_{a}(t,\mathbf{k})}{dt}\\!dt$

where ${\gamma_{a}(t,\mathbf{k})}$ is a first Brillouin zone accumulated phase
given by

$\gamma_{a}(t,\mathbf{k})=i\\!\int\limits_{-\pi/a}^{+\pi/a}\\!\\!\\!A_{x}(t,\mathbf{k})\,dk_{x},$
(77)

where ${\displaystyle
A_{x}(t,\mathbf{k})=i\left\langle\\!\Phi_{a}(t,\mathbf{k})\,|\,\frac{\partial\Phi_{a}(t,\mathbf{k})}{\partial
k_{x}}\\!\right\rangle}$ is the Berry connection of the Floquet modes. The
Floquet modes ${\left|\Phi_{a}(t,\mathbf{k})\right\rangle}$ are periodic with
respect to the time period $T$,
${\left|\Phi_{a}(t+T,\mathbf{k})\right\rangle=\left|\Phi_{a}(t,\mathbf{k})\right\rangle}$.
If we further assume that the momentum derivative of the Floquet modes
${\displaystyle\left|\frac{\partial\Phi_{a}(t,\mathbf{k})}{\partial
k_{x}}\\!\right\rangle}$ is also periodic with respect to the time period $T$,
that is, ${\displaystyle\left|\frac{\partial\Phi_{a}(t+T,\mathbf{k})}{\partial
k_{x}}\\!\right\rangle=\left|\frac{\partial\Phi_{a}(t,\mathbf{k})}{\partial
k_{x}}\\!\right\rangle}$, which is true whenever the Floquet modes are
analytic functions of ${k_{x}}$ in the ${t\\!\times\\!k_{x}}$ space, then, it
is evident that
$\displaystyle\left\langle\\!\Phi_{a}(T,\mathbf{k})\,|\,\frac{\partial\Phi_{a}(T,\mathbf{k})}{\partial
k_{x}}\\!\right\rangle=\left\langle\\!\Phi_{a}(0,\mathbf{k})\,|\,\frac{\partial\Phi_{a}(0,\mathbf{k})}{\partial
k_{x}}\\!\right\rangle$ which results into
${\displaystyle\int\limits_{0}^{T}\frac{\gamma_{a}(t,\mathbf{k})}{dt}\\!dt=0}$,
thus, the pumped charge is then given only in terms of the winding of the
Floquet band. This is the formal case that is encountered in the
literatureShih and Niu (1994); Privitera  _et al._ (2018), recovered here as a
special instance of our more general theory.

##  VI Electronic polarization

By using the DHFT and the extended velocity operator definition, we now
calculate the adiabatic electronic polarization change of a fully occupied
valence band ${\Delta\mathbf{P}_{n}}$. In this framework, we calculate the
collective (adiabatic) electrons’ displacement
${\sum_{i}^{N}\Delta\left\langle\mathbf{r}_{i}\right\rangle}$ over the
boundaries (while a crystal is modified or deformed) that gives the induced
dipole moment change of the material, which in turn provides the quantum
electronic induced polarization change when divided by the volume of the
material.

One of the main differences between polarization change and charge pumping
processes, is the voltage-bias (non-zero mean force) that is present in the
polarization phenomena while is absent in pumping counterparts. In this
manner, the presence of voltage-bias enforces a breaking of periodic Born-von
Kármán boundary conditions between the edges of the material, hence
${\left|\Psi_{n}(\mathbf{r}+\mathbf{L},t)\right|^{2}\neq\left|\Psi_{n}(\mathbf{r},t)\right|^{2}}$,
although the wavefunction may retain its Bloch form in the bulk. According to
Ref. Resta and Vanderbilt, 2007, similar breaking of the periodic Born-von
Kármán boundary conditions emerges in piezoelectric phenomena due to strain,
which is not merely a perturbing term in the Hamiltonian.

For an applied electric filed that is now spatially uniform, each spinless
electron’s quantum state evolves in time according to the TDSE

$i\hbar\frac{d}{dt}\Psi(\mathbf{r},t)=\left(\frac{1}{2m}\mathbf{p}^{2}-e\mathbf{r}\\!\cdot\\!\mathbf{E}(t)+V_{crys}(\mathbf{r},t)\right)\Psi(\mathbf{r},t),$
(78)

where ${\mathbf{E}(t)}$ is the applied homogenous electric field that is also
assumed to be time-periodic and to satisfy ${\mathbf{E}(0)=\mathbf{E}(T)=0}$.
The crystal potential ${V_{crys}(\mathbf{r},t)}$ is assumed to have time-
dependence, either as a result of the externally applied electric field, or
due to applied strain, both deforming the crystal’s ions orderly arrangement.
For example, when computing the spontaneous polarization $\mathbf{P}_{s}$ in a
ferroelectric material such as ${\text{PbTiO}_{3}}$, one can let
${V_{crys}(\mathbf{r},0)}$ refer to the centrosymmetric cubic crystal
structure and ${V_{crys}(\mathbf{r},T)}$ to the noncentrosymetric,
ferroelectric crystal structure.

Due to the lack of translational invariance of the Hamiltonian, the wave
vector ${\mathbf{k}}$ is not a conserved quantity of the dynamics. We assume
the ansatz

$\psi(\mathbf{r},t,\mathbf{k})=e^{\displaystyle
i\mathbf{k}(t)\\!\cdot\\!\mathbf{r}}u(\mathbf{r},t,\mathbf{k}),$ (79)

where the wave vector ${\mathbf{k}(t)}$ is a continuous time-dependent vector
parameter. We define the time-dependent parameter ${\mathbf{R}(t)}$ as equal
to the time-dependent wave vector ${\mathbf{R}(t)\equiv\mathbf{k}(t)}$, and
then replace Eq. (79) into Eq. (78). The time-dependent wave-vector
(parameter) is defined to satisfy the equation of motion

$\frac{\partial\mathbf{k}(t)}{\partial t}=\frac{e}{\hbar}\mathbf{E}(t)$ (80)

which guarantees that, if the electron is initially in a Bloch state, it will
continue to be in a Bloch state after the electric field is turned on.
Therefore, the parameter evolves in time according to

$\mathbf{k}(t)=\mathbf{k}_{o}+\frac{e}{\hbar}\int_{0}^{\displaystyle
t}\mathbf{E}(t^{\prime})dt^{\prime}$ (81)

where ${\mathbf{k}_{o}}$ is the initial value of the parameter. Substituting
Eq. (79) into Eq. (78), and by using Eq. (80) gives

$i\hbar\frac{d}{dt}u(\mathbf{r},t,\mathbf{k})=H_{k}(\mathbf{r},t,\mathbf{k})u(\mathbf{r},t,\mathbf{k}),$
(82)

where the quantum state ${u(\mathbf{r},t,\mathbf{k})}$ evolves over time by
the Hamiltonian

$H_{k}(\mathbf{r},t,\mathbf{k})=\frac{1}{2m}\left(\mathbf{p}+\hbar\mathbf{k}(t)\right)^{2}\,+\,V_{crys}(\mathbf{r},t).$
(83)

We assume that each electron is initially in a ground, Bloch state, and that
the electric field turns on and changes very slowly over time
${T\rightarrow\infty}$, resulting to adiabatic time-evolution, namely

$u(\mathbf{r},t,\mathbf{k})\equiv e^{\displaystyle
i\Theta_{n}(t,\mathbf{k})}u_{n}(\mathbf{r},t,\mathbf{k}),$

where ${\Theta_{n}(t,\mathbf{k})}$ is the total phase of the wavefunction (the
sum of the dynamic and the geometric adiabatic phase) and
${u_{n}(\mathbf{r},t,\mathbf{k})}$ is the instantaneous eigenstate of the
Hamiltonian ${H_{k}(\mathbf{r},t,\mathbf{k})}$ satisfying the eigenvalue
equation

$H_{k}(\mathbf{r},t,\mathbf{k})u_{n}(\mathbf{r},t,\mathbf{k})=E_{n}(t,\mathbf{k})u_{n}(\mathbf{r},t,\mathbf{k}).$
(84)

We apply the adiabatic form of the DHFT Eq. (34) on Eq. (82) where the assumed
Hamiltonian is ${H_{k}(\mathbf{r},t,\mathbf{k})}$. The gradient of the
Hamiltonian with respect to the wave vector gives
${\displaystyle\bm{\nabla}_{\mathbf{k}}H_{k}(\mathbf{r},t,\mathbf{k})=\hbar\left(\frac{i}{\hbar}[H_{k}(\mathbf{r},t,\mathbf{k}),\mathbf{r}]\right)=\hbar\mathbf{v}}$,
where $\mathbf{v}$ is the standard velocity operator. Therefore, we find

$\displaystyle\hbar\left\langle\mathbf{v}\right\rangle_{n}$ $\displaystyle=$
$\displaystyle\bm{\nabla_{\mathbf{R}}}E_{n}(t,\mathbf{k})\,+\,\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k})$
$\displaystyle\,-\,\hbar\,\bm{\mathcal{E}}_{n}(t,\mathbf{k})\,-\,\hbar\,\frac{\partial\mathbf{k}}{\partial
t}\\!\times\bm{\mathcal{B}}_{n}(t,\mathbf{k}),$

where all involved quantities are evaluated with respect to the time-dependent
eigenstates ${u_{n}(\mathbf{r},t,\mathbf{k})}$. Specifically, the curvatures
are given by

$\displaystyle\bm{\mathcal{E}}_{n}(t,\mathbf{k})$ $\displaystyle=$
$\displaystyle
i\left\langle\bm{\nabla_{\mathbf{k}}}u_{n}(t,\mathbf{k})\,|\,\frac{\partial
u_{n}(t,\mathbf{k})}{\partial t}\right\rangle$ (86) $\displaystyle-$
$\displaystyle i\left\langle\frac{\partial u_{n}(t,\mathbf{k})}{\partial
t}|\,\bm{\nabla_{\mathbf{k}}}u_{n}(t,\mathbf{k})\right\rangle$

and

$\bm{\mathcal{B}}_{n}(t,\mathbf{k})=i\left\langle\bm{\nabla_{\mathbf{k}}}u_{n}(t,\mathbf{k})\right|\times\left|\bm{\nabla_{\mathbf{k}}}u_{n}(t,\mathbf{k})\right\rangle,$

and the non-Hermitian contribution is given from

$\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k})=\left\langle
u_{n}(t,\mathbf{k})\,|\,\left(H_{k}(t,\mathbf{k})^{+}-H_{k}(t,\mathbf{k})\right)\bm{\nabla_{\mathbf{k}}}u_{n}(t,\mathbf{k})\right\rangle.$
(87)

We note that, when the instantaneous eigenstates
${u_{n}(\mathbf{r},t,\mathbf{k})}$ do not have explicit time-dependence the
electric type of Berry curvature ${\bm{\mathcal{E}}_{n}(t,\mathbf{k})}$ is by
definition zero.

A crucial difference with respect to the theoretical approach that was made in
the study of the adiabatic topological charge pumping (Thouless pump) is that
the states ${u_{n}(\mathbf{r},t,\mathbf{k})}$ are assumed to be cell periodic
in the bulk of the material at every instant, but need not satisfy periodic
boundary conditions over the boundaries of the material, therefore
${\displaystyle\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k})}$ is not zero,
due to the presence of voltage-bias or strain. We can then use Eq. (47) that
relates the non-Hermitian contributions, that is,
${\displaystyle\bm{\mathcal{S}}_{n}(t,\mathbf{k})=\hbar\left\langle\,\mathbf{v}_{b}(t,\mathbf{k})\,\right\rangle_{n}\,+\,\widetilde{\bm{\mathcal{S}}_{n}}(t,\mathbf{k}),}$
where ${\left\langle\,\mathbf{v}_{b}(t,\mathbf{k})\,\right\rangle_{n}}$ is the
non-Hermitian boundary velocity expectation value (that is non-zero when the
position operator $\mathbf{r}$ becomes anomalous) and
${\bm{\mathcal{S}}_{n}(t,\mathbf{k})}$ is evaluated with respect to the
instantaneous Bloch state, namely,

$\displaystyle\bm{\mathcal{S}}_{n}(t,\mathbf{k})$ $\displaystyle=$
$\displaystyle\left\langle\psi_{n}(t,\mathbf{k})\,|\,\left(H(\mathbf{r},t)^{+}-H(\mathbf{r},t)\right)\bm{\nabla_{\mathbf{k}}}\psi_{n}(t,\mathbf{k})\right\rangle$
(88) $\displaystyle=$
$\displaystyle\frac{i\hbar}{2}\oiint_{S}\mathbf{n}\\!\cdot\\!\left(\,(\mathbf{v}\,\psi_{n})^{\dagger}+\psi_{n}^{\dagger}\,\mathbf{v}\,\right)\\!\,\bm{\nabla_{\mathbf{k}}}\psi_{n}\,dS,$

where ${\psi_{n}\equiv\psi_{n}(\mathbf{r},t,\mathbf{k})}$, $\mathbf{v}$ is the
standard velocity operator and $\mathbf{n}$ is the unit vector that is locally
normal to the surface $S$ that encloses the material. The non-Hermitian
contribution ${\bm{\mathcal{S}}_{n}(t,\mathbf{k})}$ is always a real quantity
according to Eq. (32). In this framework Eq. (VI) takes the form

$\displaystyle\hbar\left\langle\,\mathbf{v}\,\right\rangle_{n}+\hbar\left\langle\,\mathbf{v}_{b}\,\right\rangle_{n}$
$\displaystyle=$
$\displaystyle\bm{\nabla_{\mathbf{k}}}E_{n}(t,\mathbf{k})\,+\,\bm{\mathcal{S}}_{n}(t,\mathbf{k})$
$\displaystyle\,-\,\hbar\,\bm{\mathcal{E}}_{n}(t,\mathbf{k})\,-\,\hbar\,\frac{\partial\mathbf{k}}{\partial
t}\\!\times\bm{\mathcal{B}}_{n}(t,\mathbf{k}).$

By then using the extended velocity operator definition
${\displaystyle\left\langle\,\mathbf{v}_{ext}\,\right\rangle_{n}=\left\langle\,\mathbf{v}\,\right\rangle_{n}+\left\langle\,\mathbf{v}_{b}\,\right\rangle_{n}=\displaystyle\frac{d}{dt}\\!\left\langle\,\mathbf{r}\,\right\rangle_{n},}$
we obtain the electron’s mean position time derivative that is given by

$\displaystyle\frac{d}{dt}\\!\left\langle\,\mathbf{r}\,\right\rangle_{n}$
$\displaystyle=$
$\displaystyle\bm{\upsilon}_{n}(t,\mathbf{k})\,+\,\frac{1}{\hbar}\bm{\mathcal{S}}_{n}(t,\mathbf{k})$
$\displaystyle\,-\,\bm{\mathcal{E}}_{n}(t,\mathbf{k})\,-\,\frac{\partial\mathbf{k}}{\partial
t}\\!\times\bm{\mathcal{B}}_{n}(t,\mathbf{k}).$

where
${\displaystyle\bm{\upsilon}_{n}(t,\mathbf{k})=\frac{1}{\hbar}\bm{\nabla_{\mathbf{k}}}E_{n}(t,\mathbf{k})}$
is the group velocity.

Therefore, each electron’s adiabatic displacement after the cyclic evolution
of the Hamiltonian over a period of time $T$ is given by

$\Delta\left\langle
u_{n}(t,\mathbf{k})|\,\mathbf{r}\,|u_{n}(t,\mathbf{k})\right\rangle=\int_{0}^{T}\left(\bm{\upsilon}_{n}(t,\mathbf{k})\,+\,\frac{1}{\hbar}\bm{\mathcal{S}}_{n}(t,\mathbf{k})\,-\,\bm{\mathcal{E}}_{n}(t,\mathbf{k})\,-\,\frac{\partial\mathbf{k}}{\partial
t}\\!\times\bm{\mathcal{B}}_{n}(t,\mathbf{k})\right)dt.$ (91)

By then using Eq.(3) in the thermodynamic limit, the induced electronic
polarization change for a single fully occupied valence band is given by

$\Delta\mathbf{P}_{n}=\frac{e}{(2\pi)^{3}}\int\limits_{0}^{T}dt\iiint\limits_{BZ}\left(\bm{\upsilon}_{n}(t,\mathbf{k})\,+\,\frac{1}{\hbar}\bm{\mathcal{S}}_{n}(t,\mathbf{k})\,-\,\bm{\mathcal{E}}_{n}(t,\mathbf{k})\,-\,\frac{\partial\mathbf{k}}{\partial
t}\\!\times\bm{\mathcal{B}}_{n}(t,\mathbf{k})\right)d^{3}k.$ (92)

Some remarks are then in order regarding the way we have transformed the sum
of Eq. (3) into an integration with respect to the time-dependent parameter
${\mathbf{k}(t)}$ in Eq. (92). Mathematically speaking, the equation of motion
of the time-dependent parameter in the DHFT can always be defined arbitrarily.
Each elementary volume in parameter space (assumed to be occupied by a single
spinless electron) expands or dilates according to the equation
${\displaystyle d^{3}k(t)=\text{Det}J(t,t_{o})\ d^{3}k_{o}}$ where
${\text{Det}J(t,t_{o})}$ denotes the determinant of the Jacobian matrix of the
transformation that satisfies ${\text{Det}J(t_{o},t_{o})=1}$. The Jacobian
matrix, captures the transformation of the initial value of the parameter to a
subsequent one. The time evolution of the Jacobian matrix is solely determined
by the way that one defines the velocity of the parameter and is given by
${\displaystyle\frac{d}{dt}\text{Det}J(t,t_{o})=\left(\bm{\nabla_{\mathbf{k}}}\\!\cdot\\!\frac{\partial\mathbf{k}}{\partial
t}\right)\\!\text{Det}J(t,t_{o}).}$ In this framework, the equation of motion
of the elementary volume of the parameter is determined by
${\displaystyle\frac{d}{dt}\left(d^{3}k(t)\right)=\left(\bm{\nabla_{\mathbf{k}}}\\!\cdot\\!\frac{\partial\mathbf{k}}{\partial
t}\right)d^{3}k(t)}$. Owing to our definition of the parameter Eq. (80),
${\displaystyle\bm{\nabla_{\mathbf{k}}}\\!\cdot\\!\frac{\partial\mathbf{k}}{\partial
t}=0}$, which indicates that the elementary volume per particle in parameter
space remains constant ${d^{3}k_{o}=d^{3}k(t)}$, so, no modification is needed
when we transform the sum into an integration over time-dependent wavevectors.

Comparing Eq. (92) to the polarization change formula of the Modern’s Theory
of Polarization King-Smith and Vanderbilt (1993); Vanderbilt and King-Smith
(1993), we have rigorously found three extra contributions: one due to the
group velocity ${\bm{\upsilon}_{n}(t,\mathbf{k})}$, another one due to the
non-Hermitian boundary contribution ${\bm{\mathcal{S}}_{n}(t,\mathbf{k})}$,
and the last one due to the Berry curvature
${\bm{\mathcal{B}}_{n}(t,\mathbf{k})}$.

A question now emerges: is each one of these extra contributions non-zero
separately? To answer the question based on simple criteria, we must rely on
the localization of the assumed electrons’ wavefunctions (in general the
realistic boundary conditions of the Bloch eigenstates) and the discrete
symmetries of the assumed states. Due to localization, for example, in the
simplest model, when encountering only bulk localized states, in the sense
that each electron completely avoids the boundaries of the material, the non-
Hermitian contribution is zero ${\bm{\mathcal{S}}_{n}(t,\mathbf{k})=0}$
according to Eq. (88). On the other hand, in order to take into account the
consequences of the involved discrete symmetries, we present in Table 1 the
transformations of the quantities entering into Eq. (92) as a result of, time-
reversal symmetry (T.R.S.), space-inversion symmetry (S.I.S.) and mirror
symmetry along $x$ (M.S.); each symmetry is derived in Appendix C.

Table 1: Discrete symmetries for quantities involved in the polarization change. Quantity | T.R.S. | S.I.S. | M.S.
---|---|---|---
${\bm{\upsilon}_{n}(t,\mathbf{k})}$ | ${-\bm{\upsilon}_{n}(-t,-\mathbf{k})}$ | ${-\bm{\upsilon}_{n}(t,-\mathbf{k})}$ | ${-\upsilon_{n,x}(t,-k_{x},k_{y},k_{z})}$
${\bm{\mathcal{S}}_{n}(t,\mathbf{k})}$ | ${-\bm{\mathcal{S}}_{n}(-t,-\mathbf{k})}$ | ${-\bm{\mathcal{S}}_{n}(t,-\mathbf{k})}$ | ${-\mathcal{S}_{n,x}(t,-k_{x},k_{y},k_{z})}$
${\bm{\mathcal{E}}_{n}(t,\mathbf{k})}$ | ${-\bm{\mathcal{E}}_{n}(-t,-\mathbf{k})}$ | ${-\bm{\mathcal{E}}_{n}(t,-\mathbf{k})}$ | ${-\mathcal{E}_{n,x}(t,-k_{x},k_{y},k_{z})}$
${\bm{\mathcal{B}}_{n}(t,\mathbf{k})}$ | ${-\bm{\mathcal{B}}_{n}(-t,-\mathbf{k})}$ | ${\bm{\mathcal{B}}_{n}(t,-\mathbf{k})}$ | ${\mathcal{B}_{n,x}(t,-k_{x},k_{y},k_{z})}$

For example, in a simple model where no applied external electric field is
present, thus ${\displaystyle\frac{\partial\mathbf{k}}{\partial t}=0}$, but
the crystal potential is time-dependent ${V_{crys}(\mathbf{r},t)}$ due to
applied time-periodic strain, if the crystal potential remains space-inversion
symmetric at any instant, then, the change of electronic polarization over the
time interval $T$ is zero ${\Delta\mathbf{P}_{n}=0}$ as a result of Eq (92)
and the discrete symmetries of the involved quantities presented in Table 1.
Further polarization change arguments, based on discrete-symmetries, are
discussed below.

####  VI.0.1 Time-reversal symmetry

When the Hamiltonian ${H_{k}(\mathbf{r},t,\mathbf{k})}$ given by Eq. (83) is
time-reversal symmetric, namely, is invariant under the inversions
${(i,\mathbf{r},t,\mathbf{k}(t))\longmapsto(-i,\mathbf{r},-t,-\mathbf{k}(t))}$,
and provided that the Hamiltonian ${H(\mathbf{r},t)}$ of Eq. (78) has an
intermediate time period symmetry
${H(\mathbf{r},0)=H(\mathbf{r},T)=H(\mathbf{r},T/2)}$, (although
${\displaystyle\frac{dH(\mathbf{r},0)}{dt}\neq\frac{dH(\mathbf{r},T/2)}{dt}}$),
implying that the wavefunctions differ by a phase
${\psi_{n}(\mathbf{r},0)=e^{\displaystyle
i\Theta_{1k}}\psi_{n}(\mathbf{r},T)=e^{\displaystyle
i\Theta_{2k}}\psi_{n}(\mathbf{r},T/2)}$, leads to zero polarization change
${\Delta\mathbf{P}_{n}=0}$ for the time interval of $T$. We point out that,
the inversions actually assumed are
${(i,\mathbf{r},t,\mathbf{k}_{o})\longmapsto(-i,\mathbf{r},-t,-\mathbf{k}_{o})}$
where ${\mathbf{k}_{o}}$ is the initial value of the wavevector. The time
evolution of the wavevector
${\displaystyle\mathbf{k}(t)=\mathbf{k}_{o}+\frac{e}{\hbar}\int_{0}^{\displaystyle
t}\mathbf{E}(t^{\prime})dt^{\prime}}$, indicates that
${\mathbf{k}(t)\longmapsto-\mathbf{k}(t)}$ when
${\mathbf{E}(-t)=\mathbf{E}(t)}$, which is true due to the assumed time-
reversal symmetry. This, for example, gives zero polarization change for the
whole hysteresis loop for an assumed time-reversal electric field

$\displaystyle\displaystyle\mathbf{E}(t)=\left\\{\begin{array}[]{c c
r}\mathbf{E}_{o}\,\left|sin\left(\frac{2\pi}{T}t\right)\right|&\mbox{for}&0\leq
t\leq\frac{T}{2}\\\\[8.0pt]
-\mathbf{E}_{o}\,\left|sin\left(\frac{2\pi}{T}t\right)\right|&\mbox{for}&\frac{T}{2}\leq
t\leq{T},\end{array}\right.$

while for the half hysteresis loop ${0\leq t\leq T/2}$, this gives non-zero
polarization change along the direction of the electric field. Time-reversal
symmetry, can for example be broken when the material is magnetically ordered;
in this framework, each now spinful electron’s Hamiltonian Eq. (83) may assume
to have an extra term ${-J_{ex}\,\bm{\sigma}\\!\cdot\\!\hat{\mathbf{M}}}$,
where ${J_{ex}}$ is the exchange coupling strength between the electron spin
and the magnetization and ${\hat{\mathbf{M}}}$ is the magnetization unit
vector.

####  VI.0.2 Space-inversion and mirror symmetries

Space-inversion symmetry can for example be broken by an externally applied
electric field or by a time-dependent strain that deforms the crystal. The
time-dependent electric field breaks the space-inversion symmetry of the
Hamiltonian: under the inversions
${(i,\mathbf{r},t,\mathbf{k}_{o})\longmapsto(i,-\mathbf{r},t,-\mathbf{k}_{o})}$,
the time-dependent wavevector
${\displaystyle\mathbf{k}(t)=\mathbf{k}_{o}+\frac{e}{\hbar}\int_{0}^{\displaystyle
t}\mathbf{E}(t^{\prime})dt^{\prime}}$ does not transform like
${\mathbf{k}(t)\longmapsto-\mathbf{k}(t)}$ due to the electric field
${\mathbf{E}(t)}$ that is not inverted, therefore, the Hamiltonian
${H_{k}(\mathbf{r},t,\mathbf{k})}$ does not remain invariant. Similarly, an
assumed time-dependent strain that deforms the crystal breaks the space-
inversion symmetry of the Hamiltonian for the time interval that
${V_{crys}(\mathbf{r},t)\neq V_{crys}(-\mathbf{r},t)}$.

Likewise, mirror inversion symmetry is broken in the direction of the
externally applied field and in the directions that the crystal potential is
modified either due to lower equilibrium symmetry of the crystal or as a
result of deformation by an externally applied strain.

###  VI.1 Longitudinal and transverse polarization change

Assuming that the externally applied electric field is in the $x$ direction
${\displaystyle\frac{\partial\mathbf{k}(t)}{\partial
t}=\frac{e}{\hbar}E(t)\,\mathbf{e}_{x}}$, the polarization change Eq. (92) can
be separated into two normal contributions, the longitudinal one in the $x$
direction ${\Delta{P}_{n,x}=\Delta\mathbf{P}_{n}\\!\cdot\mathbf{e}_{x}}$ and
the transverse one normal to the electric field
${\Delta{P}_{n,y}=\Delta\mathbf{P}_{n}\\!\cdot\mathbf{e}_{y}}$.

The longitudinal polarization change in the $x$ direction, where the mirror
symmetry is broken, is given by

$\Delta{P}_{n,x}=\frac{e}{(2\pi)^{3}}\int\limits_{0}^{T}dt\iiint\limits_{BZ}\upsilon_{n,x}(t,\mathbf{k})d^{3}k\
\ +\
\frac{e}{\hbar(2\pi)^{3}}\int\limits_{0}^{T}dt\iiint\limits_{BZ}\mathcal{S}_{n,x}(t,\mathbf{k})d^{3}k\
-\
\frac{e}{(2\pi)^{3}}\int\limits_{0}^{T}dt\iiint\limits_{BZ}\mathcal{E}_{n,\,x}(t,\mathbf{k})d^{3}k,$
(95)

where the last term on the right side is the Modern’s Theory of Polarization
well-known contribution King-Smith and Vanderbilt (1993); Vanderbilt and King-
Smith (1993). On the other hand, the transverse polarization change is given
by

$\displaystyle\Delta{P}_{n,y}$ $\displaystyle=$
$\displaystyle\frac{e}{(2\pi)^{3}}\int\limits_{0}^{T}dt\iiint\limits_{BZ}\upsilon_{n,y}(t,\mathbf{k})d^{3}k\
\ +\
\frac{e}{\hbar(2\pi)^{3}}\int\limits_{0}^{T}dt\iiint\limits_{BZ}\mathcal{S}_{n,y}(t,\mathbf{k})d^{3}k\
-\
\frac{e}{(2\pi)^{3}}\int\limits_{0}^{T}dt\iiint\limits_{BZ}\mathcal{E}_{n,\,y}(t,\mathbf{k})d^{3}k$
(96)
$\displaystyle+\frac{e^{2}}{\hbar(2\pi)^{3}}\int\limits_{0}^{T}dt\,\left(E(t)\iiint\limits_{BZ}\mathcal{B}_{n,\,z}(t,\mathbf{k})d^{3}k\right)$

and has an extra bulk contribution, that is
${\frac{e^{2}}{\hbar(2\pi)^{3}}\int\limits_{0}^{T}dt\,\left(E(t)\iiint\limits_{BZ}\mathcal{B}_{n,\,z}(t,\mathbf{k})\right),}$
to be compared to the longitudinal counterpart Eq. (95). Assuming that, during
the polarization change process, the mirror symmetry remains unbroken along
the transverse $y$ direction, the first three contributions on the right side
of Eq. (96) integrate to zero according to Table 1. On the other hand, the
last contribution does not turn into zero due to mirror symmetry; instead, it
integrates to zero due to time-reversal symmetry if it is present. Therefore,
in magnetically ordered materials, where time-reversal symmetry is broken, the
last contribution is not zero (although mirror symmetry may be present along
$y$ direction). The above transverse polarization change formula, is plausible
to be compared to the one found by Ref. Goryo and Kohmoto, 2002 where the
electrons were subject to static and homogeneous electric and magnetic fields.

###  VI.2 Experimental setups and the non-Hermitian contribution

It is helpful to investigate the existence of the non-Hermitian contributions
in two important polarization change experiments. Namely, in the polarization
change experiment performed with a ferroelectric material and in the induced
polarization change experiment performed with a piezoelectric material.

####  VI.2.1 Spontaneous polarization change in a ferroelectric material

In an insulating ferroelectric material, the spontaneous polarization change
is experimentally evaluated by applying a time-periodic bias
${V(\mathbf{r},t)}$ at the opposite boundaries of the material (which is in a
capacitor configuration) that results to a time-varying current ${I(t)}$ that
flows through a shorten ammeter. Taking into account the hysteresis-loop of
the experiment, the measured switching current during the reversal of the
polarization is used in order to calculate the spontaneous polarization change
at zero electric field ${\mathbf{E}=0}$. Therefore, in the spontaneous
polarization change experiments, we have both time-varying externally applied
electric field and time-varying realistic boundary conditions (non-
periodicones) over the boundaries of the ferroelectric material. In this
framework (assuming non-magnetically ordered material), both integrands
containing ${\bm{\mathcal{S}}_{n}(t,\mathbf{k})}$ and
${\bm{\mathcal{E}}_{n}(t,\mathbf{k})}$ on the right side of Eq. (92) are not
zero and have to be taken into account. We note that, the space-coordinate
integration for calculating ${\bm{\mathcal{E}}_{n}(t,\mathbf{k})}$ is taken