# Microscopic formulas for thermoelectric transport coefficients in lattice
systems

Anton Kapustin<EMAIL_ADDRESS>Lev Spodyneiko
<EMAIL_ADDRESS>California Institute of Technology, Pasadena, CA 91125,
United States

###### Abstract

A macroscopic description of thermoelectric phenomena involves several
tensorial transport coefficients. Textbook microscopic Kubo formulas for them
are plagued with ambiguities in the definitions of the current operators and
the magnetization. We derive a version of these formulas for lattice systems
which is free from ambiguities but contains additional terms compared to the
textbook results. For symmetric components of thermoelectric tensors, we
identify a large class of lattice systems for which the additional terms
vanish with a natural choice of the energy current. To eliminate ambiguities
in the skew-symmetric components, one needs to interpret them as relative
quantities: only their differences for pairs of materials are well-defined.

##  I Introduction

Hydrodynamics in a general sense is a coarse-grained description of matter in
local thermal equilibrium which incorporates all local conservation laws. It
depends on the microscopic details only through transport coefficients.
Dissipative transport coefficients such as viscosity, thermal conductivity and
electric conductivity are particularly important for applications and their
theory is better known. More recently a lot of attention has been devoted to
non-dissipative time-reversal-odd transport coefficients such as electric Hall
conductivity and Hall viscosity. This interest is explained in part by the
contribution that gapless edge modes make to some of these coefficients at low
temperatures.

Theoretical studies of transport phenomena often take Kubo formulas as their
starting point. These are microscopic formulas for transport coefficients in
terms of correlators of currents of conserved quantities. Although Kubo
formulas go back to Kubo  _et al._ (2012); Luttinger (1964) and can be found
in many textbooks, there are a number of subtleties in their derivation. It is
well appreciated by now that naive Kubo formulas for skew-symmetric parts of
the transport tensors must be supplemented with additional terms involving
magnetization and ”energy magnetization” Cooper  _et al._ (1997). Such terms
affect thermal Hall conductivity and the skew-symmetric parts of
thermoelectric coefficients. Since magnetizations are intrinsically ambiguous,
it is not obvious how to evaluate such terms, see Qin  _et al._ (2011) for a
thorough discussion.

The theory of transport coefficients developed in Luttinger (1964); Cooper
_et al._ (1997); Qin  _et al._ (2011) applies to continuum systems. It cannot
be directly applied to lattice systems because it assumes certain scaling
relations for electric and energy currents which do not hold on a lattice. (In
fact, they do not hold for interacting continuum systems either, except after
some spatial averaging Cooper  _et al._ (1997)). Since tight-binding models
and other lattice Hamiltonians are ubiquitous in theoretical condensed matter
physics, it is important to develop an analogous formalism for them. In
Kapustin and Spodyneiko (2020) we described such a formalism and used it to
derive microscopic formulas for the Hall conductivity and thermal Hall
conductivity of lattice systems. In this paper we extend it to other transport
coefficients. For simplicity, we focus on the case of 2d lattice systems and
only briefly discuss the 3d case.

The main results of the paper are as follows. Our formulas for the symmetric
parts of conductivity and thermal conductivity tensors Kapustin and Spodyneiko
(2020) are completely analogous to continuum formulas. Surprisingly, this not
the case for the symmetric parts of thermoelectric tensors. In general,
microscopic formulas for them contain local terms as well as the expected Kubo
term. We show that these extra terms are in fact required to ensure that
transport coefficients are unaffected by the ambiguities in the definition of
the microscopic energy density. We also show that in special cases, such as
systems of free particles or systems with only density-dependent interactions,
the additional terms vanish with a natural definition of currents.

The skew-symmetric parts of all transport tensors except conductivity contain
contributions from magnetizations. Since magnetizations are defined only up to
additive constants, this leads to ambiguities. In the case of thermal
conductivity, a way to resolve the ambiguities was described in Kapustin and
Spodyneiko (2020) (building on the results of Cooper  _et al._ (1997)), and
the same approach works for thermoelectric tensors. Namely, although skew-
symmetric parts of these tensors are ”contaminated” with edge effects,
ambiguities cancel when one considers differences of transport tensors for two
materials. We express this by saying that skew-symmetric tensors are relative
transport coefficients. Microscopic formulas for relative transport
coefficients take a more complicated form: they are integrals of differential
1-forms along a path in the space of parameters. These issues do not affect
the skew-symmetric part of the conductivity tensor because one can, in
principle, measure it in a torus geometry, where no boundaries are present.
This is not possible to do even in principle for the skew-symmetric parts of
other transport coefficients.

The content of the paper is as follows. In Section II, we explain how
ambiguities in definition of transport and magnetization currents leads to a
natural separation of transport coefficients into absolute and relative ones.
In Section III, we derive microscopic formulas for thermoelectric transport
coefficients. We end with a discussion of possible generalizations of our
results in Section IV. In one of the appendices, we specialize our formulas to
the case of non-interacting fermions and express thermoelectric coefficients
in terms of zero-temperature 1-particle Green’s functions in coordinate space.

The work was supported in part by the U.S. Department of Energy, Office of
Science, Office of High Energy Physics, under Award Number DE-SC0011632. A. K.
was also supported by the Simons Investigator Award.

##  II Relative and absolute transport coefficients

The total current densities are usually divided into two parts:

$\displaystyle\begin{split}{\bf j}^{N}_{\rm tot}&={\bf j}^{N}+{\bf j}^{N}_{\rm
mag},\\\ {\bf j}^{E}_{\rm tot}&={\bf j}^{E}+{\bf j}^{E}_{\rm mag},\end{split}$
(1)

where magnetization currents are by definition divergence-free vector fields
which do not contribute to net currents across any section of the system.
Therefore they must have the form

$\displaystyle{\bf j}^{N}_{\rm mag}$ $\displaystyle=\nabla\times{\bf M}^{N},$
(2) $\displaystyle{\bf j}^{E}_{\rm mag}$ $\displaystyle=\nabla\times{\bf
M}^{E},$ (3)

where ${\bf M}^{N,E}$ are defined by these equations and are usually called
magnetization density and energy magnetization density, respectively. In the
following we will use the same term magnetization for both magnetization and
magnetization density which should not lead to a confusion.

The magnetization currents can be present even in an equilibrium state. The
transport currents ${\bf j}^{N,E}$, on the other hand, can be present only in
a non-equilibrium steady state created by slowly varying gradients of external
electric potential and temperature (for simplicity of presentation we assume
the chemical potential to be constant). This follows from the Bloch theorem
Watanabe (2019); Bohm (1949) and its energy analogue Kapustin and Spodyneiko
(2019). This constrains the form of transport currents.

Further constraints arise from gauge-invariance. It requires the transport
electric current ${\bf j}^{N}$ and the transport “heat current” ${\bf
j}^{E}-\varphi{\bf j}^{N}$ to be invariant under shifts of the electrical
potential $\varphi$ by a constant. This follows from the equations defining
the current operators ${\bf J}^{N,E}({\bf r})$

$i[H,h({\bf r})]=-\nabla\cdot{\bf J}^{E}({\bf r}),\quad i[H,\rho({\bf
r})]=-\nabla\cdot{\bf J}^{N}({\bf r}),$ (4)

and fact that under a constant gauge transformation $\varphi({\bf
r})\mapsto\varphi({\bf r})+c$ the energy density also transforms as $h({\bf
r})\mapsto h({\bf r})+c\rho({\bf r})$, where $\rho$ is the electric charge
density.

Taking all these considerations into account, one finds that to leading order
in the derivative expansion the transport electric current is given by

${j}^{N}_{k}=-\sigma_{km}\partial_{m}\varphi-\nu_{km}\partial_{m}T,$ (5)

where the conductivity tensor $\sigma_{km}$ and the thermoelectric tensor
$\nu_{km}$ are functions of temperature only. For the energy current the
expansion is

${j}^{E}_{k}=\varphi{j}^{N}_{k}-\eta_{km}\partial_{m}\varphi-\kappa_{km}\partial_{m}T.$
(6)

A crucial point for this paper is that the separation of the current densities
in (1) is ambiguous. One can always remove a curl of a vector field from ${\bf
j}^{N,E}$ and add it to ${\bf j}^{N,E}_{\rm mag}$ without affecting the
conservation equation and the form of the transport equations (5,6). While
this should have no effect on physically observable quantities, it can affect
the transport coefficients. Let us specialize to the 2d case and decompose all
two-index tensors into symmetric and anti-symmetric parts:
$\sigma_{km}=\sigma^{S}_{km}+\epsilon_{km}\sigma^{A}$ and similarly for the
tensors $\nu,\eta,$ and $\kappa$. Taking into account the requirement of
gauge-invariance, the allowed redefinitions of the transport currents have the
form

$\displaystyle j^{N}_{k}$ $\displaystyle\mapsto
j^{N}_{k}+\epsilon_{km}\partial_{m}\left(\sigma_{0}\varphi+f(T)\right),$ (7)
$\displaystyle j^{E}_{k}$ $\displaystyle\mapsto
j^{E}_{k}+\epsilon_{km}\partial_{m}\left(\frac{1}{2}\sigma_{0}\varphi^{2}+f(T)\varphi+g(T)\right).$
(8)

Here $\sigma_{0}$ is a constant and $f(T),g(T)$ are arbitrary functions of
$T$. Simultaneously magnetizations are redefined as follows:

$\displaystyle M^{N}$ $\displaystyle\mapsto M^{N}-\sigma_{0}\varphi-f(T),$ (9)
$\displaystyle M^{E}$ $\displaystyle\mapsto
M^{E}-\frac{1}{2}\sigma_{0}\varphi^{2}-f(T)\varphi-g(T).$ (10)

After the redefinition transport coefficients change:

$\displaystyle\begin{split}\sigma^{A}&\mapsto\sigma^{A}-\sigma_{0},\\\
\nu^{A}&\mapsto\nu^{A}-\frac{df(T)}{dT},\\\ \eta^{A}&\mapsto\eta^{A}-f(T),\\\
\kappa^{A}&\mapsto\kappa^{A}-\frac{dg(T)}{dT}.\end{split}$ (11)

Using such a redefinition we can always set $\sigma^{A}(T)$ to vanish at $T=0$
and make $\kappa^{A}$ and $\eta^{A}$ vanish identically for any homogeneous
material. Instead of setting $\eta^{A}(T)=0$, one can choose to set
$\nu^{A}(T)=0$ and use the remaining freedom to set $\eta^{A}(0)=0$. Note also
that $\dfrac{d{\sigma}^{A}}{dT}$ and $\nu^{A}-\dfrac{d\eta^{A}}{dT}$ are
invariant under such redefinitions.

So far we have ignored the vector potential, or equivalently gauge
transformations which depend on the spatial coordinates. Allowing such gauge
transformations changes the analysis as follows. Transport electric current
${\bf j}^{N}$ and transport heat current ${\bf j}^{E}-\varphi{\bf j}^{N}$ are
now required to depend on $\varphi$ only through the electric field
$E_{k}=-\partial_{k}\varphi-\frac{\partial A_{k}}{\partial t}$. The only
difference this makes is that only transformations (7) with $\sigma_{0}=0$ are
allowed. As a result, the Hall conductivity $\sigma^{A}$ is now free from
ambiguities.

