††thanks: These authors contributed equally.††thanks: These authors
contributed equally.

# Kramers’ degenerate magnetism and superconductivity

Adil Amin  Hao Wu  Tatsuya Shishidou  Daniel F. Agterberg  1Department of
Physics, University of Wisconsin–Milwaukee, Milwaukee, Wisconsin 53201, USA

###### Abstract

Motivated by the recent discovery of odd-parity multipolar antiferromagnetic
order in $\text{CeRh}_{2}\text{As}_{2}$, we examine the coexistence of such
translation invariant Kramers’ degenerate magnetic states and
superconductivity. We show that the presence of such magnetic states
generically suppresses superconductivity, whether it be spin-singlet or spin-
triplet, unless the magnetic state drives a symmetry-required pair density
wave (PDW) superconducting order. We apply our results to
$\text{CeRh}_{2}\text{As}_{2}$, where no pair density wave order appears; and
to the loop current order in the cuprates, where such pair density wave
superconductivity must appear together with Bogoliubov Fermi surfaces. In the
former case, we explain why superconductivity is not suppressed.

Introduction The discovery of likely field-induced odd-parity
superconductivity in CeRh2As2 [1] has generated a great deal of interest [2,
3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16]. Subsequent to this, in
addition to quadrupolar order [8], CeRh2As2 has been found to host odd-parity
multipolar antiferromagnetic (AFM) order [17, 18]. In particular, nuclear
quadrupolar resonance (NQR) and nuclear magnetic resonance (NMR) have observed
that this order develops at a temperature below the superconducting transition
temperature ($T_{c}$). The symmetry of such a state is unusual since it
independently breaks parity ($\mathcal{P}$) and time reversal ($\mathcal{T}$)
symmetries but is invariant under the action of their product $\mathcal{PT}$.
Since $\mathcal{PT}$ is still preserved, electronic states retain the well-
known two-fold Kramers’ degeneracy, hence we shall refer to such magnetic
states as Kramers’ degenerate magnets. Such order has been referred to as odd
parity multipole order [19] and magnetic toroidal order [20, 21] in different
contexts. Since neither $\mathcal{P}$ nor $\mathcal{T}$ symmetries are
present, the electronic dispersion $\xi_{\boldsymbol{k}}$ does not satisfy the
usual relation $\xi_{\boldsymbol{k}}=\xi_{-\boldsymbol{k}}$. The asymmetry
between $\xi_{\boldsymbol{k}}$ and $\xi_{-\boldsymbol{k}}$
($\xi_{\boldsymbol{k}}\neq\xi_{-\boldsymbol{k}}$) naturally opens up the
question of what effect Kramers’ degenerate magnetic order has on the
superconducting state. This question has been addressed previously for multi-
band systems [19, 22] and also for loop current order in cuprates [23, 24,
25]. Superconductivity is usually a low energy scale phenomenon, thus it is
desirable to develop a general single band theory to capture the interplay of
Kramers’ degenerate order and superconductivity. Here we provide such a
general theoretical framework which allows us to encompass arbitrary energy
dispersions and gap functions, allowing us to draw broad conclusions.

Our theory shows that both spin-singlet and spin-triplet superconductivity are
strongly suppressed by such an order (such states can still survive in one
dimensional systems [19]). However, if the magnetic state belongs to a special
symmetry class then superconductivity can still appear. In particular, if the
magnetic state belongs to a vector representation of the point group, then a
symmetry dictated pair density wave (PDW) superconducting state is stabilized.
Furthermore, we find this magnetic order often implies the existence of Fermi
surfaces in the superconducting states. Such Fermi surfaces are called
Bogoliubov Fermi surfaces (BFS). When the magnetic state leads to symmetry
required PDW superconductivity, we find that the PDW order reduces the size of
the Bogoliubov Fermi surfaces, providing a natural microscopic mechanism for
the appearance of such PDW order. We apply this framework to superconductivity
coexisting with odd-parity multipolar AFM order in CeRh2As2 and coexisting
with the loop current order in the cuprates. In our companion paper, we also
apply the general framework developed to obtain a nematic BFS in the
tetragonal phase of FeSe1-xSx [21]. Such a BFS has recently been observed in
ARPES measurements for FeSe1-xSx [26].

Odd parity energy dispersion As discussed above, Kramers’ degenerate magnetic
order gives rise to asymmetry between the electronic dispersions
$\xi_{\boldsymbol{k}}$ and $\xi_{-\boldsymbol{k}}$. This asymmetry is the key
ingredient in our theory and we define
$\xi_{-,\boldsymbol{k}}=(\xi_{\boldsymbol{k}}-\xi_{-\boldsymbol{k}})/2$ to
quantify this. The ${\boldsymbol{k}}$ dependence of $\xi_{-,\boldsymbol{k}}$
is governed by the symmetry of the Kramers’ degenerate magnetic order [20].
Since these magnetic orders are translation invariant, they correspond to odd-
parity irreducible representations of the crystallographic point group.
Therefore, $\xi_{-,\boldsymbol{k}}$ can be constructed by forming products of
momentum components $\boldsymbol{k}_{i}$ that belong to the corresponding
irreducible representation of the crystallographic point group. Some examples
of the predicted form of $\xi_{-,\boldsymbol{k}}$ can be found in [27, 28, 22,
19, 29, 30, 31, 32, 33] and the recent experimental observation of this has
been reported in MnAu2 [27]. Here we discuss what form
$\xi_{-,\boldsymbol{k}}$ takes for the observed odd-parity AFM order in
CeRh2As2 and for the loop current order in the cuprates.

In CeRh2As2, with space group P4/nmm, there are two Ce atoms per unit cell.
Including a single Kramers’ doublet on each Ce site leads to the tight-binding
Hamiltonian [1],

$\displaystyle H_{0}=\ $
$\displaystyle\epsilon_{00,{\bm{k}}}\tau_{0}\sigma_{0}+\epsilon_{x0,{\bm{k}}}\tau_{x}\sigma_{0}+\epsilon_{y0,{\bm{k}}}\tau_{y}\sigma_{0}+\epsilon_{zx,{\bm{k}}}\tau_{z}\sigma_{x}+\epsilon_{zy,{\bm{k}}}\tau_{z}\sigma_{y}+\epsilon_{zz,{\bm{k}}}\tau_{z}\sigma_{z}$
$\displaystyle=\ $ $\displaystyle[t_{1}(\cos k_{x}+\cos
k_{y})-\mu]\tau_{0}\sigma_{0}+t_{c,1}\cos\frac{k_{z}}{2}\cos\frac{k_{x}}{2}\cos\frac{k_{y}}{2}\tau_{x}\sigma_{0}+t_{c,2}\sin\frac{k_{z}}{2}\cos\frac{k_{x}}{2}\cos\frac{k_{y}}{2}\tau_{y}\sigma_{0}$
$\displaystyle+\alpha_{R}(\sin k_{y}\tau_{z}\sigma_{x}-\sin
k_{x}\tau_{z}\sigma_{y})+\lambda\sin k_{z}(\cos k_{x}-\cos k_{y})\sin
k_{x}\sin k_{y}\tau_{z}\sigma_{z}.$ (1)

The $\tau_{i}$ Pauli matrices encode the Ce site basis, and the $\sigma_{i}$
Pauli matrices encode the spin basis. Here $\alpha_{R}$ and $\lambda$ are the
Rashba and Ising spin-orbit couplings. To this, we add the minimal coupling to
the odd-parity AFM order. This order is of opposite sign on the two Ce sites
in the unit cell and the Ce moments are oriented along the $z$ axis [17, 18].
Consequently, we write this coupling as $H_{c}=M_{z}\tau_{z}\sigma_{z}$. The
resulting Hamiltonian gives rise to two pseudospin (or doubly-degenerate)
bands. Treating $H_{c}$ as a perturbation to $H_{0}$ yields

$\xi_{-,\boldsymbol{k}}=\pm\frac{M_{z}\lambda\sin k_{z}(\cos k_{x}-\cos
k_{y})\sin k_{x}\sin
k_{y}}{\sqrt{\epsilon_{x0}^{2}+\epsilon_{y0}^{2}+\epsilon_{zx}^{2}+\epsilon_{zy}^{2}+\epsilon_{zz}^{2}}},$
(2)

where the $\pm$ denotes that $\xi_{-,\boldsymbol{k}}$ has a different sign on
the two bands. A key feature of the $\xi_{-,\boldsymbol{k}}$ is that even
though it has the same form for both pseudospin partners, it is generated
through the spin-orbit coupling term $\lambda$. Formally,
$\xi_{-,\boldsymbol{k}}$ has a momentum structure that belongs to the $A_{1u}$
representation of the $D_{4h}$ point group.

The form of $\xi_{-,\boldsymbol{k}}$ for loop current order in the cuprates
has been derived earlier and takes the form [23, 34, 35, 36, 37]

$\displaystyle\xi_{-,\boldsymbol{k}}=2J\left[\sin k_{x}-\sin
k_{y}-\sin\left(k_{x}-k_{y}\right)\right].$ (3)

where $J$ is the strength of the order. Physically, this term originates from
orbital currents that form closed loops within the CuO2 unit cell. In this
case, unlike for CeRh2As2 above, spin-orbit coupling does not play a role in
the origin of $\xi_{-,\boldsymbol{k}}$. Formally, for the loop current order,
the momentum structure of $\xi_{-,\boldsymbol{k}}$ belongs to the $E_{u}$
representation of the $D_{4h}$ point group. In contrast to the odd-parity
multipole AFM order in CeRh2As2, here the momentum structure of
$\xi_{-,\boldsymbol{k}}$ belongs to a vector representation of the $D_{4h}$
point group. This, as we shall see later, leads to the stabilization of PDW
superconductivity for the cuprates.

General Theory The Kramers’ degeneracy that is retained with the presence of
Kramers’ degenerate magnetic order allows a single band theory for
superconductivity [38, 39]. In the weak coupling limit, this can be formulated
for arbitrary Fermi surfaces, gap functions, and $\xi_{-,\boldsymbol{k}}$ and
this further will allow this physics to be treated within the framework of
quasi-classical theory, a powerful framework within which to examine
superconductivity [40]. We consider the following single-band Hamiltonian

$H=\sum_{\boldsymbol{ks}}\varepsilon_{\boldsymbol{k}}c_{\boldsymbol{ks}}^{\dagger}c_{\boldsymbol{ks}}+\frac{1}{2}\sum_{\boldsymbol{k}\boldsymbol{k}^{\prime}\alpha\beta\alpha^{\prime}\beta^{\prime}}V_{\alpha\beta\alpha^{\prime}\beta^{\prime}}\left(\boldsymbol{k},\boldsymbol{k}^{\prime}\right)c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\alpha}^{\dagger}c_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\beta}^{\dagger}c_{-\boldsymbol{k}^{\prime}+\frac{\boldsymbol{q}}{2},\beta^{\prime}}c_{\boldsymbol{k}^{\prime}+\frac{\boldsymbol{q}}{2},\alpha^{\prime}}$
(4)

where the operator
$c_{\boldsymbol{ks}}^{\dagger}\left(c_{\boldsymbol{ks}}\right)$ creates
(annihilates) electrons with momentum $\boldsymbol{k}$ and spin $s$;
$s,\alpha,\beta,\alpha^{\prime},\beta^{\prime}=\uparrow,\downarrow$;
$\boldsymbol{q}$ is the center-of-mass momentum of the Cooper pairs;
$\varepsilon_{\boldsymbol{k}}\equiv\xi_{\boldsymbol{k}}-\mu=\xi_{+,\boldsymbol{k}}+\xi_{-,\boldsymbol{k}}-\mu$,
and $\xi_{+,\boldsymbol{k}}=(\xi_{\boldsymbol{k}}+\xi_{-\boldsymbol{k}})/2$ is
the usual $\boldsymbol{k}$ symmetric dispersion included for single-band
theories, $\mu$ is the chemical potential. The role of Kramers’ degenerate
magnetic order is included through the addition of $\xi_{-,\boldsymbol{k}}$.
Here we explicitly include the Cooper center of mass momentum since later we
find that symmetry requires this must be nonzero for Kramers’ degenerate
magnetic order belonging to a vector representation of the point group.

Prior to presenting the results of our analysis, we note that generally
Kramers’ degenerate magnetic order leads to mixing between even parity,
pseudospin-singlet, and odd parity, pseudospin-triplet, superconductivity.
This can be understood from Ginzburg-Landau (GL) theory. In particular,
consider an even parity order parameter $\psi$, an odd parity order parameter
$\eta$, and a Kramers’ degenerate magnetic order parameter $M$. Then symmetry
allows the following term $iM(\psi\eta^{*}-\eta\psi^{*})$ in the GL free
energy, which ensures a singlet-triplet mixed order parameter of the form
$\psi+i\eta$. However, in the single band limit we consider here, it is
possible to show that this pseudospin singlet-triplet mixing vanishes if we
assume that the pairing interactions are not changed by the Kramers’
degenerate magnetic order. This is justified if $\xi_{-,\boldsymbol{k}}$ has a
much smaller energy scale than the electronic bandwidth. Hence we consider
each pairing channel independently. For pseudospin singlet order, we take
$V_{\alpha\beta\alpha^{\prime}\beta^{\prime}}\left(\boldsymbol{k},\boldsymbol{k}^{\prime}\right)=\frac{V_{s}}{2}\left[f_{\boldsymbol{k}}i{\sigma}_{\mathrm{2}}\right]_{\alpha\beta}\left[f_{\boldsymbol{k}^{\prime}}i{\sigma}_{\mathrm{2}}\right]_{\alpha^{\prime}\beta^{\prime}}^{\dagger}$
while for pseudospin triplet order we take
$V_{\alpha\beta\alpha^{\prime}\beta^{\prime}}\left(\boldsymbol{k},\boldsymbol{k}^{\prime}\right)=-\frac{V_{t}}{2}\left[\boldsymbol{d}_{\boldsymbol{k}}\cdot{\boldsymbol{\sigma}}i{\sigma_{2}}\right]_{\alpha\beta}\left[\boldsymbol{d}_{\boldsymbol{k}^{\prime}}\cdot{\boldsymbol{\sigma}}i{\sigma_{2}}\right]_{\alpha^{\prime}\beta^{\prime}}^{\dagger}$
and assume that both $f_{\boldsymbol{k}}$ and
$\boldsymbol{d}_{\boldsymbol{k}}$ can be chosen real - in practice this
applies to singly degenerate superconducting irreducible representations. It
is not difficult to generalize this to higher-dimensional irreducible
representations.

Since we are interested in the evolution of the superconducting state both at
low temperatures and near $T_{c}$, we compute the free energy within mean-
field theory

$\Omega_{\mathrm{s}}-\Omega_{\mathrm{n}}=\left\\{-2k_{B}\mathrm{~{}T}\sum_{\boldsymbol{k}}\ln\left[\cosh\left(\frac{\beta
E_{\boldsymbol{k},\boldsymbol{q}}}{2}\right)\right]+\frac{\left|\Delta_{\boldsymbol{q}}\right|^{2}}{V_{s}}\right\\}-\left\\{-2k_{B}\mathrm{~{}T}\sum_{\boldsymbol{k}}\ln\left[\cosh\left(\frac{\beta
E_{\boldsymbol{k},\boldsymbol{q},0}}{2}\right)\right]\right\\}$ (5)

