# Pair-breaking scattering interference as a mechanism for superconducting gap
modulation

Zhi-Qiang Gao  Department of Physics, University of California, Berkeley,
California 94720, USA  Yu-Ping Lin  Department of Physics, University of
California, Berkeley, California 94720, USA  Dung-Hai Lee  Department of
Physics, University of California, Berkeley, California 94720, USA  Materials
Sciences Division, Lawrence Berkeley National Laboratory, Berkeley, California
94720, USA

###### Abstract

We propose the “pair-breaking scattering interference” as a general source of
coherence peak modulations in superconductors. Assuming this mechanism, we
present a simple physical picture for the coherence peak modulations in
overdoped cuprate Bi2Sr2Ca2Cu3O10+δ (Bi-2223), ferromagnetic iron pnictide
EuRbFe4As4 (Eu-1144), and kagome metals $A$V3Sb5 ($A=$ K, Rb, and Cs).
Specifically, we explain the wavevectors, the particle-hole symmetry, and the
dependence on the internal or external Zeeman-field of the coherence peak
modulations. This work is intended as a cautious reminder to the scientific
community when asserting the existence of a pair density wave phenomenon in
the absence of tunneling conductance modulations in the normal state.

Introduction. The pair density wave (PDW), a superconducting (SC) state
characterized by the condensation of Cooper pairs possessing a non-zero
center-of-mass momentum [1, 2], has been a sought-after realm of inquiry.
Owing to the advance of scanning tunneling microscope (STM) techniques,
experimental observations relevant to PDW can be made by detecting the spatial
modulations of either the height or the energy positions (the SC gap) of the
coherence peaks. For example, uni-directional modulations in the height and
the energy positions of the coherence peaks have been observed in the halo of
vortex cores for Bi2Sr2CaCu2O8 (Bi-2212) [3]. Bi-directional coherence peak
modulations have been recently observed in overdoped Bi2Sr2CaCu2O8+δ (Bi-2212)
and Bi2Sr2Ca2Cu3O10+δ (Bi-2223) (where the modulation wavevector
$\approx(\pi,\pi)$) [4], heavy-fermion superconductor UTe2 [5], and kagome
metal CsV3Sb5 [6]. Adding to this list, uni-directional modulations of the SC
gap have been reported for the ferromagnetic (FM) iron pnictide EuRbFe4As4
(Eu-1144) [7], and along the defect domain walls formed when the monolayer
Fe(Se,Te) is grown on SrTiO3 [8].

This paper excludes the discussion of coherence peak modulations that could be
attributed to a secondary order arising from the mere coexistence of a normal-
state charge density wave (CDW) and a uniform SC order. The observed “PDW” in
UTe2 falls within this category. In contrast, the coherence peak modulations
in Bi-2212, Bi-2223, Eu-1144, and $A$V3Sb5 with $A=$ K, Rb, and Cs belong to
the category we shall discuss.

It is very important to note that the concept to be expounded below was
initially introduced in Ref. [4] for Bi-2212 and Bi-2223, and in this paper,
we extend its applicability to encompass the cases of Eu-1144 and $A$V3Sb5. We
have coined this concept the “pair-breaking scattering interference” (PBSI),
to differentiate it from the PDW phenomenon. This work serves as a cautious
reminder to the scientific community when asserting the existence of a PDW
under the situation that a similar modulation in the local density of states
is absent in the normal state.

The experimental phenomenology. To begin, we briefly summarize the phenomenon
of the coherence peak modulation arsing from the spatial modulation of the
superconducting gap. As a function of spatial coordinate $\mathbf{r}$, the
positive and negative coherence peaks detected at the SC gap edges are
situated at energies

$\displaystyle\Delta_{+}(\mathbf{r})$
$\displaystyle=\Delta_{0}+\Delta_{1}\cos(\mathbf{Q}\cdot\mathbf{r}),$ (1)
$\displaystyle\Delta_{-}(\mathbf{r})$
$\displaystyle=-[\Delta_{0}+\rho\Delta_{1}\cos(\mathbf{Q}\cdot\mathbf{r})].$

Here $2\Delta_{0}>0$ denotes the average energy distance between the coherence
peaks, and $\Delta_{1}\cos(\mathbf{Q}\cdot\mathbf{r})$ represents the
modulation component with amplitude $\Delta_{1}>0$ (usually small compared
with $\Delta_{0}$) at wavevector $\mathbf{Q}$. The coefficient $\rho=\pm 1$
controls the particle-hole symmetry/asymmetry nature of the modulation. In
Eu-1144, the modulation is uni-directional and particle-hole symmetric, hence
$\rho=1$ [7]. Meanwhile, Bi-2212 and Bi-2223 exhibits bi-directional, ${\bf
Q}\approx(\pi,\pi)$, particle-hole asymmetric modulations with $\rho=-1$ [4].
To understand how the coherence peak modulations arise, we look into the
electronic structure and the pairing symmetry associated with the uniform SC
order in these materials (Fig. 1).

Figure 1: Fermi surfaces, gap function signs, and impurity-scattering
wavevectors of (a) Eu-1144 and (b) Bi-2223. For clarity, the BZ of Eu-1144 is
drawn with the center at $M$, and the BZ edge is indicated by the black dashed
lines. The Fermi surfaces are colored, where the red and blue segments have
positive and negative gap functions (labeled with $\pm$), respectively. The
red dashed curves in panel (a) are equal-energy contours of the SC states. The
gray dashed lines in (b) indicate the gap-vanishing axes of $d$-wave pairing.
The arrows $\mathbf{Q}$ represent the wavevectors of the coherence peak
modulations.

Eu-1144: The Fermi surfaces are composed of hole pockets and electron pockets
at the Brillouin zone (BZ) center $\Gamma$ and corner $M$, respectively [9].
The SC state is nodeless, and the widely believed pairing symmetry is
$s_{\pm}$ (meaning the signs of the gap function are opposite between the
$\Gamma$\- and $M$-point Fermi pockets). The observed SC gap modulation (with
wavevector $\sim$ 8 lattice constant) roughly agrees with $2\pi$ over the
momentum distance between the tips of the ellipsoidal electron pockets
centered at M [see Fig. 1 (a)]. We attribute the lack of the $C_{4}$ rotation
symmetry to the nematic order pervasive in iron pnictides. In the subsequent
discussions we assume the SC gap modulation is due to the PBSI of the
Bogoliubov quasi-particles (QPs) at the tips of the electron pockets. This
assumption is reasonable because the area enclosed by SC state equal-energy
contours [red dashed curves in Fig. 1(a)] expands the fastest at the tips of
ellipsoid, and the density of states is highest there. We refer to the tips of
the ellipsoid as “hot spots”.

Bi-2212 and Bi-2223: The Fermi surfaces of the overdoped Bi-2212 and Bi-2223
are schematically shown in Fig. 1 (b). The different colors indicate that the
gap function has opposite signs on the respective segments of the Fermi
surface. The superconducting gap is maximal on the Fermi surfaces closest to X
and Y. Moreover, due to the proximity of the van Hove singularity, the density
of states is the largest, and hence the hot spots are situated there. The bi-
directional coherence peak-modulation wavevectors $\mathbf{Q}\approx(\pi,\pi)$
[4] again approximately fit the wavevectors connecting the inequivalent hot
spots.

The PBSI mechanism for the periodic SC gap modulation. Consider the SC mean-
field Hamiltonian of a spin-singlet superconductor

$H=\sum_{\mathbf{k}}\Psi_{\mathbf{k}}^{\dagger}(\epsilon_{\mathbf{k}}ZI-
hZZ-\Delta^{R}_{\mathbf{k}}YY-\Delta^{I}_{\mathbf{k}}XY)\Psi_{\mathbf{k}}.$
(2)

The Nambu spinor at momentum ${\bf k}$ is defined as
$\Psi_{\mathbf{k}}=(c_{\mathbf{k}\uparrow},c_{\mathbf{k}\downarrow},c_{-\mathbf{k}\uparrow}^{\dagger},c_{-\mathbf{k}\downarrow}^{\dagger})^{T}$.
$I$, $X$, $Y$, $Z$ are the Pauli matrices. Two Pauli matrices standing next to
each other, e.g., XY, represents the tensor product, with the first Pauli
matrix acting in the particle-hole space and the second one acting in the spin
space.

The dispersion relation $\epsilon_{\mathbf{k}}$ is defined with respect to the
Fermi level $\epsilon_{\mathbf{k}_{F}}=E_{F}=0$, where the Fermi momenta
$\mathbf{k}_{F}$ constitute the Fermi surface. Meanwhile, the gap function
$\Delta_{\mathbf{k}}=\Delta^{R}_{\mathbf{k}}+i\Delta^{I}_{\mathbf{k}}=|\Delta_{\mathbf{k}}|e^{i\theta_{\mathbf{k}}}$
is included to capture the Cooper pairing. A Zeeman term with $h>0$ is also
introduced in Eq.(2), which accounts for the coupling to a spontaneous
generated internal field due to ferromagnetism or an external Zeeman field.
Since the Zeeman term $\sim ZZ$ commutes with the Hamiltonian $H$, the spins
$\uparrow$ and $\downarrow$ are the good quantum numbers with $ZZ=+1$ and
$-1$, respectively. Thus, the Hamiltonian in Eq.(2) decomposes into two spin
parts

$\displaystyle
H_{\uparrow}=\Psi_{\mathbf{k}\uparrow}^{\dagger}[\epsilon_{\mathbf{k}}Z-hI+\Delta^{R}_{\mathbf{k}}X-\Delta^{I}_{\mathbf{k}}Y)]\Psi_{\mathbf{k}\uparrow},$
$\displaystyle
H_{\downarrow}=\Psi_{\mathbf{k}\downarrow}^{\dagger}[\epsilon_{\mathbf{k}}Z+hI-\Delta^{R}_{\mathbf{k}}X+\Delta^{I}_{\mathbf{k}}Y)]\Psi_{\mathbf{k}\downarrow},$
(3)

with
$\Psi_{\mathbf{k}\uparrow}=(c_{\mathbf{k}\uparrow},c_{-\mathbf{k}\downarrow}^{\dagger})^{T}$
and
$\Psi_{\mathbf{k}\downarrow}=(c_{\mathbf{k}\downarrow},c_{-\mathbf{k}\uparrow}^{\dagger})^{T}$.

Figure 2: Evolution of quasiparticle energy levels under Zeeman field and the
pair-breaking scattering interference (PBSI). In the ideal superconducting
state the quasiparticles exhibit particle-hole symmetric energy levels
$\pm\Delta_{0}$ at the two hot spots $\mathbf{k}_{1,2}$. The Zeeman field $h$
breaks the spin degeneracy. The two brown energy levels $\pm(-h+\Delta_{0})$
are most important to the physics at the gap edge. In the presence of PBSI
(the crossed dashed lines), the energy levels at the two hot spots get
hybridized and split. The original two-fold degenerate black and brown energy
levels are split into the (green, orange) and the (red, blue) energy levels.

The eigenvalues of the BdG Hamiltonians $H_{\uparrow,\downarrow}$ are
symmetrically centered about energies $\epsilon=\mp h$. At the Fermi momenta
$\mathbf{k}_{F}$, the QP operators are

$\displaystyle\chi_{\mathbf{k}_{F}\uparrow,\downarrow}$
$\displaystyle=\frac{1}{\sqrt{2}}(c_{\mathbf{k}_{F}\uparrow,\downarrow}\pm
e^{i\theta_{\mathbf{k}_{F}}}c_{-\mathbf{k}_{F}\downarrow,\uparrow}^{\dagger}),$
(4) $\displaystyle\eta_{\mathbf{k}_{F}\uparrow,\downarrow}$
$\displaystyle=\frac{1}{\sqrt{2}}(c_{\mathbf{k}_{F}\uparrow,\downarrow}\mp
e^{i\theta_{\mathbf{k}_{F}}}c_{-\mathbf{k}_{F}\downarrow,\uparrow}^{\dagger}).$

Assuming a small Zeeman splitting $h<|\Delta_{\mathbf{k}_{F}}|$ [10], the QP
energy levels associated with $\chi_{\mathbf{k}_{F}\uparrow,\downarrow}$, with
the energies $|\Delta_{\mathbf{k}_{F}}|\mp h>0$, always sit above the Fermi
level (Fig. 2). Conversely, the QP energy levels associated with
$\eta_{\mathbf{k}_{F}\uparrow,\downarrow}$ carrying the energies
$-|\Delta_{\mathbf{k}_{F}}|\mp h<0$ stay below the Fermi level. Note that the
two QP levels associated with $\chi_{\mathbf{k}_{F}\uparrow}$ and
$\eta_{\mathbf{k}_{F}\downarrow}$ have energies closest to the Fermi level,
namely, $\pm(|\Delta_{\mathbf{k}_{F}}|-h)$. They are the ones most relevant to
the gap-edge physics.

Now we study the effect of impurity scattering on the Bolgoiubov QP levels
associated with the two hot spots $\mathbf{k}_{1}$ and $\mathbf{k}_{2}$. In
the following we define $\mathbf{Q}=\mathbf{k}_{1}-\mathbf{k}_{2}$. The gap
functions are
$\Delta_{\mathbf{k}_{1}}=s\Delta_{\mathbf{k}_{2}}=\Delta_{0}e^{i\theta}$,
where $\Delta_{0}>0$ and $s=\pm 1$ represents the relative sign of the gap
function at the hot spots.

Associated with $\chi_{\mathbf{k}_{1,2}\uparrow\downarrow}$ there are four QP
levels above the Fermi energy $E_{F}$, and associated with
$\eta_{\mathbf{k}_{1,2}\uparrow\downarrow}$ there are four QP levels below the
Fermi level. In the following we shall ignore the impurity induced mixing
between the $\chi$ levels and the $\eta$ levels, while only focus on the
mixing among the levels below or above $E_{F}$. The reason for doing so is
that such mixing has the effect of reducing the superconducting gap, hence is
pair-breaking. Without affecting the nature of magnetic scattering we shall
assume it preserves the $z$-component of the spin, under which the impurity-
induced pair-breaking (PB) Hamiltonians for spin-up ($\uparrow$) and spin-down
($\downarrow$) read

$\displaystyle H_{\text{PB},\uparrow\downarrow}$
$\displaystyle=\bigg{[}V^{\text{s}}_{\mathbf{Q}}\bigg{(}\frac{1-s}{2}\bigg{)}\pm
V^{\text{m}}_{\mathbf{Q}}\bigg{(}\frac{1+s}{2}\bigg{)}\bigg{]}$ (5)
$\displaystyle\quad\times(\chi_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}\chi_{\mathbf{k}_{2}\uparrow\downarrow}+\eta_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}\eta_{\mathbf{k}_{2}\uparrow\downarrow})+\text{h.c.}.$