There is a natural way to fix ambiguities in $M^{N}$ and $M^{E}$ and therefore
also in $\nu^{A},\eta^{A}$ and $\kappa^{A}$ Cooper  _et al._ (1997); Kapustin
and Spodyneiko (2020). If we consider a material with a boundary, the
magnetizations as well as all transport coefficients can be set to zero
outside. This removes all ambiguities from transport coefficients, but
obscures the fact that some transport coefficients are defined relative to
vacuum, while others do no depend on any choices and can be measured in the
bulk. We will call them relative and absolute transport coefficients,
respectively. According to the above analysis, all symmetric transport
coefficients as well as $\sigma^{A}$ are absolute, while $\nu^{A},\eta^{A},$
and $\kappa^{A}$ are relative. The combination $\nu^{A}-\dfrac{d\eta^{A}}{dT}$
is also absolute.

This distinction has consequences for the microscopic formulas that can be
derived for transport coefficients which are usually called Kubo formulas
Kapustin and Spodyneiko (2020). As we just explained, determination of
relative transport coefficients require considering a system with boundaries.
On the other hand, as we show in the paper, derivatives of relative transport
coefficients with respect to temperature or the parameters of the Hamiltonian
involves only correlation functions of a system without boundary. The values
of relative transport coefficients for any particular material can be found by
integrating this differential over the parameters and/or temperature. The non-
uniqueness in the choice of the base point of the integral reflects the
ambiguity in the definition of the magnetization currents and can be fixed by
choosing the base point to be a trivial insulator. The resulting microscopic
formula for a relative transport coefficient is manifestly independent of the
choice of boundary conditions at the cost of depending on the correlation
functions of a whole family of systems which interpolates between the system
of interest and a trivial insulator. On the other hand, microscopic formulas
for absolute transport coefficients depend only on the correlation function of
the system at a fixed temperature and values of all parameters.

As a consistency check, let us verify that physical bulk quantities depend
only on absolute transport coefficients. For the time derivatives of charge
and energy densities we get

$\displaystyle\frac{\partial\rho^{N}}{\partial t}=-\nabla\cdot{\bf j}^{N}$
$\displaystyle=-\sigma_{km}^{S}\partial_{k}{\bf
E}_{m}-\frac{{d\sigma}_{km}}{dT}{\bf
E}_{m}\partial_{k}T-\nu^{S}_{km}\partial_{k}\partial_{m}T-\frac{{d\nu}^{S}_{km}}{dT}\partial_{k}T\partial_{m}T,$
(12) $\displaystyle\begin{split}\frac{\partial\rho^{E}}{\partial
t}=-\nabla\cdot{\bf j}^{E}&=-\varphi\nabla\cdot{\bf j}^{N}+\sigma^{S}_{km}{\bf
E}_{k}{\bf E}_{m}-\left(\nu_{km}+\frac{d{\eta}_{mk}}{dT}\right){\bf
E}_{k}\partial_{m}T-\eta^{S}_{km}\partial_{k}{\bf E}_{m}\\\
&+\kappa^{S}_{km}\partial_{k}\partial_{m}T+\frac{d{\kappa}^{S}_{km}}{dT}\partial_{k}T\partial_{m}T,\end{split}$
(13)

where the external fields are assumed to be time-independent and thus
$\nabla\times{\bf E}=0$. As expected, these time derivatives are unaffected by
the transformations (11).

The time-derivative of entropy density is usually written as follows Landau
and Lifshitz (1984):

$\frac{\partial s}{\partial t}=\frac{1}{T}\sigma^{S}_{km}{\bf E}_{k}{\bf
E}_{m}+\frac{1}{T^{2}}\kappa^{S}_{km}\partial_{k}T\partial_{m}T-\frac{1}{T^{2}}(T\nu_{km}+\eta_{mk}){\bf
E}_{k}\partial_{m}T-\nabla\cdot{\tilde{\bf j}}^{S}_{k},$ (14)

where the entropy current density is

${\tilde{\bf j}}^{S}_{k}=\frac{1}{T}\eta_{km}{\bf
E}_{m}-\frac{1}{T}\kappa^{S}_{km}\partial_{m}T.$ (15)

The r.h.s. of eq. (14) seems to depend on some relative transport
coefficients. However, if one redefines the entropy current as follows:

${\bf j}^{S}_{k}=\frac{1}{T}\eta^{S}_{km}{\bf
E}_{m}-\frac{1}{T}\kappa^{S}_{km}\partial_{m}T,$ (16)

then one can write

$\frac{\partial s}{\partial t}=\frac{1}{T}\sigma^{S}_{km}{\bf E}_{k}{\bf
E}_{m}+\frac{1}{T^{2}}\kappa^{S}_{km}\partial_{k}T\partial_{m}T-\frac{1}{T^{2}}\left(T\nu^{S}_{km}+\eta^{S}_{mk}+\epsilon_{km}T\left(\nu^{A}-\frac{d{\eta}^{A}}{dT}\right)\right){\bf
E}_{k}\partial_{m}T-\nabla\cdot{\bf j}^{S}_{k}.$ (17)

It is manifest now that both the entropy production rate and the entropy
current depend only on the absolute transport coefficients.

Since only absolute transport coefficients enter the expressions for the
divergences of currents, measuring net currents through closed curves (or
surfaces, if we are discussing a 3d material) does not allow to determine
relative transport coefficients. This applies even to infinite curves with a
boundary at infinity, provided $\varphi$ and $T$ tend to fixed values at
infinity. The latter condition must be imposed to eliminate the contribution
of magnetization currents. For example, if we compute the electric current
$I^{N}_{x}$ through a vertical line $x=0$, then the contribution of $\nu^{A}$
drops out because

$\int_{-\infty}^{\infty}dy\,\nu^{A}\partial_{y}Tdy=0.$ (18)

The above considerations apply to a homogeneous material whose transport
coefficients are constants. If one considers a heterogeneous material, such as
an interface between two homogeneous ones, then the expressions for net
currents will involve differences between relative transport coefficients. For
example, consider a sample such that $\nu^{A}$ interpolates between
$\nu^{A}_{1}$ for $y\ll 0$ and $\nu^{A}_{2}$ for $y\gg 0$. Suppose that the
temperature is a function of $y$ only which is equal to $T_{b}$ throughout the
interface region and approaches $T_{\infty}$ at $y\rightarrow\pm\infty$. Then
the $\nu^{A}$-dependent contribution to the net electric current in the $x$
direction is

$\displaystyle-\int_{-\infty}^{\infty}dy\,\nu^{A}\partial_{y}T$
$\displaystyle=(T_{b}-T_{\infty})(\nu_{2}^{A}-\nu_{1}^{A}),$ (19)

Similarly, the contribution of $\kappa^{A}$ to the net heat current is
$(T_{b}-T_{\infty})(\kappa_{2}^{A}-\kappa_{1}^{A}).$ By creating an electric
potential $\varphi$ which is equal to $\varphi_{b}$ in the interface region
and approaches $\varphi_{\infty}$ at $y\rightarrow\pm\infty$, one can also
measure $\sigma^{A}_{2}-\sigma^{A}_{1}$ and $\eta_{2}^{A}-\eta^{A}_{1}$.

In the case of the electric Hall conductivity one can do better by utilizing a
time-dependent vector potential rather than a scalar potential and working in
a cylinder geometry or a torus geometry. Then in principle one can determine
$\sigma^{A}$ for a single material by measuring the net flow of electric
charge across a section of a cylinder or a torus as one inserts a unit of
magnetic flux through this section. This does not work for $\eta^{A}$ because
the physical quantity that needs to be measured is the net amount of heat
transferred to the heat bath as one inserts a unit of magnetic flux (see Fig.
1). Therefore heat transport will receive a contribution from the work of the
electromotive force $\mathcal{E}$ on the net electric edge currents
$I^{N}_{\rm edge}$. The edge currents are proportional to the jump in the
magnetization along the boundary and they make $\eta^{A}$ relative even in the
cylinder geometry.

Figure 1: Insertion of a flux $\Phi$ into a cylinder creates electric field
around it. If $\eta_{xy}\neq 0$, this electric field $\bf E$ drives a heat
current $I^{\mathcal{Q}}$ along the axis of the cylinder. The jump of the
magnetization gives rise to an edge current $I^{N}_{\rm edge}$. Work done on
the edge current by the electric field contributes to the net heat
$\mathcal{Q}$ transferred to the heat bath.

##  III Microscopic formulas for thermoelectric coefficients

###  III.1 Currents on a lattice

We follow the conventions of Kapustin and Spodyneiko (2020). We consider a
lattice system with a Hamiltonian $H=\sum_{p\in\Lambda}H_{p}$, where
$\Lambda\subset{\mathbb{R}}^{2}$ is a not necessarily regular lattice. The
operators $H_{p}$ have a finite range, i.e. there exist $R$ such that $H_{p}$
acts trivially on site $q$ if $|p-q|>R$. The space of states at each site of
the lattice is assumed to be finite-dimensional. The electric charge operator
$Q$ has the form $Q=\sum_{p\in\Lambda}Q_{p}$, where $Q_{p}$ has integral
eigenvalues (we set the electric charge of electron to be 1) and acts only on
site $p$. This means that the $U(1)$ symmetry is on-site. In particular,
$[Q_{p},Q_{q}]=0$ for all $p,q\in\Lambda$.

The electric current from site $q$ to site $p$ is defined as

$J^{N}_{pq}=i[H_{q},Q_{p}]-i[H_{p},Q_{q}].$ (20)

The energy current from site $q$ to site $p$ is defined as

$J^{E}_{pq}=-i[H_{p},H_{q}].$ (21)

These currents enter the charge and energy conservation equations

$\displaystyle\frac{dQ_{q}}{dt}$
$\displaystyle=i\sum_{p\in\Lambda}[H_{p},Q_{q}]=-\sum_{p\in\Lambda}J^{N}_{pq},$
(22) $\displaystyle\frac{dH_{q}}{dt}$
$\displaystyle=i\sum_{p\in\Lambda}[H_{p},H_{q}]=-\sum_{p\in\Lambda}J^{E}_{pq}.$
(23)

The net current from a subset $B\subset\Lambda$ to its complement
$A=\Lambda\backslash B$ is given by

$\displaystyle J^{N,E}(A,B)=\sum_{p\in A}\sum_{q\in B}J^{N,E}_{pq}.$ (24)

This observable measures the net current across the boundary of $A$ and $B$.
In this paper we will need its mild generalization. For a given skew-symmetric
function $\eta(p,q):\Lambda\times\Lambda\rightarrow\mathbb{R}$ satisfying

$\eta(p,q)+\eta(q,r)+\eta(r,p)=0,\quad\forall p,q,r\in\Lambda,$ (25)

define

$\displaystyle
J^{N,E}(\eta)=\frac{1}{2}\sum_{p,q\in\Lambda}\eta(p,q)J^{N,E}_{pq}.$ (26)

Current $J(A,B)$ from $B$ to $A=\Lambda\backslash B$ corresponds to the case
$\eta(p,q)=\chi_{B}(q)-\chi_{B}(p)$, where the function $\chi_{B}(p)$ is equal
1 on $B$ and 0 otherwise.

The equilibrium expectation value of the currents satisfy

$\sum_{p\in\Lambda}\langle J^{N,E}_{pq}\rangle=0.$ (27)

This equation is a lattice analog of the continuum equation

$\displaystyle\nabla\cdot\langle{\bf J}^{N,E}({\bf r})\rangle=0.$ (28)

In the continuum the general solution to this equation

$\displaystyle\langle J_{k}^{N,E}({\bf
r})\rangle=-\epsilon_{kj}\partial_{j}M^{N,E}({\bf r})$ (29)

defines the magnetization $M^{N}$ and energy magnetization $M^{E}$.
Analogously, on a lattice the solution of equation (27) is

$\langle J^{N,E}_{qr}\rangle=\sum_{p\in\Lambda}M^{N,E}_{pqr},$ (30)