where
$E_{\boldsymbol{k},\boldsymbol{q}}\equiv(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}+\sqrt{\epsilon_{\boldsymbol{k},\boldsymbol{q}}^{2}+\left|\Delta_{\boldsymbol{q}}f_{\boldsymbol{k}}\right|^{2}}$
with
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}\equiv\frac{1}{2}\left(\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}-\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}\right)$
and
$\epsilon_{\boldsymbol{k},\boldsymbol{q}}\equiv\frac{1}{2}\left(\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}\right)$,
$E_{\boldsymbol{k},\boldsymbol{q},0}\equiv(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}+\sqrt{\epsilon_{\boldsymbol{k},\boldsymbol{q}}^{2}}$,
$V_{s}$ is the pairing interaction. Eq. 5 applies to pseudospin-singlet order,
for the pseudospin-triplet order, $\left|f_{\boldsymbol{k}}\right|^{2}$ is
replaced by $\left|\boldsymbol{d}_{\boldsymbol{k}}\right|^{2}$. Minimizing
this free energy with respect to $\left|\Delta_{\boldsymbol{q}}\right|$ and
$\boldsymbol{q}$ leads to the gap equation and the condition on
$\boldsymbol{q}$ that ensures a vanishing supercurrent

$\ln\left(\frac{T}{T_{c0}}\right)=\left\langle 2\pi
k_{B}T\left|f_{\boldsymbol{k}}\right|^{2}Re\sum_{n=0}^{\infty}\left(\frac{1}{\delta}-\frac{1}{\hbar\omega_{n}}\right)\right\rangle_{k}$
(6)
$\left\langle\frac{\partial(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}}{\partial\boldsymbol{q}}\left|f_{\boldsymbol{k}}\right|^{2}\operatorname{Im}\sum_{n=0}^{\infty}\frac{1}{\delta\left(\delta+\hbar\omega_{n}+i(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}\right)}\right\rangle_{k}=0$
(7)

where
$\delta\equiv\sqrt{\left(\hbar\omega_{n}+i(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}\right)^{2}+\left|\Delta_{\boldsymbol{q}}\right|^{2}\left|f_{k}\right|^{2}}$,
and $\omega_{n}$ is the Matsubara frequency which satisfies
$\hbar\omega_{n}=(2n+1)\pi k_{B}T,n\in\mathbb{Z}$. Throughout this work, we
adopt the convention that the average over the Fermi surface
$\left\langle\left|f_{k}\right|^{2}\right\rangle_{k}=1$. This average over the
Fermi surface is weighted by the momentum-dependent density of states. By
numerically solving this set of equations, we can find the gap and the optimal
$\boldsymbol{q}$ as functions of $T$. As shown in the Appendix, the Green’s
function and the anomalous Green’s function are

$\hat{G}=-\frac{\hbar\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)\sigma_{0}}{\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)\left(\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}-i\hbar\omega_{n}\right)+|\psi(\boldsymbol{k},\boldsymbol{q})|^{2}}$
(8)
$\hat{F}=\frac{\hbar\hat{\Delta}}{\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}-i\hbar\omega_{n}\right)\left(\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)+|\psi(\boldsymbol{k},\boldsymbol{q})|^{2}}$
(9)

where $\sigma_{0}$ is the identity matrix and $\hat{\Delta}$ is the gap
matrix. As can be also seen in the Appendix, for the pseudospin-triplet
$|\psi(\boldsymbol{k},\boldsymbol{q})|^{2}$ is replaced by
$|\boldsymbol{d}(\boldsymbol{k},\boldsymbol{q})|^{2}$

To gain an understanding of the role of the Kramers’ degenerate magnetic
order, we initially consider $T$ near $T_{c}$. The equation determining the
critical temperature $T_{c}^{q}$ is

$\ln\left(\frac{T_{c}^{q}}{T_{c0}}\right)=\left\langle\left|f_{k}\right|^{2}\left\\{\psi\left(\frac{1}{2}\right)-Re\left[\psi\left(\frac{1}{2}+\frac{i(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}}{2\pi
k_{B}T_{c}^{q}}\right)\right]\right\\}\right\rangle_{k}$ (10)

where $T_{c0}$ is the critical temperature for
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}=0$, the digamma function
$\psi(z)=-\gamma+\sum_{n=0}^{\infty}\left(\frac{1}{n+1}-\frac{1}{n+z}\right)$,
and $\gamma$ is the Euler-Mascheroni constant. Assuming
$|\boldsymbol{q}|\ll|\boldsymbol{k}|$, we have
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}=\xi_{-,\boldsymbol{k}}+\frac{\boldsymbol{q}}{2}\cdot\boldsymbol{\nabla}\xi_{+,\boldsymbol{k}}$.
The expression
$\psi\left(\frac{1}{2}\right)-Re\left[\psi\left(\frac{1}{2}+\frac{i(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}}{2\pi
k_{B}T_{c}^{q}}\right)\right]$ is intrinsically negative, which yields
$T_{c}^{q}<T_{c0}$. From the monotonic decreasing dependence of the Digamma
function, we find that $T_{c}^{q}$ is decreasing with increasing
$|(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}|$. Eventually,
$T_{c}^{q}$ will become zero, i.e., superconductivity is destroyed. Using the
asymptotic expansion of the Digamma function, the critical value for
$|(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}|$ corresponding to
$T_{c}=0$ is given by

$\left\langle\left|f_{\boldsymbol{k}}\right|^{2}\ln\left(\frac{|(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}|}{2\pi
k_{B}T_{c0}}\right)\right\rangle_{k}=\psi\left(\frac{1}{2}\right)$ (11)

We conclude that once the magnitude of
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}$ becomes on the order of
$k_{B}T_{c0}$ superconductivity is destroyed. We note that the usual Pauli
limiting field ,$H_{p}$, for an $s$-wave superconductor is given by
$\mu_{B}H_{P}=\Delta_{0}/\sqrt{2}$ hence, throughout this paper to we choose
to express $(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}$ in units of
$\mu_{B}H_{P}$ to provide some context. Importantly, since
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}$ depends not only on
$\xi_{-,\boldsymbol{k}}$, but also on $\boldsymbol{q}$, it is possible that
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}$ can become small even if
$\xi_{-,\boldsymbol{k}}$ is larger than the superconducting gap. To illustrate
this, it is useful to consider an example with
$\xi_{+,\boldsymbol{k}}=\frac{\hbar^{2}\boldsymbol{k}^{2}}{2m}$ and
$\xi_{-,\boldsymbol{k}}=\alpha k_{x}$. Then,
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}=\alpha
k_{x}+\frac{\hbar^{2}}{2m}\left(q_{x}k_{x}+q_{y}k_{y}\right)$. The choice
$\boldsymbol{q}=\left(q_{x},q_{y}\right)=\left(-\frac{2m\alpha}{\hbar^{2}},0\right)$
will cancel the $\xi_{-,\boldsymbol{k}}$ term making
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}=0$. Therefore, the
critical temperature remains unchanged from $T_{c0}$ in spite of a
$\xi_{-,\boldsymbol{k}}$ that can be much larger than the gap.

To understand in more detail how a PDW state with a non-zero $\boldsymbol{q}$
can resurrect a suppressed $T_{c}$ due to $\xi_{-,\boldsymbol{k}}$, we carry
out a small $(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}$ expansion.
The critical temperature is then determined by

$\ln\left(\frac{T_{c}^{q}}{T_{c0}}\right)=\frac{7\zeta(3)}{\left(2\pi
k_{B}T_{c}^{q}\right)^{2}}\\\
\times\left\langle\left|f_{\boldsymbol{k}}\right|^{2}\left[-\xi_{-,\boldsymbol{k}}^{2}-\frac{\hbar^{2}(\boldsymbol{q}\cdot\boldsymbol{v})^{2}}{4}-\hbar\xi_{-,\boldsymbol{k}}\boldsymbol{q}\cdot\boldsymbol{v}\right]\right\rangle_{k}$
(12)

where
$\boldsymbol{v}(\boldsymbol{k})\equiv\frac{1}{\hbar}\nabla\xi_{+,\boldsymbol{k}}$
is the velocity. Only the third term
($-\hbar\xi_{-,\boldsymbol{k}}\boldsymbol{q}\cdot\boldsymbol{v}$) can be
positive and hence increase $T_{c}^{q}$. Thus, the existence of a PDW state
depends on a nonzero average
$\left\langle\left|f_{\boldsymbol{k}}\right|^{2}\left[-\hbar\xi_{-,\boldsymbol{k}}\boldsymbol{q}\cdot\boldsymbol{v}\right]\right\rangle_{k}$.
Since the components of the velocity $\boldsymbol{v}$ belong to a vector
representation of the point group, the only way to get this average nonzero is
to require that the momentum structure of $\xi_{-,\boldsymbol{k}}$ also
belongs to a vector representation. This implies that the Kramers’ degenerate
magnetic order belongs to a vector representation of the point group. Indeed,
when the Kramers’ degenerate magnetic order belongs to a vector
representation, symmetry implies that the Ginzburg-Landau free energy contains
a non-vanishing Lifshitz invariant of the form
$\sum_{j}\gamma_{j}\left[\psi\left(D_{j}\psi\right)^{*}+\psi^{*}\left(D_{j}\psi\right)\right]$
with $D_{j}\equiv-i\nabla_{j}-2eA_{j}$, which guarantees the appearance of a
PDW state with non-zero $\boldsymbol{q}$. This symmetry dictated PDW order
typically takes the form $\psi=\psi_{0}e^{i\boldsymbol{q}\cdot\boldsymbol{r}}$
with a single-plane wave, unlike the multiple plane-wave solutions often
associated with PDW order [41]. The existence of Lifshitz invariants, and the
concomitant finite momentum PDW pairing, is important for the superconducting
diode effect [42, 43, 44, 45]. This suggests that superconductivity coexisting
with the Kramers’ degenerate magnetic order provides a route towards creating
the superconducting diode effect.

In addition to the generating PDW states, Kramers’ degenerate magenetic order
can also lead to novel low energy excitation spectra. In particular, by
finding the poles of the Green’s function $G$, we can obtain the quasiparticle
dispersion

$E=(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}\pm\sqrt{\epsilon_{\boldsymbol{k},\boldsymbol{q}}^{2}+\left|\Delta_{\boldsymbol{q}}f_{\boldsymbol{k}}\right|^{2}}$
(13)

In principle, because of the appearance of
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}$ in this expression,
Bogoliubov Fermi surfaces (BFS) can exist. These BFS are given by $E=0$ and
appear for $\boldsymbol{k}$ where
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}$ is larger than
$|\Delta_{\boldsymbol{q}}f_{\boldsymbol{k}}|$. Indeed, these are guaranteed to
appear when there are gap nodes present i.e. $f_{\boldsymbol{k}}=0$ (or
$\left|\boldsymbol{d}_{\boldsymbol{k}}\right|=0$) provided that
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}$ does not also vanish on
these nodes. When this occurs,
$(\Delta\varepsilon)_{\boldsymbol{k},\boldsymbol{q}}$ inflates the nodes of
the original gap functions into Bogoliubov Fermi surfaces [46, 47, 48, 49]. In
our companion paper [21] we carry out a detailed analysis of the existence of
Bogoliubov Fermi surfaces in FeSe1-xSx.

Application to CeRh2As2 In addition to Eq. 2 for $\xi_{-,\boldsymbol{k}}$, we
require a description of the normal state to apply our results to CeRh2As2.
DFT results [7, 50] yield a Fermi surface as depicted in Fig. 1. They further
reveal that the majority of the density of states (DOS) ($\sim 80\%$) is
concentrated on the Fermi surface structures near the X–M line, which we refer
to as beans. Here we include only those beans and, for simplicity, ignore the
$c$-axis dependence of $\xi_{+_{,}\boldsymbol{k}}$ (this assumption does not
qualitatively change our results). Furthermore, we assume that a
$d_{x^{2}-y^{2}}$ gap structure as depicted in Fig. 1 appears when no magnetic
field is applied. This assumption also does not qualitatively change the
results unless the gap function has nodes on the beans. We briefly comment on
this possibility later.

Figure 1:  Sketch of the normal state Fermi surface of CeRh2As2 for $k_{z}=0$
found from renormalized DFT calculations. The majority of the density of
states sit on the beans. The $d_{x^{2}-y^{2}}$ pairing symmetry considered
here is denoted by the $\Delta$ values as shown.

To describe theses beans, we develop a power series expansion in powers of
$\delta k_{x},\delta k_{y}$ centered about the minima of the band along the
X–M line. Explicitly considering a bean centered at
$(k_{x},k_{y})=\left(\frac{\pi}{a},\frac{0.4\pi}{a}\right)$, we have the
dispersion

$\xi_{+,\boldsymbol{k}}=\alpha\left(\delta k_{x}\right)^{2}+\beta\left(\delta
k_{y}\right)^{2}+\gamma_{1}\left(\delta k_{x}\right)^{2}\left(\delta
k_{y}\right)+\gamma_{2}\left(\delta k_{y}\right)^{3}.$ (14)

Mirror symmetry prevents terms odd in $\delta k_{x}$ appearing in this
expression. We also expand $\xi_{-,\boldsymbol{k}}$ in Eq. 2 around the
centers of the beans. For the pocket centered at
$\left(\frac{\pi}{a},\frac{0.4\pi}{a}\right),\xi_{-,\boldsymbol{k}}=$
$\tilde{\lambda}\delta k_{x}a\sin k_{z}$.

The numerical solution of the linearized gap equation is shown in Fig. 2.
Because the Kramers’ degenerate magnetic order does not belong to a vector
representation, a PDW state is not required by symmetry in this case. Thus
when $\xi_{-,\boldsymbol{k}}$ becomes sufficiently large, superconductivity is
suppressed as seen in Fig. 2. However, this is not observed in experiment [17,
18], which raises the question, how does superconductivity still survive. Here
we suggest that the explanation for this persistence hinges on the value of
$\tilde{\lambda}$ defined above. $\tilde{\lambda}$ on the beans is determined
by the ratio of the Ising coupling ($\lambda$) to the Rashba spin-orbit
coupling ($\alpha_{R}$) which occur in Eq. 1. In particular, using effective
$g$ factors from DFT based renormalized band structure calculations we
constrain $\tilde{\lambda}<0.2M_{z}$ [51]. The critical temperature of $M_{z}$
is of the order of the critical temperature of superconductivity, hence
$M_{z}\approx\mu_{B}H_{P}$ thus $\tilde{\lambda}<0.2\mu_{B}H_{P}$. We further
estimate $a\text{k}_{F}\approx\frac{1}{2}$. Hence
$\frac{a\text{k}_{F}\tilde{\lambda}}{\mu_{B}H_{P}}<0.1$. This is shown as the
shaded region in Fig. 2. Thus superconductivity still persists despite the
presence of Kramers’ degenerate magnetic order.

While the pairing state we consider above has $d_{x^{2}-y^{2}}$ symmetry, it
is nodeless on the Fermi surface beans which do not intersect the
$d_{x^{2}-y^{2}}$ nodes. However, the $d_{x^{2}-y^{2}}$ gap will have nodes
along on the $\Gamma$-centered Fermi surfaces, possibly allowing for
Bogoliubov Fermi surfaces to appear in the superconducting state. However,
these nodes will not be inflated to Bogoliubov Fermi surfaces since
$(\Delta\varepsilon)_{\boldsymbol{k}}$ vanishes along the diagonal the nodes.
We can further consider other gap functions with nodes to see if Bogoliubov
Fermi surfaces are expected. For example, nodes can appear on the beans for a
$d_{xy}$ like gap. The beans will then have nodes where they intersect the
Brillouin Zone face. However, these nodes will also not be inflated to form
Bogoliubov Fermi surfaces because again $(\Delta\varepsilon)_{\boldsymbol{k}}$
vanishes on this face. Consequently, we do not expect Bogoliubov Fermi
surfaces to appear as a consequence of the Kramer’s degenerate magnetic order
in CeRh2As2 for any pairing symmetry.