Here $V^{\text{s}}_{\mathbf{Q}}$ and $V^{\text{m}}_{\mathbf{Q}}$ denote the
scalar and magnetic impurity potentials, respectively. See Supplemental
Material (S. M.) [11] Section I for a detailed derivation of this Hamiltonian.
Importantly, (i) if the signs of the gap function at the hot spots are the
same ($s=1$), the scalar impurity does not break pairs. This is consistent
with Anderson’s theorem. In contrast, the magnetic impurity breaks Cooper
pairs. On the other hand, if (ii) the gap sign at the hot spots are opposite
($s=-1$), the roles of the scalar and magnetic impurities are switched.
Eu-1144, and overdoped Bi-2212 and Bi-2223, precisely correspond to the two
cases (i) and (ii) above, respectively.

In the presence of the pair-breaking scattering, the QP states at the hot
spots are mixed. This mixing induces the splitting of QP energy levels (Fig.
2). For weak scattering such splitting is small compared with the other energy
scales such as $h$ and $\Delta_{0}$. Defining the pair-breaking impurity-
scattering potential as $V^{\text{s/m}}_{\mathbf{Q}}=Ve^{i\phi}$ with
magnitude $V>0$ and phase $\phi$, the new QP operators are

$\displaystyle\chi_{\uparrow,S}$
$\displaystyle=\chi_{\mathbf{k}_{1}\uparrow}-e^{i\phi}\chi_{\mathbf{k}_{2}\uparrow},\quad\chi_{\uparrow,L}=\chi_{\mathbf{k}_{1}\uparrow}+e^{i\phi}\chi_{\mathbf{k}_{2}\uparrow},$
(6) $\displaystyle\chi_{\downarrow,S}$
$\displaystyle=\chi_{\mathbf{k}_{1}\downarrow}+se^{i\phi}\chi_{\mathbf{k}_{2}\downarrow},\quad\chi_{\downarrow,L}=\chi_{\mathbf{k}_{1}\downarrow}-se^{i\phi}\chi_{\mathbf{k}_{2}\downarrow},$
$\displaystyle\eta_{\uparrow,S}$
$\displaystyle=\eta_{\mathbf{k}_{1}\uparrow}+e^{i\phi}\eta_{\mathbf{k}_{2}\uparrow},\quad\eta_{\uparrow,L}=\eta_{\mathbf{k}_{1}\uparrow}-e^{i\phi}\eta_{\mathbf{k}_{2}\uparrow},$
$\displaystyle\eta_{\downarrow,S}$
$\displaystyle=\eta_{\mathbf{k}_{1}\downarrow}-se^{i\phi}\eta_{\mathbf{k}_{2}\downarrow},\quad\eta_{\downarrow,L}=\eta_{\mathbf{k}_{1}\downarrow}+se^{i\phi}\eta_{\mathbf{k}_{2}\downarrow}.$

Here $S$ and $L$ label the eigenstates with smaller and larger energy gaps
from the Fermi level, respectively. The $1/\sqrt{2}$ factors are neglected for
clarity. Among these eigenstates, the group of four corresponding to
$\chi_{\uparrow,S,L}$ and $\eta_{\downarrow,S,L}$ (see Fig. 2) are the most
relevant to physics at the gap edge.

Since the eigenstates associated with the QP operators in Eq.(6) are linear
combinations of those at the two different hot-spot momenta, the absolute
value square of them, which appears in the expression of the local density of
states (LDOS), will have a cross term that modulate at the momentum ${\bf
Q=k_{1}-k_{2}}$. This is the PBSI we refer to earlier.

Let us first analyze the limit where we can completely ignore the four QP
levels farther from the Fermi energy. At spatial position $\mathbf{r}$ and
energy bias $\omega$, we derive the LDOS from the four low-energy eigenstates
corresponding to $\chi_{\uparrow,S,L}$ and $\eta_{\downarrow,S,L}$. See S. M.
[11] Section II for a detailed derivation. The LDOS at positive bias and
negative bias are

$\displaystyle D_{+}(\omega,\mathbf{r})$
$\displaystyle\sim\left[1-\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)\right]\delta(\omega-\Delta_{S})$
(7)
$\displaystyle+\left[1+\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)\right]\delta(\omega-\Delta_{L}),$
$\displaystyle D_{-}(\omega,\mathbf{r})$
$\displaystyle\sim\left[1-s\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)\right]\delta(\omega+\Delta_{S})$
$\displaystyle+\left[1+s\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)\right]\delta(\omega+\Delta_{L}),$

where $\Delta_{L,S}=\Delta_{0}-h\pm V$. In Eq.(7) the proportionality constant
is the tunneling matrix element. For $s=1$ (the case of Eu-1144), the
coherence peaks exhibit particle-hole symmetric modulations [Fig. 3(a)],
consistent with the experimental observation [7]. On the other hand, for
$s=-1$ (the case of Bi-2223 and Bi-2212) the coherence peaks exhibit particle-
hole asymmetric modulation [Fig. 3(b)] as seen in Ref. [4].

Figure 3: The coherence peak modulations derived from the energy levels at the
hot spots, shown by the LDOS (label omitted for clarity) as functions of
spatial position $\mathbf{r}$ and energy bias $\omega$. Here we use the
Lorentzian function representing the contribution of each energy level to the
local density of states. (a1) When the gap function has the same sign at the
two hot spots ($s=1$), such as in Eu-1144, the modulations are particle-hole
symmetric. (b1) If the gap function has opposite signs at the two hot spots
($s=-1$), as in the case of Bi-2223, the modulations are particle-hole
asymmetric. From (a2) and (b2), it can be seen that the height of coherence
peaks also modulates for both $s=1$ and $s=-1$ cases, and shows a maximum when
the modulation amplitude is reached.

Now we take all QP levels into account. In Eu-1144, the coherence peak
modulations only occur below a FM transition temperature which is lower than
the SC $T_{c}$ [7]. In contrast, in Bi-2212 and Bi-2223, the coherence peak
appears in the absence of the magnetic field or spontaneous magnetization [4].
To understand the difference in these two cases, we inspect how the QP states
corresponding to $\chi_{\downarrow,S,L}$ and $\eta_{\uparrow,S,L}$ enter the
coherence peak modulations. The LDOS contributed by these states are

$\displaystyle D_{+}^{\prime}(\omega,\mathbf{r})$
$\displaystyle\sim[1+s\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega-\Delta_{S}^{\prime})$
(8)
$\displaystyle+[1-s\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega-\Delta_{L}^{\prime}),$
$\displaystyle D_{-}^{\prime}(\omega,\mathbf{r})$
$\displaystyle\sim[1+\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega+\Delta_{S}^{\prime})$
$\displaystyle+[1-\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega+\Delta_{L}^{\prime}),$

where $\Delta^{\prime}_{L,S}=\Delta_{0}+h\pm V$. Again, the proportionality
constant is the tunneling matrix element. Assuming that $D_{\pm}$ and
$D^{\prime}_{\pm}$ contribute with equal weight, the total LDOS is given by
the

$D_{{\rm
tot},\pm}(\omega,\mathbf{r})=D_{\pm}(\omega,\mathbf{r})+D^{\prime}_{\pm}(\omega,\mathbf{r}).$

When the Zeeman field is absent $h=0$, $\Delta_{S,L}^{\prime}=\Delta_{S,L}$,
yielding the total LDOS given by

$\displaystyle D_{\rm
tot,+}(\omega,\mathbf{r})\sim\Big{\\{}\left[1+\frac{s-1}{2}\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)\right]\delta(\omega-\Delta_{S})$
$\displaystyle+\left[1-\frac{s-1}{2}\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)\right]\delta(\omega-\Delta_{L})\Big{\\}},$
$\displaystyle D_{\rm
tot,-}(\omega,\mathbf{r})\sim\Big{\\{}\left[1-\frac{s-1}{2}\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)\right]\delta(\omega+\Delta_{S})$
$\displaystyle+\left[1+\frac{s-1}{2}\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)\right]\delta(\omega+\Delta_{L})\Big{\\}}.$
(9)

The modulations vanish in the case $s=1$ (Eu-1144), due to the destructive
interference between $D_{\pm}$ and $D^{\prime}_{\pm}$. This explains the
absence of coherence peak modulations above the FM transition temperature in
Eu-1144 [7], where the spontaneously generated Zeeman field is zero. The spin-
splitting Zeeman field, wherever it originates from, serves as a minimal
requirement for the modulations to occur in the case $s=1$. On the other hand,
$D_{\pm}$ and $D^{\prime}_{\pm}$ constructively interfere in the case $s=-1$.
Therefore, the coherence peak modulations exist without Zeeman field, which is
consistent with the experiments of Bi-2212 and Bi-2223 [4].

So far, we have studied two extreme cases: (i) the high-energy QP eigenstates
are completely ignored, and (ii) the high and the low energy QP eigenstates
contribute with equal weight. In reality, the total LDOS under Zeeman field
$h\neq 0$ should be somewhere in-between. In the $s=1$ case (Eu-1144), the
total LDOSs are

$\displaystyle D_{{\rm tot},\pm}^{\text{Eu-1144}}(\omega,\mathbf{r})$
$\displaystyle\sim[1-W(h)\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega\mp\Delta_{S})$
(10)
$\displaystyle+[1+W(h)\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega\mp\Delta_{L}).$

Here the weight function $W(h)\in[0,1]$ is a smooth monotonic function with
$W(0)=0$. Eq.(10) can explain the emergence of coherence peak modulation below
the ferromagnetic transition temperature [7] in Eu-1144. If an external
magnetic field is applied, the modulation should be visible above the FM
transition temperature. Such a possibility can be verified in future
experiments.

On the other hand, the LDOS in the case where $s=-1$ (Bi-2212 and Bi-2223)
read

$\displaystyle D_{\pm}^{\text{Bi-2212/2223}}(\omega,\mathbf{r})$
$\displaystyle\sim[1\mp\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega\mp\Delta_{S})$
(11)
$\displaystyle+[1\pm\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega\mp\Delta_{L}).$

Thus, the Zeeman field does not affect the particle-hole asymmetric coherence
peak modulations, which can be checked experimentally by applying an in-plane
field.

A clarification is in order, namely, where the magnetic scattering of QP comes
from in Eu-1144. Below the ferromagnetic transition temperature, where a
spontaneous magnetization breaking the time reversal symmetry exists, any
impurity scattering should acquire some magnetic nature.

Kagome metals $A$V3Sb5. Recently, the coherence peak modulations at three
inequivalent wavevectors are observed in the kagome metal CsV3Sb5 [6]. The
observed wavevectors coincide with the $3/4$ of the Bragg wavevectors.
Interestingly, these are equivalent to the wavevectors connecting the three
inequivalent small Fermi pockets after CDW orders reconstruct the Fermi
surfaces [12, 13, 14, 15]. Importantly, the coherence peak modulation
wavevectors are not observed in the normal phase above SC $T_{c}$. This
suggests our PBSI scenario as a promising candidate for the origin of
coherence peak modulations. Identifying the centers of the small Fermi pockets
as the hot spots, the PBSI mechanism is again applicable. Note that the
particle-hole symmetry/asymmetry nature of the modulation is determined solely
by the symmetry of the gap function and the impurity type. Since the nodeless
$s$-wave pairing is proposed experimentally [16], we think the same-sign gap
function at the hot spots is relevant. Meanwhile, the $\mu$SR experiments
suggest the presence of local internal magnetic field [17, 18], which can
serve as the intrinsic Zeeman field and the source of magnetic scattering.

Based on these properties, we derive the particle-hole symmetric coherence
peak modulations at three inequivalent wavevectors in S. M. Section III [11].
This result is consistent with the experiments [6, 15].

It should be noted that for CsV3Sb5, one CDW wavevector, which is 1/4 of a
reciprocal lattice vector, occurs above SC $T_{c}$ [19]. Since the 1/4 and the
3/4 wavevectors differed by a reciprocal lattice vector are equivalent, it is
possible that the coherence peak modulation is induced by the CDW and and a
uniform SC. However, we deem this as unlikely because there are three
inequivalent 3/4 wavevectors appearing in the coherence peak modulations. More
recently, a similar STM experiment has been done on KV3Sb5 [15].
Interestingly, the 1/4 wavevector is absent in the normal phase. Hence if
future experiments still show 3/4 coherence peak modulation wavevectors, it
will serve as a better candidate that the PBSI scenario might apply.

Conclusion. We have proposed the pair-breaking scattering interference (PBSI)
as an alternative explanation of coherence peak modulations in
superconductors. This explanation involves (i) the momentum locations of the
hot spots, which can be obtained experimentally from, say, angle-resolved
photoemission spectroscopy (ARPES); (ii) the sign structure of the gap
function at the hot spots, which requires phase-sensitive measurements such as
phase-sensitive quasiparticle interference STM spectroscopy; and (iii) the
nature of impurities, which can be obtained through, say, NMR measurement and
the phase-sensitive quasiparticle interference in STM. The predictions of the
pair-breaking scattering interference can be examined experimentally in STM
measurements. This work serves as a cautious reminder to the scientific
community when asserting the existence of a PDW phenomenon.

###### Acknowledgements.

D.-H.L. thanks Yayu Wang for collaboration on Ref. [4] which motivates the
current manuscript. He also thanks J. C. Séamus Davis for comments and the
question concerning the modulation in the coherence peak height. D.-H.L. is
grateful to Steve Kivelson for bringing the Eu-1144 experiment to his
attention. Z.-Q.G. acknowledges Steve Kivelson, Srinivas Raghu, and Zhaoyu Han
for fruitful discussions. This work is primarily funded by the U.S. Department
of Energy, Office of Science, Office of Basic Energy Sciences, Materials
Sciences and Engineering Division under Contract No. DE-AC02-05-CH11231
(Theory of Materials program KC2301). Y.-P.L. acknowledges the fellowship
support from the Gordon and Betty Moore Foundation through the Emergent
Phenomena in Quantum Systems (EPiQS) program.