where $M^{N,E}_{pqr}$ are skew-symmetric functions of the lattice points
$p,q,r\in\Lambda$. These are lattice analogs of the magnetization and the
energy magnetisation. Physically, in continuum case $M^{N,E}(r)$ represent the
circulating currents of the system in equilibrium. Similarly, $M^{N,E}_{pqr}$
physically can be thought as quantity which measures the circulating current
around a triangle formed by $p,q,r$.

Unfortunately, $M^{N,E}_{pqr}$ is not unique: one can always redefine

$M^{N,E}_{pqr}\mapsto M^{N,E}_{pqr}+\sum_{s\in\Lambda}N_{pqrs},$ (31)

where $N_{pqrs}$ is skewsymmetric function of its subscripts which decays
whenever any two of them are far apart. This corresponds to ambiguity in
splitting of the circulating currents into contributions of magnetization from
the different triangles and it is absent in continuum case.

There is an additional ambiguity corresponding to existence of solutions to
the equation $\sum_{r\in\Lambda}M^{N,E}_{pqr}=0$ which are not of the form
$\sum_{s\in\Lambda}N_{pqrs}$. It corresponds to an ambiguity of addition of a
constant to the magnetization in continuum case. A standard method to deal
with the later is to consider a system with a boundary and fix the
magnetization to be zero outside of the system. In this paper, we want to
think about all transport coefficients as manifestly bulk quantities and avoid
dealing with boundaries. While magnetization itself suffers from ambiguities
and depends non-locally on the boundary conditions, the variation of
magnetization with respect to parameters of the Hamiltonian is local. Indeed,
consider the variation of the equation (30) with respect to a parameter
$\lambda^{\ell}$ of the Hamiltonian

$\frac{\partial}{\partial\lambda^{\ell}}\langle
J^{N,E}_{pq}\rangle=\sum_{r\in\Lambda}\mu^{N,E}_{pqr,\ell},$ (32)

where $\mu^{N,E}_{pqr,\ell}=\dfrac{\partial
M^{N,E}_{pqr}}{\partial\lambda^{\ell}}$ and is given by Kitaev (2006)

$\mu^{N,E}_{pqr,\ell}=-\beta\langle\langle\frac{\partial
H_{p}}{\partial\lambda^{\ell}};J^{N,E}_{qr}\rangle\rangle-\beta\langle\langle\frac{\partial
H_{r}}{\partial\lambda^{\ell}};J^{N,E}_{pq}\rangle\rangle-\beta\langle\langle\frac{\partial
H_{q}}{\partial\lambda^{\ell}};J^{N,E}_{rp}\rangle\rangle,$ (33)

where $\langle\langle A;B\rangle\rangle$ denotes the Kubo canonical pairing
Kubo  _et al._ (2012). Using the properties of the Kubo pairing (see Appendix
A), one can easily verify the identity (32). In the following we will combine
derivatives of magnetizations with respect to different parameters into
1-forms on the parameter space
$\mu_{pqr}^{N,E}=\sum_{\ell}\mu^{N,E}_{pqr,\ell}d\lambda^{\ell}$.

###  III.2 Equilibrium conditions and driving forces

In the following sections we will follow Luttinger Luttinger (1964) and study
the behavior of the system coupled to external potentials

$H_{p}^{\psi,\varphi}=(1+\psi(p))(H_{p}+\varphi(p)Q_{p}),$ (34)

where $\varphi(p)$ is external electric potential and $\psi(p)$ can be thought
of as gravitational potential. The potentials are assumed to infinitesimally
small slowly-varying functions of $p$ which vanish at infinity. After coupling
to external potentials the system will eventually relax into a state with
density matrix

$\rho\sim\exp\left(-\frac{H^{\psi,\varphi}-\mu_{0}Q}{T_{0}}\right),$ (35)

where $T_{0}$ and $\mu_{0}$ are the temperature and local chemical potential
of the system at infinity. On the other hand, on physical grounds we expect
local observables supported in some small but macroscopic region around site
$p$ to be described by a thermal density matrix

$\rho(p)\sim\left(-\frac{H-\mu(p)Q}{T(p)}\right),$ (36)

where the local temperature $T(p)$ and chemical potential $\mu(p)$ are slowly
varying functions of $p$. The equilibrium conditions can be found to be
Luttinger (1964); Cooper  _et al._ (1997); Qin  _et al._ (2011)

$\displaystyle(1+\psi(p))(\mu(p)+\varphi(p))$ $\displaystyle=\mu_{0},$ (37)
$\displaystyle(1+\psi(p))T(p)$ $\displaystyle=T_{0}.$ (38)

These relations together with the absence of transport currents in equilibrium
can be used to derive Einstein relations between transport coefficients
Luttinger (1964). The latter also leads to transport currents being
proportional to the gradients of the left hand sides of eqs. (37) and (38).
The fact that the driving forces depend only on specific combinations of
$\psi,\varphi,T,\mu$ will be used to relate the response to variations of the
thermodynamic parameter $T$ to the response to variations of the external
field $\psi$.

Note that, in section II for simplicity we assumed that chemical potential is
constant. We will continue doing so in the rest of the paper and we will
assume $\mu(p)$ to be constant. Generalization to a more general case of non-
constant local chemical potential is straightforward, but clumsy due to
ambiguities in the separation of electrochemical potential $\mu+\varphi$ into
two parts.

###  III.3 Nernst effect

In order to find the thermoelectric coefficient coefficient $\nu_{xy}$ we
deform the Hamiltonian density by

$\displaystyle\Delta H_{p}=\epsilon e^{st}{\overset{\frown}{g}}(p)H_{p},$ (39)

where ${\overset{\frown}{g}}(p)$ is a hat-shaped function as in Fig. 2,
$\epsilon$ is an infinitesimal parameter, and $s$ is a small positive number
which controls how fast the perturbation is turned on. We will consider the
so-called fast regime Luttinger (1964) in which the characteristic time $1/s$
is large but not large enough in order for the two slopes of the hat
${\overset{\frown}{g}}(p)$ to come into equilibrium.

(a) (a) (b) (c)

Figure 2:  (a) Heat map of the function $g(p)$ corresponding to a temperature
gradient in a horizontal strip. (b) Heat map of the function
${\overset{\frown}{g}}(p)$ corresponding to two horizontal strips with the
opposite signs of the temperature gradient. (c) The function $g(p)$ restricted
to the dashed line in (a). (d) The function ${\overset{\frown}{g}}(p)$
restricted to the dashed line in (b). Dotted line in (d) represents the
dependence of some parameter $\lambda$ of the Hamiltonian on $y(p)$.

The change of the state of the system can be found as follows. The density
matrix

$\displaystyle\rho(t)=\rho_{0}+\Delta\rho(t)$ (40)

satisfies the quantum Liouville equation

$\frac{d\Delta\rho}{dt}=-i[\Delta H,\rho_{0}]+\dots,$ (41)

where $\rho_{0}$ is the equilibrium density matrix at $t=-\infty$ and dots
represent higher order terms in $\epsilon$. The solution to this equation is

$\displaystyle\Delta\rho(t=0)=-\rho_{0}\int_{0}^{\infty}dt\int_{0}^{\beta}d\tau\Delta\dot{H}(-t-i\tau),$
(42)

where the dot denotes the time derivative. The change of the observable $A$
can be found to be

$\displaystyle\Delta\langle A\rangle=\langle\Delta
A\rangle-\beta\int_{0}^{\infty}dt\langle\langle
A;\Delta\dot{H}(-t)\rangle\rangle,$ (43)

where we used the Kubo pairing notation (see Appendix A) and $\Delta A$ is the
variation of the operator arising from explicit dependence of $A$ on the
Hamiltonian. Using explicit form of the perturbation (39), energy conservation
law and properties of the Kubo pairing we can rewrite this formula as

$\displaystyle\Delta\langle A\rangle=\langle\Delta
A\rangle+\epsilon\beta\int_{0}^{\infty}dte^{-st}\langle\langle
A(t);J^{E}(\delta{\overset{\frown}{g}})\rangle\rangle,$ (44)

where we neglected term proportional to small $s$.

The change in the electric current across a vertical line $x=a$ is

$\Delta\langle J^{N}(\delta f)\rangle=\langle\Delta J^{N}(\delta
f)\rangle+\epsilon\beta\int_{0}^{\infty}dte^{-st}\langle\langle J^{N}(\delta
f,t);J^{E}(\delta{\overset{\frown}{g}})\rangle\rangle,$ (45)

where $f(p)=\theta(a-x(p))$ is a step function. The explicit variation of the
current is

$\displaystyle\begin{split}&\langle\Delta J^{N}(\delta
f)\rangle=\frac{i\epsilon}{2}\sum_{p,q\in\Lambda}\big{\langle}{\overset{\frown}{g}}(q)[H_{q},Q_{p}]-{\overset{\frown}{g}}(p)[H_{p},Q_{q}]\big{\rangle}(f(q)-f(p))\\\
&=\frac{\epsilon}{4}\sum_{p,q\in\Lambda}\langle
J^{N}_{pq}\rangle({\overset{\frown}{g}}(p)+{\overset{\frown}{g}}(q))(f(q)-f(p))\\\
&\qquad+\frac{i\epsilon}{4}\sum_{p,q\in\Lambda}\langle[H_{q},Q_{p}]+[H_{p},Q_{q}](f(q)-f(p))({\overset{\frown}{g}}(q)-{\overset{\frown}{g}}(p))\rangle,\end{split}$
(46)

where in the last line we separated the result into two contribution formally
skew-symmetric and symmetric in $f,{\overset{\frown}{g}}$. The second term
depends only on difference of values of ${\overset{\frown}{g}}$ at different
sites as expected for a transport current. On the other hand, the first term
depends on the value of ${\overset{\frown}{g}}$ and does not seem to be of the
form expected for a transport current. As was explained in Cooper  _et al._
(1997) this contribution is related to magnetization currents. Indeed we can
rewrite

$\displaystyle\frac{\epsilon}{4}\sum_{p,q\in\Lambda}\langle
J^{N}_{pq}\rangle({\overset{\frown}{g}}(p)+{\overset{\frown}{g}}(q))(f(q)-f(p))=\sum_{p,q,r\in\Lambda}M^{N}_{pqr}(f(q)-f(p))({\overset{\frown}{g}}(q)-{\overset{\frown}{g}}(r)),$
(47)

where we used skew-symmetry of $M_{pqr}^{N}$.

Written in this form, the response is proportional to the differences of
${\overset{\frown}{g}}$ at different points and therefore receives appreciable
contributions only from regions $I$ and $II$ in Fig. 2. However, the
transformation (47) contains an important subtlety. The right-hand side hand
side contains a magnetization contribution which is ambiguously defined while
left-hand side is unambiguous. There is no contradiction because the ambiguity
in the region $I$ will compensate the one in the region $II$. Moreover, in a
homogeneous system the response will be zero, because these two regions
compensate each other exactly. One way to deal with it is to introduce a
boundary somewhere in between the two regions and enforce the magnetization
$M$ to be 0 outside of the sample. This approach is used in Cooper  _et al._
(1997) in the continuum case. However, an explicit boundary introduces
additional computational challenges and makes the bulk nature of the Nernst
effect obscure.

In this paper, we will use an alternative approach proposed in Kapustin and
Spodyneiko (2020). Instead of introducing a sharp boundary, we will make one
parameter of the Hamiltonian $\lambda$ to have slightly different values in
regions $I$ and $II$ (see Fig. 2). We can write the hat-shaped function
${\overset{\frown}{g}}$ as a difference of two functions $g_{I}$ and $g_{II}$,
${\overset{\frown}{g}}(p)=g_{II}(p)-g_{I}(p)$, each of which is a translate of
the function $g(p)$ which depends only on $y(p)$ and is shown in Fig. 2. The
functions $g_{I,II}$ are non-constant only in regions $I$ and $II$
respectively. Then the magnetization contribution can be rewritten as