Figure 2:  Transition temperature as a function of the antisymmetric coupling
parameter $\lambda$ in CeRh2As2. The parameter $a$ is the lattice constant,
$k_{F}\equiv\sqrt{\frac{\mu}{\sqrt{\alpha\beta}}}$ is the natural unit for
$\delta k_{x}$ and $\delta k_{y}$. We choose the values
$\alpha=\frac{1}{(1.25)^{2}}$, $\beta=\frac{1}{(0.55)^{2}}$, $\mu=0.7$, and
$\gamma_{1}=0.75$, $\gamma_{2}=-0.9$. The Pauli limiting field energy scale
$\mu_{B}H_{P}=\Delta_{0}/\sqrt{2}$, where $\Delta_{0}=1.764k_{B}T_{c0}$.

Application to loop current order in cuprates The normal state is given by the
dispersion $\xi_{+,\boldsymbol{k}}=-2t\left(\cos k_{x}+\cos
k_{y}\right)-4t^{\prime}\cos\left(k_{x}\right)\cos\left(k_{y}\right)$, where
$t$ and $t^{\prime}$ are the nearest and next nearest neighbor hoppings on the
square 2D lattice [23]. Eq. 3 above gives $\xi_{-,\boldsymbol{k}}$. The gap
function is taken to be
$\Delta_{\boldsymbol{k},\boldsymbol{q}}=\Delta_{\boldsymbol{q}}f_{\boldsymbol{k}}=\Delta_{\boldsymbol{q}}C\left(\cos
k_{x}-\cos k_{y}\right)$ with normalization factor $C$ (determined by
$\left\langle\left|f_{k}\right|^{2}\right\rangle_{k}=1$). Since the Kramers’
degenerate magnetic order belongs to a vector representation, our general
results imply that a PDW state is stabilized by the loop current order. In
Fig. 3, we compare the calculated values for $T_{c}$ with and without the PDW
order.

Figure 3:  Transition temperature as a function of the loop current order
parameter for both $\boldsymbol{q}=0$ and finite $\boldsymbol{q}$ theories.
The Pauli limiting field energy scale $\mu_{B}H_{P}=\Delta_{0}/\sqrt{2}$,
where $\Delta_{0}=1.764k_{B}T_{c0}$. Using hopping parameter $t$ as the energy
scale, we choose $t^{\prime}=-0.25$, $\mu=-0.9$. We can observe a nonmonotonic
behavior in $\boldsymbol{q}=0$ plot. This kind of behavior typically suggests
the existence of a PDW state.

Remarkably, the maximum $J$ value for the PDW state is about 15 times larger
than the $\boldsymbol{q}=0$ state. In addition, since the $d$-wave gap has
nodes, the loop current order will give rise to BFSs. This has been examined
previously when there is no PDW order [23, 52, 25]. As shown in Fig. 4, we
find that the BFS are strongly shrunk due to the PDW order. The gain in
condensation energy due to the shrinking of the BFS provides an energetic
mechanism to stabilize the PDW state.

Figure 4:  Bogoliubov Fermi surfaces. The upper plot is for $q=0$, the lower
plot is for $q\neq 0$.

Conclusions We considered the interplay of $\mathcal{PT}$ translation
invariant Kramers’ degenerate magnetic order and superconductivity. The
absence of $\mathcal{P}$ and $\mathcal{T}$ symmetry leads to an asymmetry in
the electronic dispersion relation ($\xi_{\textbf{k}}\neq\xi_{-\textbf{k}}$).
We develop a general single-band framework to capture the effect of the
Kramers’ degenerate magnetic order on superconductivity. We find that when the
magnetic states belong to a vector representation of the point group, pair
density wave states are stabilized. In the absence of such symmetries,
superconductivity is generically suppressed. We apply our framework to
$\text{CeRh}_{2}\text{As}_{2}$, where we explain why superconductivity still
persists in the presence of Kramers’ degenerate magnetic order, as well as the
cuprates where pair density wave superconductivity appears accompanied with
Bogoliubov Fermi surfaces.

Acknowledgments A.A., H.W., and D.F.A. acknowledge support through DE-
SC0021971. We would like to thank Erez Berg, Elena Hassinger, Seunghyun Khim,
Steve Kivelson, Christian Parsons, and Michael Weinert, for useful
discussions.

## References

  * Khim  _et al._ [2021] S. Khim, J. Landaeta, J. Banda, N. Bannor, M. Brando, P. M. R. Brydon, D. Hafner, R. Küchler, R. Cardoso-Gil, U. Stockert, A. P. Mackenzie, D. F. Agterberg, C. Geibel, and E. Hassinger, Field-induced transition within the superconducting state of $\mathrm{Ce}{\mathrm{Rh}}_{2}{\mathrm{As}}_{2}$, Science 373, 1012 (2021).
  * Schertenleib  _et al._ [2021] E. G. Schertenleib, M. H. Fischer, and M. Sigrist, Unusual ${H}\text{$-$}{T}$ phase diagram of $\mathrm{Ce}{\mathrm{Rh}}_{2}{\mathrm{As}}_{2}$: The role of staggered noncentrosymmetricity, Phys. Rev. Research 3, 023179 (2021).
  * Ptok  _et al._ [2021] A. Ptok, K. J. Kapcia, P. T. Jochym, J. Łażewski, A. M. Oleś, and P. Piekarz, Electronic and dynamical properties of ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$: Role of ${\mathrm{Rh}}_{2}{\mathrm{As}}_{2}$ layers and expected orbital order, Phys. Rev. B 104, L041109 (2021).
  * Möckli and Ramires [2021] D. Möckli and A. Ramires, Two scenarios for superconductivity in ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, Phys. Rev. Research 3, 023204 (2021).
  * Skurativska  _et al._ [2021] A. Skurativska, M. Sigrist, and M. H. Fischer, Spin response and topology of a staggered-Rashba superconductor, Phys. Rev. Research 3, 033133 (2021).
  * Nogaki  _et al._ [2021] K. Nogaki, A. Daido, J. Ishizuka, and Y. Yanase, Topological crystalline superconductivity in locally noncentrosymmetric ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, Physical Review Research 3, L032071 (2021).
  * Cavanagh  _et al._ [2022] D. Cavanagh, T. Shishidou, M. Weinert, P. M. R. Brydon, and D. F. Agterberg, Nonsymmorphic symmetry and field-driven odd-parity pairing in $\mathrm{Ce}{\mathrm{Rh}}_{2}{\mathrm{As}}_{2}$, Physical Review B 105, L020505 (2022).
  * Landaeta  _et al._ [2022] J. Landaeta, P. Khanenko, D. Cavanagh, C. Geibel, S. Khim, S. Mishra, I. Sheikin, P. M. R. Brydon, D. F. Agterberg, M. Brando, _et al._ , Field-angle dependence reveals odd-parity superconductivity in ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, Physical Review X 12, 031001 (2022).
  * Nogaki and Yanase [2022] K. Nogaki and Y. Yanase, Even-odd parity transition in strongly correlated locally noncentrosymmetric superconductors: Application to ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, Physical Review B 106, L100504 (2022).
  * Möckli and Ramires [2021] D. Möckli and A. Ramires, Superconductivity in disordered locally noncentrosymmetric materials: An application to ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, Physical Review B 104, 134517 (2021).
  * Hazra and Coleman [2022] T. Hazra and P. Coleman, Triplet pairing mechanisms from Hund’s-Kondo models: applications to ${\mathrm{UTe}}_{3}$ and ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, arXiv preprint arXiv:2205.13529 (2022).
  * Mishra  _et al._ [2022] S. Mishra, Y. Liu, E. D. Bauer, F. Ronning, and S. M. Thomas, Anisotropic magnetotransport properties of the heavy-fermion superconductor ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, Physical Review B 106, L140502 (2022).
  * Kimura  _et al._ [2021] S.-i. Kimura, J. Sichelschmidt, and S. Khim, Optical study of the electronic structure of locally noncentrosymmetric ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, Physical Review B 104, 245116 (2021).
  * Machida [2022] K. Machida, Violation of Pauli-Clogston limit in the heavy-fermion superconductor ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$: Duality of itinerant and localized 4 f electrons, Physical Review B 106, 184509 (2022).
  * Semeniuk  _et al._ [2023] K. Semeniuk, D. Hafner, P. Khanenko, T. Lühmann, J. Banda, J. F. Landaeta, C. Geibel, S. Khim, E. Hassinger, and M. Brando, Superconductivity versus quadrupole density wave in ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, arXiv preprint arXiv:2301.09151 (2023).
  * Onishi  _et al._ [2022] S. Onishi, U. Stockert, S. Khim, J. Banda, M. Brando, and E. Hassinger, Low-temperature thermal conductivity of the two-phase superconductor ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, Frontiers in Electronic Materials 2 (2022).
  * Kibune  _et al._ [2022] M. Kibune, S. Kitagawa, K. Kinjo, S. Ogata, M. Manago, T. Taniguchi, K. Ishida, M. Brando, E. Hassinger, H. Rosner, _et al._ , Observation of antiferromagnetic order as odd-parity multipoles inside the superconducting phase in ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, Physical Review Letters 128, 057002 (2022).
  * Kitagawa  _et al._ [2022] S. Kitagawa, M. Kibune, K. Kinjo, M. Manago, T. Taniguchi, K. Ishida, M. Brando, E. Hassinger, C. Geibel, and S. Khim, Two-dimensional xy-type magnetic properties of locally noncentrosymmetric superconductor ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, Journal of the Physical Society of Japan 91, 043702 (2022).
  * Sumita and Yanase [2016] S. Sumita and Y. Yanase, Superconductivity in magnetic multipole states, Physical Review B 93, 224507 (2016).
  * Hayami  _et al._ [2018] S. Hayami, M. Yatsushiro, Y. Yanagi, and H. Kusunose, Classification of atomic-scale multipoles under crystallographic point groups and application to linear response tensors, Physical Review B 98, 165110 (2018).
  * wu [2] H. Wu, A. Amin, Y. Yu, and D. F. Agterberg, submitted to Physical Review Letters (2023).
  * Sumita  _et al._ [2017] S. Sumita, T. Nomoto, and Y. Yanase, Multipole superconductivity in nonsymmorphic Sr2IrO4, Physical Review Letters 119, 027001 (2017).
  * Berg  _et al._ [2008] E. Berg, C.-C. Chen, and S. A. Kivelson, Stability of nodal quasiparticles in superconductors with coexisting orders, Physical Review Letters 100, 027003 (2008).
  * Allais and Senthil [2012] A. Allais and T. Senthil, Loop current order and d-wave superconductivity: Some observable consequences, Physical Review B 86, 045118 (2012).
  * Wang and Vafek [2013] L. Wang and O. Vafek, Quantum oscillations of the specific heat in d-wave superconductors with loop current order, Physical Review B 88, 024506 (2013).
  * Nagashima  _et al._ [2022] T. Nagashima, T. Hashimoto, S. Najafzadeh, S.-i. Ouchi, T. Suzuki, A. Fukushima, S. Kasahara, K. Matsuura, M. Qiu, Y. Mizukami, K. Hashimoto, Y. Matsuda, T. Shibauchi, S. Shin, and K. Okazaki, Discovery of nematic Bogoliubov Fermi surface in an iron-chalcogenide superconductor, Research Square preprint https://doi.org/10.21203/rs.3.rs-2224728/v1 (2022).
  * Fedchenko  _et al._ [2022] O. Fedchenko, L. Šmejkal, M. Kallmayer, Y. Lytvynenko, K. Medjanik, S. Babenkov, D. Vasilyev, M. Klaeui, J. Demsar, G. Schönhense, _et al._ , Direct observation of antiferromagnetic parity violation in the electronic structure of Mn2Au, Journal of Physics: Condensed Matter 34, 425501 (2022).
  * Šmejkal  _et al._ [2017] L. Šmejkal, J. Železnỳ, J. Sinova, and T. Jungwirth, Electric control of dirac quasiparticles by spin-orbit torque in an antiferromagnet, Physical Review Letters 118, 106402 (2017).
  * Urru  _et al._ [2022] A. Urru, J.-R. Soh, N. Qureshi, A. Stunault, B. Roessli, H. M. Rønnow, and N. A. Spaldin, Neutron scattering from local magnetoelectric multipoles: a combined theoretical, computational, and experimental perspective, arXiv preprint arXiv:2212.06779 (2022).
  * Bhowal and Spaldin [2021] S. Bhowal and N. A. Spaldin, Revealing hidden magnetoelectric multipoles using compton scattering, Physical Review Research 3, 033185 (2021).
  * Spaldin and Ramesh [2019] N. Spaldin and R. Ramesh, Advances in magnetoelectric multiferroics, Nature Materials 18, 203 (2019).
  * Yatsushiro  _et al._ [2022] M. Yatsushiro, R. Oiwa, H. Kusunose, and S. Hayami, Analysis of model-parameter dependences on the second-order nonlinear conductivity in $\mathcal{PT}$-symmetric collinear antiferromagnetic metals with magnetic toroidal moment on zigzag chains, Phys. Rev. B 105, 155157 (2022).
  * Hayami  _et al._ [2022] S. Hayami, M. Yatsushiro, and H. Kusunose, Nonlinear spin hall effect in $\mathcal{PT}$-symmetric collinear magnets, Phys. Rev. B 106, 024405 (2022).
  * Varma [1997] C. Varma, Non-Fermi-liquid states and pairing instability of a general model of copper oxide metals, Physical Review B 55, 14554 (1997).
  * Varma [1999] C. Varma, Pseudogap phase and the quantum-critical point in copper-oxide metals, Physical Review Letters 83, 3538 (1999).
  * Simon and Varma [2002] M. Simon and C. Varma, Detection and implications of a time-reversal breaking state in underdoped cuprates, Physical Review Letters 89, 247003 (2002).
  * Varma [2006] C. Varma, Theory of the pseudogap state of the cuprates, Physical Review B 73, 155113 (2006).
  * Sigrist and Ueda [1991] M. Sigrist and K. Ueda, Phenomenological theory of unconventional superconductivity, Reviews of Modern physics 63, 239 (1991).
  * Gor’kov and Rashba [2001] L. P. Gor’kov and E. I. Rashba, Superconducting 2d system with lifted spin degeneracy: mixed singlet-triplet state, Physical Review Letters 87, 037004 (2001).
  * Serene and Rainer [1983] J. Serene and D. Rainer, The quasiclassical approach to superfluid 3He , Physics Reports 101, 221 (1983).
  * Agterberg  _et al._ [2020] D. F. Agterberg, J. S. Davis, S. D. Edkins, E. Fradkin, D. J. V. Harlingen, S. A. Kivelson, P. A. Lee, L. Radzihovsky, J. M. Tranquada, and Y. Wang, The physics of pair-density waves: Cuprate superconductors and beyond, Annual Review of Condensed Matter Physics 11, 231 (2020).
  * Yuan and Fu [2022] N. F. Yuan and L. Fu, Supercurrent diode effect and finite-momentum superconductors, Proceedings of the National Academy of Sciences 119, e2119548119 (2022).
  * He  _et al._ [2022] J. J. He, Y. Tanaka, and N. Nagaosa, A phenomenological theory of superconductor diodes, New Journal of Physics 24, 053014 (2022).
  * Daido  _et al._ [2022] A. Daido, Y. Ikeda, and Y. Yanase, Intrinsic superconducting diode effect, Physical Review Letters 128, 037001 (2022).
  * Pal  _et al._ [2022] B. Pal, A. Chakraborty, P. K. Sivakumar, M. Davydova, A. K. Gopi, A. K. Pandeya, J. A. Krieger, Y. Zhang, M. Date, S. Ju, _et al._ , Josephson diode effect from cooper pair momentum in a topological semimetal, Nature Physics 18, 1228 (2022).
  * Agterberg  _et al._ [2017] D. F. Agterberg, P. M. R. Brydon, and C. Timm, Bogoliubov Fermi surfaces in superconductors with broken time-reversal symmetry, Physical review letters 118, 127001 (2017).
  * Brydon  _et al._ [2018] P. M. R. Brydon, D. F. Agterberg, H. Menke, and C. Timm, Bogoliubov Fermi surfaces: General theory, magnetic order, and topology, Physical Review B 98, 224509 (2018).
  * Link and Herbut [2020] J. M. Link and I. F. Herbut, Bogoliubov-Fermi surfaces in noncentrosymmetric multicomponent superconductors, Physical Review Letters 125, 237004 (2020).
  * Timm  _et al._ [2017] C. Timm, A. Schnyder, D. F. Agterberg, and P. M. R. Brydon, Inflated nodes and surface states in superconducting half-heusler compounds, Physical Review B 96, 094526 (2017).
  * Hafner  _et al._ [2022] D. Hafner, P. Khanenko, E.-O. Eljaouhari, R. Küchler, J. Banda, N. Bannor, T. Lühmann, J. Landaeta, S. Mishra, I. Sheikin, _et al._ , Possible quadrupole density wave in the superconducting kondo lattice ${\mathrm{CeRh}}_{2}{\mathrm{As}}_{2}$, Physical Review X 12, 011023 (2022).
  * Shishidou [2023] T. Shishidou (2023), in preparation.
  * Kivelson and Varma [2012] S. Kivelson and C. Varma, Fermi pockets in a d-wave superconductor with coexisting loop-current order, arXiv preprint  arXiv:1208.6498 (2012).