## References

  * Fulde and Ferrell [1964] P. Fulde and R. A. Ferrell, Superconductivity in a strong spin-exchange field, Phys. Rev. 135, A550 (1964).
  * Larkin and Ovchinnikov [1964] A. I. Larkin and Y. N. Ovchinnikov, Nonuniform state of superconductors, Zh. Eksp. Teor. Fiz. 47, 1136 (1964).
  * Edkins  _et al._ [2019] S. D. Edkins, A. Kostin, K. Fujita, A. P. Mackenzie, H. Eisaki, S. Uchida, S. Sachdev, M. J. Lawler, E.-A. Kim, J. C. S. Davis, and M. H. Hamidian, Magnetic field–induced pair density wave state in the cuprate vortex halo, Science 364, 976 (2019).
  * Zou  _et al._ [2022] C. Zou, Z. Hao, X. Luo, S. Ye, Q. Gao, M. Xu, X. Li, P. Cai, C. Lin, X. Zhou, D.-H. Lee, and Y. Wang, Particle–hole asymmetric superconducting coherence peaks in overdoped cuprates, Nat. Phys. 18, 551 (2022).
  * Gu  _et al._ [2023] Q. Gu, J. P. Carroll, S. Wang, S. Ran, C. Broyles, H. Siddiquee, N. P. Butch, S. R. Saha, J. Paglione, J. C. Séamus Davis, and X. Liu, Detection of a pair density wave state in UTe2, Nature 618, 921 (2023).
  * Chen  _et al._ [2021] H. Chen, H. Yang, B. Hu, Z. Zhao, J. Yuan, Y. Xing, G. Qian, Z. Huang, G. Li, Y. Ye, S. Ma, S. Ni, H. Zhang, Q. Yin, C. Gong, Z. Tu, H. Lei, H. Tan, S. Zhou, C. Shen, X. Dong, B. Yan, Z. Wang, and H.-J. Gao, Roton pair density wave in a strong-coupling kagome superconductor, Nature 599, 222 (2021).
  * Zhao  _et al._ [2023] H. Zhao, R. Blackwell, M. Thinel, T. Handa, S. Ishida, X. Zhu, A. Iyo, H. Eisaki, A. N. Pasupathy, and K. Fujita, Smectic pair-density-wave order in EuRbFe4As4, Nature 618, 940 (2023).
  * Liu  _et al._ [2023] Y. Liu, T. Wei, G. He, Y. Zhang, Z. Wang, and J. Wang, Pair density wave state in a monolayer high-tc iron-based superconductor, Nature 618, 934 (2023).
  * Hemmida  _et al._ [2021] M. Hemmida, N. Winterhalter-Stocker, D. Ehlers, H.-A. K. von Nidda, M. Yao, J. Bannies, E. D. L. Rienks, R. Kurleto, C. Felser, B. Büchner, J. Fink, S. Gorol, T. Förster, S. Arsenijevic, V. Fritsch, and P. Gegenwart, Topological magnetic order and superconductivity in $\mathrm{Eu}\mathrm{Rb}{\mathrm{Fe}}_{4}{\mathrm{As}}_{4}$, Phys. Rev. B 103, 195112 (2021).
  * [10] This is indeed the case for Eu-1144 [4], since the Zeeman field originates from the weaker ferromagnetism in the stronger superconductivity. However, the results remain the same even if the reverse is true.
  * [11] See Supplemental Material at [URL] for the pair-breaking Hamiltonian, derivation of LDOS at zero temperature, and superconducting gap modulations in $A\text{V}_{3}\text{Sb}_{5}$.
  * Lin and Nandkishore [2021] Y.-P. Lin and R. M. Nandkishore, Complex charge density waves at Van Hove singularity on hexagonal lattices: Haldane-model phase diagram and potential realization in the kagome metals $A{V}_{3}{\mathrm{Sb}}_{5}$ ($A$=K, Rb, Cs), Phys. Rev. B 104, 045122 (2021).
  * Lin and Nandkishore [2022] Y.-P. Lin and R. M. Nandkishore, Multidome superconductivity in charge density wave kagome metals, Phys. Rev. B 106, L060507 (2022).
  * Zhou and Wang [2022] S. Zhou and Z. Wang, Chern fermi pocket, topological pair density wave, and charge-4e and charge-6e superconductivity in kagomé superconductors, Nat. Commun. 13, 7288 (2022).
  * Li  _et al._ [2023] H. Li, D. Oh, M. Kang, H. Zhao, B. R. Ortiz, Y. Oey, S. Fang, Z. Ren, C. Jozwiak, A. Bostwick, E. Rotenberg, J. G. Checkelsky, Z. Wang, S. D. Wilson, R. Comin, and I. Zeljkovic, Small fermi pockets intertwined with charge stripes and pair density wave order in a kagome superconductor, Phys. Rev. X 13, 031030 (2023).
  * Xu  _et al._ [2021] H.-S. Xu, Y.-J. Yan, R. Yin, W. Xia, S. Fang, Z. Chen, Y. Li, W. Yang, Y. Guo, and D.-L. Feng, Multiband superconductivity with sign-preserving order parameter in kagome superconductor ${\mathrm{csv}}_{3}{\mathrm{sb}}_{5}$, Phys. Rev. Lett. 127, 187004 (2021).
  * Yu  _et al._ [2021] L. Yu, C. Wang, Y. Zhang, M. Sander, S. Ni, Z. Lu, S. Ma, Z. Wang, Z. Zhao, H. Chen, K. Jiang, Y. Zhang, H. Yang, F. Zhou, X. Dong, S. L. Johnson, M. J. Graf, J. Hu, H.-J. Gao, and Z. Zhao, Evidence of a hidden flux phase in the topological kagome metal CsV3Sb5, arXiv e-prints , arXiv:2107.10714 (2021), arXiv:2107.10714 [cond-mat.supr-con] .
  * Guguchia  _et al._ [2023] Z. Guguchia, C. Mielke III, D. Das, R. Gupta, J.-X. Yin, H. Liu, Q. Yin, M. H. Christensen, Z. Tu, C. Gong, N. Shumiya, M. S. Hossain, T. Gamsakhurdashvili, M. Elender, P. Dai, A. Amato, Y. Shi, H. C. Lei, R. M. Fernandes, M. Z. Hasan, H. Luetkens, and R. Khasanov, Tunable unconventional kagome superconductivity in charge ordered RbV3Sb5 and KV3Sb5, Nat. Commun. 14, 153 (2023).
  * Zhao  _et al._ [2021] H. Zhao, H. Li, B. R. Ortiz, S. M. L. Teicher, T. Park, M. Ye, Z. Wang, L. Balents, S. D. Wilson, and I. Zeljkovic, Cascade of correlated electron states in the kagome superconductor CsV3Sb5, Nature 599, 216 (2021).

##  I Supplemental Material for “Pair-breaking scattering interference as a
mechanism for superconducting gap modulation”

##  II The pair-breaking Hamiltonian

In the solids, the electrons
$c_{\mathbf{r}}=(c_{\mathbf{r}\uparrow},c_{\mathbf{r}\downarrow})^{T}$ can
scatter by the scalar and magnetic impurities (for simplicity we assume the
magnetic scattering to be $S_{z}$-conserving without losing generality)

$H_{\text{imp}}=\sum_{\mathbf{r}}c_{\mathbf{r}}^{\dagger}[V^{\text{s}}(\mathbf{r})I+V^{\text{m}}(\mathbf{r})Z]c_{\mathbf{r}},$
(S1)

where $V^{\text{s,m}}(\mathbf{r})$ are the scalar and the magnetic impurity-
scattering potentials, respectively. A Fourier transform
$c_{\mathbf{r}}=(1/\sqrt{N})\sum_{\mathbf{k}}c_{\mathbf{k}}e^{i\mathbf{k}\cdot\mathbf{r}}$,
where $N$ is the number of sites, yields the momentum-space form

$H_{\text{imp}}=\sum_{\mathbf{k}\mathbf{q}}c_{\mathbf{k}}^{\dagger}[V^{\text{s}}_{\mathbf{q}}I+V^{\text{m}}_{\mathbf{q}}Z]c_{\mathbf{k}-\mathbf{q}},$
(S2)

with
$V^{\text{s,m}}_{\mathbf{q}}=(1/N)\sum_{\mathbf{r}}V^{\text{s,m}}(\mathbf{r})e^{-i\mathbf{q}\cdot\mathbf{r}}$.
The expression in terms of the Nambu operator $\psi_{\bf k}=(c_{\bf
k\uparrow},c_{\bf
k\downarrow},c^{\dagger}_{\bf-k\uparrow},c^{\dagger}_{\bf-k\downarrow})^{T}$
reads

$H_{\text{imp}}=\sum_{\mathbf{k}\mathbf{q}}\Psi_{\mathbf{k}}^{\dagger}[V^{\text{s}}_{\mathbf{q}}ZI+V^{\text{m}}_{\mathbf{q}}ZZ]\Psi_{\mathbf{k}-\mathbf{q}}.$
(S3)

Since $H_{\text{imp}}$ commutes with $ZZ$, we can decompose it into the
$\uparrow$ ($ZZ=+1$) and $\downarrow$ ($ZZ=-1$) sectors with

$H_{\text{imp},\uparrow\downarrow}=\sum_{\mathbf{k}\mathbf{q}}\Psi_{\mathbf{k}\uparrow\downarrow}^{\dagger}[V^{\text{s}}_{\mathbf{q}}Z\pm
V^{\text{m}}_{\mathbf{q}}I]\Psi_{\mathbf{k}-\mathbf{q}\uparrow\downarrow},$
(S4)

where $\psi_{\bf k\uparrow}=(c_{\bf
k\uparrow},c^{\dagger}_{\bf-k\uparrow})^{T}$ and $\psi_{\bf
k\downarrow}=(c_{\bf k\downarrow},c^{\dagger}_{\bf-k\downarrow})^{T}$.

Now we derive the pair-breaking Hamiltonian. Here we focus on the hot spots
$\mathbf{k}_{1}$ and $\mathbf{k}_{2}$ on the Fermi surface, where the
quasiparticle operators $\chi_{\mathbf{k}_{1,2}\uparrow\downarrow}$ and
$\eta_{\mathbf{k}_{1,2}\uparrow\downarrow}$ are given by

$\Psi_{\mathbf{k}_{1,2}\uparrow\downarrow}=\frac{1}{\sqrt{2}}\left(\begin{array}[]{c}\chi_{\mathbf{k}_{1,2}\uparrow\downarrow}+\eta_{\mathbf{k}_{1,2}\uparrow\downarrow}\\\
\pm
e^{-i\theta_{\mathbf{k}_{1,2}}}[\chi_{\mathbf{k}_{1,2}\uparrow\downarrow}-\eta_{\mathbf{k}_{1,2}\uparrow\downarrow}]\end{array}\right).$
(S5)

The impurity-scattering potential can be transformed into the quasiparticle
eigenbasis

$H_{\text{imp},\uparrow\downarrow}=\frac{1}{\sqrt{2}}\left(\begin{array}[]{cc}[\chi_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}+\eta_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}],&\pm
e^{i\theta_{\mathbf{k}_{1}}}[\chi_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}-\eta_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}]\end{array}\right)[V^{\text{s}}_{\mathbf{Q}}Z\pm
V^{\text{m}}_{\mathbf{Q}}I]\frac{1}{\sqrt{2}}\left(\begin{array}[]{c}\chi_{\mathbf{k}_{2}\uparrow\downarrow}+\eta_{\mathbf{k}_{2}\uparrow\downarrow}\\\
\pm
e^{-i\theta_{\mathbf{k}_{2}}}[\chi_{\mathbf{k}_{2}\uparrow\downarrow}-\eta_{\mathbf{k}_{2}\uparrow\downarrow}]\end{array}\right)+\text{h.c.}.$
(S6)

The pair-breaking (PB) part of the Hamiltonian is obtained by omitting the
crossing terms between $\chi_{\mathbf{k}_{F}\uparrow\downarrow}$ and
$\eta_{\mathbf{k}_{F}\uparrow\downarrow}$, yielding

$\displaystyle H_{\text{PB},\uparrow\downarrow}$
$\displaystyle=\chi_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}\left(\begin{array}[]{cc}1,&\pm
e^{i\theta_{\mathbf{k}_{1}}}\end{array}\right)\frac{1}{2}\left(\begin{array}[]{cc}V^{\text{s}}_{\mathbf{Q}}\pm
V^{\text{m}}_{\mathbf{Q}}&0\\\ 0&-V^{\text{s}}_{\mathbf{Q}}\pm
V^{\text{m}}_{\mathbf{Q}}\end{array}\right)\left(\begin{array}[]{c}1\\\ \pm
e^{-i\theta_{\mathbf{k}_{2}}}\end{array}\right)\chi_{\mathbf{k}_{2}\uparrow\downarrow}$
(S7)
$\displaystyle\quad+\eta_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}\left(\begin{array}[]{cc}1,&\mp
e^{i\theta_{\mathbf{k}_{1}}}\end{array}\right)\frac{1}{2}\left(\begin{array}[]{cc}V^{\text{s}}_{\mathbf{Q}}\pm
V^{\text{m}}_{\mathbf{Q}}&0\\\ 0&-V^{\text{s}}_{\mathbf{Q}}\pm
V^{\text{m}}_{\mathbf{Q}}\end{array}\right)\left(\begin{array}[]{c}1\\\ \mp
e^{-i\theta_{\mathbf{k}_{2}}}\end{array}\right)\eta_{\mathbf{k}_{2}\uparrow\downarrow}+\text{h.c.}$
$\displaystyle=\chi_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}\frac{1}{2}[(V^{\text{s}}_{\mathbf{Q}}\pm
V^{\text{m}}_{\mathbf{Q}})+(-V^{\text{s}}_{\mathbf{Q}}\pm
V^{\text{m}}_{\mathbf{Q}})e^{i(\theta_{\mathbf{k}_{1}}-\theta_{\mathbf{k}_{2}})}]\chi_{\mathbf{k}_{2}\uparrow\downarrow}$
$\displaystyle\quad+\eta_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}\frac{1}{2}[(V^{\text{s}}_{\mathbf{Q}}\pm
V^{\text{m}}_{\mathbf{Q}})+(-V^{\text{s}}_{\mathbf{Q}}\pm
V^{\text{m}}_{\mathbf{Q}})e^{i(\theta_{\mathbf{k}_{1}}-\theta_{\mathbf{k}_{2}})}]\eta_{\mathbf{k}_{2}\uparrow\downarrow}+\text{h.c.}$
$\displaystyle=\bigg{\\{}V^{\text{s}}_{\mathbf{Q}}\bigg{[}\frac{1-e^{i(\theta_{\mathbf{k}_{1}}-\theta_{\mathbf{k}_{2}})}}{2}\bigg{]}\pm
V^{\text{m}}_{\mathbf{Q}}\bigg{[}\frac{1+e^{i(\theta_{\mathbf{k}_{1}}-\theta_{\mathbf{k}_{2}})}}{2}\bigg{]}\bigg{\\}}(\chi_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}\chi_{\mathbf{k}_{2}\uparrow\downarrow}+\eta_{\mathbf{k}_{1}\uparrow\downarrow}^{\dagger}\eta_{\mathbf{k}_{2}\uparrow\downarrow})+\text{h.c.}.$

This is the pair-breaking Hamiltonian we use in the main text.

##  III Derivation of LDOS at zero temperature

In the STM measurements, the LDOS at spatial position $\mathbf{r}$ and bias
$\omega$ is determined by the local differential conductance

$D(\omega,\mathbf{r})\sim\frac{dI}{dV}(\omega,\mathbf{r}),$ (S8)

where $I$ and $V$ are the local current and voltage, respectively. At positive
and negative biases, the differential conductance corresponding to the
electron and the hole tunnelings are