$\epsilon(\lambda_{II}-\lambda_{I})\mu^{N}_{\lambda}(\delta g\cup\delta
f)+O((\lambda_{II}-\lambda_{I})^{2}),$ (48)

where $\mu^{N}_{pqr,\lambda}=\frac{\partial M^{N}_{pqr}}{\partial\lambda}$,
and we introduced a notation

$\displaystyle\mu^{N}_{\lambda}(\delta f_{1}\cup\delta
f_{2})=\frac{1}{6}\sum_{p,q,r\in\Lambda}\mu^{N}_{pqr,\lambda}(f_{1}(q)-f_{2}(p))(f_{2}(r)-f_{2}(q)).$
(49)

Combining this with other contributions we find

$\displaystyle\begin{split}\delta\langle J^{N}(\delta
f)\rangle&\approx\epsilon(\lambda_{II}-\lambda_{I})\Bigg{\\{}\frac{\partial}{\partial\lambda}\Big{[}\beta\int_{0}^{\infty}dte^{-st}\langle\langle
J^{N}(\delta f,t);J^{E}(\delta g)\rangle\rangle+U(\delta f,\delta
g)\Big{]}+\mu^{N}_{\lambda}(\delta g\cup\delta f)\Bigg{\\}},\end{split}$ (50)

where

$\displaystyle U(\delta f,\delta
g)=\frac{i}{4}\sum_{p,q\in\Lambda}\langle[H_{q},Q_{p}]+[H_{p},Q_{q}](f(q)-f(p))(g(q)-g(p))\rangle.$
(51)

The function $g$ as in Fig. 2 is not compactly supported and thus it takes an
infinite time for the system to equilibrate. Therefore, one can take the limit
$s\rightarrow 0$ while staying in the “fast” regime. From eqs. (37,38) one
finds that electric current after perturbation by gravitational potential
$\epsilon g(p)$ is equal to the current generated by

$\displaystyle T(p)$ $\displaystyle=\epsilon g(p)T_{0},$ (52)
$\displaystyle\varphi(p)$ $\displaystyle=\epsilon g(p)\mu_{0}.$ (53)

From continuum phenomenological equation (5) one finds the current across
$x=a$ line to be

$-\epsilon
T_{0}\int_{-\infty}^{\infty}\nu_{xy}\partial_{y}gdy-\epsilon\mu_{0}\int_{-\infty}^{\infty}\sigma_{xy}\partial_{y}gdy=\epsilon
T_{0}\int_{\lambda_{I}}^{\lambda_{II}}\frac{\partial\nu_{xy}}{\partial\lambda}d\lambda+\epsilon\mu_{0}\int_{\lambda_{I}}^{\lambda_{II}}\frac{\partial\sigma_{xy}}{\partial\lambda}d\lambda.$
(54)

Comparing it to (50) we find

$\displaystyle\begin{split}d\nu_{xy}&=d\Big{[}\beta^{2}\lim_{s\rightarrow
0}\int_{0}^{\infty}dte^{-st}\langle\langle J^{N}(\delta
f,t);J^{\mathcal{Q}}(\delta g)\rangle\rangle+\beta U(\delta f,\delta
g)\Big{]}-\beta\mu^{N}(\delta f\cup\delta g),\end{split}$ (55)

where we introduced the notation $J^{\mathcal{Q}}=J^{E}-\mu J^{N}$ for the
heat current, inverse temperature $\beta=\frac{1}{T}$, we dropped the
subscript 0 from $T_{0}$ and $\mu_{0}$ since this formula contains correlation
functions of the unperturbed system in equilibrium. We combined differential
with respect to parameter into the differential form $\mu^{N}(\delta
f\cup\delta g)=\mu^{N}_{\lambda}(\delta f\cup\delta g)d\lambda$ and the
derivation can be straightforwardly extended to involve several parameters.
The exterior derivative
$d=\sum_{\ell}d\lambda\frac{\partial}{\partial\lambda^{\ell}}$ acts on the
parameter space.

Since rescaling the temperature is equivalent to rescaling the Hamiltonian, we
can extend this 1-form to the enlarged parameter space which includes $T$.
Then we can define the difference of coefficients $\eta_{xy}$ for any two 2d
materials, regardless of the temperature. Explicitly, let us define the
rescaled Hamiltonian $H_{\lambda_{0}}=\lambda_{0}H$, where we introduced an
additional scaling parameter $\lambda_{0}$. Then

$\left(T\frac{\partial}{\partial
T}+\left.{\lambda_{0}}\frac{\partial}{\partial{\lambda_{0}}}\right|_{{\lambda_{0}}=1}\right)\nu_{xy}({\lambda_{0}},T)=0.$
(56)

Therefore we can define the $T$-component of the 1-form on the enlarged
parameter space as follows:

$\frac{d\nu_{xy}}{dT}=\frac{\partial}{\partial
T}\left[\beta^{2}\int_{0}^{\infty}dte^{-st}\langle\langle J^{N}(\delta
f,t);J^{\mathcal{Q}}(\delta g)\rangle\rangle+\beta U(\delta f,\delta
g)\right]-\beta^{2}\tau^{N}(\delta f\cup\delta g),$ (57)

where $\tau^{N}(\delta f\cup\delta g)$ is given by eq. (49) with $\mu_{pqr}$
replaced with

$\tau^{N}_{pqr}=\beta\langle\langle
H_{p};J^{N}_{qr}\rangle\rangle+\beta\langle\langle
H_{r};J^{N}_{pq}\rangle\rangle+\beta\langle\langle
H_{q};J^{N}_{rp}\rangle\rangle,$ (58)

which is obtained from $\mu^{N}$ by replacing $dH_{p}$ with $-H_{p}$.

###  III.4 Ettingshausen effect

One can derive a formula for the coefficient $\eta_{xy}$ in a similar way. In
this section, we will display only the key steps, since all arguments are the
same.

In order to find the coefficient we deform the Hamiltonian density by

$\displaystyle\Delta H_{p}=\epsilon e^{st}{\overset{\frown}{g}}(p)Q_{p},$ (59)

where ${\overset{\frown}{g}}(p)$ is a hat-shaped function of $y(p)$ as in Fig.
2.

The change in the energy current across a vertical line $x=a$ is

$\Delta\langle J^{E}(\delta f)\rangle=\langle\Delta J^{E}(\delta
f)\rangle+\epsilon\beta\int_{0}^{\infty}dte^{-st}\langle\langle J^{E}(\delta
f,t);J^{N}(\delta{\overset{\frown}{g}})\rangle\rangle,$ (60)

where $f=\theta(a-x(p))$. The explicit variation of the current is

$\displaystyle\begin{split}\langle\Delta J^{E}(\delta
f)\rangle&=\frac{i\epsilon}{2}\sum_{p,q\in\Lambda}\big{\langle}{\overset{\frown}{g}}(p)[H_{q},Q_{p}]-{\overset{\frown}{g}}(q)[H_{p},Q_{q}]\big{\rangle}(f(q)-f(p))\\\
&\qquad=\frac{\epsilon}{4}\sum_{p,q\in\Lambda}\langle
J^{N}_{pq}\rangle({\overset{\frown}{g}}(p)+{\overset{\frown}{g}}(q))(f(q)-f(p))-U(\delta
f,\delta{\overset{\frown}{g}}),\end{split}$ (61)

The first term can be expressed in terms of magnetization as in (47). We write
${\overset{\frown}{g}}(p)$ as a difference
${\overset{\frown}{g}}(p)=g_{II}(p)-g_{I}(p)$, where $g_{I,II}$ are translates
of a smeared step-function $g(p)$. Then we rewrite the response as a
difference of conductivities of different materials:

$\displaystyle\begin{split}\delta\langle J^{E}(\delta
f)\rangle&\approx\epsilon(\lambda_{II}-\lambda_{I})\Bigg{\\{}\frac{\partial}{\partial\lambda}\Big{[}\beta\int_{0}^{\infty}dte^{-st}\langle\langle
J^{E}(\delta f,t);J^{N}(\delta g)\rangle\rangle-U(\delta f,\delta
g)\Big{]}+\mu^{N}_{\lambda}(\delta g\cup\delta f)\Bigg{\\}}.\end{split}$ (62)

On the other hand, from the continuum phenomenological equation (6) one finds
the current across the line $x=a$ to be

$-\epsilon\int_{-\infty}^{\infty}\eta_{xy}\partial_{y}gdy-\epsilon\mu\int_{-\infty}^{\infty}\sigma_{xy}\partial_{y}gdy=\epsilon\int_{\lambda_{I}}^{\lambda_{II}}\frac{\partial\eta_{xy}}{\partial\lambda}d\lambda+\epsilon\mu\int_{\lambda_{I}}^{\lambda_{II}}\frac{\partial\sigma_{xy}}{\partial\lambda}d\lambda,$
(63)

where the second term originates from the contribution of $\varphi{\bf j}^{N}$
to ${\bf j}^{E}$. Comparing it to (62) we find

$\displaystyle\begin{split}d\eta_{xy}&=d\Big{[}\beta\lim_{s\rightarrow
0}\int_{0}^{\infty}dte^{-st}\langle\langle J^{\mathcal{Q}}(\delta
f,t);J^{N}(\delta g)\rangle\rangle-U(\delta f,\delta
g)\Big{]}-\beta\mu^{N}(\delta f\cup\delta g).\end{split}$ (64)

This 1-form can be extended to include temperature as a parameter using the
scaling relation

$\left(T\frac{\partial}{\partial
T}+\left.\lambda_{0}\frac{\partial}{\partial\lambda_{0}}\right|_{\lambda_{0}=1}\right)\frac{\eta_{xy}(\lambda_{0},T)}{T}=0.$
(65)

Therefore we can define the $T$-component of the 1-form on the enlarged
parameter space as follows:

$\frac{d}{dT}\frac{\eta_{xy}}{T}=\frac{\partial}{\partial
T}\left[\beta^{2}\int_{0}^{\infty}dte^{-st}\langle\langle
J^{\mathcal{Q}}(\delta f,t);J^{N}(\delta g)\rangle\rangle-\beta U(\delta
f,\delta g)\right]-\beta^{2}\tau^{N}(\delta f\cup\delta g),$ (66)

where $\tau^{N}$ is given by (58).

###  III.5 Symmetric parts of transport coefficients

Note that the 1-forms $\mu^{N}(\delta f\cup\delta g)$ and $\tau^{N}(\delta
f\cup\delta g)$ are formally skew-symmetric under the exchange of $f$ and $g$.
To make use of this symmetry, we need to argue that $f$ can be replaced with a
smeared step-function in the $x$ direction. This can be argued using matching
between microscopic theory and hydrodynamics. Namely, we expect that the
microscopic linear response can be used to compute the properties of a Non-
Equilibrium Steady State (NESS). Replacing $f$ with a smeared step-function
changes the operators $J^{N}(\delta f)$ and $J^{E}(\delta f)$ by
$J^{N}(\delta{\overset{\frown}{f}})$ and $J^{E}(\delta{\overset{\frown}{f}})$.
The latter operators can be equivalently written as
$\frac{dQ({\overset{\frown}{f}})}{dt}$ and
$\frac{dH({\overset{\frown}{f}})}{dt}$, where

$Q({\overset{\frown}{f}})=\sum_{p\in\Lambda}{\overset{\frown}{f}}(p)Q_{p},\quad
H({\overset{\frown}{f}})=\sum_{p\in\Lambda}{\overset{\frown}{f}}(p)H_{p},$
(67)