## Appendices

##  Appendix A Appendix A: Derivation of the Green’s functions

Starting from the single-band Hamiltonian

$H=\sum_{\boldsymbol{ks}}\varepsilon_{\boldsymbol{k}}c_{\boldsymbol{ks}}^{\dagger}c_{\boldsymbol{ks}}+\frac{1}{2}\sum_{\boldsymbol{k}\boldsymbol{k}^{\prime}\alpha\beta\alpha^{\prime}\beta^{\prime}}V_{\alpha\beta\alpha^{\prime}\beta^{\prime}}\left(\boldsymbol{k},\boldsymbol{k}^{\prime}\right)c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\alpha}^{\dagger}c_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\beta}^{\dagger}c_{-\boldsymbol{k}^{\prime}+\frac{\boldsymbol{q}}{2},\beta^{\prime}}c_{\boldsymbol{k}^{\prime}+\frac{\boldsymbol{q}}{2},\alpha^{\prime}}$
(A.1)

within mean-filed approximation, we get the mean-field Hamiltonian

$H=\sum_{\boldsymbol{ks}}\varepsilon_{\boldsymbol{k}}c_{\boldsymbol{ks}}^{\dagger}c_{\boldsymbol{ks}}-\frac{1}{2}\sum_{\boldsymbol{k}\alpha\beta}\Delta_{\alpha\beta}(\boldsymbol{k},\boldsymbol{q})c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\alpha}^{\dagger}c_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\beta}^{\dagger}-\frac{1}{2}\sum_{\boldsymbol{k}\alpha\beta}\Delta_{\alpha\beta}^{*}(\boldsymbol{k},\boldsymbol{q})c_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\beta}c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\alpha}$
(A.2)

where we have defined the gap function as

$\Delta_{\alpha\beta}(\boldsymbol{k},\boldsymbol{q})\equiv-\sum_{\boldsymbol{k}^{\prime}\alpha^{\prime}\beta^{\prime}}V_{\alpha\beta\alpha^{\prime}\beta^{\prime}}\left(\boldsymbol{k},\boldsymbol{k}^{\prime}\right)\left\langle
c_{-\boldsymbol{k}^{\prime}+\frac{\boldsymbol{q}}{2},\beta^{\prime}}c_{\boldsymbol{k}^{\prime}+\frac{\boldsymbol{q}}{2},\alpha^{\prime}}\right\rangle$
(A.3)

and then

$\Delta_{\alpha\beta}^{*}(\boldsymbol{k},\boldsymbol{q})=-\sum_{\boldsymbol{k}^{\prime}\alpha^{\prime}\beta^{\prime}}\left[V_{\alpha\beta\alpha^{\prime}\beta^{\prime}}\left(\boldsymbol{k},\boldsymbol{k}^{\prime}\right)\right]^{*}\left\langle
c_{-\boldsymbol{k}^{\prime}+\frac{\boldsymbol{q}}{2},\beta^{\prime}}c_{\boldsymbol{k}^{\prime}+\frac{\boldsymbol{q}}{2},\alpha^{\prime}}\right\rangle^{*}=-\sum_{\boldsymbol{k}^{\prime}\alpha^{\prime}\beta^{\prime}}V_{\alpha^{\prime}\beta^{\prime}\alpha\beta}\left(\boldsymbol{k}^{\prime},\boldsymbol{k}\right)\left\langle
c_{\boldsymbol{k}^{\prime}+\frac{\boldsymbol{q}}{2},\alpha^{\prime}}^{\dagger}c_{-\boldsymbol{k}^{\prime}+\frac{\boldsymbol{q}}{2},\beta^{\prime}}^{\dagger}\right\rangle$
(A.4)

where $\left\langle...\right\rangle$ stands for a grand canonical ensemble
average. Introduce the finite-temperature Green’s function

$G_{\alpha\beta}(\boldsymbol{k},\boldsymbol{q},\tau)\equiv-\left\langle
T_{\tau}c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\alpha}(\tau)c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\beta}^{\dagger}(0)\right\rangle$
(A.5)

The operator $A(\tau)$ denotes $e^{H\tau/\hbar}Ae^{-H\tau/\hbar}$, $A(0)$ is
for $\tau=0$, where $H$ is the Hamiltonian, $\tau$ is the ‘imaginary time’
which is treated as a real quantity. $T_{\tau}$ is the time ordering operator
with respect to the imaginary time. Using the Heaviside step function

$\theta(x)=\begin{cases}0&x<0\\\ 1&x>0\end{cases}$ (A.6)

we can express the Green’s function as

$G_{\alpha\beta}(\boldsymbol{k},\boldsymbol{q},\tau)=-\theta(\tau)\left\langle
c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\alpha}(\tau)c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\beta}^{\dagger}(0)\right\rangle+\theta(-\tau)\left\langle
c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\beta}^{\dagger}(0)c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\alpha}(\tau)\right\rangle$
(A.7)

Using
$\frac{\partial}{\partial\tau}G_{\alpha\beta}(\boldsymbol{k},\boldsymbol{q},\tau)$,
we get the equation of motion for the Green’s function

$\left(\frac{\partial}{\partial\tau}+\frac{\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}}{\hbar}\right)G_{\alpha\beta}(\boldsymbol{k},\boldsymbol{q},\tau)+\frac{1}{\hbar}\sum_{\alpha^{\prime\prime}}\Delta_{\alpha\alpha^{\prime\prime}}(\boldsymbol{k},\boldsymbol{q})\left\langle
T_{\tau}c_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\alpha^{\prime\prime}}^{\dagger}(\tau)c_{\boldsymbol{k}+\frac{\mathfrak{\boldsymbol{q}}}{2},\beta}^{\dagger}(0)\right\rangle=-\delta(\tau)\delta_{\alpha\beta}$
(A.8)

We can define the anomalous Green’s function

$F_{\alpha\beta}^{\dagger}(\boldsymbol{k},\boldsymbol{q},\tau)\equiv-\left\langle
T_{\tau}c_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\alpha}^{\dagger}(\tau)c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\beta}^{\dagger}(0)\right\rangle$
(A.9)

Similarly, using
$\frac{\partial}{\partial\tau}F_{\alpha\beta}^{\dagger}(\boldsymbol{k},\boldsymbol{q},\tau)$,
we get the equation of motion for the anomalous Green’s function

$\left(\frac{\partial}{\partial\tau}-\frac{\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}}{\hbar}\right)F_{\alpha\beta}^{\dagger}(\boldsymbol{k},\boldsymbol{q},\tau)-\frac{1}{\hbar}\sum_{\alpha^{\prime\prime}}\Delta_{\alpha^{\prime\prime}\alpha}^{*}(\boldsymbol{k},\boldsymbol{q})G_{\alpha^{\prime\prime}\beta}(\boldsymbol{k},\boldsymbol{q},\tau)=0$
(A.10)

Transforming the Green’s function and the anomalous Green’s function from the
imaginary-time-momentum space to the frequency-momentum space, we get the
Gor’kov equations

$\left(-i\omega_{n}+\frac{\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}}{\hbar}\right)G_{\alpha\beta}\left(\boldsymbol{k},\boldsymbol{q},\omega_{n}\right)-\frac{1}{\hbar}\sum_{\alpha^{\prime\prime}}\Delta_{\alpha\alpha^{\prime\prime}}(\boldsymbol{k},\boldsymbol{q})F_{\alpha^{\prime\prime}\beta}^{\dagger}\left(\boldsymbol{k},\boldsymbol{q},\omega_{n}\right)=-\delta_{\alpha\beta}$
(A.11)
$\left(-i\omega_{n}-\frac{\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}}{\hbar}\right)F_{\alpha\beta}^{\dagger}\left(\boldsymbol{k},\boldsymbol{q},\omega_{n}\right)-\frac{1}{\hbar}\sum_{\alpha^{\prime\prime}}\Delta_{\alpha^{\prime\prime}\alpha}^{*}(\boldsymbol{k},\boldsymbol{q})G_{\alpha^{\prime\prime}\beta}\left(\boldsymbol{k},\boldsymbol{q},\omega_{n}\right)=0$
(A.12)

where $\delta_{\alpha\beta}$ is the Kronecker delta. We can define a second
anomalous Green’s function

$F_{\alpha\beta}(\boldsymbol{k},\boldsymbol{q},\tau)\equiv-\left\langle
T_{\tau}c_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\alpha}(\tau)c_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2},\beta}(0)\right\rangle$
(A.13)

Combining the definition of the gap function, the expression of the complex
conjugate of the gap function, and the definitions of the anomalous Green’s
function and the second anomalous Green’s function in imaginary-time space and
their transformations in frequency space, we can get a relationship

$F_{\beta^{\prime}\alpha^{\prime}}^{\dagger}\left(\boldsymbol{k}^{\prime},\boldsymbol{q},\omega_{n}\right)=\left[F_{\alpha^{\prime}\beta^{\prime}}\left(\boldsymbol{k}^{\prime},\boldsymbol{q},\omega_{n}\right)\right]^{*}$
(A.14)

Then, we can replace the anomalous Green’s functions by the second anomalous
Green’s functions in the Gor’kov equations to get

$\left(-i\omega_{n}+\frac{\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}}{\hbar}\right)G_{\alpha\beta}\left(\boldsymbol{k},\boldsymbol{q},\omega_{n}\right)-\frac{1}{\hbar}\sum_{\alpha^{\prime\prime}}\Delta_{\alpha\alpha^{\prime\prime}}(\boldsymbol{k},\boldsymbol{q})\left[F_{\beta\alpha^{\prime\prime}}\left(\boldsymbol{k},\boldsymbol{q},\omega_{n}\right)\right]^{*}=-\delta_{\alpha\beta}$
(A.15)
$\left(-i\omega_{n}-\frac{\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}}{\hbar}\right)\left[F_{\beta\alpha}\left(\boldsymbol{k},\boldsymbol{q},\omega_{n}\right)\right]^{*}-\frac{1}{\hbar}\sum_{\alpha^{\prime\prime}}\Delta_{\alpha^{\prime\prime}\alpha}^{*}(\boldsymbol{k},\boldsymbol{q})G_{\alpha^{\prime\prime}\beta}\left(\boldsymbol{k},\boldsymbol{q},\omega_{n}\right)=0$
(A.16)

The Gor’kov equations can be written in a matrix form

$\left(-i\hbar\omega_{n}+\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}\right)\hat{G}-\hat{\Delta}\hat{F}^{\dagger}=-\hbar\sigma_{0}$
(A.17)
$\left(-i\hbar\omega_{n}-\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}\right)\hat{F}^{\dagger}-\hat{\Delta}^{\dagger}\hat{G}=0$
(A.18)

where $\hat{\Delta}$ is the gap matrix and $\sigma_{0}$ is the identity
matrix. Solving this system of matrix equations, we obtain the Green’s
function and the anomalous Green’s function as

$\hat{G}=-\hbar\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)\left[\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)\left(\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}-i\hbar\omega_{n}\right)\sigma_{0}+\hat{\Delta}\hat{\Delta}^{\dagger}\right]^{-1}$
(A.19)
$\hat{F}=\hbar\left[\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}-i\hbar\omega_{n}\right)\left(\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)\sigma_{0}+\hat{\Delta}\hat{\Delta}^{\dagger}\right]^{-1}\hat{\Delta}$
(A.20)

For singlet pairing, the gap matrix
$\hat{\Delta}(\boldsymbol{k},\boldsymbol{q})$ is anti-symmetric which can be
characterized by a single even function $\psi(\boldsymbol{k},\boldsymbol{q})$,

$\hat{\Delta}(\boldsymbol{k},\boldsymbol{q})=\psi(\boldsymbol{k},\boldsymbol{q})i\sigma_{2}=\left[\begin{array}[]{cc}0&\psi(\boldsymbol{k},\boldsymbol{q})\\\
-\psi(\boldsymbol{k},\boldsymbol{q})&0\end{array}\right]$ (A.21)

Substituting $\hat{\Delta}\hat{\Delta}^{\dagger}=|\psi|^{2}\sigma_{0}$ into
$\hat{G}$ and $\hat{F}$, we get the explicit forms of the Green’s function and
the anomalous Green’s function

$\hat{G}=-\frac{\hbar\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)\sigma_{0}}{\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)\left(\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}-i\hbar\omega_{n}\right)+|\psi(\boldsymbol{k},\boldsymbol{q})|^{2}}$
(A.22)
$\hat{F}=\frac{\hbar\hat{\Delta}}{\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}-i\hbar\omega_{n}\right)\left(\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)+|\psi(\boldsymbol{k},\boldsymbol{q})|^{2}}$
(A.23)

For triplet pairing, the gap matrix
$\hat{\Delta}(\boldsymbol{k},\boldsymbol{q})$ is symmetric which can be
parametrized by an odd vector function
$\boldsymbol{d}(\boldsymbol{k},\boldsymbol{q})$,

$\hat{\Delta}(\boldsymbol{k},\boldsymbol{q})=\boldsymbol{d}(\boldsymbol{k},\boldsymbol{q})\cdot\boldsymbol{\sigma}i\sigma_{2}=\left[\begin{array}[]{cc}-d_{x}(\boldsymbol{k},\boldsymbol{q})+\boldsymbol{i}d_{y}(\boldsymbol{k},\boldsymbol{q})&d_{z}(\boldsymbol{k},\boldsymbol{q})\\\
d_{z}(\boldsymbol{k},\boldsymbol{q})&d_{x}(\boldsymbol{k},\boldsymbol{q})+\boldsymbol{i}d_{y}(\boldsymbol{k},\boldsymbol{q})\end{array}\right]$
(A.24)