$\displaystyle D_{+}(\omega,\mathbf{r})$
$\displaystyle\sim\sum_{n}\sum_{\sigma}|\langle
n|c_{\mathbf{r}\sigma}^{\dagger}|\Omega\rangle|^{2}\delta(\omega-E_{n}),$ (S9)
$\displaystyle D_{-}(\omega,\mathbf{r})$
$\displaystyle\sim\sum_{n}\sum_{\sigma}|\langle
n|c_{\mathbf{r}\sigma}|\Omega\rangle|^{2}\delta(\omega+E_{n}).$

Here $|{\Omega}\rangle$ is the SC ground state and $|{n}\rangle$ are
eigenstates with energies $E_{n}$. After Fourier transform, the tunneling
amplitudes read

$\displaystyle\langle n|c_{\mathbf{r}\sigma}^{\dagger}|\Omega\rangle$
$\displaystyle=\frac{1}{\sqrt{N}}\sum_{\mathbf{k}}\langle
n|c_{\mathbf{k}\sigma}^{\dagger}|\Omega\rangle
e^{-i\mathbf{k}\cdot\mathbf{r}},$ (S10) $\displaystyle\langle
n|c_{\mathbf{r}\sigma}|\Omega\rangle$
$\displaystyle=\frac{1}{\sqrt{N}}\sum_{\mathbf{k}}\langle
n|c_{\mathbf{k}\sigma}|\Omega\rangle e^{i\mathbf{k}\cdot\mathbf{r}}.$

On the Fermi surface, we express the $c_{\mathbf{k}_{F}\sigma}$ in terms of
the quasiparticle operators as

$c_{\mathbf{k}_{F}\sigma}=\frac{1}{\sqrt{2}}(\chi_{\mathbf{k}_{F}\sigma}+\eta_{\mathbf{k}_{F}\sigma}).$
(S11)

Note that the ground state obeys the condition
$\chi_{\mathbf{k}_{F}\sigma}|{\Omega}\rangle=\eta_{\mathbf{k}_{F}\sigma}^{\dagger}|{\Omega}\rangle=0$
(at zero temperature, the quasiparticle states above and below the Fermi level
are fully filled and empty, respectively). The tunneling amplitudes then
become

$\displaystyle\langle n|c_{\mathbf{k}_{F}\sigma}^{\dagger}|\Omega\rangle$
$\displaystyle=\frac{1}{\sqrt{2}}\langle
n|\chi_{\mathbf{k}_{F}\sigma}^{\dagger}|\Omega\rangle,$ (S12)
$\displaystyle\langle n|c_{\mathbf{k}_{F}\sigma}|\Omega\rangle$
$\displaystyle=\frac{1}{\sqrt{2}}\langle
n|\eta_{\mathbf{k}_{F}\sigma}|\Omega\rangle.$

Within the QP approximation, the relevant eigenstates $|{n}\rangle$ are the
one-quasiparticle and one-quasihole excitations

$|{n}\rangle=\sum_{\mathbf{k}_{F}}\sum_{\sigma}(\alpha_{n\mathbf{k}_{F}\sigma}^{*}\chi_{\mathbf{k}_{F}\sigma}^{\dagger}+\beta_{n\mathbf{k}_{F}\sigma}\eta_{\mathbf{k}_{F}\sigma})|{\Omega}\rangle.$
(S13)

According to the fact
$\langle\Omega|\chi_{\mathbf{k}\sigma}\chi_{\mathbf{k}^{\prime}\sigma^{\prime}}^{\dagger}|\Omega\rangle=\langle\Omega|\eta_{\mathbf{k}\sigma}^{\dagger}\eta_{\mathbf{k}^{\prime}\sigma^{\prime}}|\Omega\rangle=\delta_{\mathbf{k}\mathbf{k}^{\prime}}\delta_{\sigma\sigma^{\prime}}$,
the scattering amplitudes take the form

$\displaystyle\langle n|\chi_{\mathbf{k}_{F}\sigma}^{\dagger}|\Omega\rangle$
$\displaystyle=\alpha_{n\mathbf{k}_{F}\sigma},$ (S14) $\displaystyle\langle
n|\eta_{\mathbf{k}_{F}\sigma}|\Omega\rangle$
$\displaystyle=\beta_{n\mathbf{k}_{F}\sigma}^{*}.$

Inserting this result back to the real-space tunneling amplitudes, we get

$\displaystyle\langle n|c_{\mathbf{r}\sigma}^{\dagger}|\Omega\rangle$
$\displaystyle=\frac{1}{\sqrt{2N}}\sum_{\mathbf{k}_{F}}\alpha_{n\mathbf{k}_{F}\sigma}e^{-i\mathbf{k}_{F}\cdot\mathbf{r}},$
(S15) $\displaystyle\langle n|c_{\mathbf{r}\sigma}|\Omega\rangle$
$\displaystyle=\frac{1}{\sqrt{2N}}\sum_{\mathbf{k}_{F}}\beta_{n\mathbf{k}_{F}\sigma}^{*}e^{i\mathbf{k}_{F}\cdot\mathbf{r}}.$

Now we derive the LDOS associated with the hot spots $\mathbf{k}_{1}$ and
$\mathbf{k}_{2}$. The relevant tunneling amplitudes read

$\displaystyle\langle n|c_{\mathbf{r}\sigma}^{\dagger}|\Omega\rangle$
$\displaystyle=\frac{1}{\sqrt{2N}}\left[\alpha_{n\mathbf{k}_{1}\sigma}e^{-i\mathbf{Q}\cdot\mathbf{r}}+\alpha_{n\mathbf{k}_{2}\sigma}\right]e^{-i\mathbf{k}_{2}\cdot\mathbf{r}},$
(S16) $\displaystyle\langle n|c_{\mathbf{r}\sigma}|\Omega\rangle$
$\displaystyle=\frac{1}{\sqrt{2N}}\left[\beta_{n\mathbf{k}_{1}\sigma}^{*}e^{i\mathbf{Q}\cdot\mathbf{r}}+\beta_{n\mathbf{k}_{2}\sigma}^{*}\right]e^{i\mathbf{k}_{2}\cdot\mathbf{r}},$

where ${\bf Q=k_{1}-k_{2}}$. The LDOS then take the form

$\displaystyle D_{+}(\omega,\mathbf{r})$
$\displaystyle\sim\frac{1}{2N}\sum_{n}\sum_{\sigma}\left|\alpha_{n\mathbf{k}_{1}\sigma}e^{-i\mathbf{Q}\cdot\mathbf{r}}+\alpha_{n\mathbf{k}_{2}\sigma}\right|^{2}\delta(\omega-
E_{n}),$ (S17) $\displaystyle D_{-}(\omega,\mathbf{r})$
$\displaystyle\sim\frac{1}{2N}\sum_{n}\sum_{\sigma}\left|\beta_{n\mathbf{k}_{1}\sigma}^{*}e^{i\mathbf{Q}\cdot\mathbf{r}}+\beta_{n\mathbf{k}_{2}\sigma}^{*}\right|^{2}\delta(\omega+E_{n}).$

The contributions from the four low-energy quasiparticle eigenstates
$\chi_{\uparrow,S,L}$ and $\eta_{\downarrow,S,L}$ are

$\displaystyle D_{+}(\omega,\mathbf{r})$
$\displaystyle\sim\left|e^{-i\mathbf{Q}\cdot\mathbf{r}}-e^{i\phi}\right|^{2}\delta(\omega-\Delta_{S})+\left|e^{-i\mathbf{Q}\cdot\mathbf{r}}+e^{i\phi}\right|^{2}\delta(\omega-\Delta_{L})$
(S18)
$\displaystyle\sim[1-\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega-\Delta_{S})+[1+\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega-\Delta_{L}),$
$\displaystyle D_{-}(\omega,\mathbf{r})$
$\displaystyle\sim\left|e^{i\mathbf{Q}\cdot\mathbf{r}}-se^{-i\phi}\right|^{2}\delta(\omega+\Delta_{S})+\left|e^{i\mathbf{Q}\cdot\mathbf{r}}+se^{-i\phi}\right|^{2}\delta(\omega+\Delta_{L})$
$\displaystyle\sim[1-s\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega+\Delta_{S})+[1+s\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega+\Delta_{L}).$

Similarly, we can do the calculation for the four high-energy quasiparticle
eigenstates $\chi_{\downarrow,S,L}$ and $\eta_{\uparrow,S,L}$

$\displaystyle D_{+}^{\prime}(\omega,\mathbf{r})$
$\displaystyle\sim\left|e^{-i\mathbf{Q}\cdot\mathbf{r}}+se^{i\phi}\right|^{2}\delta(\omega-\Delta_{S}^{\prime})+\left|e^{-i\mathbf{Q}\cdot\mathbf{r}}-se^{i\phi}\right|^{2}\delta(\omega-\Delta_{L}^{\prime})$
(S19)
$\displaystyle\sim[1+s\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega-\Delta_{S}^{\prime})+[1-s\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega-\Delta_{L}^{\prime}),$
$\displaystyle D_{-}^{\prime}(\omega,\mathbf{r})$
$\displaystyle\sim\left|e^{i\mathbf{Q}\cdot\mathbf{r}}+e^{-i\phi}\right|^{2}\delta(\omega+\Delta_{S}^{\prime})+\left|e^{i\mathbf{Q}\cdot\mathbf{r}}-e^{-i\phi}\right|^{2}\delta(\omega+\Delta_{L}^{\prime})$
$\displaystyle\sim[1+\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega+\Delta_{S}^{\prime})+[1-\cos(\mathbf{Q}\cdot\mathbf{r}+\phi)]\delta(\omega+\Delta_{L}^{\prime}).$

These are the LDOS formula we use in the main text.

##  IV Superconducting gap modulations in $A\text{V}_{3}\text{Sb}_{5}$

Figure S1: Small Fermi pockets (yellow) of kagome metals $A$V3Sb5 at
$M_{1,2,3}$ points. The gap function has the same sign on all these pockets.
The pockets at opposite momenta are equivalent, since they are connected by a
reciprocal lattice vector of the folded BZ. The three inequivalent Fermi
pockets are connected by the wavevectors $\mathbf{Q}_{1,2,3}$. The Bragg peaks
of the CDW phase are marked by black crosses. For clarity, the wavevectors of
the CDW order connecting these peaks are omitted.

In the kagome metals $A$V3Sb5, the $2\times 2$ CDW order folds the BZ into the
$1/2\times 1/2$ one (the area enclosed by the dashed lines in Fig. S1). The
reciprocal lattice vectors of the folded BZ connect the origin and the black
crosses in Fig. S1, with the CDW ordering wavevectors connecting the black
crosses. After the CDW-induced band reconstruction, small Fermi pockets appear
at three inequivalent momenta $M_{1,2,3}$ on the edges of the folded BZ. The
inequivalent Fermi pockets are connected by three wavevectors
$\mathbf{Q}_{1,2,3}$. These wavevectors are $1/4$ of the original reciprocal
lattice vectors. They are related to the modulation wavevectors of the
coherence peak, modulo the receiprocal lattice vectors. It is important to
note that $\mathbf{Q}_{1,2,3}$ are distinct from the CDW ordering wavevectors.
In the following we identify the centers of the small Fermi pockets as the hot
spots. They give rise to the coherence peak modulations wavevectors in the
superconducting state via the PBSI mechanism.

Our analysis assumes the same-sign gap functions at all hot spots. According
to the main text, the relevant pair-breaking impurity is the magnetic one.
Similar to Eq.(5) in the main text and Eq.(S7) in the previous section, the
pair-breaking part of the impurity scattering Hamiltonian is

$H_{\text{PB},\uparrow,\downarrow}=\chi_{\uparrow\downarrow}^{\dagger}\mathcal{H}_{\text{PB},\uparrow,\downarrow}\chi_{\uparrow\downarrow}+\eta_{\uparrow\downarrow}^{\dagger}\mathcal{H}_{\text{PB},\uparrow,\downarrow}\eta_{\uparrow\downarrow}.$
(S20)

Here
$\chi_{\uparrow\downarrow}=(\chi_{M_{1}\uparrow\downarrow},\chi_{M_{2}\uparrow\downarrow},\chi_{M_{3}\uparrow\downarrow})^{T}$
and
$\eta_{\uparrow\downarrow}=(\eta_{M_{1}\uparrow\downarrow},\eta_{M_{2}\uparrow\downarrow},\eta_{M_{3}\uparrow\downarrow})^{T}$
represent the 3-component quasiparticle operators at the hot spots
$M_{1,2,3}$, and

$\mathcal{H}_{\text{PB},\uparrow,\downarrow}=\pm\begin{pmatrix}0&V^{\text{m}}_{\mathbf{Q}_{3}}&(V^{\text{m}}_{\mathbf{Q}_{2}})^{\dagger}\\\
(V^{\text{m}}_{\mathbf{Q}_{3}})^{\dagger}&0&V^{\text{m}}_{\mathbf{Q}_{1}}\\\
V^{\text{m}}_{\mathbf{Q}_{2}}&(V^{\text{m}}_{\mathbf{Q}_{1}})^{\dagger}&0\end{pmatrix}.$
(S21)

Due to the involvement of three hot spots, the energy splitting is more
complicated than the 2-hot-spot results. Here we consider a simplified case
where the $C_{3}$ rotation symmetry is preserved. The pair-breaking
Hamiltonian can be chosen, for example, as

$\mathcal{H}_{\text{PB},\uparrow,\downarrow}=\pm
V\begin{pmatrix}0&e^{i\phi}&e^{-i\phi}\\\ e^{-i\phi}&0&e^{i\phi}\\\
e^{i\phi}&e^{-i\phi}&0\end{pmatrix}.$ (S22)

As in the 2-hot-spot case, we focus on the lowest-energy branches
$\chi_{\uparrow}$ and $\eta_{\downarrow}$. The eigenstates and corresponding
energies are

$\displaystyle E_{+n}(\phi)$ $\displaystyle=(\Delta_{0}-h)+2V\cos(\phi+2\pi
n/3),\quad\chi_{\uparrow n}=\frac{1}{\sqrt{3}}(e^{i2\pi
n/3}\chi_{\mathbf{M}_{1}\uparrow}+e^{-i2\pi
n/3}\chi_{\mathbf{M}_{2}\uparrow}+\chi_{\mathbf{M}_{3}\uparrow}),$ (S23)
$\displaystyle E_{-n}(\phi)$ $\displaystyle=-(\Delta_{0}-h)-2V\cos(\phi+2\pi
n/3),\quad\eta_{\downarrow n}=\frac{1}{\sqrt{3}}(e^{i2\pi
n/3}\eta_{\mathbf{M}_{1}\downarrow}+e^{-i2\pi
n/3}\eta_{\mathbf{M}_{2}\downarrow}+\eta_{\mathbf{M}_{3}\downarrow}),$