and ${\overset{\frown}{f}}(p)$ is a hat-shaped function which depends only on
$x(p)$. On the other hand, if the microscopic linear response is to reproduce
the expected properties of a NESS, the expectation value of these observables
must be zero. Thus we can replace $f$ with a smeared step-function in the
$x$-direction without affecting $\nu_{xy}$ or $\eta_{xy}$.

After this has been done, exchanging $f$ and $g$ is equivalent to exchanging
$x$ and $y$. Thus $\mu^{N}(\delta f\cup\delta g)$ does not enter the
microscopic formulas for the symmetrized thermoelectric coefficients. These
formulas then can be integrated, giving

$\displaystyle\nu^{S}_{xy}$
$\displaystyle=\frac{\beta^{2}}{2}\lim_{s\rightarrow
0}\int_{0}^{\infty}dte^{-st}\left[\langle\langle J^{N}(\delta
f,t);J^{\mathcal{Q}}(\delta g)\rangle\rangle+\langle\langle J^{N}(\delta
g,t);J^{\mathcal{Q}}(\delta f)\rangle\rangle\right]+\beta U(\delta f,\delta
g),$ (68) $\displaystyle\eta^{S}_{xy}$
$\displaystyle=\frac{\beta}{2}\lim_{s\rightarrow
0}\int_{0}^{\infty}dte^{-st}\left[\langle\langle J^{\mathcal{Q}}(\delta
f,t);J^{N}(\delta g)\rangle\rangle+\langle\langle J^{\mathcal{Q}}(\delta
g,t);J^{N}(\delta f)\rangle\rangle\right]-U(\delta f,\delta g).$ (69)

As we show in Appendix C, the two terms on the right-hand side of these
equation are not separately invariant under Hamiltonian density redefinition,
but the full transport coefficients are invariant. The correction term
$U(\delta f,\delta g)$ is zero for fermionic systems with only density-
dependent interactions (see Appendix D for more details).

###  III.6 Skew-symmetric parts of transport coefficients

In a similar way one can find skew-symmetric parts of transport coefficients