where the Pauli matrix vector
$\boldsymbol{\sigma}=\left(\sigma_{1},\sigma_{2},\sigma_{3}\right)$. Assuming
real $\boldsymbol{d}(\boldsymbol{k},\boldsymbol{q})$, we have
$\hat{\Delta}\hat{\Delta}^{\dagger}=|\boldsymbol{d}|^{2}\sigma_{0}$.
Substitute this into the solutions for $\hat{G}$ and $\hat{F}$, we get the
explicit forms for the Green’s functions

$\hat{G}=-\frac{\hbar\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)\sigma_{0}}{\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)\left(\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}-i\hbar\omega_{n}\right)+|\boldsymbol{d}(\boldsymbol{k},\boldsymbol{q})|^{2}}$
(A.25)
$\hat{F}=\frac{\hbar\hat{\Delta}}{\left(\varepsilon_{-\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}-i\hbar\omega_{n}\right)\left(\varepsilon_{\boldsymbol{k}+\frac{\boldsymbol{q}}{2}}+i\hbar\omega_{n}\right)+|\boldsymbol{d}(\boldsymbol{k},\boldsymbol{q})|^{2}}$
(A.26)


C_{1}\epsilon^{1/p},$

where $C_{1}$ is a constant that depends only on the domain $\Omega$.
Therefore we have that for $|z-z_{0}|=2\epsilon$, then $|F(z)|\geq
C_{0}C_{1}\epsilon^{1+1/p}/\operatorname{diam}(\Omega)=:d$. This implies
$F(Q)$ does not intersect a disc of radius $d$ centered at 0. On the other
hand, note that $|F(z)|\leq 2e^{2\pi/m}$ for $z\in\Omega$. Putting this all
together, we find that

$\displaystyle\mathcal{E}_{1}^{\bullet}(\theta)$ $\displaystyle\leq
10\sum_{Q:\text{dist}(Q,z_{0})\geq 2\epsilon}\int_{F(Q)}\frac{1}{|z|^{2}}\,dA$
$\displaystyle\leq
10\int_{0}^{2\pi}\int_{d}^{2e^{2\pi/m}}\frac{1}{s^{2}}s\,ds\,dt$
$\displaystyle\leq
C\log\left(\frac{\operatorname{diam}(\Omega)}{\epsilon}\right),$

for $\epsilon\leq\epsilon_{0}$, where $C$ and $\epsilon_{0}$ only depend on
$\Omega$.

∎

We need modified versions of the three results in [11, Section 7]. The logical
relationship between these results is the following:

$\text{Lemma 7.3}\implies\text{Lemma 7.2}\implies\text{Proposition 7.1}.$

Since these results make use of the fact that the boundary value problem, in
their case, is defined on the whole topological boundary, we will need to
modify them. This will have a ripple effect, causing us to update all three
results, although the proofs connecting the results will follow the logic of
[11] and we will omit some of the details.

We begin with stating the modification of [11, Proposition 7.1]:

######  Proposition 32.

Assume we are in the setting of Theorem 3, and $G_{1}=(V^{\bullet}\sqcup
V^{\circ},E)$ is a sub-orthodiagonal map of $G$ which is either equal to $G$
or as in Lemma 29. Let $x,y\in V^{\bullet}$. Assume we can choose $r,R$ so
that $r\geq\frac{1}{2}|x-y|$, $R\geq 2r+3\epsilon$, and that there is a path
from $x$ to $y$ in $G^{\bullet}\cap D_{r}(z)$, where $z=\frac{1}{2}(x+y)$.
Further, assume that either $D_{R}(z)$ does not intersect $\partial
V^{\bullet}$ or, if $G_{1}=G$, then
$R\leq\frac{1}{4}\text{{dist}}(T_{1},T_{2})$.

Let $h:V^{\bullet}\to\mathbb{R}$ be discrete harmonic on ${\rm
int}(V^{\bullet})$. Set

$\kappa=\max_{v,v^{\prime}\in\partial V^{\bullet}\cap
D_{R}(z)}|h(v)-h(v^{\prime})|,$

with $\kappa=0$, if $\partial V^{\bullet}\cap D_{R}(z)=\emptyset$. Then there
exists a universal constant $C<\infty$ such that

$|h(x)-h(y)|\leq\frac{C\,\mathcal{E}_{1}^{\bullet}(h)^{1/2}}{\log^{1/2}\left[R/(r+\epsilon)\right]}+\kappa.$

The proof of this result is the same as the proof of Proposition 7.1 in [11],
simply substituting Lemma 33 in the place of [11, Lemma 7.2]. The rough idea
is that for a lower bound on the energy, we must be able to find a thick
family of curves connecting the red dashed line to the blue dashed line in
Figure 8. This is what Lemma 33 does. For this, we introduce one more
definition from [11]. For $\rho>0$, a $\rho$-edge of $G^{\bullet}$ is an edge
$e\in E^{\bullet}$ whose dual edge $e^{\circ}\in E^{\circ}$ connects two
vertices $w,w^{\prime}\in V^{\circ}$ with $|w|<\rho\leq|w^{\prime}|$.

######  Lemma 33.

Under the assumptions of Proposition 32, suppose that the origin is located at
$\frac{1}{2}(x+y)$ and that $h(x)<h(y)$. Define

$\displaystyle X$ $\displaystyle=\\{v\in V^{\bullet}\,:\,h(v)\leq h(x)\\}$
$\displaystyle Y$ $\displaystyle=\\{v\in V^{\bullet}\,:\,h(v)\geq h(y)\\}$

and set $Y^{\prime}=Y\cup(\partial V^{\bullet}\cap D_{R}(0))$. Assume that
$X\cap(\partial V^{\bullet}\cap D_{R}(0))=\emptyset$. Then for every
$\rho\in(r+\epsilon,R-\epsilon)$, there is a path in $G^{\bullet}$ from $X$ to
$Y^{\prime}$ consisting entirely of $\rho$-edges.

This lemma will follow from the following result about plane graphs. The proof
of Lemma 33 from Lemma 34 follows the same logic that was used in [11] to
prove Lemma 7.2 from Lemma 7.3, and we omit the proof. See Figure 8 for an
illustration of the setting of Lemma 34.

######  Lemma 34.

Let $H$ be a finite connected plane graph with planar dual graph
$H^{\dagger}$. Let $W,W^{\dagger}$ be the vertex sets of $H,H^{\dagger}$. For
each $w\in W^{\dagger}$, let $f(w)$ denote the face of $H$ that contains $w$.
Assume that each inner face $f(w)$ is contained in the interior of
$D_{\epsilon}(w)$. Fix $0<r<R$ with $r+\epsilon<R-\epsilon$. Let $\zeta\in
W^{\dagger}$ be the dual vertex contained in the outer face of $H$, and assume
that $\zeta\notin D_{R}=D_{R}(0)$.

Assume that $\partial H$, the topological boundary of the outer face of $H$,
contains a cycle formed by 2 or 4 nonempty paths $S_{1},T_{1}$ or
$S_{1},T_{1},S_{2},T_{2}$. Let $x,y\in W$ satisfy $|x|,|y|\leq r$. Let $X,Y$
be subsets of $W$ such that there is a path in $H$ from $x$ to $\cup S_{k}$
whose vertices are all in $X$, and there is a path in $H$ from $y$ to $\cup
S_{k}$ whose vertices are all in $Y$. Assume there is a path in $H$ from $x$
to $y$ that is contained in $D_{r}=D_{r}(0)$. Further assume that $S_{k}\cap
D_{R}^{c}\neq\emptyset$.

Set $Y^{\prime}=Y\cup[D_{R}\cap(\cup S_{k})]$. If $X\cap
Y^{\prime}=\emptyset$, then for each $\rho\in(r+\epsilon,R-\epsilon)$, there
is a path in $H\setminus(\partial H)$ from $X$ to $Y^{\prime}$ consisting
entirely of $\rho$-edges.

$S_{1}$$S_{2}$$T_{1}$$T_{2}$$D_{R}$$x$$y$ Figure 8: A sketch of the setting
for Lemma 34: the topological boundary of $H$ subdivided into paths
$S_{1},T_{1},S_{2},T_{2}$ (solid black), the path $\gamma_{x}$ from $x$ to
$\cup S_{k}$ with vertices in $X$ (dashed red), the path $\gamma_{y}$ from $y$
to $\cup S_{k}$ with vertices in $Y$ (dashed blue), and the path $\gamma_{xy}$
from $x$ to $y$ in $D_{r}$ (dashed black).

###### Proof.

We begin by following the proof of Lemma 7.3 in [11] to obtain a cycle $K$ in
$H$ such that all of its edges are $\rho$-edges. For convenience, we briefly
repeat that argument here. Let $u\in W^{\dagger}$ be a dual vertex contained
in an inner face of $H$ with $x$ incident to this face. Then $u\neq\zeta$ and
$|u|\leq r+\epsilon<\rho$. Let $\sigma_{\rho}$ be the set of vertices $w\in
W^{\dagger}$ so that there is a path from $u$ to $w$ in $H^{\dagger}$ with all
vertices of this path in the open disk of radius $\rho$ centered at 0. Let
$H[\sigma_{\rho}]$ be the subgraph of $H$ whose vertices and edges make up the
boundaries of the faces containing $w$ for each $w\in\sigma_{\rho}$. The graph
$H[\sigma_{\rho}]$ is connected and contained in $D_{R}$. Let
$f_{\rho}^{\text{out}}$ be the outer face of $H[\sigma_{\rho}]$, which
contains the complement of $D_{R}$. Let $K$ be the boundary of
$f_{\rho}^{\text{out}}$, which is a connected subgraph. One can check that all
of the edges in $K$ are $\rho$-edges. Further, we note that the edges in $K$
are either in $\partial H$, or within distance $\epsilon$ of $C_{\rho}$, the
circle of radius $\rho$ centered at 0. The former situation can happen when
$C_{\rho}$ extends into the outer face of $H$. However, when an arc of
$C_{\rho}$ is contained entirely within the interior of $H$, then there will
be edges of $K$ within distance $\epsilon$ of this arc.

Let $\gamma_{x}$ denote the path in $H$ from $x$ to $\cup S_{k}$ whose
vertices are all in $X$, let $\gamma_{y}$ denote the path in $H$ from $y$ to
$\cup S_{k}$ whose vertices are all in $Y$, and let $\gamma_{xy}$ denote the
path in $H$ from $x$ to $y$ that is contained in $D_{r}$. Since
$\gamma_{x}\subset X$ is disjoint from $(\cup S_{k})\cap D_{R}\subset
Y^{\prime}$, we must have that $\gamma_{x}$ connects $x$ to a point in
$\partial H$ that is outside $D_{R}$. Hence $\gamma_{x}$ must cross the
annulus $D_{R}\setminus D_{r}$. For the sake of simplicity, assume that
$\gamma_{y}$ connects to a vertex in $S_{2}$ (otherwise renumber the boundary
paths.) Then $\gamma_{y}\cup S_{2}$ must also cross the annulus
$D_{R}\setminus D_{r}$, since $S_{k}\cap D_{R}^{c}\neq\emptyset$. Further
$\gamma_{x}$ and $\gamma_{y}\cup(\cup S_{k}\cap D_{R})$ must be disjoint
within $D_{R}$, since the vertices of the former are in $X$ and the vertices
of the later are in $Y^{\prime}$. Hence we can find a cycle $J$ composed of
three disjoint arcs: (1) an arc contained in $\gamma_{xy}$ within $D_{r}$, (2)
an arc is contained in $\gamma_{x}$ that crosses $D_{R}\setminus D_{r}$, and
(3) an arc contained in $\gamma_{y}\cup(\cup S_{k})\cup(\cup T_{k})$ that
crosses $D_{R}\setminus D_{r}$ and whose intersection with $D_{R}$ is
contained in $\gamma_{y}\cup(\cup S_{k})$.

Now the interior of $J$ is contained in the interior of $H$ by construction.
Further it contains an arc of $C_{\rho}$. Hence there must be a path contained
in $K$ that crosses the interior of $J$ from a vertex in $\gamma_{x}$ to a
vertex in $\gamma_{y}\cup(\cup S_{k})$. This path could potentially contain
some edges in $\partial H$, but the only way for this to happen is for the
path to contain edges from $\gamma_{x}$ or $\gamma_{y}\cup(\cup S_{k})$.
Hence, we can shorten the path by starting at the last vertex it encounters in
$\gamma_{x}$ and ending at the first vertex it encounters in
$\gamma_{y}\cup(\cup S_{k})$. Then all the edges in the path are strictly
within the interior of $J$ and hence in the interior of $H$. This gives us the
desired path.

∎

Lastly, we need the following result that allows us to convert a curve in
$\hat{G}$ to a nearby path in $G^{\bullet}$. This result follows exactly the
proof of Lemma 3.2 in [11], and so we omit the proof.

######  Lemma 35.

Let $G$ be a finite orthodiagonal map with edge lengths bounded by $\epsilon$,
let $z,w\in V^{\bullet}$, and let $\gamma_{zw}$ be a curve from $z$ to $w$ in
$\hat{G}$. Then there is a path $\gamma$ in $G^{\bullet}$ so that
$\operatorname{diam}(\gamma)\leq\operatorname{diam}(\gamma_{zw})+2\epsilon.$

###  4.2 Proof of Theorem 3

This section is devoted to the proof of Theorem 3. Following the same
structure as in [11], we begin with an outline of the overall logic and then
conclude by establishing Claims 36-40.

###### Proof of Theorem 3..

Recall that our goal is to show that (22) holds for $\lambda<C^{\prime}$ for
some $C^{\prime}>0$ which depends only on $\Omega$. The first step is to show
that the desired bound holds near the boundary arcs $A$ and $B$.

######  Claim 36.

For $\lambda$ small enough (depending only on $\Omega$), the bound (22) holds
for all $z\in V^{\bullet}\cap\overline{\Omega}$ such that dist$(z,A\cup B)\leq
2\lambda^{1/4}$.

For the rest of the proof, we will assume that there are $z\in
V^{\bullet}\cap\overline{\Omega}$ with dist$(z,A\cup B)>2\lambda^{1/4}$. Since
diam$(\Omega)=1$, this will imply that $2\lambda^{1/4}\leq 1$. This is
equivalent to $\lambda\leq\frac{1}{16}$ or
$\lambda\leq\frac{1}{8}\lambda^{1/4}$, which we will use at times. We will
also assume at times that $\lambda$ is small enough, where this depends only
on the analytic quadrilateral $\Omega$.

Next we need to construct a sub-orthodiagonal map which has distance at least
$\lambda^{1/4}$ from $A\cup B$, putting us in the setting for many of the
tools developed in the previous section. This task is accomplished in the
following construction and claim. Let $P$ be the set of all inner faces $Q$ of
$G$ such that dist$(Q,A\cup B)\geq\lambda^{1/4}$, and let $G[P]$ be the
subgraph of $G$ that is the union of the boundaries of the faces in $P$. Note
that $G[P]$ may be composed of more than one block (i.e. maximal 2-connected
component).

######  Claim 37.

The set $P$ is nonempty. The inner faces of $G[P]$ are precisely the elements
of $P$, and every point $p$ on the boundary of the outer face of $G[P]$ is
either in $T_{1}\cup T_{2}$ or satisfies

$\mathrm{dist}(p,A\cup B)\leq\frac{5}{4}\lambda^{1/4}.$ (27)