where $n=0,1,2$ labels the eigenstates, $\Delta_{0}$ is the average
superconducting gap, and $h$ is the Zeeman field. Note that the three-fold
degeneracy is lifted by the pair-breaking scattering.

The LDOS can be derived similarly to the 2-hot-spot case. From the tunneling
amplitudes

$\displaystyle|\langle n|c_{\mathbf{r}\sigma}^{\dagger}|\Omega\rangle|^{2}$
$\displaystyle=\frac{1}{2N}|\alpha_{n\mathbf{M}_{1}\sigma}e^{-i\mathbf{M}_{1}\cdot\mathbf{r}}+\alpha_{n\mathbf{M}_{2}\sigma}e^{-i\mathbf{M}_{2}\cdot\mathbf{r}}+\alpha_{n\mathbf{M}_{3}\sigma}e^{-i\mathbf{M}_{3}\cdot\mathbf{r}}|^{2}$
(S24) $\displaystyle=\frac{1}{2N}\frac{1}{3}|e^{i2\pi
n/3}e^{-i\mathbf{M}_{1}\cdot\mathbf{r}}+e^{-i2\pi
n/3}e^{-i\mathbf{M}_{2}\cdot\mathbf{r}}+e^{-i\mathbf{M}_{3}\cdot\mathbf{r}}|^{2}$
$\displaystyle\sim
3+2\sum_{\alpha=1}^{3}\cos(\mathbf{Q}_{\alpha}\cdot\mathbf{r}+2\pi n/3),$
$\displaystyle|\langle n|c_{\mathbf{r}\sigma}|\Omega\rangle|^{2}$
$\displaystyle=\frac{1}{2N}|\beta_{n\mathbf{M}_{1}\sigma}^{*}e^{i\mathbf{M}_{1}\cdot\mathbf{r}}+\beta_{n\mathbf{M}_{2}\sigma}^{*}e^{i\mathbf{M}_{2}\cdot\mathbf{r}}+\beta_{n\mathbf{M}_{3}\sigma}^{*}e^{i\mathbf{M}_{3}\cdot\mathbf{r}}|^{2}$
$\displaystyle=\frac{1}{2N}\frac{1}{3}|e^{-i2\pi
n/3}e^{i\mathbf{M}_{1}\cdot\mathbf{r}}+e^{i2\pi
n/3}e^{i\mathbf{M}_{2}\cdot\mathbf{r}}+e^{i\mathbf{M}_{3}\cdot\mathbf{r}}|^{2}$
$\displaystyle\sim
3+2\sum_{\alpha=1}^{3}\cos(\mathbf{Q}_{\alpha}\cdot\mathbf{r}+2\pi n/3),$

we obtain the corresponding LDOS

$\displaystyle D_{\pm}(\omega,\mathbf{r},\phi)$
$\displaystyle\sim\sum_{n=0}^{2}\left[3+2\sum_{\alpha=1}^{3}\cos(\mathbf{Q}_{\alpha}\cdot\mathbf{r}+2\pi
n/3)\right]\delta\left(\omega\mp[(\Delta_{0}-h)+2V\cos(\phi+2\pi
n/3)]\right).$ (S25)

Note that the result is particle-hole symmetric, as expected for the situation
of same-sign gap functions under magnetic impurity scattering. The spatially
varying superconducting gap, i.e. the energy position of coherence peaks,
$\Delta_{\pm}(\mathbf{r},\phi)$ is determined by the first local maxima of
LDOS (See Figs. S2 for illustration)

$\left.\frac{\partial}{\partial\omega}D_{\pm}(\omega,\mathbf{r},\phi)\right|_{\omega=\Delta_{\pm}(\mathbf{r},\phi)}=0.$
(S26)

Without loss of generality, we restrict $\phi\in[0,\frac{\pi}{6})$. The
profile of positive- and negative-bias superconducting gap functions
$\Delta_{\pm}(\mathbf{r},\phi)$ is plotted for $\phi=0.0$, $0.2$, and $0.4$ in
Figs. S2, where the modulations show the $C_{3}$ symmetry explicitly. Notably,
the wavevectors $\mathbf{Q}_{1,2,3}$ are not affected by the phase $\phi$,
i.e. the details of the impurity scattering. Although the modulation profile
changes for different $\phi$, the particle-hole symmetry and the wavevectors
remain unchanged.

Figure S2: The energy separation of the coherence peaks as a function of
spatial coordinates at the positive bias $\Delta_{+}(\mathbf{r},\phi)$ with
(a1) $\phi=0.0$, (b1) $\phi=0.2$, and (c1) $\phi=0.4$, and at the negative
bias $\Delta_{-}(\mathbf{r},\phi)$ with (a2) $\phi=0.0$, (b2) $\phi=0.2$, and
(c2) $\phi=0.4$. These plots are calculated in the presence of a Zeeman field,
since the kagome materials fall into the $s=1$ category in the main text and
hence coherence-peak modulations only occur when Zeeman field is non-zero. It
can be observed that the change of $\phi$ does not affect the wavevectors of
the modulations which is particle-hole symmetric. Illustration of the
superconducting gaps as the first local maxima of $D_{+}(\omega,\mathbf{r})$
(red) and $D_{-}(\omega,\mathbf{r})$ (blue) for a fixed $\mathbf{r}$ with (a3)
$\phi=0.0$, (b3) $\phi=0.2$ and (c3) $\phi=0.4$, respectively, where the first
maxima are pointed by arrows.


1–20

# Impulsive impact of a submerged body

Yuriy A. Semenov  Guo Xiong Wu  Author’s affiliation  Y. A.SEMENOV1 Email
address for correspondence<EMAIL_ADDRESS>Y. N. SAVCHENKO1
G.Y.SAVCHENKO1 1National Academy od Science of Ukraine, Institute of
Hydromechanics, Kyiv, 03057, Ukraine

(??; ?? and in revised form ??)

###### Abstract

An analytical solution of the impulsive impact of a cylindrical body submerged
below a calm water surface is obtained by solving a free boundary problem. The
shape of the cross section of the body is arbitrary. The integral hodograph
method is applied to derive the complex velocity potential defined in a
parameter plane. The boundary-value problem is reduced to a Fredholm integral
equation of the first kind in the velocity magnitude on the free surface. The
velocity field, the impulsive pressure on the body surface, and the added mass
are determined in a wide range of depths of submergence for various cross-
sectional shapes, such as a flat plate, a circular cylinder, and a rectangle.

††volume: ???

##  1 Introduction

The concept of impulsive fluid/structure interaction is widely used to study
the initial stage of violent water impact flows. Among the earliest work is
von Karman’s (1929) analytical solution of the impulsive impact of a flat
plate floating on a free surface. Havelock (1949) studied an impulsively
starting motion of a cylinder, with a constant velocity and a constant
acceleration, respectively. He applied a linearized free surface boundary
condition and investigated the full time evolution of the free surface. Tyvand
and Miloh (1995) studied an unsteady nonlinear free-surface flow using the
method of small-time series expansion taking into consideration orders high
enough to account for the leading gravitational effects on the surface
elevation and predict the hydrodynamic force acting on the cylinder. The
related problem of water-entry (water-exit) of a circular cylinder was studied
theoretically and experimentally by Greenhow and Lin (1983) and Greenhow
(1987).

The impulsive concept is used to predict wave impacts on marine and coastal
structures (Cooker (1995)), ship slamming (Faltinsen (2005)), the impulsive
vertical motion of a body initially floating on a flat free surface (Iafrati
and Korobkin  (2005)), dam-break flows (Korobkin and Yilmaz  (2009)),
impulsive sloshing in containers and tanks (Tyvand and Miloh (2012)), drops
that hit a solid or liquid surface in impulsive impact (Hjelmervik and Tyvand
(2017)). A solution based on the impulsive concept may contain a singularity
on the three-phase contact line. In such cases, the solution is used as an
outer solution which has to be matched with an inner solution in the vicinity
of the singularity using the method of matched asymptotic expansion (King and
Needham  (1994); Needham et al.  (2007)).

In contrast to the previous studies, we consider the impulsive motion of a
body fully submerged in a liquid. The motivation for this research comes from
the naval hydrodynamics of high-speed hydrofoil crafts, whose foil system may
experience sudden vertical impacts caused by waves hitting the main body of
the craft. The impulsive pressure on the body, the velocity on the free
surface, and the associated added mass are determined in a wide range of
depths of submergence and for various shapes including a flat plate, a
circular cylinder, and a rectangle.

##  2 Boundary-value problem

A sketch of the physical domain is shown in figure 1($a$). The body submerged
below a calm free surface is symmetric about the $Y-$axis, and thus only a
half of the flow region is considered. Before the impact, $t=0$, the body and
the liquid are at rest. At time $t=0^{+}$ the body is suddenly set into motion
with acceleration $a$ directed downward such that during infinitesimal time
interval $\Delta t\rightarrow 0$, the speed of the body reaches the value $U$.
The problem of a rigid body moving in a fluid body is kinematically equivalent
to the problem of a fluid body moving around a fixed rigid body with
acceleration $a$ at infinity. We define a noninertial Cartesian system of
coordinates $XY$ attached to the body at point $A$, and an inertial system of
coordinates $X^{\prime}Y^{\prime}$ in which the velocity of the liquid at
infinity equals zero. The body is assumed to have an arbitrary shape, which
can be defined by the slope of the body as a function of the arc length
coordinate $S$, $\beta_{b}=\beta_{b}(S)$. The liquid is assumed to be ideal
and incompressible. The flow starting from the rest remains irrotational at
subsequent times. The gravity and surface tension are ignored.

Figure 1: ($a$) Sketch of the physical plane, ($b$) the parameter, or
$\zeta-$plane.

We can introduce complex potentials $W(Z)=\Phi(X,Y)+i\Psi(X,Y)$ and
$W^{\prime}(Z)=\Phi^{\prime}(X,Y)+i\Psi^{\prime}(X,Y)$ with $Z=X+iY$. In these
systems, the velocity fields are related as follows

$\frac{dW}{dZ}=\frac{dW^{\prime}}{dZ}-iat,$ (1)

where $a$ is the acceleration, $0<t<\Delta t$ and $\Delta t\rightarrow 0$.
Integrating equation (1), we can find

$W=W^{\prime}-iatZ,\quad\frac{\partial W}{\partial t}=\frac{\partial
W^{\prime}}{\partial t}-iaZ,\quad\frac{\partial\Phi^{\prime}}{\partial
t}=\frac{\partial\Phi}{\partial t}-aY.$ (2)

By substituting the last equation into Bernoulli’s equation

$\frac{\partial\Phi^{\prime}}{\partial
t}+\frac{p}{\rho}+\frac{|V^{\prime}|^{2}}{2}=\frac{p_{a}}{\rho}$ (3)

integrating it over infinitesimal time interval $\Delta t\rightarrow 0$, and
taking into account that the integral of the third term tends to zero, one can
obtain

$P=\int_{0}^{\Delta t}pdt=-\rho\Phi+\rho UY,$ (4)

where $P$ is the impulsive pressure and $U=a\Delta t$. Here,
$|V^{\prime}|<\infty$ is the velocity magnitude, $p$ and $p_{a}$ are the
hydrodynamic pressure and the pressure on the free surface, respectively.

We introduce dimensionless quantities normalized to $U$, $L$, and $\rho$:Then,
$v=|V|/U$, $x=X/L$, $y=Y/L$, $h=H/L$, $\phi(s)=\Phi(S)/(LU)$. The vertical
impulse force $F_{y}$ is obtained by integrating the impulse pressure over the
body surface,

$F_{y}=-2\rho L^{2}U\int_{s_{A}}^{s_{C}}\phi(s)\cos(n,y)ds-2\rho UA=\rho
mL^{2}U,$ (5)

where $s$ is the arc length coordinate along the body surface, $s_{A}$ and
$s_{C}$ are the arc lengths of points $A$ and $C$, $m$ is the coefficient of
the added mass and $A$ is the cross-sectional area . The factor $"2"$ accounts
for the force acting on the part of the body symmetric about the $y-$axis. Due
to the symmetry of the body about the $y$-axis, the horizontal impulse force
equals zero.

It is desired to determine the velocity potential $\phi(s)$ immediately after
the impact.

##  3 Conformal mapping

A direct finding of the complex potential $w(z)$ is a challenge; therefore, we
introduce an auxiliary parameter plane, or $\zeta-$plane as suggested by
Michell (1890) and Joukovskii (1890). We formulate boundary-value problems for
the complex velocity function, $dw/dz$, and for the derivative of the complex
potential, $dw/d\zeta$, both defined in the $\zeta-$plane. Then the derivative
of the mapping function is obtained as $dz/d\zeta=(dw/d\zeta)/(dw/dz)$, and
its integration provides the mapping function $z=z(\zeta)$ relating the
coordinates in the parameter and physical planes.

We choose the first quadrant of the $\zeta-$plane in figure 1$b$ as the region
corresponding to the fluid region in the physical plane (figure 1$a$). The
conformal mapping theorem allows us to arbitrarily choose the locations of
three points, which are points $O$ at the origin ($\zeta=0$), $D$
($D^{\prime}$) at infinity, and $B$ at $\zeta=1$ (see figure 1$b$). The
position of points $A$ ($\zeta=a$) and $C$ ($\zeta=c$) has to be determined
from the solution of the problem and physical considerations.

###  3.1 Expressions for the complex velocity and the derivative of the
complex potential

The body is considered to be fixed; therefore, the velocity direction and the
slope of the body coincide. Besides, at this stage we assume that the velocity
magnitude on the free surface is a known function of the parameter variable,
$v(\eta)$. Then the boundary-value problem for the complex velocity in the
first quadrant of the parameter plane can be written as follows