$\displaystyle d\nu^{A}$
$\displaystyle=d\Big{[}\frac{\beta^{2}}{2}\lim_{s\rightarrow
0}\int_{0}^{\infty}dte^{-st}\left(\langle\langle J^{N}(\delta
f,t);J^{\mathcal{Q}}(\delta g)-\langle\langle J^{N}(\delta
g,t);J^{\mathcal{Q}}(\delta f\rangle\rangle\right)\Big{]}-\beta\mu^{N}(\delta
f\cup\delta g),$ (70) $\displaystyle d\eta^{A}$
$\displaystyle=d\Big{[}\frac{\beta}{2}\lim_{s\rightarrow
0}\int_{0}^{\infty}dte^{-st}\left(\langle\langle J^{\mathcal{Q}}(\delta
f,t);J^{N}(\delta g)\rangle\rangle-\langle\langle J^{\mathcal{Q}}(\delta
g,t);J^{N}(\delta f)\rangle\rangle\right)\Big{]}-\beta\mu^{N}(\delta
f\cup\delta g).$ (71)

These formulas give only derivatives of transport coefficients with respect to
parameters. Integration of these formulas over parameters or temperature gives
the difference of relative transport coefficients at different values of
parameters. It is natural to define relative transport coefficients of a
trivial insulator to be zero. Determination of the relative transport
coefficient in this case would correspond to integration over a path in the
parameter space from a trivial insulator to the material of interest.

There are many paths which one can use to deform a system into a trivial one.
Consistency requires the integral of $d\nu^{A}$ or $d\eta^{A}$ to depend only
on the endpoints of the path. This means that the 1-form $\mu^{N}(\delta
f\cup\delta g)$ must be exact. This can be proved using the techniques of
Kapustin and Spodyneiko (2020) where $\mu^{E}(\delta f\cup\delta g)$ was shown
to be exact.

##  IV Discussion

In this paper we have derived microscopic formulas for “transverse”
thermoelectric coefficients of general 2d lattice systems. It was convenient
to decompose them into symmetric and anti-symmetric parts, since they have
qualitatively different behavior: the former are absolute transport
coefficients, while the latter are relative. Similar formulas for electric
Hall conductivity and thermal Hall conductivity have already been derived in
Kapustin and Spodyneiko (2020).

The usual Kubo formulas for transport coefficients require averaging the
correlators of currents over the whole space. In contrast, our microscopic
formulas involve net currents through two perpendicular lines. This is because
a current on a 2d lattice is a function of a pair of points, and the natural
observable associated to it is localized on a line rather than at a point.
Despite this, after the limit $s\rightarrow 0$ has been taken, our formula
computes the same quantity as the usual continuum Kubo formula.

It is natural to ask whether longitudinal components of conductivity, thermal
conductivity, and thermoelectric tensors of a 2d lattice system can be
computed in a similar manner. This is easily achieved: one simply replaces two
perpendicular lines with two lines making a nonzero angle $\theta$. It is easy
to see determine from hydrodynamics which linear combination of components of
the transport tensors describes the corresponding linear response. For
example, if both currents involved are electric currents, one of the lines is
given by $y=0$, and the other one is $y=x\cdot\tan\theta$, then the correlator

$\beta\lim_{s\rightarrow 0+}\int_{0}^{\infty}e^{-st}\langle\langle
J^{N}(\delta f,t);J^{N}(\delta g)\rangle\rangle dt$ (72)

measures $\sigma_{xy}-\frac{1}{\tan\theta}\sigma_{yy}$. Thus by changing the
functions $f,g$ one can extract all four components of the conductivity
tensor. The same is true about other transport coefficients.

In this paper we focused on the case of 2d materials, but the 3d case can be
accommodated as well. One can simply replace a lattice in ${\mathbb{R}}^{2}$
with a lattice in ${\mathbb{R}}^{2}\times[0,L]$, impose periodic boundary
conditions in the third direction, divide all formulas by $L$, and take the
limit $L\rightarrow\infty$. The functions $f,g$ remain independent of the
third coordinate. It should not matter whether the limit $L\rightarrow\infty$
is taken before or after the limit $s\rightarrow 0$, since the problem is
translationally invariant in the third direction.

In the 2d case, the quantity $\nu^{A}(T)$ (normalized relative to the vacuum)
is dimensionless, and Bloch’s theorem implies that for gapped systems
$\nu^{A}(0)$ does not change under the variations of the Hamiltonian which do
not close the gap. Thus if $\nu^{A}(0)$ were nonzero, it would represent a new
topological invariant of gapped 2d phases of matter. However, one can show
that on very general grounds $\nu^{A}(0)$ vanishes for all gapped systems
Levin  _et al._ . By Onsager reciprocity, the $T\rightarrow 0$ limit of
$\eta^{A}(T)/T$ also vanishes. Thus topological invariants of gapped 2d
systems arise only from the Hall conductivity and the thermal Hall
conductivity.

##  Appendix A Kubo canonical pairing

Kubo canonical pairing of two operators $A,B$ is defined as follows:

$\langle\langle A;B\rangle\rangle=\frac{1}{\beta}\int_{0}^{\beta}\langle
A(-i\tau)B\rangle d\tau-\langle A\rangle\langle B\rangle.$ (73)

Here $\langle\ldots\rangle$ denotes average over a Gibbs state at temperature
$T=1/\beta$ (or more generally, over a state satisfying the Kubo-Martin-
Schwinger condition), and $A(-i\tau)=e^{H\tau}Ae^{-H\tau}$. Kubo paring
determines static linear response: if the Hamiltonian is perturbed by $\lambda
B$, where $\lambda$ is infinitesimal, then the change in the expectation value
of $A$ due to the change in the equilibrium density matrix is

$\Delta\langle A\rangle=-\beta\lambda\langle\langle A;B\rangle\rangle.$ (74)

Kubo pairing is symmetric, $\langle\langle A;B\rangle\rangle=\langle\langle
B;A\rangle\rangle$, and satisfies

$\beta\langle\langle i[H,A];B\rangle\rangle=\langle i[B,A]\rangle.$ (75)

In finite volume, one can write it in terms of the energy eigenstates as
follows:

$\langle\langle A;B\rangle\rangle=Z^{-1}\sum_{n,m}\langle
n|\bar{A}|m\rangle\langle m|\bar{B}|n\rangle\frac{e^{-\beta E_{m}}-e^{-\beta
E_{n}}}{\beta(E_{n}-E_{m})},$ (76)

where $\bar{A}=A-\langle A\rangle,$ and $\bar{B}=B-\langle B\rangle.$

##  Appendix B Onsager reciprocity revisited

Derivations of Onsager relations are based on the analysis of hydrodynamic
fluctuations, so it might seem that they should put constraints only on those
transport coefficients which enter the hydrodynamic equations of motion. On
closer inspection, one finds Casimir (1945) that the derivation involves net
currents which measure the rate of change of conserved quantities in a finite
volume and thus require understanding boundary contributions. As a result,
Onsager reciprocity constrains both absolute transport coefficients and
relative transport coefficients defined relative to the vacuum. Equivalently,
it imposes conditions on the derivatives of relative transport coefficients
with respect to parameters. To illustrate how this works, let us discuss the
constraints imposed by Onsager reciprocity on relative transport coefficients
of time-reversal-invariant 2d systems. For the skew-symmetric thermal
conductivity $\kappa^{A}$, time-reversal-invariance implies

$\frac{\partial}{\partial\lambda}\kappa^{A}=0,$ (77)

where $\lambda$ is a parameter of the Hamiltonian. Thus $\kappa^{A}$ can be a
function of temperature only. Further, if we treat $T$ as a parameter, then
scaling analysis gives

$\frac{\partial}{\partial T}\frac{\kappa^{A}}{T}=0.$ (78)

Hence $\kappa^{A}(T)=aT$, where $a$ does not depend on parameters. The
parameter $a$ has no physical significance, but it is natural to set it to
zero, so that the vacuum has zero thermal Hall conductivity. Thus we reach the
standard conclusion that for a system with time-reversal invariance
$\kappa^{A}=0$.

The case of thermoelectric coefficients is slightly different. Usually one
says that Onsager reciprocity requires $\nu_{km}=T^{-1}\eta_{mk}$, which
implies $\nu^{A}+T^{-1}\eta^{A}=0$ Landau and Lifshitz (1984). Since both
$\nu^{A}$ and $\eta^{A}$ are relative transport coefficients, one should
interpret this as

$\frac{\partial}{\partial\lambda}\left(\nu^{A}+T^{-1}\eta^{A}\right)=0.$ (79)

Hence $\nu^{A}+T^{-1}\eta^{A}$ can depend only on the temperature. If we treat
temperature as a parameter, then the scaling analysis gives

$\frac{\partial}{\partial T}\left(\nu^{A}+T^{-1}\eta^{A}\right)=0.$ (80)

Hence $\nu^{A}+T^{-1}\eta^{A}=a$, where $a$ is a constant which is independent
of any parameters or temperature and has no physical significance. One can
choose it to be zero. Then $\nu^{A}=-T^{-1}\eta^{A}$. So for a time-reversal-
invariant 2d system there is only one independent skew-symmetric
thermoelectric transport coefficient, namely $\nu^{A}$.

##  Appendix C Invariance under Hamiltonian density redefinition

For a given Hamiltonian, there are many ways to define a Hamiltonian density.
A typical example of this is the ambiguity in splitting an interaction term
between two sites $p$ and $q$ into $H_{p}$ and/or $H_{q}$. In this appendix,
we will show that our microscopic formulas for physically observable transport
coefficients are independent of the choice of the Hamiltonian density, even
though individual terms in the microscopic formulas are not invariant. For
some systems this can be used to simplify the microscopic formulas.

###  C.1 Invariance of the electric current

Consider the following change of the Hamiltonian density

$\displaystyle H_{p}\rightarrow H_{p}+\sum_{r\in\Lambda}A_{rp},$ (81)

where $A_{rp}$ is skew-symmetric in $r,p$. We want the final Hamiltonian to be
$U(1)$-invariant. Therefore, we have to impose

$[Q,\sum_{r\in\Lambda}A_{rp}]=0.$ (82)

For a general choice of $A_{pq}$ a stronger condition

$[Q,A_{pq}]=0,$ (83)

will not hold. However, one can always redefine $A_{pq}$ (by subtracting the
$U(1)$-non-invariant part) in such a way that (83) holds without affecting
$H_{p}$. In the following we will assume this was done and (83) is true.

Under the transformation (81) the electric current changes as

$\displaystyle J^{N}_{pq}\rightarrow
J^{N}_{pq}+i\sum_{r\in\Lambda}\left([A_{rq},Q_{p}]-[A_{rp},Q_{q}]\right).$
(84)

Even though the current density changes, the net current through any section
is invariant. Indeed,

$J^{N}(\delta f)\rightarrow J^{N}(\delta
f)+\frac{i}{2}\sum_{p,q,r\in\Lambda}\left([A_{rq},Q_{p}]-[A_{rp},Q_{q}]\right)(f(q)-f(p)),$
(85)

and the last term is zero since

$\sum_{p,q,r\in\Lambda}\left([A_{rq},Q_{p}]-[A_{rp},Q_{q}]\right)(f(q)-f(p))\\\
=\sum_{p,q,r\in\Lambda}\left([A_{rq},Q_{p}]+[A_{pr},Q_{q}]+[A_{qp},Q_{r}]\right)(f(q)-f(p))\\\
=\frac{1}{3}\sum_{p,q,r\in\Lambda}\left([A_{rq},Q_{p}]+[A_{pr},Q_{q}]+[A_{qp},Q_{r}]\right)(f(q)-f(p)+f(p)-f(r)+f(r)-f(q))=0,$
(86)

where we have used (83) and the skew-symmetry of
$[A_{rq},Q_{p}]+[A_{pr},Q_{q}]+[A_{qp},Q_{r}]$.

###  C.2 Covariance of the energy current

Let us now consider the effect of the redefinition of the Hamiltonian density
on the energy current. Imposing an energy analog of (82) or (83)

$[H,\sum_{r\in\Lambda}A_{rp}]\stackrel{{\scriptstyle?}}{{=}}0,\qquad\textrm{or}\qquad[H,A_{pq}]\stackrel{{\scriptstyle?}}{{=}}0,$
(87)

is far too restrictive, since it would only allow changes of the Hamiltoniain
density by conserved quantities. For example, the difference between putting
the interaction term between the two sites $p$ and $q$ either into $H_{p}$ or
into $H_{q}$ corresponds to $A_{pq}$ equal to the interaction term. Obviously,
interaction terms are not integrals of motion in general. Because of this we
will not impose either of the equations in (87).

Under the redefinition of the Hamiltonian density (81) the energy current
changes as

$\displaystyle J^{E}_{pq}\rightarrow
J^{E}_{pq}+i\sum_{r\in\Lambda}\left([A_{rq},H_{p}]+[H_{q},A_{rp}]\right),$
(88)

while the net current transforms as

$\displaystyle J^{E}(\delta f)\rightarrow J^{E}(\delta
f)+\frac{i}{2}\sum_{p,q,r\in\Lambda}\left([A_{rq},H_{p}]+[H_{q},A_{rp}]\right)(f(q)-f(p)).$
(89)

The last term can be rewritten as

$\frac{i}{2}\sum_{p,q,r\in\Lambda}\left([A_{rq},H_{p}]+[H_{q},A_{rp}]\right)(f(q)-f(p))\\\
=\frac{i}{2}\sum_{p,q,r\in\Lambda}\left([A_{rq},H_{p}]+[A_{pr},H_{q}]+[A_{qp},H_{r}]\right)(f(q)-f(p))-\frac{i}{2}\sum_{p,q,r\in\Lambda}[A_{qp},H_{r}](f(q)-f(p))\\\
=\frac{i}{6}\sum_{p,q,r\in\Lambda}\left([A_{rq},H_{p}]+[A_{pr},H_{q}]+[A_{qp},H_{r}]\right)(f(q)-f(p)+f(p)-f(r)+f(r)-f(q))\\\
-\frac{i}{2}\sum_{p,q\in\Lambda}[H,A_{pq}](f(q)-f(p))=-\dot{A}(\delta f),$

where we have defined

$\displaystyle A(\delta f)=\frac{1}{2}\sum_{p,q\in\Lambda}A_{pq}(f(q)-f(p)).$
(90)

We find that the net energy current transforms as follows under a redefinition
of the Hamiltonian density:

$\displaystyle J^{E}(\delta f)\rightarrow J^{E}(\delta f)-\dot{A}(\delta f).$
(91)

But this should be expected since a redefinition of the energy density changes
how we define the energy of sub-regions and therefore should affect the net
energy current. Indeed, one can see that (89) is exactly the transformation
needed in order to satisfy the energy conservation law

$\displaystyle\dot{H}_{p}=-\sum_{q\in\Lambda}J^{E}_{pq}\quad\rightarrow\quad\dot{H}_{p}+\sum_{q\in\Lambda}\dot{A}_{qp}=-\sum_{q\in\Lambda}J^{E}_{pq}+\sum_{q\in\Lambda}\dot{A}_{qp}$
(92)

for the new energy density $H_{p}+\sum_{q\in\Lambda}A_{qp}$. By summing this
transformation law over $p$ weighted by a function $f(p)$ with a compact
support we find that

$\displaystyle\dot{H}(f)=-J^{E}(\delta
f)\quad\rightarrow\quad\dot{H}_{p}+\dot{A}(\delta f)=-J^{E}(\delta
f)+\dot{A}(\delta f),$ (93)

which reproduces (91). Here we used an identity

$\displaystyle\sum_{p,q\in\Lambda}A_{pq}f(q)=\frac{1}{2}\sum_{p,q\in\Lambda}A_{pq}(f(q)-f(p))=A(\delta
f)$ (94)

which is true for any $f$ with a compact support.

From the above discussion, one can see that energy current is not invariant
but covariant under energy density redefinitions. If we choose $f(p)$ to be 1
when $p$ is in some compact set $B$ and zero otherwise, the physical meaning
of (91) is very clear. It corresponds to ambiguities in the energy currents
due to interaction terms along the boundary of $B$. Depending on how we
distribute the interaction terms among $H_{p}$ we can change the energy stored
in the region $B$ as well as energy current through its boundary.

###  C.3 Invariance of the microscopic formulas for thermoelectic coefficients

In this section we will show that the coefficients $\nu_{xy}$ and $\eta_{xy}$
are invariant under a redefinition of the Hamiltonian density. We will start
with skew-symmetric coefficients

$\displaystyle d\nu^{A}$ $\displaystyle=\frac{1}{2}d\Big{(}\nu^{\rm
Kubo}(\delta f,\delta g)-\nu^{\rm Kubo}(\delta g,\delta
f)\Big{)}-\beta^{2}\mu^{N}(\delta f\cup\delta g),$ (95) $\displaystyle
d\eta^{A}$ $\displaystyle=\frac{1}{2}d\Big{(}\eta^{\rm Kubo}(\delta f,\delta
g)-\eta^{\rm Kubo}(\delta g,\delta f)\Big{)}-\beta\mu^{N}(\delta f\cup\delta
g).$ (96)

Here we defined the Kubo parts as

$\displaystyle\nu^{\rm Kubo}(\delta f,\delta g)$
$\displaystyle=\beta^{2}\lim_{s\rightarrow
0}\int_{0}^{\infty}dte^{-st}\langle\langle J^{N}(\delta
f,t);J^{\mathcal{Q}}(\delta g)\rangle\rangle,$ (97) $\displaystyle\eta^{\rm
Kubo}(\delta f,\delta g)$ $\displaystyle=\beta\lim_{s\rightarrow
0}\int_{0}^{\infty}dte^{-st}\langle\langle J^{\mathcal{Q}}(\delta
f,t);J^{N}(\delta g)\rangle\rangle.$ (98)

Under Hamiltonian density redefinition the Kubo parts transform as

$\displaystyle\begin{split}\nu^{\rm Kubo}(\delta f,\delta
g)&\rightarrow\nu^{\rm Kubo}(\delta f,\delta g)-\beta^{2}\lim_{s\rightarrow
0}\int_{0}^{\infty}dte^{-st}\langle\langle J^{N}(\delta f,t);\dot{A}(\delta
g)\rangle\rangle\\\ &=\nu^{\rm Kubo}(\delta f,\delta
g)-\beta^{2}\langle\langle J^{N}(\delta f);A(\delta
g)\rangle\rangle,\end{split}$ (99) $\displaystyle\begin{split}\eta^{\rm
Kubo}(\delta f,\delta g)&\rightarrow\eta^{\rm Kubo}(\delta f,\delta
g)-\beta\lim_{s\rightarrow
0}\int_{0}^{\infty}dte^{-st}\langle\langle\dot{A}(\delta f,t);J^{N}(\delta
g);\rangle\rangle\\\ &=\nu^{\rm Kubo}(\delta f,\delta g)+\beta\langle\langle
A(\delta f);J^{N}(\delta g)\rangle\rangle,\end{split}$ (100)

where we used properties of the Kubo pairing.

Before finding the variation of the magnetization term it is useful to rewrite
it slightly:

$\displaystyle\begin{split}\mu^{N}(\delta f\cup\delta
g)=\frac{1}{2}\sum_{p,q\in\Lambda}\left[\frac{1}{3}\sum_{r\in\Lambda}\mu_{pqr}(g_{p}+g_{q}+g_{r})-\frac{1}{2}\sum_{r\in\Lambda}\mu_{pqr}(g(p)+g(q))\right](f(q)-f(p))\\\
=\frac{1}{2}\sum_{p,q\in\Lambda}\left[\frac{1}{3}\sum_{r\in\Lambda}\mu_{pqr}(g_{p}+g_{q}+g_{r})-\frac{1}{2}d\langle
J^{N}_{pq}\rangle(g(p)+g(q))\right](f(q)-f(p)).\end{split}$ (101)

Note that one cannot expand the square brackets, since the two resulting sums
over $p,q$ will not converge separately.

Let us find the variation of $\dfrac{1}{2}\langle
J^{N}_{pq}\rangle(g(p)+g(q))$ under a Hamiltonian density redefinition. It
reads

$\displaystyle\begin{split}\frac{1}{2}\langle
J^{N}_{pq}\rangle(g(p)+g(q))&\rightarrow\frac{1}{2}\langle
J^{N}_{pq}\rangle(g(p)+g(q))\\\
&+\frac{i}{2}\sum_{r\in\Lambda}\langle[A_{rq},Q_{p}]-[A_{rp},Q_{q}]\rangle(g(p)+g(q))\end{split}$
(102)

The last term can be rewritten as follows:

$\displaystyle\begin{split}\frac{i}{2}\sum_{r\in\Lambda}\langle[A_{rq},Q_{p}]&-[A_{rp},Q_{q}]\rangle(g(p)+g(q))\\\
&=\frac{\beta}{2}\sum_{r\in\Lambda}\langle\langle
g(p)\dot{Q}_{p};A_{rq}\rangle\rangle+\frac{\beta}{2}\sum_{r\in\Lambda}\langle\langle\dot{Q}_{p};g(q)A_{rq}\rangle\rangle-(p\leftrightarrow
q),\end{split}$ (103)

where we used the properties of the Kubo pairing. The first term in this
expression can be rewritten as

$\sum_{r\in\Lambda}\langle\langle
g(p)\dot{Q}_{p};A_{rq}\rangle\rangle-(p\leftrightarrow
q)=-\frac{1}{2}\sum_{s,r\in\Lambda}\langle\langle
J_{sp}^{N}(g(s)+g(p));A_{rq}\rangle\rangle-(p\leftrightarrow q)\\\
=-\frac{1}{2}\sum_{s,r\in\Lambda}\langle\langle
J^{N}_{sp}(g(s)+g(p))+J^{N}_{ps}(g(s)-g(p));A_{rq}\rangle\rangle-(p\leftrightarrow
q)=\langle\langle J^{N}(\delta g);A_{qp}\rangle\rangle\\\
-\frac{1}{2}\sum_{r,s\in\Lambda}\Bigg{[}\langle\langle
J^{N}_{rp}(g(r)+g(p));A_{sq}\rangle\rangle+\langle\langle
J^{N}_{ps}(g(s)-g(p));A_{rq}\rangle\rangle+\text{2 perms}\Bigg{]},$ (104)

where ”2 perms” means the two cyclic permutations in $p,q,r$. Note that the
term in square brackets is skew-symmetric in $p,q,r$. The second term can be
rewritten as

$\sum_{r\in\Lambda}\langle\langle\dot{Q}_{p};g(q)A_{rq}\rangle\rangle-(p\leftrightarrow
q)=-\sum_{s,r\in\Lambda}\langle\langle
J^{N}_{sp};g(q)A_{rq}\rangle\rangle-(p\leftrightarrow q)\\\
=-\frac{1}{2}\sum_{s,r\in\Lambda}\langle\langle
J^{N}_{sp};A_{rq}(g(q)+g(r))+A_{qr}(g(r)-g(q))\rangle\rangle-(p\leftrightarrow
q)=\langle\langle J^{N}_{pq};A(\delta g)\rangle\rangle\\\
-\frac{1}{2}\sum_{s,r\in\Lambda}\Big{[}\langle\langle
J^{N}_{sp};A_{rq}(g(q)+g(r))\rangle\rangle+\langle\langle
J^{N}_{rp};A_{qr}(g(r)-g(q))\rangle\rangle+\text{2 perms}\Big{]}.$ (105)

Note that term in square brackets is skew-symmetric in $p,q,r$

By combining equations (101-105) we find that the magnetization contribution
changes under a redefinition of the Hamiltonian density as follows:

$\displaystyle\begin{split}\mu^{N}(\delta f\cup\delta
g)\rightarrow\mu^{N}(\delta f\cup\delta g)-\frac{\beta}{2}d\langle\langle
J^{N}(\delta f);A(\delta g)\rangle\rangle&+\frac{\beta}{2}d\langle\langle
J^{N}(\delta g);A(\delta f)\rangle\rangle\\\
&+\frac{1}{2}\sum_{p,q,r\in\Lambda}C_{pqr}(f(q)-f(p)),\end{split}$ (106)

where $C_{pqr}$ is a skew-symmetric function of $p,q,r$ which is combination
of skew-symmetric parts (and their derivatives) in the right-hand sides of
Eqs. (101-105). Due to its skew-symmetry we find that

$\frac{1}{2}\sum_{p,q,r\in\Lambda}C_{pqr}(f(q)-f(p))=\frac{1}{6}\sum_{p,q,r\in\Lambda}C_{pqr}(f(q)-f(p)+f(p)-f(s)+f(s)-f(q))=0.$
(107)

We see that the variation of the magnetization exactly compensates the
variation of the Kubo parts. Thus the skew-symmetric parts of the
thermoelectric tensors are invariant under a redefinition of the Hamiltonian
density.

Now let us consider the symmetric parts

$\displaystyle\nu^{S}_{xy}$ $\displaystyle=\frac{1}{2}\Big{(}\nu^{\rm
Kubo}(\delta f,\delta g)+\nu^{\rm Kubo}(\delta g,\delta f)\Big{)}+\beta
U(\delta f,\delta g),$ (108) $\displaystyle\eta^{S}_{xy}$
$\displaystyle=\frac{1}{2}\Big{(}\eta^{\rm Kubo}(\delta f,\delta g)+\eta^{\rm
Kubo}(\delta g,\delta f)\Big{)}-U(\delta f,\delta g).$ (109)

The variation of Kubo parts were already determined before, so we focus on the
transformation of $U$. Under (81) it transforms as follows:

$\displaystyle U(\delta f,\delta g)\rightarrow U(\delta f,\delta
g)+\frac{i}{4}\sum_{p,q\in\Lambda}\langle[\partial A_{q},Q_{p}]+[\partial
A_{p},Q_{q}]\rangle(f(q)-f(p))(g(q)-g(p)).$ (110)

We can rewrite this equation by noticing that

$\displaystyle\frac{i}{2}\langle[\partial A_{q},Q_{p}]+[\partial
A_{p},Q_{q}]\rangle(g(q)-g(p))=-\frac{\beta}{2}\langle\langle
g(p)\dot{Q}_{p};\partial
A_{q}\rangle\rangle+\frac{\beta}{2}\langle\langle\dot{Q}_{p};\partial
g(q)A_{q}\rangle\rangle-(p\leftrightarrow q).$ (111)

Then using eqs. (104, 105) we find

$U(\delta f,\delta g)\rightarrow U(\delta f,\delta
g)+\frac{\beta}{2}\langle\langle J^{N}(\delta f);A(\delta
g)\rangle\rangle+\frac{\beta}{2}\langle\langle J^{N}(\delta g);A(\delta
f)\rangle\rangle$ (112)

We see that the variation of this term cancels the varitions of the Kubo
parts.

One can do the same checks for the thermal Hall conductivity and verify that
the microscopic formula derived in Kapustin and Spodyneiko (2020) is in
invariant under a redefinition of the Hamiltonian density. To linear order in
$A_{pq}$ all the manipulations are almost the same except for the replacement
$Q_{p}\rightarrow H_{p}$ and $J^{N}\rightarrow J^{E}$.

##  Appendix D Thermoelectric coefficients for free fermions

###  D.1 Definitions and correlation functions

In this appendix we will specialize our microscopic formulas for coefficients
$\nu$ and $\eta$ to free fermionic systems. The Hamiltonian is taken to be

$H=\sum_{p,q\in\Lambda}a^{\dagger}_{p}h(p,q)a_{q},$ (113)

where an infinite matrix $h(p,q)$ is Hermitian $h(p,q)^{*}=h(q,p)$, and
$a^{\dagger}_{p},a_{p}$ are fermionic creation-annihilation operators
satisfying the standard anti-commutation relations

$a_{p}a^{\dagger}_{q}+a^{\dagger}_{q}a_{p}=\delta_{p,q},\qquad
a_{p}a_{q}+a_{q}a_{p}=a_{p}^{\dagger}a^{\dagger}_{q}+a^{\dagger}_{q}a_{p}^{\dagger}=0.$
(114)

We define the Hamiltonian density on site $p$ to be

$H_{p}=\frac{1}{2}\sum_{m\in\Lambda}\left(a^{\dagger}_{p}h(p,m)a_{m}+a^{\dagger}_{m}h(m,p)a_{p}\right).$
(115)

The charge operator on site $p$ is defined as

$Q_{p}=a^{\dagger}_{p}a_{p}.$ (116)

The electric current can be found from the conservation equation:

$J_{pq}^{N}=i(a^{\dagger}_{q}h(q,p)a_{p}-a^{\dagger}_{p}h(p,q)a_{q}).$ (117)

The net current through a section defined by $\delta f(p,q)=f(q)-f(p)$ is

$J(\delta f)=-ia^{\dagger}[h,f]a,$

where a bounded function $f\in\ell^{2}(\Lambda)$ is understood as an operator
acting on the one-particle Hilbert space $\ell^{2}(\Lambda)$ by
multiplication. Summation over sites is implicit.

The energy current operator is

$J^{E}_{pq}=\frac{-i}{4}\sum_{m\in\Lambda}\left(a^{\dagger}_{p}h(p,q)h(q,m)a_{m}-a^{\dagger}_{q}h(q,p)h(p,m)a_{m}\right.\\\
\left.-a^{\dagger}_{m}h(m,q)h(q,p)a_{q}+a^{\dagger}_{m}h(m,p)h(p,q)a_{q}\right.\\\
\left.+a^{\dagger}_{p}h(p,m)h(m,q)a_{q}-a^{\dagger}_{q}h(q,m)h(m,p)a_{p}\right).$
(118)

The net energy current is

$J^{E}(\delta f)=-\frac{i}{2}a^{\dagger}[h^{2},f]a.$

The state of the system at a temperature $T=1/\beta$ is defined via Wick’s
theorem and Gibbs distribution

$\displaystyle\langle a_{p}(t)a^{\dagger}_{q}\rangle$ $\displaystyle=$
$\displaystyle\left\langle p\left|\frac{e^{-iht}}{1+e^{-\beta
h}}\right|q\right\rangle,$ (119) $\displaystyle\langle
a_{p}(t)^{\dagger}a_{q}\rangle$ $\displaystyle=$ $\displaystyle\left\langle
q\left|\frac{e^{iht}}{1+e^{\beta h}}\right|p\right\rangle,$ (120)

where $a_{p}(t)$ are operators in the Heisenberg picture.

Using these formulas we find

$\langle J^{N}(\delta f,t)J^{N}(\delta g)\rangle=-{\rm
Tr}\left([h,f]\frac{e^{-iht}}{1+e^{-\beta h}}[h,g]\frac{e^{iht}}{1+e^{\beta
h}}\right),$

where the trace is over the 1-particle Hilbert space $\ell^{2}(\Lambda)$, and
the functions $f:\Lambda\rightarrow{\mathbb{R}}$ and
$g:\Lambda\rightarrow{\mathbb{R}}$ are operators on this Hilbert space. The
operators $[h,f]$ and $[h,g]$ have support on a vertical strip and a
horizontal strip, respectively.

Switching to the energy basis, substituting $t\rightarrow t-i\tau$, and
integrating from $0$ to $\beta$ over $\tau$ we find

$\langle\langle J^{N}(\delta f,t);J^{N}(\delta
g)\rangle\rangle=\frac{-1}{\beta}\sum_{n,m}\langle n|[h,f]|m\rangle\langle
m|[h,g]|n\rangle
e^{i(\varepsilon_{n}-\varepsilon_{m})t}\frac{e^{\beta\varepsilon_{n}}-e^{\beta\varepsilon_{m}}}{(1+e^{\beta\varepsilon_{n}})(1+e^{\beta\varepsilon_{m}})(\varepsilon_{n}-\varepsilon_{m})},$

where $\varepsilon_{n}$ are 1-particle Hamiltonian energy eigenvalues.

Multiplying this by $e^{-st}$ and integrating over $t$, we arrive at

$\sigma_{xy}=i\lim_{s\rightarrow 0}\sum_{n,m}\frac{\langle
n|[h,f]|m\rangle\langle
m|[h,g]|n\rangle}{\varepsilon_{n}-\varepsilon_{m}+is}\frac{\mathfrak{f}(\varepsilon_{n})-\mathfrak{f}(\varepsilon_{m})}{\varepsilon_{n}-\varepsilon_{m}},$

where $\mathfrak{f}(\varepsilon)=\frac{1}{1+e^{\beta(\varepsilon)}}$ is the
Fermi-Dirac distribution. We absorb the chemical potential into a shift of the
Hamiltonian.

The above expressions assume a discrete energy spectrum and thus can only be
used for finite-volume systems. To get an expression applicable to infinite-
volume systems, let us rewrite it in terms of the one-particle Green’s
functions $G_{\pm}(z)=1/(z-h\pm i0)$. Some of the useful formulas are

$\langle a^{\dagger}Aa\rangle=-\frac{1}{2\pi
i}\int_{-\infty}^{\infty}dz\,\mathfrak{f}(z){\rm
Tr}\Big{(}\big{[}G_{+}-G_{-}\big{]}A\Big{)},$
$\displaystyle\begin{split}-\beta\langle\langle
a^{\dagger}Aa;a^{\dagger}Ba\rangle\rangle=-\frac{1}{2\pi
i}\int_{-\infty}^{\infty}dz\,\mathfrak{f}(z){\rm
Tr}\Big{(}\big{[}G_{+}-G_{-}\big{]}AG_{+}B+G_{-}A\big{[}G_{+}-G_{-}\big{]}B\Big{)}\\\
=-\frac{1}{2\pi i}\int_{-\infty}^{\infty}dz\,\mathfrak{f}(z){\rm
Tr}\Big{(}G_{+}AG_{+}B-G_{-}AG_{-}B\Big{)},\end{split}$ (121)

where we have dropped $z$ for $G_{\pm}(z)$. Here $A$ and $B$ are operators
acting on the one-particle Hilbert space, and in the second formula we assumed
in addition that their average is zero: $\langle a^{\dagger}Aa\rangle=\langle
a^{\dagger}Ba\rangle=0$. Note also that

$hG_{\pm}=G_{\pm}h=zG_{\pm}-1,\quad[G_{\pm},A]=G_{\pm}[h,A]G_{\pm}.$

Using this notation the formula for the electric conductivity takes the form

$\sigma_{xy}=-\frac{1}{2\pi}\int_{-\infty}^{\infty}dz\,\mathfrak{f}(z){\rm
Tr}\big{\\{}[h,f]G_{+}^{2}[h,g](G_{+}-G_{-})-[h,f](G_{+}-G_{-})[h,g]G_{-}^{2}\big{\\}},$
(122)

where the integration is over the real axis in the $z$-plane.

###  D.2 Magnetization term

The magnetization differential for an arbitrary deformation $dh$ of the
1-particle Hamiltonian is given by

$\mu^{N}(\delta f\cup\delta
g)=\frac{1}{4\pi}\int_{-\infty}^{\infty}dz\,\mathfrak{f}(z){\rm
Tr}\Big{(}G_{+}dhG_{+}\Big{\\{}[h,f]G_{+}[h,g]-[h,g]G_{+}[h,f]\Big{\\}}\Big{)}-(G_{+}\rightarrow
G_{-}).$ (123)

For temperature variations this expression can be simplified to

$\tau^{N}(\delta f\cup\delta g)=\frac{1}{4\pi}\int_{-\infty}^{\infty}dz{\rm
Tr}\Big{(}\mathfrak{f}(z)(G_{+}-G_{-})[h^{2},f]G_{-}^{2}[h,g]-\mathfrak{f}(z)(G_{+}-G_{-})[h^{2},g]G_{+}^{2}[h,f]\\\
+\mathfrak{f}^{\prime}(z)(G_{+}-G_{-})h[h,g]G_{+}[h,f]+\mathfrak{f}^{\prime}(z)(G_{+}-G_{-})h[h,g]G_{-}[h,f]\Big{)}-(f\leftrightarrow
g).$ (124)

These expressions are needed only for the evaluation of skew-symmetric parts
of the thermoelectric coefficients.

###  D.3 $U$-term

Let us study the term (51) for free fermionic system. In this case the
relevant many-body operators become

$\left([Q_{p},H_{q}]+[Q_{q},H_{p}]\right)=a^{\dagger}[h,\delta_{p}\delta_{q}]a,$
(125)

where $\delta_{p}$ is a Kronecker delta function equal 1 on site $p$ and 0 on
all other sites. A product of two delta functions enforces $q=p$ in the
summation over $p$ and $q$. Since $U$ also involves a factor of
$(g(p)-g(q))(f(p)-f(q))$, $U(\delta f,\delta g)$ vanishes for systems of free
fermions.

More generally, one can consider a system of fermions with only density-
dependent interactions. Namely, suppose we allow the following interaction
term in the Hamitonian (113):

$H^{\rm
int}=\sum_{p_{1},\dots,p_{n}\in\Lambda}V(p_{1},\dots,p_{n})Q_{p_{1}}\dots
Q_{p_{n}},$ (126)

where $V(p_{1},\dots,p_{n})$ is a function of $n$ sites which describes the
potential energy of many-body interaction and decays rapidly when the points
$p_{1},\ldots,p_{n}$ are far from each other. One can see that this term will
leave eq. (125) unaffected since $Q_{p}$ commute with each other. We conclude
that for fermionic system with only density-dependent interactions there is no
correction originating from $U$ to the symmetric thermoelectric coefficients
provided $H_{p}$ is chosen in the manner explained above.

###  D.4 Skew-symmetric part

Consider the variation of the Kubo parts (97,98) of the skew-symmetric
thermoelectic coefficients under a rescaling of the Hamiltonian:
$dh=h\,d\lambda_{0}$. We get

$d\nu^{A}_{\rm Kubo}=\beta d\eta^{A}_{\rm Kubo}(\delta f\cup\delta
g)=\frac{d\lambda_{0}}{4\pi}\int_{-\infty}^{\infty}dz{\rm
Tr}\Big{(}\mathfrak{f}(z)(G_{+}-G_{-})[h^{2},f]G_{-}^{2}[h,g]\\\
-\mathfrak{f}(z)(G_{+}-G_{-})[h^{2},g]G_{+}^{2}[h,f]+\mathfrak{f}^{\prime}(z)(G_{+}-G_{-})h[h,g]G_{+}[h,f]+\mathfrak{f}^{\prime}(z)(G_{+}-G_{-})h[h,g]G_{-}[h,f]\\\
-2\mathfrak{f}^{\prime}(z)h^{2}(G_{+}-G_{-})[h,f]G_{+}^{2}[h,g]+2\mathfrak{f}^{\prime}(z)h^{2}(G_{+}-G_{-})[h,g]G_{-}^{2}[h,f]\Big{)}-(f\leftrightarrow
g).$ (127)

Summing up this contributions with the magnetization contribution gives

$\frac{d\nu^{A}}{dT}=\frac{d}{dT}\left(\frac{\eta^{A}}{T}\right)=\frac{1}{2\pi
T^{2}}\int_{-\infty}^{\infty}dz\mathfrak{f}^{\prime}(z)z^{2}{\rm
Tr}\Big{(}(G_{+}-G_{-})[h,f]G_{+}^{2}[h,g]\\\
-(G_{+}-G_{-})[h,g]G_{-}^{2}[h,f]-(f\leftrightarrow g)\Big{)}.$ (128)

Integrating over the temperature and using the formula

$\displaystyle\int_{T}^{\infty}\frac{dT}{T^{2}}\mathfrak{f}^{\prime}(z)=\frac{c_{1}(\mathfrak{f}(z))-\log
2}{z^{2}},$ (129)

where

$\displaystyle c_{1}(x)=\int_{0}^{x}dt\log\left(\frac{1-t}{t}\right)=-x\log
x-(1-x)\log(1-x),$ (130)

gives

$\displaystyle\begin{split}\nu^{A}=\frac{\eta^{A}}{T}=\frac{1}{2\pi}\int_{-\infty}^{\infty}dz\,c_{1}\big{(}\mathfrak{f}(z)\big{)}{\rm
Tr}\Big{(}(G_{+}&-G_{-})[h,f]G_{+}^{2}[h,g]\\\
&-(G_{+}-G_{-})[h,g]G_{-}^{2}[h,f]\Big{)}-(f\leftrightarrow g).\end{split}$
(131)

Here we normalized the thermoelectric coefficients to be 0 in the infinite-
temperature state. Note that since in the limit $T\rightarrow 0$ the Fermi-
Dirac distribution $\mathfrak{f}(z)$ becomes a step-function, and since
$c_{1}(0)=c_{1}(1)=0$, both $\nu^{A}(T)$ and $\eta^{A}(T)/T$ vanish at $T=0$
regardless of the choice of the Hamiltonian $h$.

###  D.5 Symmetric part

Symmetric parts of transverse thermoelectric coefficients are

$\displaystyle\nu^{S}_{xy}$
$\displaystyle=-\frac{\beta}{8\pi}\int_{-\infty}^{\infty}dz\,\mathfrak{f}(z){\rm
Tr}\big{\\{}[h,f]G_{+}^{2}[h^{2},g](G_{+}-G_{-})-[h,f](G_{+}-G_{-})[h^{2},g]G_{-}^{2}\big{\\}}+(f\leftrightarrow
g),$ $\displaystyle\eta^{S}_{xy}$
$\displaystyle=-\frac{1}{8\pi}\int_{-\infty}^{\infty}dz\,\mathfrak{f}(z){\rm
Tr}\big{\\{}[h^{2},f]G_{+}^{2}[h,g](G_{+}-G_{-})-[h^{2},f](G_{+}-G_{-})[h,g]G_{-}^{2}\big{\\}}+(f\leftrightarrow
g).$

As explained in the body of the paper, longitudinal parts are given by the
same formulas with a more general choice of the functions $f,g$.

## References

  * Kubo  _et al._ (2012) R. Kubo, M. Toda, and N. Hashitsume,  _Statistical Physics II: Nonequilibrium Statistical Mechanics_, Springer Series in Solid-State Sciences (Springer Berlin Heidelberg, 2012).
  * Luttinger (1964) J. M. Luttinger, “Theory of Thermal Transport Coefficients,” Phys. Rev. 135, A1505 (1964).
  * Cooper  _et al._ (1997) N. R. Cooper, B. I. Halperin, and I. M. Ruzin, “Thermoelectric response of an interacting two-dimensional electron gas in a quantizing magnetic field,” Phys. Rev. B 55, 2344 (1997).
  * Qin  _et al._ (2011) T. Qin, Q. Niu, and J. Shi, “Energy Magnetization and the Thermal Hall Effect,” Phys. Rev. Lett. 107, 236601 (2011).
  * Kapustin and Spodyneiko (2020) A. Kapustin and L. Spodyneiko, “Thermal Hall conductance and a relative topological invariant of gapped two-dimensional systems,” Phys. Rev. B101, 045137 (2020), arXiv:1905.06488 [cond-mat.str-el] .
  * Watanabe (2019) H. Watanabe, “A Proof of the Bloch Theorem for Lattice Models,” Journal of Statistical Physics 177 (2019), 10.1007/s10955-019-02386-1.
  * Bohm (1949) D. Bohm, “Note on a Theorem of Bloch Concerning Possible Causes of Superconductivity,” Phys. Rev. 75, 502 (1949).
  * Kapustin and Spodyneiko (2019) A. Kapustin and L. Spodyneiko, “Absence of energy currents in an equilibrium state and chiral anomalies,” Phys. Rev. Lett. 123, 060601 (2019).
  * Landau and Lifshitz (1984) L. D. Landau and E. M. Lifshitz,  _Electrodynamics of continuous media; 2nd ed._, Course of theoretical physics (Butterworth, Oxford, 1984).
  * Kitaev (2006) A. Kitaev, “Anyons in an exactly solved model and beyond,” Annals of Physics 321, 2 (2006).
  * (11) M. Levin, A. Kapustin, and L. Spodyneiko, “Nernst and Ettingshausen effects in gapped quantum materials,” arXiv:2103.02628 [cond-mat.str-el] .
  * Casimir (1945) H. B. G. Casimir, “On Onsager’s Principle of Microscopic Reversibility,” Rev. Mod. Phys. 17, 343 (1945).