Each $z\in V^{\bullet}\cap\overline{\Omega}$ with $\mathrm{dist}(z,A\cup
B)>2\lambda^{1/4}$ is a vertex of a unique block of $G[P]$. If $H$ is a block
of $G[P]$, then the inner faces of $H$ are all elements of $P$ and every point
$p$ on the boundary of the outer face of $H$ is either in $T_{1}\cup T_{2}$ or
satisfies (27).

Claim 37 implies that each block $H$ of $G[P]$ is a sub-orthodiagonal map. Fix
$z_{0}\in V^{\bullet}$ with dist$(z_{0},A\cup B)>2\lambda^{1/4}$, and let
$G_{1}$ be the block of $G[P]$ containing $z_{0}$. Note that if the
topological boundary $\partial G_{1}$ does not intersect $T_{1}\cup T_{2}$,
then we are in the setting addressed in [11] of no free boundary, and their
work takes care of this case. Hence we will assume that $\partial G_{1}$ does
intersect $T_{1}\cup T_{2}$. We take $\partial V_{1}^{\bullet}$ to be the
points from $V_{1}^{\bullet}$ that are in $\partial
G_{1}\setminus\left(T_{1}\cup T_{2}\right)$ and ${\rm
int}(V_{1}^{\bullet})=V_{1}^{\bullet}\setminus\partial V_{1}^{\bullet}$.
Lastly, by Lemma 29, we analytically extend $\psi$ (and hence $h_{c}$) to be
defined on a domain $\mathcal{O}$ containing $\hat{G}_{1}$ so that
$\mathrm{dist}(\hat{G}_{1},\partial\mathcal{O})\geq\frac{5}{16}\lambda^{1/4}$.

Our goal is to show that

$|h_{d}(z_{0})-h_{c}(z_{0})|\leq\frac{C}{\log^{1/2}\left(1/\lambda\right)}$
(28)

for $\lambda$ small enough (depending only on $\Omega$). To accomplish this,
we need derivative bounds on the extended $h_{c}=\operatorname{Re}\psi$ on
$\hat{G}_{1}$.

######  Claim 38.

Let $\tilde{G}_{1}$ denote the union of the convex hulls of the closures of
the inner faces of $G_{1}$. Then for $\lambda$ small enough (depending on
$\Omega$),

$||\nabla h_{c}||_{\infty,\tilde{G}_{1}}\leq C\lambda^{-1/4}\;\;\;\text{ and
}\;\;\;||Hh_{c}||_{\infty,\tilde{G}_{1}}\leq C\lambda^{-1/2}.$

Define $h_{d}^{(1)}:V_{1}^{\bullet}\to\mathbb{R}$ to be the solution to the
discrete Dirichlet problem on ${\rm int}(V_{1}^{\bullet})$ with boundary data
$h_{c}$ on $\partial V_{1}^{\bullet}$. Then

$|h_{d}(z_{0})-h_{c}(z_{0})|\leq|h_{d}(z_{0})-h_{d}^{(1)}(z_{0})|+|h_{d}^{(1)}(z_{0})-h_{c}(z_{0})|.$
(29)

We will bound each of the terms in (29) separately, and we begin with the
first term.

######  Claim 39.

The first term of (29) is bounded as follows:

$|h_{d}(z_{0})-h_{d}^{(1)}(z_{0})|\leq\frac{C}{\log^{1/2}(1/\lambda)}.$

We now turn to the second term of (29) and set $f=h_{d}^{(1)}-h_{c}$, which
satisfies that $f=0$ on $\partial V^{\bullet}_{1}$. The following claim is a
key part of the proof, and it will be used to bound $|f(z_{0})|$, which is the
second term in (29).

######  Claim 40.

Let $D_{4\epsilon}$ be the disk of radius $4\epsilon$ centered at $z_{0}$. For
$\lambda$ small enough (depending only on $\Omega$), the set $D_{4\epsilon}$
is disjoint from $\partial V^{\bullet}_{1}$, and for $y\in V^{\bullet}_{1}\cap
D_{4\epsilon}$

$|f(z_{0})-f(y)|\leq\frac{C}{\log^{1/2}(1/\lambda)}.$

Proposition 30 and Claim 38 imply that

$\mathcal{E}^{\bullet}_{1}(f)\leq
C\left[M^{2}\epsilon^{2}+L^{2}\lambda\right]\leq C\lambda^{1/2},$

where $L=||\nabla h_{c}||_{\infty,\hat{G}_{1}}$ and
$M=||Hh_{c}||_{\infty,\tilde{G}_{1}}.$ Applying Proposition 31 yields a unit
flow $\theta$ in $G^{\bullet}_{1}$ from $V_{1}^{\bullet}\cap D_{4\epsilon}$ to
$\partial V_{1}^{\bullet}$ such that

$\mathcal{E}_{1}^{\bullet}(\theta)\leq C\log\left(1/\lambda\right).$

Then $-\theta$ is a unit flow from $\partial V_{1}^{\bullet}$ to
$V_{1}^{\bullet}\cap D_{4\epsilon}$. Applying Proposition 4.11 in [11] to $f$
and $-\theta$ shows that

$\text{gap}_{\partial V_{1}^{\bullet},V_{1}^{\bullet}\cap
D_{4\epsilon}}(f)\leq\mathcal{E}_{1}^{\bullet}(-\theta)^{1/2}\mathcal{E}_{1}^{\bullet}(f)^{1/2}\leq
C\lambda^{1/4}\log^{1/2}\left(1/\lambda\right),$

where

$\displaystyle\text{gap}_{X,Y}(f):=\min_{y\in Y}f(y)-\max_{x\in X}f(x).$ (30)

Note that Proposition 4.11 only applies in the case that the gap is
nonnegative, but the inequality is immediate in the other case. Since $f=0$ on
$\partial V_{1}^{\bullet}$, there is $y\in V_{1}^{\bullet}\cap D_{4\epsilon}$
so that

$f(y)\leq C\lambda^{1/4}\log^{1/2}(1/\lambda).$

Then combining this with Claim 40, we obtain

$f(z_{0})=f(y)+[f(z_{0})-f(y)]\leq
C\lambda^{1/4}\log^{1/2}(1/\lambda)+\frac{C}{\log^{1/2}(1/\lambda)}\leq\frac{C}{\log^{1/2}(1/\lambda)}$

Applying the same argument to $-f$ yields that

$|f(z_{0})|\leq\frac{C}{\log^{1/2}(1/\lambda)}.$

All that remains is to establish Claims 36-40. Two of these, Claims 37 and 39,
follow the arguments of the corresponding results, Claims 8.2 and 8.4, in [11]
with only minor changes of constants. As a result, we omit those proofs here.

###### Proof of Claim 36..

Let $z\in V^{\bullet}\cap\overline{\Omega}$ such that dist$(z,A\cup B)\leq
2\lambda^{1/4}$. Choose $u\in A\cup B$ such that $|z-u|\leq 2\lambda^{1/4}$.
Then there is $v\in S_{1}\cup S_{2}$ such that there is a path $\gamma_{uv}$
from $u$ to $v$ in $\hat{G}$ with length at most $\delta$. Further by
following the edge containing $v$ to its endpoint in $V^{\bullet}$, there is
$w\in\partial V^{\bullet}$ such that $|v-w|\leq\epsilon$ (and let
$\gamma_{vw}$ be the path from $v$ to $w$ along the edge). Then $|z-w|\leq
2\lambda^{1/4}+\delta+\epsilon\leq 4\lambda^{1/4}$.

Since $h_{d}(w)=h_{c}(u)$, we have

$|h_{d}(z)-h_{c}(z)|\leq|h_{d}(z)-h_{d}(w)|+|h_{c}(u)-h_{c}(z)|.$

From Proposition 18,

$|h_{c}(u)-h_{c}(z)|\leq C|u-z|^{1/4}\leq
C\lambda^{1/16}\leq\frac{C}{\log^{1/2}\left(1/\lambda\right)}.$

Thus it remains to bound $|h_{d}(z)-h_{d}(w)|$.

Our goal will be to use Proposition 32, and so we must find a path from $z$ to
$w$ in $G^{\bullet}$. Since $\Omega$ is a linearly connected domain (with
constant $K$), there is a path $\gamma_{zu}$ in
$\overline{\Omega}\subset\hat{G}$ from $z$ to $u$ so that

$\operatorname{diam}(\gamma_{zu})\leq K|z-u|\leq 2K\lambda^{1/4}.$

Therefore $\gamma_{zw}=\gamma_{zu}\cup\gamma_{uv}\cup\gamma_{vw}$ is a path
from $z$ to $w$ in $\hat{G}$ with

$\operatorname{diam}(\gamma_{zw})\leq
2K\lambda^{1/4}+\delta+\epsilon\leq(2K+2)\lambda^{1/4}.$

Now by Lemma 35 the path $\gamma_{zw}$ can be converted into a path $\gamma$
in $G^{\bullet}$ so that

$\operatorname{diam}(\gamma)\leq\operatorname{diam}(\gamma_{zu})+2\epsilon\leq(2K+4)\lambda^{1/4}.$

Now a disc centered at $\frac{1}{2}(z+w)$ will contain $\gamma$ if its radius
is at least

$\frac{1}{2}|z-w|+\operatorname{diam}(\gamma)\leq
2\lambda^{1/4}+(2K+4)\lambda^{1/4}=(2K+6)\lambda^{1/4}.$

Set $r=(2K+6)\lambda^{1/4}$ and $R=(4K+15)\lambda^{1/8}$. Note that, by taking
$\lambda$ small enough,

$2r+3\epsilon\leq R\leq\frac{1}{4}\hat{d},$

where
$\hat{d}:=\min\\{\text{{dist}}(T_{1},T_{2}),\text{{dist}}(S_{1},S_{2})\\}$ and

$\frac{R}{r+\epsilon}=\frac{(4K+15)\lambda^{1/8}}{(2K+6)\lambda^{1/4}+\epsilon}\geq\lambda^{-1/8}.$

The fact that $R\leq\frac{1}{4}\hat{d}$ means that $D_{R}$, the disc of radius
$R$ centered at $\frac{1}{2}|z+w|$ can intersect at most one of $S_{1}$ or
$S_{2}$. Because $h_{d}$ is constant equal to 0 on $S_{1}\cap V^{\bullet}$ and
constant equal to 1 on $S_{2}\cap V^{\bullet}$, the term $\kappa$ from
Proposition 32 will be zero. Therefore Proposition 32 gives

$|h_{d}(z)-h_{d}(w)|\leq\frac{C\,\mathcal{E}^{\bullet}(h_{d})^{1/2}}{\log^{1/2}\left(1/\lambda\right)}.$

To complete the proof, it remains to bound $\mathcal{E}^{\bullet}(h_{d})$. We
reflect $\psi$ over $\tau_{1}$ and $\tau_{2}$, so that $\psi$ is conformal in
an open set $\mathcal{O}$ containing $\Omega$, as before. Define
$\mathcal{O}_{\text{mid}}=\\{z\in\mathcal{O}\,:\,1/3\leq\psi(z)\leq 2/3\\}$.
Consider the distance between $\tau_{k}\cap\mathcal{O}_{\text{mid}}$ and
$\partial\mathcal{O}$. By taking $\lambda$ small enough, we may assume that
$3\lambda<\text{dist}\left(\tau_{k}\cap\mathcal{O}_{\text{mid}},\,\partial\mathcal{O}\right)$.
Then the boundary of any connected component of
$\hat{G}\setminus\mathcal{O}_{\text{mid}}$ will intersect exactly one of
$\psi^{-1}(\\{\text{Re}(z)=1/3\\})$ or $\psi^{-1}(\\{\text{Re}(z)=2/3\\})$.
Let $\hat{G}_{\text{left}}$ be the union of all components whose boundary
intersects $\psi^{-1}(\\{\text{Re}(z)=1/3\\})$, and similarly define
$\hat{G}_{\text{right}}$ as the union of components whose boundary intersects
$\psi^{-1}(\\{\text{Re}(z)=2/3\\})$. Define the continuous function $f$ on
$\hat{G}$ as follows:

$f(z)=\begin{cases}0&\quad\text{if }z\in\hat{G}_{\text{left}}\\\
1&\quad\text{if }z\in\hat{G}_{\text{right}}\\\ 3h_{c}(z)-1&\quad\text{else}\\\
\end{cases}.$

Then by Proposition 4.8 in [11],

$\mathcal{E}^{\bullet}(h_{d})\leq\mathcal{E}^{\bullet}(f).$

Form $G_{\text{mid}}$ from all the quadrilaterals that intersect
$\mathcal{O}_{\text{mid}}$. Note that the only quadrilaterals in $G$ that
contribute to $\mathcal{E}^{\bullet}(f)$ are those in $G_{\text{mid}}$.
Further, for these $Q=[v_{1},w_{1},v_{2},w_{2}]$ in $G_{\text{mid}}$, we have
that

$|f(v_{1})-f(v_{2})|\leq 3\left|h_{c}(v_{1})-h_{c}(v_{2})\right|.$

Therefore

$\mathcal{E}^{\bullet}(f)\leq 9\mathcal{E}_{\text{mid}}^{\bullet}(h_{c}),$

where the energy on the right is with respect to the graph $G_{\text{mid}}$.
Since $G_{\text{mid}}$ is contained within $\mathcal{O}$, Proposition 5.2 in
[11] implies that

$\mathcal{E}_{\text{mid}}^{\bullet}(h_{c})\leq
2\int_{\hat{G}_{\text{mid}}}|\nabla
h_{c}|^{2}\,dA+2\,\text{area}(\hat{G}_{\text{mid}})(10LM\epsilon+8M^{2}\epsilon^{2}),$

where $L=||\nabla h_{c}||_{\infty,\hat{G}_{\text{mid}}}$ and
$M=||Hh_{c}||_{\infty,\tilde{G}_{\text{mid}}}$. The constants $L$ and $M$
depend on $\Omega$ but are independent of $\epsilon$ and $\delta$. Since
$\mathcal{O}$ and the extended $\psi$ are defined by reflection,
$\psi(\mathcal{O})$ is contained in the rectangle $[0,1]\times[-m,2m]$ and
$\int_{\hat{G}_{\text{mid}}}|\nabla h_{c}|^{2}\,dA\leq 3\int_{\Omega}|\nabla
h_{c}|^{2}\,dA$. However, $\int_{\Omega}|\nabla h_{c}|^{2}\,dA=m$ the modulus
of curves in $R$ from $\\{\text{Re}(z)=0\\}$ to $\\{\text{Re}(z)=1\\}$. Thus
we have shown that $\mathcal{E}^{\bullet}(h_{d})$ is bounded by a constant
(which depends only on $\Omega$), and hence

$|h_{d}(z)-h_{d}(w)|\leq\frac{C}{\log^{1/2}\left(1/\lambda\right)}.$ (31)

∎

In the proof of Claim 36, we have established the following, which we will use
in proving Theorem 4

######  Lemma 41.

Let $z\in V^{\bullet}\cap\overline{\Omega}$ such that dist$(z,A\cup B)\leq
2\lambda^{1/4}$. Then there exists $w\in\partial V^{\bullet}$ so that
$|h_{d}(z)-h_{d}(w)|\leq\frac{C}{\log^{1/2}\left(1/\lambda\right)}$.

###### Proof of Claim 38..

Note that $\sup_{y\in\Omega}|h_{c}(y)|\leq 1$, since $h_{c}$ is the real part
of the conformal map $\psi$ onto $[0,1]\times[0,m]$. For all
$q\in\tilde{G}_{1}$, Lemma 29 implies that for $\lambda$ small enough

$\text{dist}(q,\partial\mathcal{O})\geq\frac{5}{16}\lambda^{1/4}-\epsilon\geq\frac{3}{16}\lambda^{1/4}.$
(32)