$\chi(\xi)=\arg\left(\frac{dw}{dz}\right)=\left\\{{\begin{array}[]{l}-\beta_{b}(a)+\beta_{0},\quad
0\leq\xi\leq a,\\\ -\beta_{b}(\xi),\qquad\quad\;a\leq\xi\leq c,\\\
-\beta_{b}(c)-\beta_{0},\quad c\leq\xi<\infty.\end{array}}\right.$ (6)
$v(\eta)=\left|{\frac{dw}{dz}}\right|_{\zeta=i\eta},\qquad 0\leq\eta<\infty.$
(7)

where $\beta_{0}=\pi/2$, $\beta_{b}(a)=\pi$ and $\beta_{b}(c)=0$. Equation (6)
satisfies the conditions: $\chi(\xi)=-\pi/2$ along the intervals $OA$ and $CD$
on the symmetry line, and $\chi(\xi)=-\beta_{b}(\xi)$ along the body. The
argument of the complex velocity exhibits jumps $\Delta=-\pi/2$ at points $A$
and $C$ when we move along the boundary in the physical plane from point $O$
to point $D$. This boundary-value problem can be solved by applying the
following integral formula (Semenov & Yoon, 2009):

$\frac{dw}{dz}=v_{\infty}\exp\left[{\frac{1}{\pi}\int\limits_{0}^{\infty}{\frac{d\chi}{d{\xi}}\ln\left({\frac{\zeta+{\xi}}{\zeta-{\xi}}}\right)d{\xi}-\frac{i}{\pi}\int\limits_{0}^{\infty}{\frac{d\ln
v}{d{\eta}}\ln\left({\frac{\zeta-i{\eta}}{\zeta+i{\eta}}}\right)d{\eta}+i\chi_{\infty}}}}\right],$
(8)

where $v_{\infty}=v(\eta)_{\eta\rightarrow\infty}$ and
$\chi_{\infty}=\chi(\xi)_{\xi\rightarrow\infty}$. Substituting equations (6)
and (7) into (8) and evaluating the first integral over the step change at
points $\zeta=a$ and $\zeta=c$, we obtain

$\displaystyle\frac{dw}{dz}$ $\displaystyle=$ $\displaystyle
v_{\infty}\left(\frac{\zeta-a}{\zeta+a}\right)^{\frac{1}{2}}\left(\frac{\zeta-c}{\zeta+c}\right)^{\frac{1}{2}}$
$\displaystyle\times$
$\displaystyle\exp\left[\frac{1}{\pi}\int\limits_{a}^{c}{\frac{d\beta_{b}}{d{\xi}}\ln\left(\frac{\zeta-\xi}{\zeta+\xi}\right)d{\xi}}-\frac{i}{\pi}\int\limits_{0}^{\infty}{\frac{d\ln
v}{d{\eta}}\ln\left({\frac{\zeta-i{\eta}}{\zeta+i{\eta}}}\right)d{\eta}}-i\beta_{0}\right].$

It can be easily verified that for $\zeta=\xi$ the argument of the right-hand
side of (3.1) is the function $-\beta_{b}(\xi)$, while for $\zeta=i\eta$ the
modulus of (3.1) is the function $v(\eta)$, i.e. the boundary conditions (6)
and (7) are satisfied. It can also be seen that the complex velocity function
has zeros of order $1/2$, which correspond to the flow around a corner of
angle $\pi/2$ at points $A$ and $C$.

In order to derive the derivative of the complex potential, we have to analyse
its behaviour. Before the impact, the free surface is flat, and it coincides
with the $x-$axis (see figure 1$a$). The pressure along the free surface is
constant. Then, as follows from the Euler equation, the velocity generated by
the impact is perpendicular to the free surface where the pressure is
constant, or the velocity is directed in the $y-$direction. Therefore, the
$x-$component of the velocity is zero, and
$dw=(dw/dz)dz=(v_{x}-iv_{y})dx=-ivdx$. Thus, the free surface corresponds to
the interval $(-\infty,0)$ on the imaginary axis $\psi$ in the $w-$plane.
$\Im(w)=0$ along the line $OABCD$ due to the impermeability condition on the
body $ABC$ and the intervals $OA$ and $CD$ of the symmetry line. At the same
time, $\Re(w)$ changes from zero at point $O$ to $-\infty$ at infinity $D$.
Thus, the flow region in the physical plane corresponds to the third quadrant
in the $w-$plane. They are related as $w=-K\zeta$ where $K$ is a positive real
number. Then one can obtain

$\frac{dw}{d\zeta}=-K.$ (10)

The simple form of the complex potential, $w(\zeta)=-K\zeta+w_{O}$, allows one
to exclude the parameter $\zeta$ and obtain the solution in the Kirchhoff’s
form, for which the complex velocity is the function of the complex potential.
Here, $w_{O}$ is the complex potential at point $O$.

The derivative of the mapping function is obtained by dividing (10) by (3.1)

$\displaystyle\frac{dz}{d\zeta}$ $\displaystyle=$
$\displaystyle-K\left(\frac{\zeta+a}{\zeta-a}\right)^{\frac{1}{2}}\left(\frac{\zeta+c}{\zeta-c}\right)^{\frac{1}{2}}$
$\displaystyle\times$
$\displaystyle\exp\left[\frac{1}{\pi}\int\limits_{a}^{c}{\frac{-d\beta_{b}}{d{\xi}}\ln\left(\frac{\zeta-\xi}{\zeta+\xi}\right)d{\xi}}+\frac{i}{\pi}\int\limits_{0}^{\infty}{\frac{d\ln
v}{d{\eta}}\ln\left({\frac{\zeta-i{\eta}}{\zeta+i{\eta}}}\right)d{\eta}}+i\beta_{0}\right].$

The integration of this equation yields the mapping function $z=z(\zeta)$
relating the parameter and the physical planes. Equations (3.1) and (10)
include the parameters $a$, $c$, $K$ and the functions $v(\eta)$ and
$\beta_{b}(\xi)$, all to be determined from physical considerations and the
kinematic boundary condition on the free surface and on the solid boundary
$OABCD$.

###  3.2 Body boundary conditions

The arc lengths between points $A$ and $B$, $s_{AB}$, and between points $B$
and $C$, $s_{BC}$, and the depth of submergence are determined as follows

$\int_{a}^{1}\left|\frac{dz}{d\zeta}\right|_{\zeta=\xi}=s_{AB},\qquad\int_{1}^{c}\left|\frac{dz}{d\zeta}\right|_{\zeta=\xi}=s_{BC},\qquad\int_{0}^{a}\left|\frac{dz}{d\zeta}\right|_{\zeta=\xi}=h.$
(12)

The unknown function $\beta_{b}(\xi)$ is determined from the following
integro-differential equation

$\frac{d\beta_{b}}{d\xi}=\frac{d\beta_{b}}{ds}\left|\frac{dz}{d\zeta}\right|_{\zeta=\xi},$
(13)

where $\beta_{b}(s)$ is the given function. Equation (13) is solved by
iteration using $d\beta_{b}/d\xi$ in (3.1) known at the previous iteration.

###  3.3 Free surface boundary conditions

An impulsive impact is characterized by infinitesimally small time interval
$\Delta t\rightarrow 0$ such that the position of the free surface does not
change during the impact. From the Euler equations it follows that the
velocity generated by the impact is perpendicular to the free surface
($p=p_{a}$)

$\arg\left(\left.\frac{dw}{dz}\right|_{\zeta=i\eta}\right)=-\beta_{0},\qquad
0\leq\eta\leq\infty.$ (14)

Taking the argument of the complex velocity from (3.1), we obtain the
following integral equation in the function $d\ln v/d\eta$

$\int\limits_{0}^{\infty}\frac{d\ln
v}{d\eta}\ln\left|\frac{\eta^{\prime}-\eta}{\eta^{\prime}+\eta}\right|d\eta^{\prime}+\tan^{-1}\left(\frac{\eta}{a}\right)+\tan^{-1}\left(\frac{\eta}{c}\right)+\frac{2}{\pi}\int_{a}^{c}\frac{d\beta_{b}}{d\xi}\tan^{-1}\left(\frac{\eta}{\xi}\right)d\xi=0$
(15)

Equation (15) is a Fredholm integral equation of the first kind with a
logarithmic kernel. Its solution takes the form (see Appendix)

$v(\eta)=\sqrt{\eta^{2}+a^{2}}\sqrt{\eta^{2}+c^{2}}\exp\left(\frac{1}{\pi}\int_{a}^{c}\frac{d\beta_{b}}{d\xi}\ln(\eta^{2}+\xi^{2})d\xi\right).$
(16)

The system of equations (12), (13) and (15) forms a closed system of equations
in the parameters $a$, $c$, $K$ and functions $\beta_{b}(\xi)$ and $\ln
v(\eta)$.

##  4 Results

Figure 2: Velocity magnitude on the free surface (left axis) for various
depths of submergence (right axis): h=0.25 (red, solid line); 0.5 (blue,
dashed line); 1.0 (olive, dot-dashed line)

The results presented below are shown in the system of coordinates related to
the liquid at infinity. The velocity distribution on the free surface
generated by the impulsive impact of a flat plate is shown in figure 2 for
various depths of submergence. The position of the plate (right axis) is shown
as a thick horizontal line. The colours of the plate and the associated
velocity distribution are the same. After the impact, the plate gets a
velocity equal to $-1$. At a relatively small depth of submergence, $h<0.25$,
the velocity of the liquid above the plate, $|x|<1$, is close to $-1$, i.e.
the plate entrains the liquid, and they move almost together. Outside of the
plate, $|x|>1$, the plate displaces the liquid, and it moves upward providing
the balance of the liquid coming into and out of the flow region. The larger
the depth of submergence of the plate, the smaller the response of the free
surface.

Figure 3: ($a$) velocity along the free surface (left axis) ($b$) impulsive
pressure along the cylinder for various depths of submergence: $h=0.02$ (red,
solid line), $0.1$ (blue, dashed line), $1.0$ (magenta, dotted line), $3.0$
(green, dot dashed) and $10.0$ (purple, dash-dot-dot)

The velocity distribution along the free surface and the impulsive pressure
versus the arc length coordinate $s$ are shown for a circular cylinder in
figures 3$a$ and 3$b$, respectively. The coordinates $s=0$ and $s=pi$
correspond to the top (point $A$) and the bottom (point $C$) of the cylinder.
The bottom pressure is positive, and it depends only weakly on the depth of
submergence. At the top of the cylinder, the pressure is negative, but it
increases as the cylinder approaches the free surface.

The response of the free surface caused by the impulsive impact of a submerged
rectangle is shown in figure 4 for depth $h=0.25$ and various values of the
side length $b$ of the rectangle. For a small length, $b=0.1$, the velocity
distribution is close to that shown in figure 2 for a plate at depth $h=0.25$.
For larger lengths $b$, the velocity decreases, while the length of the free
surface in $x-$direction affected by the impact increases. Such behaviour
provides a balance between the liquid moving into and out of the flow region.

Figure 4: Velocity along the free surface generated by the impact of a
rectangle submerged at depth $h=0.25$ (left axis) and the position of the
rectangle (right axis): side length $b=0.1$ (red, solid line), $m=2.272$;
$b=1.0$ (blue, dashed line) , $m=2.786$; $b=3$ (magenta, dotted line),
$m=3.241$.

The added mass coefficients are shown in table 1 for various shapes of the
body. For a flat plate, $m\rightarrow\pi/2$ as $h\rightarrow 0$, which agrees
with the added mass for a flat plate floating on a free surface (von Karman
impact solution). For a large depth of submergence, $h=50$, the effect of the
free surface becomes negligible, and the coefficient of the added mass
approaches the value corresponding to the added mass in an unbounded fluid
domain.

$h$ | $0$ | $0.02$ | $0.05$ | $0.1$ | $0.3$ | $0.5$ | $1$ | $5$ | $50$ | $\infty^{\ast}$
---|---|---|---|---|---|---|---|---|---|---
$Plate$ | 1.571 |  001.6240 | 1.735 | 1.876 | 2.265 | 2.516 | 2.835 | 3.108 | 3.137 | $\pi$
$Circle$ | - |  002.090 | 2.111 | 2.162 | 2.370 | 2.531 | 2.777 | 3.104 | 3.141 | $\pi$
$Squire$ | - |  002.993 | 3.011 | 3.048 | 3.292 | 3.552 | 4.024 | 4.667 | 4.754 | 4.754

Table 1: Added mass coefficient for various body shapes and depths of
submergence. ∗ unbounded fluid domain

.

##  5 Upward impulsive impact

Here, we analyse how the direction of an impact on a submerged body affects
the velocity field. In the system of coordinates attached to the body, a
change of the impact direction results in a reversal of the velocity direction
of the liquid at infinity and on the whole solid boundary including the body
and the symmetry line. Equation (6) will keep its form if we set values
$\beta_{0}=-\pi/2$, $\beta_{b}(a)=0$ and $\beta_{b}(c)=-\pi$, or subtract
$\pi$ from $\chi(\xi)$ in (6). In this case, expression for the complex
velocity (3.1) keeps its form. The derivative of the complex potential is
obtained from (10) if we change the sign of the expression. Thus, all the
equations of the problem keep their form. Therefore, the velocity magnitude on
the free surface (16) remains the same for both upward and downward impact
directions.

In the system of coordinates related to the liquid, a change of the impact
direction results in a reversal of the velocity direction on the free surface.
Figure 5 illustrates the velocity distribution on the free surface
corresponding to an upward impact of the plate for the same depths of
submergence as shown in figure 2 for a downward impact.

Figure 5: Same as in figure 2 but opposite (upward) direction of the impact.

##  6 Conclusions

An analytical solution of the impulsive impact of a cylindrical body below an
undisturbed free surface is found using the integral hodograph method. The
cross-sectional shape of the body may be a polygon or an arbitrary smooth
shape. The results are shown for a flat plate, a circular cylinder, and a
rectangle. The equation for the complex velocity includes the velocity
magnitude on the free surface, whose analytic form has been determined.

The associated added masses are found as a function of the depth of
submergence. As the depth of submergence tends to infinity, the added mass
tends to the value corresponding to the added mass in an unbounded fluid
domain.

It is shown that upward and downward impacts generate an identical magnitude
of the velocity on the free surface and identical added mass coefficients.
However, the velocity direction is opposite. The obtained solution can be
considered as a first-order solution in solving the problem using the method
of small-time series.

## Declaration of interests

The authors report no conflict of interests.

## Appendix A

Let’s consider a polygon inscribed in a given smooth boundary of the body. Let
the $N_{p}$ be the number of the sides of length $l_{i}$ the slope to which is
$\beta_{pi}$. Let $\xi_{i}$, $i=1,2,\ldots,N_{p}$, be the points in the
parameter region corresponding to the vertexes of the polygon. Then the
function $\beta_{b}(\xi)$ can be written explicitly

$\beta_{b}(\xi)=\beta_{bi},\qquad\xi_{i-1}<\xi\leq\xi_{i},\quad i=1,\ldots
N_{p},$ (17)

where $\xi_{0}=a$, $\xi_{N_{p}}=c$, $\beta_{p1}=\pi$ and $\beta_{Np}=0$.

By substituting (17) into (3.1) and evaluating the integral over the step
change in the function $\beta_{b}(\xi)$ at points $\xi=\xi_{i}$, we obtain the
complex velocity for the polygon-shaped body

$\displaystyle\frac{dw}{dz}$ $\displaystyle=$
$\displaystyle\left(\frac{\zeta-a}{\zeta+a}\right)^{\frac{1}{2}}\left(\frac{\zeta-c}{\zeta+c}\right)^{\frac{1}{2}}\prod_{i=1}^{N_{p}}\left(\frac{\zeta-\xi_{i}}{\zeta+\xi_{i}}\right)^{\frac{\Delta\beta_{pi}}{\pi}}$
$\displaystyle\times$
$\displaystyle\exp\left[-\frac{i}{\pi}\int\limits_{0}^{\infty}{\frac{d\ln
v}{d{\eta}}\ln\left({\frac{\zeta-i{\eta}}{\zeta+i{\eta}}}\right)d{\eta}}-i\beta_{0}\right].$

The free-surface boundary condition (14) with the complex velocity (A) leads
to the integral equation

$\int\limits_{0}^{\infty}\frac{d\ln
v}{d\eta}\ln\left|\frac{\eta^{\prime}+\eta}{\eta^{\prime}-\eta}\right|d\eta^{\prime}=f(\eta),$
(19)

where

$f(\eta)=\tan^{-1}\left(\frac{\eta}{a}\right)+\tan^{-1}\left(\frac{\eta}{c}\right)+\frac{2}{\pi}\sum_{i=1}^{N_{b}}\Delta\beta_{bi}\tan^{-1}\left(\frac{\eta}{\xi_{i}}\right).$

By applying the following transformations (Polyanin and Manzhirov  (1998))

$\frac{d\ln
v}{d\eta}=-\frac{2}{\pi^{2}}\frac{d}{d\eta}\int_{\eta}^{\infty}\frac{F(u)du}{\sqrt{u^{2}-\eta^{2}}},\quad
F(u)=\frac{d}{du}\int_{0}^{u}\frac{pf(p)}{\sqrt{u^{2}-p^{2}}}dp,$ (20)

the solution of the integral equation (19) is obtained

$v(\eta)=\sqrt{\eta^{2}+a^{2}}\sqrt{\eta^{2}+c^{2}}\prod_{i=1}^{N_{p}}\left(\eta^{2}+{\xi_{i}}^{2}\right)^{\frac{\Delta\beta_{bi}}{\pi}}.$
(21)

By taking the limit of (21) for $N_{b}\rightarrow\infty$ and using
$\Delta\beta_{bi}=\left(\frac{d\beta_{b}}{d\xi}\right)_{i}{\Delta\xi}_{i}$, we
obtain equation (16).

## References

  * Cooker (1995) Cooker, M. J., Peregrine, D.H. 1995 Pressure-impulse theory for liquid impact problems, J. Fluid Mech. 297, pp. 193 – 214.
  * Faltinsen (2005) Faltinsen, O. M. 2005 Hydrodynamics of High-speed Marine Vehicles. Cambridge University Press, 454 pp.
  * Greenhow and Lin (1983) Greenhow, M. 1983 Nonlinear free surface effects: experiments and theory. Rep. 83-19. MIT, Dept. of Ocean Engineering.
  * Greenhow (1987) Greenhow, M. and Yanbao, L. 1987 Added masses for circular cylinders near or penetrating fluid boundaries-review, extension and application to water-entry, -exit and slamming. Ocean Engineering, 14 (4). pp. 325 – 348.
  * Havelock (1949) Havelock, T. H. 1949a The wave resistance of a cylinder started from rest. Q. J. Mech. Appl. Maths, 2, pp. 325 – 334.
  * Hjelmervik and Tyvand (2017) Hjelmervik, K.B., Tyvand, P.A. 2017 Incompressible impulsive wall impact of liquid cylinders. J. Eng Math 103, 159 – 171. https://doi.org/10.1007/s10665-016-9866-6
  * Joukovskii (1890) Joukovskii, N. E.  1890 Modification of Kirchhoff’s method for determination of a fluid motion in two directions at a fixed velocity given on the unknown streamline. Math. Sbornik. 15 (1), pp. 121 – 278.
  * von Karman (1929) von Karman, T. 1929 The impact of seaplane floats during landing. Washington, DC:NACA Tech. Note 321.
  * King and Needham (1994) King, A., Needham, D. 1994 The initial development of a jet caused by fluid, body and free-surface interaction. Part 1. A uniformly accelerating plate. J. of Fluid Mech. 268 pp. 89 – 101.
  * Michell (1890) Michell, J. H. 1890 On the theory of free stream lines, Phil. Trans. R. Soc. Lond. A 181, 389-431.
  * Needham et al. (2007) Needham, D., Billingham, J., King, A. 2007 The initial development of a jet caused by fluid, body and free-surface interaction. Part 2. An impulsively moved plate. J. Fluid Mech. 578 pp. 67 – 84.
  * Korobkin and Yilmaz (2009) Korobkin, A., Yilmaz, O. 2009 The initial stage of dam-break flow. J Eng Math 63, 293-308.
  * Iafrati and Korobkin (2005) Iafrati, A., Korobkin, A.A. 2005 Starting flow generated by the impulsive start of a floating wedge. J Eng Math 51, 99 – 125.
  * Korobkin and Scolan (2006) Korobkin, A. A. and Scolan, Y.-M. 2006 Three-dimensional theory of water impact. Part 2. Linearized Wagner problem. J. of Fluid Mech, 549, pp. 343 - 373.
  * Newman (1977) Newman, J.N. 1977 Marine Hydrodynamics. The MIT Press, Cambridge, Massuchesetts. 389 p.
  * Oliver (2007) Oliver J.M. 2007 Second-order Wagner theory for two-dimensional water-entry problems at small deadrise angles. J. Fluid Mech. 572, pp. 59 - 85.
  * Polyanin and Manzhirov (1998) Polyanin, A. D. and Manzhirov, A. V. 1998 Hand book of Integral Equations. CRC Press Boca Raton London New York Washington, D.C.
  * Semenov & Yoon (2009) Semenov, Y. A. & Yoon, B-S. 2009 Onset of flow separation for the oblique water impact of a wedge. _Phys. of Fluids_ 21, 112103.
  * Tyvand and Miloh (1995) Tyvand, P. A. and Miloh, T. 1995 Free-surface flow due to impulsive motion of a submerged circular cylinder. J. Fluid Mech. 286, pp. 67-101.
  * Tyvand and Miloh (2012) Tyvand, P.A. and Miloh, T. 2012 Incompressible impulsive sloshing J. Fluid Mech. 708, pp. 279-302.


Naoto Kan

# A Microscopic study of Magnetic monopoles in Topological Insulators

Shoto Aoki  Hidenori Fukaya  Mikito Koshino  Yoshiyuki Matsuki

###### Abstract

In this article, we analyze a magnetic monopole in topological insulators. The
monopole obtain a fractional electric charge because of the Witten effect. We
consider this system with a microscopic view by adding the Wilson term to the
ordinary Dirac Hamiltonian. The Wilson term yields the positive mass shift to
the effective mass of the electrons, then the curved domain-wall is
dynamically generated around the monopole. The zero-modes of the electrons are
localized on the domain-wall, which can be identified as the source of the
electric charge.

Preprint number: OU-HET-1212

##  1 Introduction

A magnetic monopole has been intensively studied in particle physics. In
Maxwell theory, the monopole can be seen as a heavy point-like object with
magnetic charge, i.e., the Dirac monopole [1], and the world line of the
monopole is referred as the ’t Hooft loop. In the modern sense, the ’t Hooft
loop operator is the charged object of the one-form magnetic symmetry in
Maxwell theory. In grand unified theories, the monopole can be appeared as the
’t Hooft–Polyakov monopole which induces the proton decay [2, 3].

We explore a behavior of the magnetic monopole in the so-called $\theta$
vacuum. It is known that the monopole dresses electric charge in the non-zero
$\theta$ vacuum, which is called the Witten effect [4]111See also related
recent studies [5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16]. The dressed
electric charge, in general, becomes fractional. The non-zero $\theta$ angle
can be realized by introducing the axion field, but the following issues we
study in this article are unclear in the effective field theoretical
description: (1) what is the origin of the electric charge (it must be
electrons, but unclear), (2) why the electron is bound to the monopole, (3)
why the electric charge is fractional.

In this article, we focus on the $\theta=\pi$ vacuum, which is known that the
vacuum describes the time-reversal symmetric topological insulator. From the
view point of a microscopic description222We find Refs. [17, 18, 19, 20] on
the microscopic approach in condensed matter physics., we address the three
issues mentioned above. As a result of introducing the Wilson term to the
Hamiltonian, we will find an interpretation of our issues. Because of the
positive correction to the effective mass from the Wilson term, the curved
domain-wall333The curved domain-wall fermion is also studied in Refs. [21, 22,
23] is dynamically created around the monopole [24]. The chiral zero-modes of
the electron are localized on the created domain-wall and bound to the
monopole. Due to the cobordism invariant nature of the Atiyah–Singer index,
another domain-wall is required. Since the 50% of the wavefunction is located
on the domain-wall around the monopole, we observe the fractional electric
charge.

##  2 A monopole in three-dimensions

###  2.1 A naive Dirac equation with a monopole

We start with solving the naive Dirac equation in the magnetic monopole
background. The Dirac Hamiltonian is

$\displaystyle
H=\gamma_{0}\left(\gamma^{i}\left(\partial_{i}-iA_{i}\right)+m\right),$ (1)

where the gamma matrices are defined by $\gamma_{0}=\sigma_{3}\otimes{\bf 1}$
and $\gamma_{i}=\sigma_{1}\otimes\sigma_{i}$. Let us put the Dirac monopole
with magnetic charge $q_{m}$ at the origin. Then the background gauge
connection is given by

$\displaystyle A_{1}=-{q_{m}y\over r(r+z)},\quad A_{2}={q_{m}x\over
r(r+z)},\quad A_{3}=0.$ (2)

Note that the magnetic charge is quantized, $q_{m}=\mathbb{Z}/2$, due to the
Dirac quantization condition. The Hamiltonian commutes with the total angular
momentum $J_{i}$,

$\displaystyle J_{i}=L_{i}+{\sigma_{i}\over 2},\quad
L_{i}=-i\epsilon_{ijk}x_{j}\left(\partial_{k}-iA_{k}\right)-{q_{m}x_{i}\over
r},$ (3)

with $[J_{i},J_{j}]=i\epsilon_{ijk}J_{k}$, where $L_{i}$ is the orbital
angular momentum. As usual, we use ${\bf 1}\otimes J_{3}$ to parameterize the
eigenstates. The square of the total angular momentum also commutes with the
Hamiltonian: $[{\bf 1}\otimes J^{2},H]=0$. In addition, there is another
operator which commutes with $H$ and ${\bf 1}\otimes J_{i}$. We define the
“spherical” operator as

$\displaystyle D^{S^{2}}=\sigma^{i}\left(L_{i}+{q_{m}x_{i}\over r}\right)+1.$
(4)

We will see a physical meaning of this operator later. We can easily check
$[\sigma_{3}\otimes D^{S^{2}},H]=0$ and $[J_{i},D^{S^{2}}]=0$.

Let us solve the Dirac equation, $H\psi=E\psi$. We find the normalizable zero-
mode ($E=0$) solution localized at the origin with $j=|q_{m}|-1/2$:

$\displaystyle\psi_{j,j_{3},0}(r,\theta,\phi)={C\over
r}e^{-|m|r}\begin{pmatrix}1\\\ \operatorname{\rm sign}(m)\operatorname{\rm
sign}(q_{m})\end{pmatrix}\otimes\chi_{j,j_{3},0}(\theta,\phi),$ (5)

where we introduced the eigenstate of $D^{S^{2}}$ satisfying
$D^{S^{2}}\chi_{j,j_{3},0}=0$. The solution is bound to monopole, but it still
not enough to explain the following issues: The solution Eq. (5) can exist for
both positive and negative mass. The Witten effect predicts that the monopole
becomes dyon in the topological insulator with $m<0$, but not in the normal
insulator with $m>0$. We cannot describe this prediction in the microscopic
sense unless we impose the chiral boundary condition at the origin by hand. It
also remains unclear why the electric charge becomes fractional. In Ref. [20],
summing all charges within the Dirac sea and comparing it to the charge
configuration in the absence of a monopole yields half of an electric charge,
yet this can only be achieved following a regularization process that disrupts
charge conservation.

###  2.2 A regularized Dirac equation

To explore the issues we saw at the end of the previous section, let us define
the regularized Dirac Hamiltonian,

$\displaystyle H_{\rm reg}=H+\gamma_{0}{D_{i}^{\dagger}D^{i}\over M_{\rm
PV}}=\gamma_{0}\left(\gamma^{i}D_{i}+m+{D_{i}^{\dagger}D^{i}\over M_{\rm
PV}}\right).$ (6)

The physical meaning of the new term we added is as follows. In continuum
field theory, to obtain finite results, we must regularize the theory. We now
regularize our theory using the Pauli–Villars (PV) regularization. Then the
partition function is given by

$\displaystyle Z=\det\left({\not{D}+m\over\not{D}+M_{\rm
PV}}\right)=\det\left[{1\over M_{\rm
PV}}\left(\not{D}+m+{D^{\dagger}_{\mu}D^{\mu}\over M_{PV}}+{\cal O}(1/M_{\rm
PV}^{2},m/M_{\rm PV},F_{\mu\nu}/M_{\rm PV})\right)\right],$ (7)

where $M_{\rm PV}$ is the mass of the PV field 444Precisely, in order to
completely regularize the theory, we have to introduce multiple PV fields.
However, a single PV field is sufficient for our current purpose.. Thus we can
interpret the term we added to the Hamiltonian as the leading term of the
large $M_{\rm PV}$ expansion. We refer the additional term as the Wilson term
since the Laplacian corresponds to the Wilson term in lattice gauge theory.

Before finding the zero-mode eigenstates of the regularized Hamiltonian
$H_{\rm reg}$, let us see the role of the Wilson term schematically. Since the
Laplacian $D_{i}^{\dagger}D^{i}$ is alway positive, the mass shift due to the
Wilson term is also always positive when we take $M_{\rm PV}$ positive. In the
topological insulator with $m<0$, the Wilson term can change the sign of the
effective mass,

$\displaystyle m<0\quad\to\quad m_{\rm eff}=m+{D_{i}^{\dagger}D^{i}\over
M_{\rm PV}}\sim m+{1\over M_{\rm PV}r_{1}^{2}}>0,$ (8)

when the magnetic flux is concentrated in the region $r<r_{1}$. It implies
that the inner region $r<r_{1}$ becomes a normal insulator, and the spherical
domain-wall is dynamically created and the chiral edge-modes appear on it.
Note that this does not happen in the normal insulator with $m>0$ since we
have no sign flip of the mass around the monopole.

Let us investigate the regularized Dirac Hamiltonian with the monopole
background Eq. (2). The rotational symmetry is still genuine, $[H_{\rm
reg},{\bf 1}\otimes J_{i}]$, but the spherical operator does not commute with
the Laplacian: $[D^{S^{2}},D_{i}^{\dagger}D^{i}]\neq 0$. We find solution of
the zero-mode for $r_{1}\to 0$,

$\displaystyle\psi_{j,j_{3}}^{\rm mono}(r,\theta,\phi)={Be^{-M_{\rm
PV}r/2}\over\sqrt{r}}I_{\nu}(\kappa r)\begin{pmatrix}1\\\ -\operatorname{\rm
sign}(q_{m})\end{pmatrix}\otimes\chi_{j,j_{3},0}(\theta,\phi),$ (9)

where $\nu=\sqrt{4|q_{m}|+1}/2$, $\kappa=M_{\rm PV}\sqrt{1+4m/M_{\rm PV}}/2$,
and $I_{\nu}(z)$ is the modified Bessel function of the first kind. We note
that since $I_{\nu}(\kappa r)\sim\exp(M_{\rm PV}r/2-|m|r)/\sqrt{\pi M_{\rm
PV}r}$ in the $M_{\rm PV}\to\infty$ limit, the solution coincides with the
naive solution Eq. (5). We also note that the solution has a peak at
$r=|q_{m}|/M_{\rm PV}$ only when the mass $m$ is negative. As is expected,
this implies that the the domain-wall is created by the Wilson term. We show
the the zero-mode solution with and without the Wilson term in Fig. 1.

Figure 1:  The plot of the zero-mode solution with the Wilson term in Eq. (9)
and without the Wilson term in Eq. (5). We set $q_{m}=1/2$, $m=-0.1$ and
$M_{\rm PV}=10$. Figure from Ref. [24].

We have set the parameters with $q_{m}=1/2$, $m=-0.1$ and $M_{\rm PV}=10$.

###  2.3 Half-integral charge

So far, we implicitly considered a $\mathbb{R}^{3}$ space, but to discuss
nature of topology, we need to compactify the space on a $S^{3}$. Then
topological insulator region with $m_{\rm eff}<0$ seems to have the topology
of a disk $D^{3}$ with a small $S^{2}$ boundary at $r=r_{1}$. However, due to
the cobordism invariance of the Atiyah–Singer index, the index must be zero:

$\displaystyle{\rm Ind}\,D^{S^{2}}={1\over 4\pi}\int_{S^{2}}F={1\over
4\pi}\int_{D^{3}}dF=0,$ (10)

which contradicts with ${\rm Ind}\,D^{S^{2}}=2q_{m}$ if there are the only
chiral zero modes. A resolution is to create another domain-wall at, say,
$r=r_{0}$, outside of the topological insulator. Another zero-modes are
localized at the outside domain-wall, and the index is kept trivial:

$\displaystyle\int_{M}dF=\int_{S^{2}_{\rm mono}\cup S^{2}_{\rm out}}F=0,$ (11)

where $M$ is the region of the topological insulator, $S^{2}_{\rm mono}$ and
$S^{2}_{\rm out}$ correspond to the domain-wall around the monopole and
outside domain-wall, respectively. The zero modes localized on $S^{2}_{\rm
mono}$ and on $S^{2}_{\rm out}$ are mixed by the tunneling effect for finite
$r_{0}$, i.e., $\psi=\alpha\psi^{\rm mono}_{j,j_{3}}+\beta\psi^{\rm
DW}_{j,j_{3}}$, where $\alpha$ and $\beta$ are constant coefficients which
parameterize the mixing, $\psi^{\rm mono/DW}_{j,j_{3}}$ are the localized zero
modes around the monopole/outside domain-wall. We can show that $(\psi^{\rm
mono}_{j,j_{3}})^{\dagger}H\psi^{\rm mono}_{j,j_{3}}=(\psi^{\rm
DW}_{j,j_{3}})^{\dagger}H\psi^{\rm DW}_{j,j_{3}}=0$ for the diagonal parts and
$(\psi^{\rm mono}_{j,j_{3}})^{\dagger}H\psi^{\rm DW}_{j,j_{3}}=(\psi^{\rm
DW}_{j,j_{3}})^{\dagger}H\psi^{\rm mono}_{j,j_{3}}=:\Delta\in\mathbb{R}$ for
the off-diagonal perts. Then we can obtain that $\alpha=\pm\beta$ and the
split energy is $E=\pm\Delta$. As a result, the 50% of the zero mode state is
located at the monopole, while the other is sit at the outside domain-wall,
thus the dressed electric charge of the monopole becomes half-integer.

##  3 Numerical results

###  3.1 Lattice setup

We prepare the three-dimensional hyper-cubic lattice with size $L=31$555In
Ref. [24], we studied the three cases with $L=23,31$ and $47$.. Imposing on
open boundary conditions, we put a monopole at $\bm{x}=\bm{x}_{\rm
m}=(L/2,L/2,L/2)$ and anti-monopole at $\bm{x}=\bm{x}_{\rm a}=(L/2,L/2,1/2)$.
The continuum vector potential at $\bm{x}=(x,y,z)$ is given by

$\displaystyle A_{1}(\bm{x})$ $\displaystyle=q_{m}\left[{-(y-y_{\rm
m})\over|\bm{x}-\bm{x}_{\rm m}|(|\bm{x}-\bm{x}_{\rm m}|+(z-z_{\rm
m}))}-{-(y-y_{\rm a})\over|\bm{x}-\bm{x}_{\rm a}|(|\bm{x}-\bm{x}_{\rm
a}|+(z-z_{\rm a}))}\right],$ (12) $\displaystyle A_{2}(\bm{x})$
$\displaystyle=q_{m}\left[{x-x_{\rm m}\over|\bm{x}-\bm{x}_{\rm
m}|(|\bm{x}-\bm{x}_{\rm m}|+(z-z_{\rm m}))}-{x-x_{\rm
a}\over|\bm{x}-\bm{x}_{\rm a}|(|\bm{x}-\bm{x}_{\rm a}|+(z-z_{\rm
a}))}\right],$ (13) $\displaystyle A_{3}(\bm{x})$ $\displaystyle=0,$ (14)

with $q_{m}\in\mathbb{Z}/2$. We assume a position-dependent mass $m(\bm{x})$
for fermions. We set negative mass $m(\bm{x})=-m_{0}$ for
$r=\sqrt{|\bm{x}-\bm{x}_{\rm m}|}<r_{0}=3L/8$, while positive mass
$m(\bm{x})=+m_{0}$ for $r>r_{0}$ with $m_{0}(L+1)=14$. The Wilson Dirac
Hamiltonian is given by

$\displaystyle H_{\rm
Wilson}=\gamma^{0}\left[\sum^{3}_{i=1}\left(\gamma^{i}{\nabla^{f}_{i}+\nabla^{b}_{i}\over
2}-{1\over 2}\nabla_{i}^{f}\nabla^{b}_{i}\right)+m(\bm{x})\right],$ (15)

where the forward covariant difference is
$\nabla_{i}^{f}\psi(\bm{x})=U_{i}(\bm{x})\psi(\bm{x}+\bm{e}_{i})-\psi(\bm{x})$,
and the backward difference is
$\nabla^{b}_{i}\psi(\bm{x})=\psi(\bm{x})-U^{\dagger}_{i}(\bm{x}-\bm{e}_{i})\psi(\bm{x}-\bm{e}_{i})$.
The link variable $U_{j}(\bm{x})$ is defined as

$\displaystyle
U_{j}(\bm{x})=\exp\left(i\int^{1}_{0}A_{j}(\bm{x}+\bm{e}_{j}l)dl\right).$ (16)

###  3.2 Numerical analises

We compute the energy eigenvalues and eigenfunctions of $H_{\rm Wilson}$ with
$q_{m}=1/2$ and as $E_{k}$ and $\phi_{k}(\bm{x})$, respectively. We plot the
eigenvalues $E_{k}$ in units of $r_{0}$ in Fig. 2. We show the results of the
continuum calculation without the Wilson term by cross symbols. We also show
the chirality expectation value measured by

$\displaystyle\sum_{\bm{x}}\phi^{\dagger}_{k}(\bm{x})\left(\sigma_{1}\otimes\sigma_{r}\right)\phi_{k}(\bm{x})$
(17)

using the color gradation. We find the two nearest zero-modes.

In the left panel of Fig. 3, we plot the local amplitude of the (nearest)
zero-mode,

$\displaystyle A_{1}(\bm{x})=\phi_{1}(\bm{x})^{\dagger}\phi_{1}(\bm{x})r^{2}$
(18)

normalized by $r^{2}$ at the $z=16$ slice. The amplitude has two peaks around
$r=0$ and $r=r_{0}$. The 50% of the state is located around the monopole,
while the other 50% is located at $r=r_{0}$, which implies that we observe the
half electric charge.

In the right of Fig. 3, the distribution of the local effective mass

$\displaystyle m_{\rm eff}(\bm{x})=\phi_{1}(\bm{x})^{\dagger}\left(-{1\over
2}\sum_{i}\nabla^{f}_{i}\nabla^{b}_{i}+m(\bm{x})\right)\phi_{1}(\bm{x})r^{2}/A_{1}(\bm{x})$
(19)

at the $z=16$ slice is plotted, as well. We see that the small island of the
normal insulator, i.e., the positive mass region appears around the monopole;
the domain-wall is dynamically created.

Figure 2:  The spectrum of the energy eigenvalue $E_{k}$ with magnetic charge
$q_{m}=1/2$. Figure from Ref. [24]

Figure 3:  Left panel: the amplitude of the nearest zero-mode $A_{1}(\bm{x})$
at $z=16$ with $q_{m}=1/2$. Right panel: the local effective mass of the
nearest zero-mode $m_{\rm eff}(\bm{x})$. Figure from Ref. [24]

##  4 Summary

In this article, we have microscopically studied the magnetic monopole in the
topological insulators. In the analysis of the naive Dirac equation, it is
possible to find the zero mode solution localized around a monopole. However,
the solution cannot satisfactorily explain why the Witten effect occurs only
when the fermion mass is negative, and why the electric charge dressed by the
monopole becomes fractional. To explain this, we have considered a regularized
Hamiltonian with the Wilson term. Due to the effects of the Wilson term, the
domain wall is created around the monopole only when the fermion mass is
negative. As a result, the zero modes localize on the domain wall.
Furthermore, the cobordism invariance of the Atiyah–Singer index theorem
requires the existence of the outside domain wall, and we have concluded that
the $50\%$ of the wavefunction of the zero modes are localized around the
monopole, so that the electric charge becomes fractional. we have also
numerically analyzed the creation of the domain-wall around the monopole and
the emergence of the chiral zero modes.

### Acknowledgements

We thank Mikio Furuta, Masahiro Hotta, Takuto Kawakami, Tsunehide Kuroki,
Okuto Morikawa, Tetsuya Onogi, Shigeki Sugimoto, Kazuki Yamamoto, Masahito
Yamazaki, and Ryo Yokokura for useful discussions. The work of SA was
supported by JSPS KAKENHI Grant Number JP23KJ1459. The work of HF and NK was
supported by JSPS KAKENHI Grant Number JP22H01219. The work of MK was
supported by JSPS KAKENHI Grant Number JP21H05236, JP21H05232, JP20H01840, and
JP20H00127, and by JST CREST Grant Number JPMJCR20T3, Japan.

## References

  * [1] P.A.M. Dirac, _Quantised singularities in the electromagnetic field,_ , _Proc. Roy. Soc. Lond. A_ 133 (1931) 60.
  * [2] C.G. Callan, Jr, _Disappearing dyons_ , _Phys. Rev. D_ 25 (1982) 2141.
  * [3] V.A. Rubakov, _Adler-Bell-Jackiw anomaly and fermion number breaking in the presence of a magnetic monopole_ , _Nucl. Phys. B_ 203 (1982) 311.
  * [4] E. Witten, _Dyons of charge e theta/2 pi_ , _Phys. Lett. B_ 86 (1979) 283.
  * [5] M. Pretko, _Higher-Spin witten effect and Two-Dimensional fracton phases_ , _Phys. Rev. B Condens. Matter_ 96 (2017) 125151.
  * [6] N. Yamamoto, _Magnetic monopoles and fermion number violation in chiral matter_ , 2005.05028.
  * [7] M. Abe, O. Morikawa and H. Suzuki, _Fractional topological charge in lattice abelian gauge theory_ , _PTEP_ 2023 (2022) 023B03.
  * [8] Y. Hidaka, M. Nitta and R. Yokokura, _Higher-form symmetries and 3-group in axion electrodynamics_ , _Phys. Lett. B_ 808 (2020) 135672\.
  * [9] H. Fukuda and K. Yonekura, _Witten effect, anomaly inflow, and charge teleportation_ , _JHEP_ 01 (2020) 119.
  * [10] A. Sekine and K. Nomura, _Axion electrodynamics in topological materials_ , _J. Appl. Phys._ 129 (2020) 141101.
  * [11] P. Agrawal, K.V. Berghaus, J. Fan, A. Hook, G. Marques-Tavares and T. Rudelius, _Some open questions in axion theory_ , 2203.08026.
  * [12] Y. Hamada, R. Kitano, R. Matsudo and K. Mukaida, _Understanding anomalous particle production in massless QED via time-varying $\vartheta$ angle_, _JHEP_ 09 (2022) 218.
  * [13] Y. Choi, H.T. Lam and S.-H. Shao, _Non-invertible gauss law and axions_ , _JHEP_ 09 (2022) 067.
  * [14] R. Yokokura, _Non-invertible symmetries in axion electrodynamics_ , 2212.05001.
  * [15] N. Kan, O. Morikawa, Y. Nagoya and H. Wada, _Higher-group structure in lattice abelian gauge theory under instanton-sum modification_ , _Eur. Phys. J. C_ 83 (2023) 481.
  * [16] M. Abe, N. Kan, O. Morikawa, Y. Nagoya, S. Onoda and H. Wada, _Higher-group symmetry in lattice gauge theories with restricted topological sectors_ , 2310.01787.
  * [17] X.-L. Qi, T. Hughes and S.-C. Zhang, _Topological field theory of Time-Reversal invariant insulators_ , _Phys. Rev. B Condens. Matter_ 78 (2008) 195424.
  * [18] T. Sasaki, E. Imai and I. Kanazawa, _Witten effect and fractional electric charge on the domain wall between topological insulators and spin ice compounds_ , _J. Phys. Conf. Ser._ 568 (2014) 052029.
  * [19] H.-G. Zirnstein and B. Rosenow, _Topological magnetoelectric effect: Nonlinear time-reversal-symmetric response, witten effect, and half-integer quantum hall effect_ , _Phys. Status Solidi B Basic Res._ 257 (2020) 1900698.
  * [20] H. Yamagishi, _Fermion-monopole system reexamined_ , _Phys. Rev. D_ 27 (1983) 2383.
  * [21] S. Aoki and H. Fukaya, _Curved domain-wall fermions_ , _PTEP_ 2022 (2022) 063B04.
  * [22] S. Aoki and H. Fukaya, _Curved domain-wall fermion and its anomaly inflow_ , _PTEP_ 2023 (2022) 033B05.
  * [23] S. Aoki, H. Fukaya and N. Kan, _A lattice regularization of weyl fermions in a gravitational background_ , _in these proceedings_ .
  * [24] S. Aoki, H. Fukaya, N. Kan, M. Koshino and Y. Matsuki, _Magnetic monopole becomes dyon in topological insulators_ , _Phys. Rev. B Condens. Matter_ 108 (2023) 155104.