Proposition 3.7 of [11] and (32) imply that

$||\nabla h_{c}||_{\infty,\tilde{G}_{1}}\leq C\lambda^{-1/4}.$

We obtain the Hessian bound $||Hh_{c}||_{\infty,\tilde{G}_{1}}\leq
C\lambda^{-1/2}$ from Proposition 3.7 of [11] and (32) exactly as in the proof
of Claim 8.3 in [11].

∎

###### Proof of Claim 40..

Recall that we have $z_{0}\in V_{1}^{\bullet}$ with dist$(z_{0},A\cup
B)>2\lambda^{1/4}$. From Claim 37, we know that every point $p\in\partial
G_{1}$ satisfies dist$(p,A\cup B)\leq\frac{5}{4}\lambda^{1/4}$. Thus the
distance from $z_{0}$ to $\partial V_{1}^{\bullet}$ is at least
$2\lambda^{1/4}-\frac{5}{4}\lambda^{1/4}=\frac{3}{4}\lambda^{1/4}\geq
6\lambda$, since $\lambda^{1/4}\geq 8\lambda$, and this implies that
$D_{4\epsilon}$ cannot intersect $\partial V_{1}^{\bullet}$.

Let $y\in V^{\bullet}_{1}\cap D_{4\epsilon}$. Then

$|f(z_{0})-f(y)|\leq|h_{d}^{(1)}(z_{0})-h_{d}^{(1)}(y)|+|h_{c}(z_{0})-h_{c}(y)|.$

Our next step is to get a path from $z_{0}$ to $y$ in $\hat{G}_{1}$ with
length bounded by $C\lambda$. Then by Claim 38, we will have that

$|h_{c}(z_{0})-h_{c}(y)|\leq||\nabla h_{c}||_{\infty,\hat{G}_{1}}\cdot
C\lambda\leq C\lambda^{3/4}.$

Since $y$ may not be in $\overline{\Omega}$, we will choose $w$ in
$\overline{\Omega}$ so that there is a path $\gamma_{wy}$ from $w$ to $y$ in
$\hat{G}_{1}$ with the length of $\gamma_{wy}$ bounded by $\delta$. Since
$\Omega$ is a linearly connected domain (with constant $K$), there is a path
$\gamma_{z_{0}w}$ in $\overline{\Omega}\subset\hat{G}_{1}$ from $z_{0}$ to $w$
so that

$\operatorname{diam}(\gamma_{z_{0}w})\leq K|z_{0}-w|\leq
K\left(|z_{0}-y|+\delta\right)\leq 4K\lambda.$

Therefore $\gamma_{z_{0}y}=\gamma_{z_{0}w}\cup\gamma_{wy}$ is a path from
$z_{0}$ to $y$ in $\hat{G_{1}}$ with

$\operatorname{diam}(\gamma_{z_{0}y})\leq(4K+1)\lambda.$

It remains to bound $|h_{d}^{(1)}(z_{0})-h_{d}^{(1)}(y)|$. Propositions 4.8
and 5.1 of [11] give

$\mathcal{E}^{\bullet}_{1}(h_{d}^{(1)})\leq\mathcal{E}^{\bullet}_{1}(h_{c})\leq
2\int_{\hat{G}_{1}}|\nabla
h_{c}|^{2}\,dA+2\text{area}(\hat{G}_{1})\,(10LM\epsilon+8M^{2}\epsilon^{2})$

where (using Claim 38)

$\displaystyle L$ $\displaystyle=||\nabla h_{c}||_{\infty,\hat{G}_{1}}\leq
C\lambda^{-1/4},$ $\displaystyle M$
$\displaystyle=||Hh_{c}||_{\infty,\tilde{G}_{1}}\leq C\lambda^{-1/2}.$

Since $\mathcal{O}$ and the extended $\psi$ are defined by reflection,
$\psi(\mathcal{O})$ is contained in the rectangle $[0,1]\times[-m,2m]$, and so

$\int_{\hat{G}_{1}}|\nabla h_{c}|^{2}\,dA\leq 3\int_{\Omega}|\nabla
h_{c}|^{2}\,dA=3m,$

where $m$ is the continuous modulus of the family of curves from $A$ to $B$ in
$\Omega$. Thus,

$\mathcal{E}^{\bullet}_{1}(h_{d}^{(1)})\leq C.$

We will utilize Proposition 32 to finish the proof, and so we must verify that
its hypotheses are satisfied. By Lemma 35, the path $\gamma_{z_{0}y}$ from
above can be converted into a path $\gamma$ in $G^{\bullet}_{1}$ so that

$\operatorname{diam}(\gamma)\leq\operatorname{diam}(\gamma_{z_{0}y})+2\epsilon\leq(4K+3)\lambda.$

Thus a disc centered at $\frac{1}{2}(z_{0}+y)$ will contain $\gamma$ if its
radius is larger than

$\frac{1}{2}|z_{0}-y|+\operatorname{diam}(\gamma)\leq\frac{3}{2}\epsilon+(4K+3)\lambda\leq(4K+5)\lambda.$

Set $r=(4K+5)\lambda$ and $R=\frac{1}{4}\lambda^{1/4}$. We have chosen
$\lambda$ to be small enough so that $2r+3\epsilon\leq R$. Note further that
since the distance from $z_{0}$ to $\partial V_{1}^{\bullet}$ is at least
$\frac{3}{4}\lambda^{1/4}$, the distance from $\frac{1}{2}(z_{0}+y)$ to
$\partial V_{1}^{\bullet}$ is at least
$\frac{3}{4}\lambda^{1/4}-\frac{3}{2}\epsilon\geq\frac{1}{2}\lambda^{1/4}$
(using that $\lambda\leq\frac{1}{8}\lambda^{1/4}$). Thus $D_{R}$ does not
intersect $\partial V_{1}^{\bullet}$ and so the term $\kappa$ in Proposition
32 is equal to zero. Therefore Proposition 32 implies that

$|h_{d}^{(1)}(z_{0})-h_{d}^{(1)}(y)|\leq\frac{C\,\mathcal{E}_{1}^{\bullet}(h_{d}^{(1)})^{1/2}}{\log^{1/2}\left[R/(r+\epsilon)\right]}\leq\frac{C}{\log^{1/2}\left(1/\lambda\right)}.$

∎

This completes the proof of Theorem 3. ∎

Note that by establishing (28) in the proof of Theorem 3, we have proven the
following result, which we formulate as a proposition to be able to reference
elsewhere.

######  Proposition 42.

Let $G_{1}$ be a sub-orthodiagonal map of $G$ so that
$\mathrm{dist}(G_{1},A\cup B)\geq\lambda^{1/4}$ and let $z_{0}$ be a vertex in
$G_{1}$. Then
$|h_{d}(z_{0})-h_{c}(z_{0})|\leq\frac{C}{\log^{1/2}\left(1/\lambda\right)}$
for $\lambda$ small enough (depending only on $\Omega$).

###  4.3 Proof of Theorem 4

###### Proof of Theorem 4..

Let $\Omega$ be an analytic quadrilateral, as defined in Definition 16, with
boundary arcs $A,\tau_{1},B,\tau_{2}$ listed counterclockwise, and let
$\Gamma$ be the curve family between $A$ and $B$ in $\Omega$. Let
$\epsilon_{n}>0$ with $\epsilon_{n}\to 0$. Let $G_{n}$ be an orthodiagonal map
with boundary arcs $S^{n}_{1},T^{n}_{1},S^{n}_{2},T^{n}_{2}$, with edge-length
at most $\epsilon_{n}$, and so that
$(G_{n};S^{n}_{1},T^{n}_{1},S^{n}_{2},T^{n}_{2})$ approximates
$(\Omega;A,\tau_{1},B,\tau_{2})$ within distance $\epsilon_{n}$ (in the sense
of Definition 17). Equip $G_{n}$ with canonical edge-weights $\sigma_{n}$ as
in Equation (9). Consider the family of all paths in $G^{\bullet}_{n}$ from
$A_{n}^{\bullet}=V_{n}^{\bullet}\cap S^{n}_{1}$ to
$B_{n}^{\bullet}=V_{n}^{\bullet}\cap S^{n}_{2}$, and let $\Gamma_{n}$ be the
unique minimal subfamily with non-crossing paths.

Let $\psi$ be a conformal map of $\Omega$ onto the rectangle
$R=[0,1]\times[0,m]$ taking the arcs $A$ and $B$ to the vertical sides of the
rectangle, and let $h=\text{Re }\psi$. Let
$h_{n}^{\bullet}:V_{n}^{\bullet}\to\mathbb{R}$ be discrete harmonic on ${\rm
int}(V_{n}^{\bullet})$, with the following boundary values:
$h_{n}^{\bullet}(z)=0$ for $z\in A_{n}$ and $h_{n}^{\bullet}(z)=1$ for $z\in
B_{n}$. Recall that

$\operatorname{Mod}_{2,\sigma_{n}}(\Gamma_{n})=\mathcal{E}_{n}^{\bullet}(h_{n}^{\bullet}).$

Similarly, the continuous modulus of $\Gamma$ equals the energy of $h$:

$\operatorname{Mod}_{2}(\Gamma)=\mathcal{E}(h)=\int_{\Omega}|\nabla
h|^{2}\,dA.$

Our goal is to show that $\operatorname{Mod}_{2,\sigma_{n}}(\Gamma_{n})$ tends
to $m:=\operatorname{Mod}_{2}(\Gamma)$ as $n$ tends to infinity. In order to
do this we introduce four other harmonic functions. First, in the continuous
case, we let $\tilde{h}$ be the harmonic conjugate of $h$, i.e.

$\tilde{h}:=\operatorname{Im}\psi.$

On the other hand, in the discrete case, we introduce three more discrete
harmonic functions

$h_{n}^{\circ}:V_{n}^{\circ}\to\mathbb{R},\;\;\tilde{h}_{n}^{\bullet}:V_{n}^{\bullet}\to\mathbb{R},\;\text{
and }\;\tilde{h}_{n}^{\circ}:V_{n}^{\circ}\to\mathbb{R}$

that satisfy the following boundary values:

$\displaystyle h_{n}^{\circ}(z)$ $\displaystyle=0\;\text{ for }\;z\in
V_{n}^{\circ}\cap S_{1}^{n}\;\text{ and }\;h_{n}^{\circ}(z)=1\;\text{ for
}z\in V_{n}^{\circ}\cap S_{2}^{n},$ $\displaystyle\tilde{h}_{n}^{\bullet}(z)$
$\displaystyle=0\;\text{ for }\;z\in V_{n}^{\bullet}\cap T_{1}^{n}\;\text{ and
}\;\tilde{h}_{n}^{\bullet}(z)=1\;\text{ for }\;z\in V_{n}^{\bullet}\cap
T_{2}^{n},$ $\displaystyle\tilde{h}_{n}^{\circ}(z)$ $\displaystyle=0\;\text{
for }\;z\in V_{n}^{\circ}\cap T_{1}^{n}\;\text{ and
}\;\tilde{h}_{n}^{\circ}(z)=1\;\text{ for }\;z\in V_{n}^{\circ}\cap
T_{2}^{n}.$

As in the proof of Theorem 3, let $P_{n}$ be the set of all inner faces $Q$ of
$G_{n}$ such that dist$(Q,A\cup B)\geq\epsilon_{n}^{1/4}$, and let
$G_{n,1}=G[P_{n}]$ be the subgraph of $G_{n}$ that is the union of the
boundaries of the faces in $P_{n}$. By Lemma 29, when $\epsilon_{n}$ is small
enough, we can analytically extend $\psi$, so that $\psi$ is defined on all
the vertices of $G_{n,1}$. We define $G_{n,2}$ in an analogous way by
switching the roles of $A,B$ and $\tau_{1},\tau_{2}$. For $k=1$ or $2$, the
notation $\mathcal{E}_{n,k}^{\bullet}$ (respectively
$\mathcal{E}_{n,k}^{\circ}$) refers to the energy with respect to the graph
$G_{n,k}^{\bullet}$ (respectively $G_{n,k}^{\circ}$).

######  Claim 43.

There exists $r_{n}\searrow 1$ so that

$\mathcal{E}_{n}^{\bullet}(h_{n}^{\bullet})\leq
r_{n}\,\mathcal{E}_{n,1}^{\bullet}(h).$ (33)

We assume Claim 43 for now and return to its proof later. Note that Claim 43
implies that there exists $r_{n}\searrow 1$ so that

$\mathcal{E}_{n}^{\circ}(h_{n}^{\circ})\leq
r_{n}\,\mathcal{E}_{n,1}^{\circ}(h),\;\;\mathcal{E}_{n}^{\bullet}(\tilde{h}_{n}^{\bullet})\leq
r_{n}\,\mathcal{E}_{n,2}^{\bullet}(\tilde{h}/m),\;\text{and
}\;\mathcal{E}_{n}^{\circ}(\tilde{h}_{n}^{\circ})\leq
r_{n}\,\mathcal{E}_{n,2}^{\circ}(\tilde{h}/m).$ (34)

Next we need some properties that follow from the Cauchy-Riemann equations for
$\psi$.

######  Claim 44.

Fix a face $Q=[v_{1},w_{1},v_{2},w_{2}]$ that is an inner face of $G_{n,1}$ or
$G_{n,2}$, and define

$\zeta_{Q}=\frac{(h(w_{2})-h(w_{1}))+i(\tilde{h}(w_{2})-\tilde{h}(w_{1}))}{|w_{2}-w_{1}|}\;\;\text{and}\;\;\eta_{Q}=\frac{(h(v_{2})-h(v_{1}))+i(\tilde{h}(v_{2})-\tilde{h}(v_{1}))}{|v_{2}-v_{1}|}.$

Then for $\epsilon_{n}$ small enough

  1. (i)

$|\zeta_{Q}-i\eta_{Q}|\leq C\epsilon_{n}^{1/2}$.

  2. (ii)

$\operatorname{Im}(\zeta_{Q}\bar{\eta}_{Q})\operatorname{Area}(Q)\leq\int_{Q}|\nabla
h|^{2}\,dA+C\epsilon_{n}^{1/4}\operatorname{Area}{Q}$.

Again, we will postpone the proof of this claim until the end, and we will
assume that we have taken $\epsilon_{n}$ small enough for its conclusion to
hold. By Claim 44(i), summing over all $Q$ that are inner faces of either
$G_{n,1}$ or $G_{n,2}$ gives

$\sum_{Q}2\operatorname{Area}(Q)|\zeta_{Q}-i\eta_{Q}|^{2}\leq C\epsilon_{n},$

since the area of $\Omega$ is finite. Expanding
$|\zeta_{Q}-i\eta_{Q}|^{2}=|\zeta_{Q}|^{2}+|\eta_{Q}|^{2}-2\operatorname{Im}(\zeta_{Q}\bar{\eta}_{Q})$
in this sum, we obtain three terms, which we will relate to different
energies. First,

$\displaystyle\mathcal{E}_{n,1}^{\circ}(h)+\mathcal{E}_{n,2}^{\circ}(\tilde{h})$
$\displaystyle\leq\sum_{Q}2\operatorname{Area}(Q)\left(|\operatorname{Re}\zeta_{Q}|^{2}+|\operatorname{Im}\zeta_{Q}|^{2}\right)$
$\displaystyle=\sum_{Q}2\operatorname{Area}(Q)|\zeta_{Q}|^{2},$

where the inequality arises since the lefthand side only uses the inner faces
of $G_{n,k}$ to compute $\mathcal{E}_{n,k}^{\circ}$. Similarly,

$\mathcal{E}_{n,1}^{\bullet}(h)+\mathcal{E}_{n,2}^{\bullet}(\tilde{h})\leq\sum_{Q}2\operatorname{Area}(Q)|\eta_{Q}|^{2}.$

Lastly by Claim 44(ii),

$\sum_{Q}2\operatorname{Area}(Q)\cdot
2\operatorname{Im}(\zeta_{Q}\bar{\eta}_{Q})\leq
4\int_{\hat{G}_{n,1}\cup\hat{G}_{n,2}}|\nabla h|^{2}\,dA+C\epsilon_{n}^{1/4}.$

Therefore collecting everything gives

$\mathcal{E}_{n,1}^{\bullet}(h)+\mathcal{E}_{n,2}^{\bullet}(\tilde{h})+\mathcal{E}_{n,1}^{\circ}(h)+\mathcal{E}_{n,2}^{\circ}(\tilde{h})\leq
4\int_{\hat{G}_{n,1}\cup\hat{G}_{n,2}}|\nabla h|^{2}\,dA+C\epsilon_{n}^{1/4},$

and from (33) and (34) we obtain

$\mathcal{E}_{n}^{\bullet}(h_{n}^{\bullet})+m^{2}\mathcal{E}_{n}^{\bullet}(\tilde{h}_{n}^{\bullet})+\mathcal{E}_{n}^{\circ}(h_{n}^{\circ})+m^{2}\mathcal{E}_{n}^{\circ}(\tilde{h}_{n}^{\circ})\leq
4r_{n}\int_{\hat{G}_{n,1}\cup\hat{G}_{n,2}}|\nabla
h|^{2}\,dA+C\epsilon_{n}^{1/4}.$

Now, we use the fact that Fulkerson duality (Lemma 24(i)) gives

$\mathcal{E}_{n}^{\bullet}(h_{n}^{\bullet})\mathcal{E}_{n}^{\circ}(\tilde{h}_{n}^{\circ})=1\qquad\text{and}\qquad\mathcal{E}_{n}^{\circ}(h_{n}^{\circ})\mathcal{E}_{n}^{\bullet}(\tilde{h}_{n}^{\bullet})=1.$

So, dividing by $m=\int_{\Omega}|\nabla h|^{2}\,dA$, we get

$\left(\frac{\mathcal{E}_{n}^{\bullet}(h_{n}^{\bullet})}{m}+\frac{m}{\mathcal{E}_{n}^{\bullet}(h_{n}^{\bullet})}\right)+\left(\frac{\mathcal{E}_{n}^{\circ}(h_{n}^{\circ})}{m}+\frac{m}{\mathcal{E}_{n}^{\circ}(h_{n}^{\circ})}\right)\leq
4r_{n}\frac{1}{m}\int_{\hat{G}_{n,1}\cup\hat{G}_{n,2}}|\nabla
h|^{2}\,dA+C\epsilon_{n}^{1/4}.$

Letting $\epsilon_{n}$ tend to zero, the right-hand side approaches 4, and
hence

$\min\left\\{\limsup_{n\rightarrow\infty}\left(\frac{\mathcal{E}_{n}^{\bullet}(h_{n}^{\bullet})}{m}+\frac{m}{\mathcal{E}_{n}^{\bullet}(h_{n}^{\bullet})}\right),\,\limsup_{n\rightarrow\infty}\left(\frac{\mathcal{E}_{n}^{\circ}(h_{n}^{\circ})}{m}+\frac{m}{\mathcal{E}_{n}^{\circ}(h_{n}^{\circ})}\right)\right\\}\leq
2.$

But, the function $f(x)=\frac{x}{m}+\frac{m}{x}$, for $x>0$, has an absolute
minimum at $x=m$ with value $2$. So, this shows that

$\lim_{n\rightarrow\infty}\mathcal{E}_{n}^{\bullet}(h_{n}^{\bullet})=\lim_{n\rightarrow\infty}\mathcal{E}_{n}^{\circ}(h_{n}^{\circ})=m.$

To finish, we much establish the two claims.

###### Proof of Claim 43.

It suffices to prove the statement for $\epsilon_{n}$ small enough. Let
$\partial V_{n,1}^{\bullet}$ consist of the vertices in $V_{n,1}^{\bullet}$
that are in the topological boundary of the outer face of $G_{n,1}$ but not in
$T_{1}^{n}\cup T_{2}^{n}$. If $z\in\partial V^{\bullet}_{n,1}$, then
$\text{dist}(z,A\cup B)\leq\frac{5}{4}\epsilon_{n}^{1/4}$ by Claim 37. Set
$A_{n}^{\prime}=\\{z\in\partial
V_{n,1}^{\bullet}\,:\,\text{dist}(z,A)\leq\frac{5}{4}\epsilon_{n}^{1/4}\\}$
and define $B_{n}^{\prime}$ similarly (by replacing $A$ with $B$ in the above
definition). We assume that $\epsilon_{n}$ is small enough so that
$A_{n}^{\prime}$ and $B_{n}^{\prime}$ are disjoint. Let $z\in A_{n}^{\prime}$.
Then we see from Lemma 41 that there exists $w\in A_{n}$ so that
$|h_{n}^{\bullet}(z)-h_{n}^{\bullet}(w)|\leq C\log^{-1/2}(1/\epsilon_{n})$.
Using this along with Theorem 3 shows that, for $z\in A_{n}^{\prime}$,

$|h(z)|=|h(z)-h_{n}^{\bullet}(w)|\leq|h(z)-h_{n}^{\bullet}(z)|+|h_{n}^{\bullet}(z)-h_{n}^{\bullet}(w)|\leq\frac{C}{\log^{1/2}(1/\epsilon_{n})}.$

Similarly, for $z\in B_{n}^{\prime}$,

$|h(z)-1|\leq\frac{C}{\log^{1/2}(1/\epsilon_{n})}.$

Therefore,

$\text{gap}_{A_{n}^{\prime},B_{n}^{\prime}}\left(h\right)\geq
1-\frac{C}{\log^{1/2}(1/\epsilon_{n})},$

where gap is defined in (30) and we assume that $\epsilon_{n}$ is small enough
so that the righthand side of this equation is positive.

Let $g$ be the solution to the discrete Dirichlet problem on
$G_{n,1}^{\bullet}$ that is 0 on $A_{n}^{\prime}$ and 1 on $B_{n}^{\prime}$.
Then Dirichlet’s Principle gives that

$\mathcal{E}_{n,1}^{\bullet}(g)=\inf_{f}\frac{\mathcal{E}_{n,1}^{\bullet}(f)}{\text{gap}_{A_{n}^{\prime},B_{n}^{\prime}}(f)^{2}}\leq
r_{n}\,\mathcal{E}_{n,1}^{\bullet}(h),$

where $r_{n}=\left(1-\frac{C}{\log^{1/2}(1/\epsilon_{n})}\right)^{-2}$. Since
every path in $\Gamma_{n}$ contains a path between $A_{n}^{\prime}$ and
$B_{n}^{\prime}$ in $G_{n,1}^{\bullet}$, Statement 5.1 in [5] implies that

$\operatorname{Mod}_{2,\sigma_{n}}(\Gamma_{n})\leq\operatorname{Mod}_{2,\sigma_{n}}\left(\text{paths
between }A_{n}^{\prime}\text{ and }B_{n}^{\prime}\text{ in
}G_{n,1}^{\bullet}\right)=\mathcal{E}_{n,1}^{\bullet}(g).$

Therefore

$\mathcal{E}_{n}^{\bullet}(h_{n}^{\bullet})=\operatorname{Mod}_{2,\sigma_{n}}(\Gamma_{n})\leq\mathcal{E}_{n,1}^{\bullet}(g)\leq
r_{n}\,\mathcal{E}_{n,1}^{\bullet}(h).$

∎

###### Proof of Claim 44.

Let $q$ be the intersection point of $e_{Q}^{\bullet}$ and $e^{\circ}_{Q}$,
and let ${\bf v}$ and ${\bf w}$ are the unit vectors in the directions of
$\vv{v_{1}v_{2}}$ and $\vv{w_{1}w_{2}}$, respectively. We will assume that
$\epsilon_{n}$ is small enough for Claim 38 to hold. Then since $\langle\nabla
h(q),{\bf w}\rangle=-\langle\nabla\tilde{h}(q),{\bf v}\rangle$,

$\displaystyle\left|\frac{h(w_{2})-h(w_{1})}{|w_{2}-w_{1}|}+\frac{\tilde{h}(v_{2})-\tilde{h}(v_{1})}{|v_{2}-v_{1}|}\right|$
$\displaystyle\hskip
79.49744pt\leq\left|\frac{h(w_{2})-h(w_{1})}{|w_{2}-w_{1}|}-\langle\nabla
h(q),{\bf
w}\rangle\right|+\left|\frac{\tilde{h}(v_{2})-\tilde{h}(v_{1})}{|v_{2}-v_{1}|}-\langle\nabla\tilde{h}(q),{\bf
v}\rangle\right|$ $\displaystyle\hskip 79.49744pt\leq
2||Hh||_{\tilde{Q}}\,\epsilon_{n}$ $\displaystyle\hskip 79.49744pt\leq
C\epsilon_{n}^{1/2}$

where Lemma 5.3 in [11] gives the penultimate inequality and Claim 38 gives
the last inequality. Similarly, we can use that $\langle\nabla h(q),{\bf
v}\rangle=\langle\nabla\tilde{h}(q),{\bf w}\rangle$ to show
$\left|\frac{\tilde{h}(w_{2})-\tilde{h}(w_{1})}{|w_{2}-w_{1}|}-\frac{h(v_{2})-h(v_{1})}{|v_{2}-v_{1}|}\right|\leq
C\epsilon_{n}^{1/2}$, and this establishes (i).

For (ii), note that

$\operatorname{Im}(\zeta_{Q}\bar{\eta}_{Q})=\frac{h(v_{2})-h(v_{1})}{|v_{2}-v_{1}|}\frac{\tilde{h}(w_{2})-\tilde{h}(w_{1})}{|w_{2}-w_{1}|}-\frac{h(w_{2})-h(w_{1})}{|w_{2}-w_{1}|}\frac{\tilde{h}(v_{2})-\tilde{h}(v_{1})}{|v_{2}-v_{1}|}.$

Using the same tools as above (Cauchy-Riemann equations, Lemma 5.3 in [11],
and Claim 38), one can show that

$\left|\operatorname{Im}(\zeta_{Q}\bar{\eta}_{Q})-|\nabla h(q)|^{2}\right|\leq
C\epsilon_{n}^{1/4}.$

Further, Lemma 5.3 in [11] and Claim 38 also give that $|\nabla
h(q)|^{2}\leq|\nabla h(z)|^{2}+C\epsilon_{n}^{1/4}$ for any $z\in Q$. Hence

$\operatorname{Im}(\zeta_{Q}\bar{\eta}_{Q})\operatorname{Area}(Q)\leq|\nabla
h(q)|^{2}\operatorname{Area}(Q)+C\epsilon_{n}^{1/4}\operatorname{Area}(Q)\leq\int_{Q}|\nabla
h|^{2}\,dA+C\epsilon_{n}^{1/4}\operatorname{Area}(Q).$

∎

By proving Claims 43 and 44, we have established Theorem 4. ∎

## References

  * [1] Ahlfors, L. V. Conformal invariants: topics in geometric function theory.  McGraw-Hill Book Co., New York-Düsseldorf-Johannesburg, 1973.  McGraw-Hill Series in Higher Mathematics.
  * [2] Albin, N., Brunner, M., Perez, R., Poggi-Corradini, P., and Wiens, N. Modulus on graphs as a generalization of standard graph theoretic quantities.  Conform. Geom. Dyn. 19 (2015), 298–317.  http://dx.doi.org/10.1090/ecgd/287.
  * [3] Albin, N., Clemens, J., Fernando, N., and Poggi-Corradini, P. Blocking duality for $p$-modulus on networks and applications.  Ann. Mat. Pura Appl. (4) 198, 3 (2019), 973–999.  https://doi-org.er.lib.k-state.edu/10.1007/s10231-018-0806-0.
  * [4] Albin, N., and Poggi-Corradini, P. Minimal subfamilies and the probabilistic interpretation for modulus on graphs.  The Journal of Analysis (2016), 1–26.  http://dx.doi.org/10.1007/s41478-016-0002-9.
  * [5] Albin, N., Sahneh, F. D., Goering, M., and Poggi-Corradini, P. Modulus of families of walks on graphs.  In Complex analysis and dynamical systems VII, vol. 699 of Contemp. Math. Amer. Math. Soc., Providence, RI, 2017, pp. 35–55.
  * [6] Alrayes, N. Approximation of $p$-modulus in the plane with discrete grids.  PhD thesis, Kansas State University, 2018.
  * [7] Badger, M. Beurling’s criterion and extremal metrics for Fuglede modulus.  Ann. Acad. Sci. Fenn. Math. 38, 2 (2013), 677–689.
  * [8] Brooks, R. L., Smith, C. A. B., Stone, A. H., and Tutte, W. T. The dissection of rectangles into squares.  Duke Math. J. 7 (1940), 312–340.
  * [9] Eriksson-Bique, S., and Poggi-Corradini, P. On the sharp lower bound for duality of modulus.  Unpublished (2021).  https://arxiv.org/pdf/2102.03035.pdf.
  * [10] Garnett, J. B., and Marshall, D. E. Harmonic measure, vol. 2 of New Mathematical Monographs.  Cambridge University Press, Cambridge, 2008.  Reprint of the 2005 original.
  * [11] Gurel-Gurevich, O., Jerison, D. C., and Nachmias, A. The Dirichlet problem for orthodiagonal maps.  Adv. Math. 374 (2020), 107379, 53.
  * [12] Hinkkanen, A. Modulus of continuity of harmonic functions.  J. Analyse Math 51 (1988), 1–29.
  * [13] Lovász, L. Graphs and geometry, vol. 65 of American Mathematical Society Colloquium Publications.  American Mathematical Society, Providence, RI, 2019.
  * [14] Lyons, R., and Peres, Y. Probability on trees and networks, vol. 42 of Cambridge Series in Statistical and Probabilistic Mathematics.  Cambridge University Press, New York, 2016.
  * [15] Lyons, R., and Peres, Y. Probability on Trees and Networks.  Cambridge University Press, New York, 2016.
  * [16] Pommerenke, C. Boundary behavior of conformal maps.  Springer-Verlag Berlin Heidelberg, 1992.
  * [17] Prasolov, M., and Skopenkov, M. Tiling by rectangles and alternating current.  J. Combin. Theory Ser. A 118, 3 (2011), 920–937.
  * [18] Rogers, L. C. G., and Williams, D. Diffusions, Markov Processes, and Martingales, 2 ed., vol. 1 of Cambridge Mathematical Library.  Cambridge University Press, 2000.
  * [19] Shewchuk, J. R. Triangle: Engineering a 2D Quality Mesh Generator and Delaunay Triangulator.  In Applied Computational Geometry: Towards Geometric Engineering, M. C. Lin and D. Manocha, Eds., vol. 1148 of Lecture Notes in Computer Science. Springer-Verlag, May 1996, pp. 203–222.  From the First ACM Workshop on Applied Computational Geometry.
  * [20] Skopenkov, M. The boundary value problem for discrete analytic functions.  Adv. Math. 240 (2013), 61–87.
  * [21] Skopenkov, M. The boundary value problem for discrete analytic functions.  Adv. Math. 240 (2013), 61–87.
  * [22] Werness, B. Discrete analytic functions on non-uniform lattices without global geometric control.  arXiv: Complex Variables (2015).