# Vaccinating according to the maximal endemic equilibrium achieves herd
immunity

Jean-François Delmas  Jean-François Delmas, CERMICS, École des Ponts, France
<EMAIL_ADDRESS>,  Dylan Dronnier  Dylan Dronnier, Université de
Neuchâtel, Switzerland<EMAIL_ADDRESS>and  Pierre-André Zitt
Pierre-André Zitt, LAMA, Université Gustave Eiffel, France  pierre-
<EMAIL_ADDRESS>
###### Abstract.

We consider the simple epidemiological SIS model for a general heterogeneous
population introduced by Lajmanovich and Yorke (1976) in finite dimension, and
its infinite dimensional generalization we introduced in previous works. In
this model the basic reproducing number $R_{0}$ is given by the spectral
radius of an integral operator. If $R_{0}>1$, then there exists a maximal
endemic equilibrium. In this very general heterogeneous SIS model, we prove
that vaccinating according to the profile of this maximal endemic equilibrium
ensures herd immunity. Moreover, this vaccination strategy is critical: the
resulting effective reproduction number is exactly equal to one. As an
application, we estimate that if $R_{0}=2$ in an age-structured community with
mixing rates fitted to social activity, applying this strategy would require
approximately 29% less vaccine doses than the strategy which consists in
vaccinating uniformly a proportion $1-1/R_{0}$ of the population.

From a dynamical systems point of view, we prove that the non-maximality of an
equilibrium $g$ is equivalent to its linear instability in the original
dynamics, and to the linear instability of the disease-free state in the
modified dynamics where we vaccinate according to $g$.

##  1\. Introduction

Increasing the prevalence of immunity from contagious disease in a population
limits the circulation of the infection among the individuals who lack
immunity. This so-called “herd effect” plays a fundamental role in
epidemiology, for example it has had a major impact in the eradication of
smallpox and rinderpest or the near eradication of poliomyelitis
[HerdImmunityFine2011]. Our aim is to present a targeted vaccination strategy
based on the heterogeneity of the infection spreading in the population which
allows to eradicate the epidemic. We consider for simplicity the deterministic
infinite-dimensional SIS model (with S=Susceptible and I=Infectious) and the
effect of a perfect vaccine. However, we take into account a very general
model for the heterogeneous population based on the infinite-dimensional model
introduced in [delmas_infinite-dimensional_2020], that encompasses the meta-
population SIS models developed by Lajmanovich and Yorke in their pioneer
paper [lajmanovich1976deterministic] or SIS model on graphs. More precisely,
the probability $u_{t}(x)$ of an individual of type $x\in\Omega$ to be
infected at time $t$ is the solution of the (infinite dimensional) ordinary
differential equation:

(1)
$\partial_{t}u_{t}(x)=(1-u_{t}(x))\int_{\Omega}k(x,y)\,u_{t}(y)\,\mu(\mathrm{d}y)-\gamma(x)u_{t}(x)\quad\text{for
$t\geq 0$ and $x\in\Omega$,}$

where $k$ is the transmission rate kernel of the disease, $\gamma$ is the
recovery rate function and $\mu(\mathrm{d}y)$ is the probability for an
individual taken at random to be of type $y\in\Omega$; see Equation (11).

In an homogeneous population, the _basic reproduction number_ of an infection,
denoted by $R_{0}$, is defined as the number of secondary cases one individual
generates on average over the course of its infectious period, in an otherwise
uninfected (susceptible) population. Intuitively, the disease should die out
if $R_{0}<1$ and invade the population if $R_{0}>1$. For the heterogenous
generalization of many classical models in epidemiology (including the
heterogeneous SIS model), it is still possible to define a meaningful basic
reproduction number $R_{0}$, as the number of secondary cases generated by a
typical infectious individual when all other individuals are uninfected and
the threshold phenomenon occurs [Diekmann1990]. In the setting of
[delmas_infinite-dimensional_2020], the reproduction number $R_{0}$ then
corresponds to the spectral radius of the next-generation operator defined as
the integral operator associated to the kernel $k(x,y)/\gamma(y)$.

After a vaccination campaign, let the vaccination strategy $\eta$ denote the
(non necessarily homogeneous) proportion of the non-vaccinated population, and
let the effective reproduction number $R_{e}(\eta)$ denote the corresponding
reproduction number of the non-vaccinated population. Following
[delmas_infinite-dimensional_2020, Section 5.3.], the effective reproduction
number $R_{e}(\eta)$ is given by the spectral radius of the effective next-
generation operator defined as the integral operator associated to the kernel
$k(x,y)\eta(y)/\gamma(y)$, where $\eta(y)$ is the proportion of individuals of
type $y$ which are not vaccinated.

The vaccination strategy $\eta$ is called _critical_ if $R_{e}(\eta)=1$.
Assuming $R_{0}>1$, suppose now that only a proportion $\eta^{\mathrm{uni}}$
of the population can catch the disease, the rest being perfectly immunized.
An infected individual will now only generate $\eta^{\mathrm{uni}}R_{0}$ new
cases, since a proportion $1-\eta^{\mathrm{uni}}$ of previously successful
infections will be prevented. Therefore, the new _effective reproduction
number_ is equal to $R_{e}(\eta^{\mathrm{uni}})=\eta^{\mathrm{uni}}R_{0}$.
This fact led to the recognition by Smith in 1970 [smith_prospects_1970] and
Dietz in 1975 [dietz1975transmission] of a simple threshold theorem: the
incidence of an infection declines if the proportion of non-immune individuals
is reduced below $\eta^{\mathrm{uni}}_{\mathrm{crit}}=1/R_{0}$. This effect is
called _herd immunity_ , and the corresponding proportion
$1-\eta^{\mathrm{uni}}_{\mathrm{crit}}$ of people that have to be vaccinated
is called the _herd immunity threshold_ [smith_concepts_2010,
somerville_public_2016].

##  2\. Critical vaccination given by the endemic equilibrium

However, herd immunity can also be achieved using a non-uniform vaccination
strategy when the population is heterogeneous. For example, the discussion of
vaccination control of gonorrhea in [hethcote, Section 4.5] suggests that it
may be better to prioritize the vaccination of people that have already caught
the disease: this leads us to consider a vaccination strategy guided by the
equilibrium state. For the SIS model in heterogeneous population with
$R_{0}>1$, there exists a maximal endemic equilibrium, say $\mathfrak{g}$,
where $\mathfrak{g}(x)$ represents the fraction of infected people in the
group with feature $x$. In other words, the function $\mathfrak{g}$ is the
maximal $[0,1]$-valued solution $g$ of:

(2)
$(1-g(x))\int_{\Omega}k(x,y)\,g(y)\,\mu(\mathrm{d}y)=\gamma(x)g(x)\quad\text{for
$x\in\Omega$.}$

Let us mention that if there exist isolated subpopulations, it is possible to
have other endemic equilibria, i.e., solutions to Equation (2) that are not
equal to $0$ for all $x$. Irreducibility conditions on the kernel $k$ ensure
however the uniqueness of the endemic equilibrium
[lajmanovich1976deterministic, delmas_infinite-dimensional_2020]. Consider the
vaccination strategy, denoted by $\eta^{\mathrm{equi}}$, corresponding to
vaccinating a fraction $\mathfrak{g}(x)$ of people in the group with feature
$x$, for all groups. In our mathematical framework, this amounts to setting:

(3) $\eta^{\mathrm{equi}}(x)=1-\mathfrak{g}(x)\quad\text{for $x\in\Omega$.}$

The following result ensures that this strategy reaches herd immunity; see
Theorem 4.1 in Section 4 for a precise mathematical statement.

######  Theorem.

In the heterogeneous SIS model with non-zero maximal endemic equilibrium, the
vaccination strategy $\eta^{\mathrm{equi}}$ is critical:

(4) $R_{e}(\eta^{\mathrm{equi}})=1.$

Let us stress that implementing the critical vaccination strategy
$\eta^{\mathrm{equi}}$ can be achieved without estimating the transmission
rate kernel and the recovery rate.

The proof of the theorem relies on the study of the spectral bound of the
linearized operator associated to equation (1) near an equilibrium. When
$R_{0}>1$, this spectral bound is non-positive at the maximal equilibrium and
positive at all other equilibria; see Proposition 4.4 ii. Thus, the non-
maximality of an equilibrium is equivalent to its linear instability in the
original dynamics. We also prove the linear instability of the disease-free
state in the modified dynamics where we vaccinate according to a non maximal
equilibrium; see Proposition 4.4 iv.

##  3\. Discussion

We expect the results obtained here for the SIS model to be generic, in the
sense that similar behaviours should also be observed in more realistic and
complex models in epidemiology for non-homogeneous populations: when an
endemic equilibrium exists, vaccinating the population according to the
maximal endemic profile should protect the population from the disease.

We refer to [ddz-theo] for a general framework for cost comparison of
vaccination strategies and the notions of “best” and “worst” vaccination
strategies; see also [ddz-Re, ddz-cordon, ddz-reg] for further comments and
various examples of optimal vaccinations.

Consider a general cost function $C$ which measures the cost for the society
of a vaccination strategy (production and diffusion). A simple and natural
choice is the uniform cost given by the overall proportion of vaccinated
individuals:

(5)
$C(\eta)=\int_{\Omega}(\mathbb{1}-\eta)\,\mathrm{d}\mu=1-\int_{\Omega}\eta\,\mathrm{d}\mu.$

We have that $C(\eta^{\mathrm{uni}}_{\mathrm{crit}})$ is equal to the herd
immunity threshold $1-1/R_{0}$ while $C(\eta^{\mathrm{equi}})$ is equal to the
proportion of people in the endemic state in an SIS infection
$\int_{\Omega}\mathfrak{g}\,\mathrm{d}\mu$. It is not possible to determine
which strategy is cheaper in general. However, in the following examples, we
are able to compare their costs for mixing structures that are redundant in
the epidemiologic litterature.

######  Example 3.1 (Homogeneous mixing).

If the population is homogeneous (which corresponds to the one-dimensional SIS
model where $\Omega$ is a singleton), then the maximal equilibrium is constant
equal to $1-1/R_{0}$. It follows that
$C(\eta^{\mathrm{uni}}_{\mathrm{crit}})=C(\eta^{\mathrm{equi}})$.

######  Example 3.2 (Proportionate mixing structure with two subpopulations).

The proportionate mixing is a classical mixing structure introduced by
[HeterogeneityINold1980] and used in many different epidemiological models. It
assumes that the number of adequate contacts between two subpopulations is
proportional to the relative activity levels of the two subpopulations. Thus
individuals in more active subpopulations will have more adequate contacts.
Let us consider the simple case where there are only two subpopulations. Then
the contact matrix is given by:

$K=\begin{pmatrix}a^{2}&ab\\\ ab&b^{2}\end{pmatrix}$

where $a$ and $b$ are positive constants that correspond to the activity
levels of the first and second subpopulations respectively. Denote by
$\mu_{1}$ and $\mu_{2}$ their respective relative size, suppose that the
recovery rate $\gamma$ is equal to $1$ for both subpopulations, and assume
without loss of generality that $a\geq b$. In this case, we get that:

$R_{0}=a^{2}\mu_{1}+b^{2}\mu_{2},\quad
R_{e}(\eta)=a^{2}\eta_{1}\mu_{1}+b^{2}\eta_{2}\mu_{2}\quad\text{and}\quad
C(\eta)=1-(\eta_{1}\mu_{1}+\eta_{2}\mu_{2})$

for the vaccination strategy $\eta=(\eta_{1},\eta_{2})$. If $R_{0}>1$, then
the (unique) non-zero equilibrium satisfies:

$(1-\mathfrak{g}_{i})\,\sum_{j=1}^{2}K_{i,j}\,\mathfrak{g}_{j}\,\mu_{j}=\mathfrak{g}_{i}\quad\text{for
$i=1,2$},$

and the corresponding vaccination strategy
$\eta^{\mathrm{equi}}=1-\mathfrak{g}$ is given by:

$\eta^{\mathrm{equi}}=\left(\frac{1}{1+ac},\frac{1}{1+bc}\right),$

where $c\in[(1-R_{0})/a,(1-R_{0})/b]$ is the unique positive solution of the
second-order equation $R_{e}(\eta^{\mathrm{equi}})=1$. It is elementary to
check that in this case $C(\eta^{\mathrm{equi}})<C(\eta^{\mathrm{uni}})$, with
an equality if and only if $a=b$. However, the critical vaccination strategy
with minimal cost, say $\eta^{\mathrm{opt}}$, corresponds to vaccinating in
priority the population with the highest activity rate, that is, if $a>b$:

$\eta^{\mathrm{opt}}=\left(\frac{1-\min(1,b^{2}\mu_{2})}{a^{2}\mu_{1}},\frac{1}{\max(1,b^{2}\mu_{2})}\right).$

######  Example 3.3 (Age and activity structure).

In [AMathematicalBritto2020], Britton, Ball and Trapman study an SEIR model,
where immunity can be obtained through infection. Using parameters derived
from real-world data, these authors noticed that the disease-induced herd
immunity level can, for some models, be substantially lower than the classical
herd immunity threshold $1-1/R_{0}$. This can be reformulated in term of
targeted vaccination strategies: prioritizing the individuals that are more
likely to get infected in a SEIR epidemic may be more efficient than
distributing uniformly the vaccine in the population.

We use the same age and activity structures to determine which strategy
between $\eta^{\mathrm{equi}}$ and $\eta^{\mathrm{uni}}_{\mathrm{crit}}$ is
less costly. More precisely, the community is categorized into six age groups
and contact rates between them are derived from an empirical study of social
contacts [UsingDataOnSWallin2006]. For the activity structure, individuals are
categorized into three different activity levels, which are arbitrary and
chosen for illustration purposes: 50% of each age cohort have normal activity,
25% have low activity corresponding to half as many contacts compared with
normal activity, and 25% have high activity corresponding to twice as many
contacts as normal activity. Note that when the population is only structured
by activity, the mixing is proportionate. Assuming that the recovery rate is
constant equal to $1$, we solved numerically Equation (2) and computed in
Table 1 the cost of the uniform and the equilibrium strategies for different
values of $R_{0}$ and different population structures. In Table 2, the
fractions of vaccinated individuals in the different age activity groups when
following the strategy $\eta^{\mathrm{equi}}$. This is done by assuming
$R_{0}=2$. Note than in this case, only three subpopulations need to be
vaccinated at a level higher than $1-1/R_{0}$.

Table 1. Cost of the equilibrium vaccination compared to the herd immunity level for different population structures. Numbers correspond to percentage. | $R_{0}=2$ |  | $R_{0}=2.5$ |  | $R_{0}=3$
---|---|---|---|---|---
| $C(\eta^{\mathrm{equi}})$ | $C(\eta^{\mathrm{uni}})$ |  | $C(\eta^{\mathrm{equi}})$ | $C(\eta^{\mathrm{uni}})$ |  | $C(\eta^{\mathrm{equi}})$ | $C(\eta^{\mathrm{uni}})$
Homogeneous | $50$ | $50$ |  | $60$ | $60$ |  | $66.7$ | $66.7$
Age structure | $46.6$ | $50$ |  | $56.7$ | $60$ |  | $63.9$ | $66.7$
Activity structure | $40.1$ | $50$ |  | $50$ | $60$ |  | $57$ | $66.7$
Age and activity structure | $35.7$ | $50$ |  | $45.2$ | $60$ |  | $52.2$ | $66.7$

Table 2. Fraction of vaccinated individuals in different groups for the strategy $\eta^{\mathrm{equi}}$. The population structure includes both age and activity. Numbers correspond to percentage. These values assume that $R_{0}=2$, so that the uniform critical vaccination consists in vaccinating $50\%$ of the population: only three groups require more vaccine in the targeted strategy than in the uniform strategy. Age group | Low activity | Average activity | High activity
---|---|---|---
0–5 years | $12.0$ | $21.4$ | $35.3$
6–12 years | $18.5$ | $31.2$ | $47.5$
13–19 years | $22.9$ | $37.3$ | $54.3$
20–39 years | $29.1$ | $45.1$ | $62.1$
40–59 years | $20.9$ | $34.6$ | $51.4$
$\geq$ 60 years | $12.4$ | $22.1$ | $36.2$

##  4\. Proof

The differential equations governing the epidemic dynamics in meta-population
SIS models were developed in paper [lajmanovich1976deterministic] in finite
dimension and generalized in [delmas_infinite-dimensional_2020].

###  4.1. The heterogeneous SIS model

Let $(\Omega,\mathscr{F},\mu)$ be a probability space, where $x\in\Omega$
represents a feature and the probability measure $\mu(\mathrm{d}x)$ represents
the fraction of the population with feature $x$. The parameters of the SIS
model are given by a _recovery rate function_ $\gamma$, which is a positive
bounded measurable function defined on $\Omega$, and a _transmission rate
kernel_ $k$, which is a non-negative measurable function defined on
$\Omega^{2}$.

In accordance with [delmas_infinite-dimensional_2020], we consider for a
kernel $\mathrm{k}$ on $\Omega$ and $q\in(1,+\infty)$ its norm:

$\left\lVert\,\mathrm{k}\,\right\rVert_{\infty,q}=\sup_{x\in\Omega}\,\left(\int_{\Omega}\mathrm{k}(x,y)^{q}\,\mu(\mathrm{d}y)\right)^{1/q}.$

For a kernel $\mathrm{k}$ on $\Omega$ such that
$\left\lVert\,\mathrm{k}\,\right\rVert_{\infty,q}$ is finite for some
$q\in(1,+\infty)$, we define the integral operator $\mathcal{T}_{\mathrm{k}}$
on the set $\mathscr{L}^{\infty}$ of bounded measurable real-valued function
on $\Omega$ by:

(6)
$\mathcal{T}_{\mathrm{k}}(g)(x)=\int_{\Omega}\mathrm{k}(x,y)g(y)\,\mu(\mathrm{d}y)\quad\text{for
}g\in\mathscr{L}^{\infty}\text{ and }x\in\Omega.$

By convention, for $f,g$ two non-negative measurable functions defined on
$\Omega$ and $\mathrm{k}$ a kernel on $\Omega$, we denote by $f\mathrm{k}g$
the kernel on $\Omega$ defined by:

(7) $f\mathrm{k}g:(x,y)\mapsto f(x)\,\mathrm{k}(x,y)g(y).$

We shall consider the kernel ${\bm{\mathrm{k}}}=k/\gamma$ (corresponding to
$k\gamma^{-1}$, which differs in general from $\gamma^{-1}k$), which is thus
defined by:

(8) $\boxed{{\bm{\mathrm{k}}}(x,y)=k(x,y)\,\gamma(y)^{-1}.}$

We shall assume that:

(9)
$\boxed{\left\lVert\,{\bm{\mathrm{k}}}\,\right\rVert_{\infty,q}<\infty\quad\text{for
some $q\in(1,+\infty)$.}}$

The integral operator $\mathcal{T}_{{\bm{\mathrm{k}}}}$ is the so called
_next-generation operator_.

Let $\Delta=\\{f\in\mathscr{L}^{\infty}\,\colon\,0\leq f\leq 1\\}$ be the
subset of non-negative functions bounded by $1$, and let $\mathbb{1}\in\Delta$
be the constant function equal to 1. The SIS dynamics considered in
[delmas_infinite-dimensional_2020] follows the vector field $F$ defined on
$\Delta$ by:

(10) $\boxed{F(g)=(\mathbb{1}-g)\mathcal{T}_{k}(g)-\gamma g.}$

More precisely, we consider $u=(u_{t},t\in\mathbb{R})$, where $u_{t}\in\Delta$
for all $t\in\mathbb{R}_{+}$ such that:

(11) $\boxed{\partial_{t}u_{t}=F(u_{t})}\quad\text{for }t\in\mathbb{R}_{+},$

with initial condition $u_{0}\in\Delta$. The value $u_{t}(x)$ models the
probability that an individual of feature $x$ is infected at time $t$; it is
proved in [delmas_infinite-dimensional_2020] that such a solution $u$ exists
and is unique.

An _equilibrium_ of (11) is a function $g\in\Delta$ such that $F(g)=0$.
According to [delmas_infinite-dimensional_2020], there exists a maximal
equilibrium $\mathfrak{g}$, i.e., an equilibrium such that all other
equilibria $h\in\Delta$ are dominated by $\mathfrak{g}$: $h\leq\mathfrak{g}$.
It is to this maximal equilibrium that the process stabilizes when started
from a situation where all the population is infected, that is, if
$u_{0}=\mathbb{1}$, then we have:

$\boxed{\lim_{t\rightarrow\infty}u_{t}=\mathfrak{g}.}$

For $T$ a bounded operator on $\mathscr{L}^{\infty}$ endowed with its usual
supremum norm, we define by
$\left\lVert\,T\,\right\rVert_{\mathscr{L}^{\infty}}$ its operator norm and
its spectral radius is given by:

$\rho(T)=\lim_{n\rightarrow\infty}\left\lVert\,T^{n}\,\right\rVert_{\mathscr{L}^{\infty}}^{1/n}.$

The _reproduction number_ $R_{0}$ associated to the SIS model given by (11) is
the spectral radius of the next-generation operator:

(12) $\boxed{R_{0}=\rho(\mathcal{T}_{{\bm{\mathrm{k}}}}),}$

If $R_{0}\leq 1$ (sub-critical and critical case), then $u_{t}$ converges
pointwise to $0$ when $t\to\infty$. In particular, the maximal equilibrium
$\mathfrak{g}$ is equal to $0$ everywhere. If $R_{0}>1$ (super-critical case),
then $0$ is still an equilibrium but different from the maximal equilibrium
$\mathfrak{g}$, as $\int_{\Omega}\mathfrak{g}\,\mathrm{d}\mu>0$.

###  4.2. Vaccination strategies

A _vaccination strategy_ $\eta$ of a vaccine with perfect efficiency is an
element of $\Delta$, where $\eta(x)$ represents the proportion of _non-
vaccinated_ individuals with feature $x$. Notice that $\eta\,\mathrm{d}\mu$
corresponds in a sense to the effective population. In particular, the
“strategy” that consists in vaccinating no one (resp. everybody) corresponds
to $\eta=\mathbb{1}$, the constant function equal to 1, (resp.
$\eta=\mathbb{0}$, the constant function equal to 0).

Recall the definition of the kernel $f{\bm{\mathrm{k}}}g$ from (7). For
$\eta\in\Delta$, the kernel $\mathrm{k}\eta=k\eta/\gamma$ has finite norm
$\left\lVert\,\cdot\,\right\rVert_{\infty,q}$, so we can consider the bounded
positive operators $\mathcal{T}_{{\bm{\mathrm{k}}}\eta}$ and
$\mathcal{T}_{k\eta}$ on $\mathscr{L}^{\infty}$. According to
[delmas_infinite-dimensional_2020, Section 5.3.], the SIS equation with
vaccination strategy $\eta$ is given by (11), where $F$ is replaced by
$F_{\eta}$ defined by:

(13) $\boxed{F_{\eta}(g)=(\mathbb{1}-g)\mathcal{T}_{k\eta}(g)-\gamma g.}$

We denote by $u^{\eta}=(u^{\eta}_{t},t\geq 0)$ the corresponding solution with
initial condition $u_{0}^{\eta}\in\Delta$. We recall that $u_{t}^{\eta}(x)$
represents the probability for an non-vaccinated individual of feature $x$ to
be infected at time $t$. We define the _effective reproduction number_
$R_{e}(\eta)$ associated to the vaccination strategy $\eta$ as the spectral
radius of the effective next-generation operator
$\mathcal{T}_{{\bm{\mathrm{k}}}\eta}$:

(14) $\boxed{R_{e}(\eta)=\rho(\mathcal{T}_{{\bm{\mathrm{k}}}\eta}).}$

For example, for the trivial vaccination strategies we get
$R_{e}(\mathbb{1})=R_{0}$ and $R_{e}(\mathbb{0})=0$.

We denote by $\mathfrak{g}_{\eta}$ the corresponding maximal equilibrium. In
particular, we have:

(15) $F_{\eta}(\mathfrak{g}_{\eta})=0.$

In particular, we have:

$R_{e}(\mathbb{1})=R_{0}\quad\text{and}\quad\mathfrak{g}=\mathfrak{g}_{\mathbb{1}}.$

###  4.3. Critical vaccination strategies

If $R_{0}\geq 1$, then a vaccination strategy $\eta$ is called _critical_ if
it achieves precisely herd immunity, that is $R_{e}(\eta)=1$.

As the spectral radius is positively homogeneous (that is, $\rho(\lambda
t)=\lambda\rho(T)$ for $\lambda\geq 0$), we also get, when $R_{0}\geq 1$, that
the uniform strategy that corresponds to the constant function:

$\boxed{\eta^{\mathrm{uni}}_{\mathrm{crit}}=\frac{1}{R_{0}}\mathbb{1}}$

is critical, as $R_{e}(\eta^{\mathrm{uni}}_{\mathrm{crit}})=1$. This is
consistent with results obtained in the homogeneous model.

As hinted in [hethcote, Section 4.5] for vaccination control of gonorrhea, it
is interesting to consider vaccinating people with feature $x$ with
probability $\mathfrak{g}(x)$. This corresponds to the strategy based on the
maximal equilibrium:

$\boxed{\eta^{\mathrm{equi}}=\mathbb{1}-\mathfrak{g}.}$

The following result entails that this strategy is critical and thus achieves
herd immunity. Recall that in the (infinite dimensional) SIS model (10) on the
probability space $(\Omega,\mathscr{F},\mu)$ the recovery rate function
$\gamma$ is positive and bounded, the transmission rate $k$ is non-negative
and the norm $\left\lVert\,{\bm{\mathrm{k}}}\,\right\rVert_{\infty,q}$ of the
kernel ${\bm{\mathrm{k}}}=k/\gamma$ is finite for some $q\in(1,+\infty)$.

######  Theorem 4.1 (The maximal equilibrium yields a critical vaccination).

Consider the SIS model (10) under the boundedness assumption (9). If
$R_{0}\geq 1$, then the vaccination strategy $\eta^{\mathrm{equi}}$ is
critical, that is, $R_{e}(\eta^{\mathrm{equi}})=1$.

This result will be proved below as a part of Proposition 4.4.

###  4.4. Proof of Theorem 4.1

For an operator $A$, we denote by $A^{\top}$ its adjoint. We first give a
preliminary lemma. For the convenience of the reader, we only use references
to the results recalled in [delmas_infinite-dimensional_2020] for positive
operators on Banach spaces. In particular, we shall use that if $\mathrm{k}$
and $\mathrm{k}^{\prime}$ are two (non-negative) kernels on $\Omega$ with
finite norms $\left\lVert\,\cdot\,\right\rVert_{\infty,q}$ for some
$q\in(1,+\infty)$, then we have that:

(16)
$\mathrm{k}\geq\mathrm{k}^{\prime}\,\Longrightarrow\,\rho(\mathcal{T}_{\mathrm{k}})\geq\rho(\mathcal{T}_{\mathrm{k}^{\prime}}),$

see for example [delmas_infinite-dimensional_2020, Theorem 3.5(i)] as the
operator $\mathcal{T}_{\mathrm{k}}-\mathcal{T}_{\mathrm{k}^{\prime}}$ is
positive. We shall also used that for two bounded operators $T$ and $S$ on
$\mathscr{L}^{\infty}$:

(17) $\rho(TS)=\rho(ST).$

We first state two technical lemmas.

######  Lemma 4.2.

Let $\mathrm{k}$ be a non-negative kernel on $\Omega$ such that
$\left\lVert\,\mathrm{k}\,\right\rVert_{\infty,q}$ is finite for some
$q\in(1,+\infty)$ and and consider the positive bounded linear integral
operator $\mathcal{T}_{\mathrm{k}}$ on $\mathscr{L}^{\infty}$. If there exists
$g\in\mathscr{L}^{\infty}_{+}$, with $\int_{\Omega}g\,\mathrm{d}\mu>0$ and
$\lambda>0$ satisfying:

$\mathcal{T}_{\mathrm{k}}(g)(x)>\lambda g(x),\quad\text{for all}\quad x\text{
such that }g(x)>0,$

then we have $\rho(\mathcal{T}_{\mathrm{k}})>\lambda$.

###### Proof.

We simply write $\mathcal{T}$ for $\mathcal{T}_{\mathrm{k}}$. Let
$A=\set{g>0}$ be the support of the function $g$. Let $\mathcal{T}^{\prime}$
be the bounded operator defined by
$\mathcal{T}^{\prime}(f)=\mathbb{1}_{A}\mathcal{T}(\mathbb{1}_{A}f)$. Since
$\mathcal{T}^{\prime}(g)=\mathbb{1}_{A}\mathcal{T}(\mathbb{1}_{A}g)=\mathbb{1}_{A}\mathcal{T}(g)>\lambda
g$, we deduce from the Collatz-Wielandt formula [delmas_infinite-
dimensional_2020, Proposition 3.6] that
$\rho(\mathcal{T}^{\prime})\geq\lambda>0$. According to [delmas_infinite-
dimensional_2020, Lemma 3.7 (v)], there exists $v\in
L^{q}_{+}\setminus\\{0\\}$, seen as an element of the topological dual of
$\mathscr{L}^{\infty}$, a left Perron eigenfunction of $\mathcal{T}^{\prime}$,
that is such that
$(\mathcal{T}^{\prime})^{\top}(v)=\rho(\mathcal{T}^{\prime})v$. In particular,
we have $v=\mathbb{1}_{A}\,v$ and thus $\int_{A}v\,\mathrm{d}\mu>0$ and
$\int_{\Omega}vg\,\mathrm{d}\mu>0$. We obtain:

$(\rho(\mathcal{T}^{\prime})-\lambda)\braket{v,g}=\braket{v,\mathcal{T}^{\prime}(g)-\lambda
g}>0.$

As $\mathcal{T}^{\prime}=\mathcal{T}_{\mathrm{k}^{\prime}}$ with
$\mathrm{k}^{\prime}=\mathbb{1}_{A}\mathrm{k}\mathbb{1}_{A}\leq\mathrm{k}$, we
deduce from (16) that
$\rho(\mathcal{T})\geq\rho(\mathcal{T}^{\prime})>\lambda$. ∎

######  Lemma 4.3.

Consider the SIS model (10) under the boundedness assumption (9). Let
$\eta,g\in\Delta$. If $F_{\eta}(g)\geq 0$, then we have
$g\leq\mathfrak{g}_{\eta}$.

###### Proof.

Consider the solution $u_{t}$ of the SIS model
$\partial_{t}u_{t}=F_{\eta}(u_{t})$ with vaccination $\eta$ and initial
condition $u_{0}=g$. According to [delmas_infinite-dimensional_2020,
Proposition 2.10], this solution is non-decreasing since $F_{\eta}(g)\geq 0$.
According to [delmas_infinite-dimensional_2020, Proposition 2.13], the
pointwise limit of $u_{t}$ is an equilibrium. As this limit is dominated by
the maximal equilibrium $\mathfrak{g}_{\eta}$ and since $u_{t}$ is non-
decreasing, this proves that $g\leq\mathfrak{g}_{\eta}$. ∎

The next result characterizes the maximal equilibrium $\mathfrak{g}$ among all
equilibria by various spectral properties; Theorem 4.1 may be viewed as a
corollary to this characterization. Recall that $R_{0}=R_{e}(\mathbb{1})$, and
that the vector field $F$ is defined by (10). Let $DF[h]$ denote the bounded
linear operator on $\mathscr{L}^{\infty}$ of the derivative of the map
$f\mapsto F(f)$ defined on $\mathscr{L}^{\infty}$ at point $h$:

$DF[h](g)=(1-h)\mathcal{T}_{k}(g)-(\gamma+\mathcal{T}_{k}(h))g\quad\text{for
$h,g\in\mathscr{L}^{\infty}$.}$

Let $s(A)$ denote the spectral bound of the bounded operator $A$, see Equation
(33) in [delmas_infinite-dimensional_2020].

######  Proposition 4.4 (Equivalent conditions for maximality).

Consider the SIS model (10) under the boundedness assumption (9). Let $h$ in
$\Delta$ be an equilibrium, that is, $F(h)=0$. The following properties are
equivalent:

  1. (i)

$h=\mathfrak{g}$,

  2. (ii)

$s(DF[h])\leq 0$,

  3. (iii)

$R_{e}((1-h)^{2})\leq 1$.

  4. (iv)

$s(DF_{(1-h)}[\mathbb{0}])\leq 0$.

  5. (v)

$R_{e}(1-h)\leq 1$.

Furthermore, $\mathfrak{g}=\mathbb{0}$ if and only if $R_{0}\leq 1$, and if
$\mathfrak{g}\neq\mathbb{0}$, then it is critical: $R_{e}(1-\mathfrak{g})=1$.

######  Remark 4.5 (On stability).

From a dynamical systems point of view, this proposition links together two
different stability properties. The (classically equivalent) conditions ii and
iii state that for the original dynamics given by (11) with vector field $F$,
the equilibrium $h$ is not linearly unstable. Similarly, conditions iv and v
both state that in the vaccinated dynamics given by the modified vector field
$F_{1-h}$ defined by (13), the disease-free equilibrium $0$ is not linearly
unstable.

In particular, in the original dynamics given by (11), equilibria that are not
maximal are necessarily linearly unstable.

###### Proof.

Let $h\in\Delta$ be an equilibrium, that is $F(h)=0$.

The equivalence between iv and v is a direct consequence of [delmas_infinite-
dimensional_2020, Proposition 4.2].

Let us show the equivalence between ii and iii. According to the same
[delmas_infinite-dimensional_2020, Proposition 4.2], $s(DF[h])\leq 0$ if and
only if:

$\rho\left(\mathcal{T}_{\mathrm{k}}\right)\leq
1\quad\text{with}\quad\mathrm{k}(x,y)=(1-h(x))\frac{k(x,y)}{\gamma(y)+\mathcal{T}_{k}(h)(y)}\cdot$

Since $F(h)=0$, we have $(1-h)/\gamma=1/(\gamma+\mathcal{T}_{k}(h))$. This
gives:

(18) $\mathsf{k}(x,y)=(1-h(x))\frac{k(x,y)(1-h(y))}{\gamma(y)}$

and thus $\mathcal{T}_{\mathrm{k}}=M_{1-h}\,\mathcal{T}_{k/\gamma}\,M_{1-h}$,
where $M_{f}$ is the multiplication operator by $f$. Recall the definition
(14) of $R_{e}$. We deduce from (17) that:

(19)
$\rho\left(\mathcal{T}_{\mathrm{k}}\right)=\rho\left(\mathcal{T}_{k/\gamma}M_{(1-h)^{2}}\right)=R_{e}((1-h)^{2}).$

This gives the equivalence between ii and iii.

We prove that i implies v. Suppose that $R_{e}(1-h)>1$. Thanks to (17), we
have
$\rho(M_{1-h}\mathcal{T}_{k/\gamma})=\rho(\mathcal{T}_{k/\gamma}M_{1-h})=R_{e}(1-h)>1$.
According to [delmas_infinite-dimensional_2020, Lemma 3.7 (v)], there exists
$v\in L^{q}_{+}\setminus\\{0\\}$ a left Perron eigenfunction of
$\mathcal{T}_{(1-h)k/\gamma}$, that is
$\mathcal{T}_{(1-h)k/\gamma}^{\top}(v)=R_{e}(1-h)v$. Using $F(h)=0$, and thus
$(1-h)\mathcal{T}_{k}(h)=\gamma h$, for the last equality, we have:

$R_{e}(1-h)\braket{v,\gamma h}=\braket{v,(1-h)\mathcal{T}_{k/\gamma}(\gamma
h)}=\braket{v,\gamma h}.$

We get $\braket{v,\gamma h}=0$ and thus $\braket{v,\mathbb{1}_{A}}=0$, where
$A=\set{h>0}$ denote the support of the function $h$. Since
$\mathcal{T}_{(1-h)k/\gamma}^{\top}(v)=R_{e}(1-h)v$ and setting
$v^{\prime}=(1-h)v$ (so that $v^{\prime}=v$ $\mu$-almost surely on $A^{c}$),
we deduce that:

$\mathcal{T}_{k^{\prime}/\gamma}^{\top}(v^{\prime})=R_{e}(1-h)v^{\prime},$

where $k^{\prime}=\mathbb{1}_{A^{c}}\,k\,\mathbb{1}_{A^{c}}$. This implies
that $\rho(\mathcal{T}_{k^{\prime}/\gamma})\geq R_{e}(1-h)$. Since
$k^{\prime}=(1-h)k^{\prime}$ and $k-k^{\prime}\geq 0$, we get that
$\mathcal{T}_{k/\gamma}-\mathcal{T}_{k^{\prime}/\gamma}$ is a positive
operator. Using (16) for the inequality as
$(1-h)k^{\prime}/\gamma\leq(1-h)k/\gamma$, we deduce that
$\rho(\mathcal{T}_{k^{\prime}/\gamma})=\rho(M_{1-h}\mathcal{T}_{k^{\prime}/\gamma})\leq\rho(M_{1-h}\mathcal{T}_{k/\gamma})=R_{e}(1-h)$.
Thus, the spectral radius of $\mathcal{T}_{k^{\prime}/\gamma}$ is equal to
$R_{e}(1-h)$. According to [delmas_infinite-dimensional_2020, Proposition
4.2], since $\rho(\mathcal{T}_{k^{\prime}/\gamma})>1$, there exists
$w\in\mathscr{L}^{\infty}_{+}\setminus\\{0\\}$ and $\lambda>0$ such that:

$\mathcal{T}_{k^{\prime}}(w)-\gamma w=\lambda w.$

This also implies that $w=0$ on $A=\set{h>0}$, that is $wh=0$ and thus
$w\mathcal{T}_{k}(h)=0$ as $\mathcal{T}_{k}(h)=\gamma h/(1-h)$. Using that
$F(h)=0$, $\mathcal{T}_{k}(w)=\mathcal{T}_{k^{\prime}}(w)=(\gamma+\lambda)w$
and $h\mathcal{T}_{k}(w)=0$, we obtain:

$F(h+w)=w(\lambda-\mathcal{T}_{k}(w)).$

Taking $\varepsilon>0$ small enough so that
$\varepsilon\mathcal{T}_{k}(w)\leq\lambda/2$ and $\varepsilon w\leq 1$, we get
$h+\varepsilon w\in\Delta$ and $F(h+\varepsilon w)\geq 0$. Then use Lemma 4.3
to deduce that $h+\varepsilon w\leq\mathfrak{g}$ and thus $h\neq\mathfrak{g}$.

To see that v implies iii, notice that $(1-h)\geq(1-h)^{2}$, and then use (16)
to deduce that
$\rho(\mathcal{T}_{{\bm{\mathrm{k}}}(1-h)})\geq\rho(\mathcal{T}_{{\bm{\mathrm{k}}}(1-h)^{2}})$
and thus $R_{e}(1-h)\geq R_{e}((1-h)^{2})$.

We prove that iii implies i. Notice that $F(g)=0$ and $g\in\Delta$ implies
that $g<1$. Assume that $h\neq\mathfrak{g}$. Notice that
$\gamma/(1-h)=\gamma+\mathcal{T}_{k}(h)$, so that
$\gamma(\mathfrak{g}-h)/(1-h)\in\mathscr{L}^{\infty}_{+}$. An elementary
computation, using $F(h)=F(\mathfrak{g})=0$ and $\mathrm{k}$ defined in (18),
gives:

$\mathcal{T}_{\mathrm{k}}\left(\gamma\frac{\mathfrak{g}-h}{1-h}\right)=(1-h)\mathcal{T}_{k}\left(\mathfrak{g}-h\right)=\gamma\frac{\mathfrak{g}-h}{1-\mathfrak{g}}=\frac{1-h}{1-\mathfrak{g}}\,\gamma\frac{\mathfrak{g}-h}{1-h}\cdot$

Since $h\neq\mathfrak{g}$ and $h\leq\mathfrak{g}$, we deduce that
$(1-h)/(1-\mathfrak{g})\geq 1$, with strict inequality on
$\set{\mathfrak{g}-h>0}$ which is a set of positive measure. We deduce from
Lemma 4.2 (with $k$ replaced by $\mathrm{k}\gamma$) that
$\rho\left(\mathcal{T}_{\mathsf{k}}\right)>1$. Then use (19) to conclude.

To conclude notice that $\mathfrak{g}=0\Longleftrightarrow R_{0}\leq 1$ is a
consequence of the equivalence between i and v with $h=0$ and
$R_{0}=R_{e}(\mathbb{1})$.

Using that $F(\mathfrak{g})=0$, we get
$\mathcal{T}_{k}(\mathfrak{g})=\gamma\mathfrak{g}/(1-\mathfrak{g})$. We deduce
that
$\mathcal{T}_{k(1-\mathfrak{g})/\gamma}(\mathcal{T}_{k}(\mathfrak{g}))=\mathcal{T}_{k}(\mathfrak{g})$.
If $\mathfrak{g}\neq 0$, we get $\mathcal{T}_{k}(\mathfrak{g})\neq 0$ (on a
set of positive $\mu$-measure). This implies that $R_{e}(1-\mathfrak{g})\geq
1$. Then use v to deduce that $R_{e}(1-\mathfrak{g})=1$ if $\mathfrak{g}\neq
0$. ∎

## Acknowledgments

This work is partially supported by Labex Bézout reference ANR-10-LABX-58.


[1]organization=IPPP, addressline=Department of Physics, Durham University,
city=Durham, postcode=DH1 3LE, country=UK

# Optical effects of domain walls

Valentin V. Khoze<EMAIL_ADDRESS>Daniel L. Milne
<EMAIL_ADDRESS>
###### Abstract

Domain walls arise in theories where there is spontaneous symmetry breaking of
a discrete symmetry such as $\mathbb{Z}_{N}$ and are a feature of many BSM
models. In this work we consider the possibility of detecting domain walls
through their optical effects and specify three different methods of coupling
domain walls to the photon. We consider the effects of these couplings in the
context of gravitational wave detectors, such as LIGO, and examine the
sensitivity of these experiments to domain wall effects. In cases where
gravitational wave detectors are not sensitive we examine our results in the
context of axion experiments and show how effects of passing domain walls can
be detected at interferometers searching for an axion signal.

##  1 Introduction

It was pointed out in Ref. [1] that gravitational wave detectors, such as
LIGO, can be used to search for potential signals of cosmological domain walls
[2]. The LIGO detector was built with the purpose of testing general
relativity by searching for gravitational waves [3]. However, it was shown by
Jaeckel et al. that domain walls could leave an optical signal at an
interferometer-based experiment that is not caused by gravity but is triggered
instead by an effective interaction between the domain wall and the photon.
Domain walls arise naturally in theories when a discrete symmetry is broken;
such discrete symmetries occur frequently in many beyond the Standard Model
(BSM) settings in particle phyisics and are readily seen in the literature.
Simple $\mathbb{Z}_{2}$ symmetries are popular in the context of two Higgs
doublet models [4, 5] and aspects of domain walls in these models were
specifically studied in [6, 7]. Many other discrete symmetries arise in many
different areas of physics such as flavour physics [8] and dark matter [9].

In this paper we extend on the work of [1] by considering different types of
interactions between the domain walls sector and QED and also accounting for
the effect of cavities in interferometers. We couple the domain wall to the
photon in three different ways and show that in each case the signal in an
interferometer is quite different and, in fact, an experiment sensitive to one
type of coupling is not in general sensitive to the others. It should be
emphasised that many previous papers have considered gravitational waves in
the context of domain walls, see e.g. Refs. [10, 11, 12], but in general these
earlier studies considered gravitational waves sourced by annihilations
between domain walls. In this paper we do not examine gravitational waves
themselves but the signal at a gravitational wave interferometer which is
actually generated by the interaction of the domain wall and the photon and
not by the distortion of the Minkowski metric.

For the domain wall itself we will use the generic simple model of Ref. [13],

$\mathcal{L}_{\rm
DW}=\frac{1}{2}\left(\partial_{\mu}\phi\right)^{2}-\frac{2m^{2}f^{2}}{N_{\phi}^{2}}\sin^{2}\left(\frac{N_{\phi}\phi}{2f}\right)$
(1)

which has domain wall solutions,

$\phi_{\rm cl}\left(z\right)=\frac{4f}{N_{\phi}}\tan^{-1}\left(e^{mz}\right).$
(2)

The Lagrangian (1) is invariant under the shift symmetry
$\phi\to\phi+\frac{2\pi}{N_{\phi}}\,f$ which is the consequence of the
$\mathbb{Z}_{N_{\phi}}$ symmetry of the underlying theory of a complex scalar
field $\Phi(x)=S(x)\,e^{i\phi(x)/f}$. Since the domain wall field $\phi$ is
the phase of a complex scalar $\Phi$, it is a pseudo-scalar, and hence the
model (1) and the corresponding classical solution (2) describe pseudo-scalar
(and for small $m$) light domain walls. The plot on the left in Fig. 1 shows
the domain wall solution (1) as a function of the normal distance $z$ to the
domain wall in km.

Figure 1: The profile of the domain wall field (2) (left) and the function
$\sin^{2}\left(\frac{N_{A}\phi}{f}\right)$ (right) plotted as functions of the
distance $z$ for $m=1$ neV, $f=1$ GeV and $N_{A}=20,\,N_{\phi}=4$.

We will consider three types of interactions between the domain wall field and
the photon field in QED that encompass couplings to all photon operators of
dimension $\leq 4$. These correspond to an effective photon mass-term
coupling,
$\mathcal{L}\supset\frac{1}{2}m_{0}^{2}\,\sin^{2}\left(\frac{N_{A}\phi}{f}\right)A_{\mu}A^{\mu}$,
an axion-like coupling,
$\mathcal{L}\supset\frac{1}{4}\tilde{g}_{DW}\,(\phi/f)\,F_{\mu\nu}\tilde{F}^{\mu\nu}$
and a coupling to the canonical photon kinetic term,
$\mathcal{L}\supset\frac{1}{4}g_{DW}\,(\phi/f)^{2}\,F_{\mu\nu}F^{\mu\nu}$. The
photon mass-term coupling is considered in Section 2, the domain wall axion-
like coupling is dealt with in Section 3, and in Section 4 we compute the
effects of the photon kinetic term coupling to domain walls. Conclusions are
presented in Section 5.

##  2 Signal from a photon mass term coupling

In this section, we review the results of [1] and extend this work by
accounting for the Fabry-Perot cavities used in LIGO. The interaction of the
domain wall field $\phi$ with the lowest dimension ($d=2$) photon field
operator is given by the effective mass term for the photon sourced by the
domain wall field. We thus consider the modified QED Lagrangian,

$\mathcal{L}=\mathcal{L}_{QED}+\frac{1}{2}m_{0}^{2}\sin^{2}\left(\frac{N_{A}\phi}{f}\right)A^{\mu}A_{\mu}\,.$
(3)

The mass parameter $m_{0}$ will ultimately be taken to be in the sub-neV
regime and the domain wall-induced function
$\sin^{2}\left(\frac{N_{A}\phi}{f}\right)$ is plotted in Fig. 1 (plot on the
right). Note that computed on the domain wall solution, this function vanishes
far from the centre of the domain wall as long as $N_{A}/N_{\phi}$ is integer
or half-integer, so that the photon mass is non-vanishing only in the vicinity
of the domain wall. This changes the dispersion relation of light which
becomes,

$\omega^{2}=k^{2}+m_{0}^{2}\sin^{2}\left(\frac{N_{A}\phi}{f}\right),$ (4)

where $k=|{\bf k}|$ is the magnitude of the photon 3-momentum.

Figure 2: A simplified representation of a Michelson interferometer used in
gravitational wave detectors such as LIGO and VIRGO. Note that ideal mirrors
produce an overall phase shift of $\pi$ which will cancel between the two
arms. We also show here the detector geometry and the definition of the angle
$\alpha$. We also point out here that the domain wall has a finite width $\sim
2/m\sim$ few km for the masses considered in this paper. This means that the
domain wall may cause additional effects in the detector apparatus (this is
addressed in Section 4).

To briefly explain the detector geometry, a typical Michelson interferometer,
as shown in Fig. 2, consists of two arms perpendicular to each other, an
effect then induces a phase difference between light in the two arms causing a
variation in the power output at the point where the arms join due to the
superposition of the two beams. Note, that an ideal mirror produces an overall
phase shift of $\pi$ which will cancel between the two arms. In this section
we take the angle $\alpha$ to be the angle the domain wall makes with the
x-axis (and we assume the normal is perpendicular to the z-axis) and we align
our and x- and y-axes with the detector arms. The domain wall profile is given
by

$\phi\,=\,\phi_{\rm cl}({\bf x}\cdot{\bf n}-z_{0}-vt)\,,$ (5)

where ${\bf n}=(\sin\alpha,\cos\alpha)$ is the unit normal vector to the wall,
$z_{0}$ is the distance from the wall to the origin at $t=0$, $v$ is the
velocity of the domain wall and $\phi_{\rm cl}$ is the domain wall solution
profile (2).

We can now use the dispersion relation (4) to find the phase velocity
$v_{P}={\omega}/{k}$ of the photon. For a light domain wall, $m\ll k$, we
have,

$v_{P}=\frac{\omega}{k}=1+\frac{m_{0}^{2}}{2k^{2}}\sin^{2}\left(\frac{N_{A}\phi}{f}\right).$
(6)

Recalling that $v_{P}=\frac{dx}{dt}$, we can then integrate both sides of this
equation to get

$2L=t-t_{0}+\frac{m_{0}^{2}}{2\omega^{2}}\int_{t_{0}}^{t}\,dt\,\sin^{2}\left(\frac{N_{A}\phi}{f}\right),$
(7)

where $t_{0}$ is the initial time when the light emitted by the laser arrives
at the beam splitter at the beginning of the detector arm of length $L$, and
$t$ is the time when the light returns to the beam splitter having travelled
along the detector arm to the mirror and back (thus covering the overall
distance $2L$). Note that we have replaced $k$ by $\omega$ in the second term
as to zeroth order in $m_{0}$, $\omega=k$. Now also note that to zeroth order
in $m_{0}$, $t_{0}=t-2L$ so if we only wish to account for
$\mathcal{O}\left(m_{0}^{2}\right)$ effects then in the limits of the integral
we can replace $t_{0}$ by $t-2L$ since this term is already proportional to
$m_{0}^{2}$. Hence we have,

$\displaystyle
2L=t-t_{0}+\frac{m_{0}^{2}}{2\omega^{2}}\int_{t-2L}^{t-L}\,dt^{\prime}\,\sin^{2}\left(\frac{N_{A}\,\phi_{\rm
cl}((t^{\prime}-t+2L)\sin\alpha-vt^{\prime})}{f}\right)+$ (8)
$\displaystyle\frac{m_{0}^{2}}{2\omega^{2}}\int_{t-L}^{t}\,dt^{\prime}\,\sin^{2}\left(\frac{N_{A}\,\phi_{\rm
cl}((t-t^{\prime})\sin\alpha-vt^{\prime})}{f}\right),$ (9)

where we have accounted for the fact that the photon is travelling along the
$X$-arm of the detector, hence in the domain wall profile (5) we take $x$ to
be the distance of a photon from the beam splitter at time $t$, we take $y=0$,
and we placed the wall at the origin at $t=0$, hence $z_{0}=0$. Also note the
necessity of splitting the domain of the integral into the region where the
light is travelling away from/towards the detector.

For our purposes it is more useful to say that if a beam of light is detected
at time $t$ then it must have been emitted at a time:

$\displaystyle
t_{0}=t-2L+\frac{m_{0}^{2}}{2\omega^{2}}\int_{t-2L}^{t-L}\,dt^{\prime}\,\sin^{2}\left(\frac{N_{A}\,\phi_{\rm
cl}((t^{\prime}-t+2L)\sin\alpha-vt^{\prime})}{f}\right)+$ (10)
$\displaystyle\frac{m_{0}^{2}}{2\omega^{2}}\int_{t-L}^{t}\,dt^{\prime}\,\sin^{2}\left(\frac{N_{A}\,\phi_{\rm
cl}((t-t^{\prime})\sin\alpha-vt^{\prime})}{f}\right).$ (11)

Importantly, gravitational wave interferometers make use of Fabry-Perot
cavities and as a result light makes more than one round trip down the
detector arms. The mean number of round trips in the cavity is related to the
finesse, where, for LIGO, we have $\mathcal{F}=450$ [15]. We may also work
with the mean amount of time light spends in the cavity, known as the storage
time, $\tau_{s}=\frac{L\mathcal{F}}{\pi c}$. Assuming the plane-wave
representation for the photon electric field, $\mathbf{E}\propto e^{-i\omega
t}$, we see that the phase difference induced by the domain wall for the light
travelling along the $X$-arm of the detector is,111Note that more accurately
the observed light is a superposition of beams of light which have made
different numbers of round trips but we neglect this effect here.

$\Delta\varphi_{x}(t)\,=\,\frac{m_{0}^{2}}{2\omega}\int_{t-\tau_{s}}^{t}\,dt^{\prime}\,\sin^{2}\left(\frac{N_{A}\,\phi_{\rm
cl}(x(t^{\prime})\sin\alpha-vt^{\prime})}{f}\right),$ (12)

and for the $Y$-arm it is,

$\Delta\varphi_{y}(t)\,=\,\frac{m_{0}^{2}}{2\omega}\int_{t-\tau_{s}}^{t}\,dt^{\prime}\,\sin^{2}\left(\frac{N_{A}\,\phi_{\rm
cl}(y(t^{\prime})\cos\alpha-vt^{\prime})}{f}\right),$ (13)

with the domain wall solution given by (2). Note that by $x(t)$ and $y(t)$ we
mean the position of the light in the respective cavity, which will involve an
appropriate splitting of the domain of the integral into regions where the
light is moving away from/towards the detector as in Eq. 8. The observable at
LIGO is then the photon phase difference between the X- and Y-arms of the
detector,

$\Delta\varphi=\Delta\varphi_{x}-\Delta\varphi_{y}\,.$ (14)

The difference between the two arms arises due to different distances to the
domain wall.

Finally, before presenting results we briefly discuss the value of various
parameters. We take $v\sim 10^{-3}$, as a generic galactic velocity, in all
scenarios considered in this paper $f$ cancels and therefore the signal is
independent of $f$. However it should be noted that that $f$ should be less
than $\sim 1$ TeV to avoid some astrophysical constraints and, similarly, we
require $m_{0}\lesssim$ few neV for astrophysical reasons as discussed in [1].
We plot the results in Fig. 3.

Figure 3: The photon phase difference, $\Delta\varphi$, between the X- and
Y-arms of gravitational waves detectors induced by the domain wall-photon mass
term interaction (1). The phase difference is plotted as a function of time
with the domain wall assumed to be at the origin at $t=0$. The chosen
parameters are $m_{0}=0.1\,\textrm{neV},\,N_{A}=20,\,N_{\phi}=4,\,v=3\times
10^{-3},\,\omega=1\,\textrm{eV},\,L=4\,\textrm{km},\,\alpha=\pi/2.2,\,\mathcal{F}=450$
and $m=0.1$ neV. Note that the result is independent of the value of $f$.

We also wish to examine in which areas of parameter space LIGO is sensitive to
a domain wall signal. In [1], it was shown that LIGO is sensitive to phase
shifts of order $\sim 10^{-10}$. Again, as in [1], we can make the
transformation $\hat{x}=mx$ and $\hat{y}=my$ in Eqs. (12)-(13) so that our
integration variables are dimensionless. We then see that the phase shift is
proportional to $\frac{m_{0}^{2}}{m\omega}$. So in order for a particular
domain wall model to produce an observable signal at LIGO we would require
$\frac{m_{0}^{2}}{m\omega}\gtrsim 10^{-10}$.

##  3 Signal from an axion-like coupling

In this section we consider an axion-like interaction of the domain wall with
the photon given by the $\phi$ coupling to the dimension-4 photon operator
via,

$\mathcal{L}_{\textrm{int}}=\frac{\tilde{g}_{DW}}{4f}\,\phi\,F_{\mu\nu}\tilde{F}^{\mu\nu}\,,$
(15)

where
$\tilde{F}^{\mu\nu}=\frac{1}{2}\epsilon^{\mu\nu\rho\sigma}F_{\rho\sigma}$ is
the dual field-strength tensor. The Euler-Lagrange equations for the QED
photon field modified by the interaction Lagrangian (15) are given by

$\partial_{\mu}\left(F^{\mu\nu}+\frac{\tilde{g}_{DW}\phi}{f}\tilde{F}^{\mu\nu}\right)=0.$
(16)

This is equivalent to modified Maxwell’s equations of the form,

$\displaystyle\nabla\cdot\mathbb{B}=0\,,$ (17)
$\displaystyle\nabla\times\mathbb{E}+\dot{\mathbb{B}}=0\,,$ (18)
$\displaystyle\nabla\cdot\mathbb{E}=\frac{\tilde{g}_{DW}}{f}\nabla\phi\cdot\mathbb{B}\,,$
(19)
$\displaystyle\nabla\times\mathbb{B}-\dot{\mathbb{E}}=\frac{\tilde{g}_{DW}}{f}\left(\mathbb{E}\times\nabla\phi-\mathbb{B}\dot{\phi}\right).$
(20)

It was shown in [16] that this leads to the rotation of the plane of
polarisation of linearly polarised light. However for our purposes it is more
useful to instead find the different dispersion relations for left and right
circular polarisations which cause this rotation as was done in a similar
search for axions [17]. We will consider an experiment of the form proposed in
[17] in this section where one arm of the detector consists of right
circularly polarised light and the other arm left circularly-polarised light.

With these new Maxwell’s equations the wave equations for the electric and
magnetic fields become [18],

$\displaystyle\Box
E_{i}+\frac{\tilde{g}_{DW}}{f}\left(\partial_{t}\phi\,\partial_{t}B_{i}-\partial_{j}\phi\,\partial_{j}B_{i}\right)=0\,,$
(21) $\displaystyle\Box
B_{i}+\frac{\tilde{g}_{DW}}{f}\left(\partial_{t}\phi\,\partial_{t}E_{i}-\partial_{j}\phi\,\partial_{j}E_{i}\right)=0\,.$
(22)

Note that the derivation of these equations requires the assumption that
second derivatives of the domain wall field are negligible compared to first
derivatives, and similarly first derivatives squared are negligible. We shall
examine this assumption at the end of this section. Substituting the usual
plane-wave ansatz for the electric and magnetic fields $\propto e^{-i\omega
t+i{\bf k}.{\bf x}}$ into these equations leads to a set of equations of the
form,

$M(\mathbb{k},\omega)(\mathbb{E},\mathbb{B})^{T}=0$ (23)

for a particular matrix $M$. This equation has non-trivial solutions only when
det $M$=0. Imposing this condition yields the relation [18]:

$\omega^{2}-k^{2}=\pm\frac{\tilde{g}_{DW}}{f}(\omega\partial_{t}\phi+\mathbf{k}.\nabla\phi).$
(24)

These are our two solutions for left and right circular polarised light,
giving the phase velocity (c.f. (6)),

$v_{P}\,=\,1\pm\frac{\tilde{g}_{DW}}{2fk^{2}}(\omega\,\partial_{t}\phi+\mathbf{k}.\nabla\phi)\,.$
(25)

Now we have $v_{P}=\frac{dx}{dt}$ so for a single round trip integration
around the detector arm we have,

$\int_{0}^{2L}dx\,=\,\int_{t_{0}}^{t}dt\,\left(1\pm\frac{\tilde{g}_{DW}}{2fk^{2}}(\omega\,\partial_{t}\phi+\mathbf{k}.\nabla\phi)\right).$
(26)

Noting that the second term of the right hand side is a total derivative
($\frac{d\phi}{dt}=\partial_{t}\phi+\nabla\phi.\hat{\bf{n}}$ where
$\hat{\bf{n}}$ is a unit vector in the direction of propogation,
$\mathbb{k}/k$ in this case) and performing the integration we obtain,

$2L\,=\,t-t_{0}\pm\frac{\tilde{g}_{DW}}{2fk}\left(\phi_{\rm cl}(vt)-\phi_{\rm
cl}(vt_{0})\right),$ (27)

where $L$ is the length of a detector arm, $t$ is the time light is detected
and $t_{0}$ is the time of emission. Using a similar argument as that leading
to Eq. (12) we obtain:

$\Delta\varphi_{x/y}\,=\,\pm\frac{\tilde{g}_{DW}}{2f}\left(\phi_{\rm
cl}(vt)-\phi_{\rm cl}(vt-v\tau_{s})\right),$ (28)

with the $\pm$ applying to the different arms of the detector and the
observable $\Delta\varphi$ being the difference between the two. Note that
here it is the different polarisations of light that lead to the phase
difference and not the different distances to the domain wall since the
observable depends only on the field strength at the observation/emission
points. It is worth bearing in mind here that since different polarisations of
light do not interfere, we must use a waveplate to change our circularly
polarised light to linearly polarised light before interfering the beams.
There are also some other waveplates required in this new setup to account for
the effects of reflection upon polarisation, for details see [17].

Let us also consider any constraints on $\tilde{g}_{DW}$. One might expect
that many of the constraints on a similar axion coupling would also be
applicable here, however the most stringent axion constraints come from the
CAST experiment [19] which examines the solar axion flux. Although we would
not expect domain walls to be produced in the sun, we would expect excitations
of the $\phi$ field to be produced. These excitations would then have the
possibility to convert into photons, the effect which is searched for by the
CAST experiment. Therefore in presenting results we will take into account
CAST constraints of $\tilde{g}_{DW}/f\lesssim 10^{-10}\,\rm{GeV}^{-1}$.
Results calculated using Eq. (28) are presented in Fig. 4. It is worth
pointing out here that these results are independent of the direction of
domain wall propagation since Eq. (28) only contains the value of the domain
wall field at the detector. In Fig. 5 we show a plot examining the values of
$m$ and $\tilde{g}_{DW}$ to which we are sensitive.

Figure 4: The value of the photon phase difference $\Delta\varphi$ to be
observed in the interferometer for the case of the axion-like interaction (15)
between the domain wall and photons using an experimental setup similar to
that suggested in [17]. $\Delta\varphi$ is plotted as a function of time for
$v=3\times
10^{-3},\,N_{\phi}=4,\,L=0.04\,\textrm{km},\,\tilde{g}_{DW}=10^{-10},\,\mathcal{F}=450$
and $m=0.1$ neV. Note that the result is independent of the value of $f$.
Figure 5: A plot showing the sensitivity of the proposed interferometer
($L=40$ m) to domain walls for various values of the domain wall velocity,
$v$. Note that the plot levels off at high masses; this is when the domain
wall becomes thin enough such that a photon experiences the full passage of
the domain wall in its storage time and so going to higher masses does not
increase the phase shift past this point. Note that a interferometer with
longer arms would improve the sensitivity to low mass domain walls. The
precise value of the CAST constraint depends on the value of $f$ (whereas our
phase shift doesn’t) but one should bear in mind that we also require
$\frac{\tilde{g}_{DW}}{f}\lesssim 10^{-10}$ $\textrm{GeV}^{-1}$.

Finally let us examine our requirement that higher order derivatives (and
higher powers of first derivatives) are negligible. This is equivalent to the
requirement that [18]

$\bigg{\lvert}\frac{\partial_{\mu}\phi}{\phi}\bigg{\rvert}\ll\omega.$ (29)

Taking the domain wall field in Eq. (2) and assuming propagation parallel to
the z-axis (a change in direction will just introduce factors of
cos/sin$\,\alpha$, where $\alpha$ is the angle between the direction of
propagation and the z-axis) we obtain

$\displaystyle\frac{\partial_{z}\phi_{\rm cl}}{\phi_{\rm
cl}}=\frac{m\,\textrm{sech}(m(z-vt))}{\arctan\left(e^{m(z-vt)}\right)}$ (30)
$\displaystyle\frac{\partial_{t}\phi_{\rm cl}}{\phi_{\rm
cl}}=\frac{mv\,\textrm{sech}(m(z-vt))}{\arctan\left(e^{m(z-vt)}\right)}.$ (31)

In Fig 6, we plot $\textrm{sech}(x)/\textrm{arctan}(e^{x})$ to show that it is
always less than two. Hence the requirement that
$\big{\lvert}\frac{\partial_{\mu}\phi}{\phi}\rvert\ll\omega$ is equivalent to
requiring that $2m\ll\omega$. This is easily satisfied for the range of $m$
and $\omega$ considered in this paper. Note that the derivative with respect
to $t$ simply introduces an extra factor of $v$ and, since
$v=\mathcal{O}(10^{-3})$, this does not affect the argument.

Figure 6: A plot of $\textrm{sech}(x)/\arctan(e^{x})$.

##  4 Signal from a canonical kinetic term coupling

Now we consider the signal from an interaction of the domain wall with the
photon kinetic term,

$\mathcal{L}_{\textrm{int}}=\frac{g_{DW}}{4f^{2}}\phi^{2}F_{\mu\nu}F^{\mu\nu}\,.$
(32)

Note that since the domain wall field is a pseudoscalar, the coupling to the
canonical kinetic term can involve only even powers of $\phi/f$. Here for
simplicity we concentrate on the $\phi^{2}$ interactions, but our analysis
below can be easily applied to the general case of any even function of
$\phi/f$. We also note that the contribution from $\phi^{2}$ is not
necessarily negligible in the EFT expansion since we expect
$\phi/f=\mathcal{O}(1)$ (c.f. Eq. (2)). The Euler-Lagrange equations for the
photon now read:

$\partial_{\mu}\left(\left(-1+\frac{g_{DW}}{f^{2}}\phi^{2}\right){F}^{\mu\nu}\right)=0.$
(33)

Such interactions were considered for a fundamental particle in a strong
magnetic field in e.g.[20, 21] but here we attempt to derive more general
results which may be applied to interferometers where there is no magnetic
field present. Instead we start, as in the previous section, with the modified
Maxwell’s equations resulting from the interaction in (32),

$\nabla.\bf{B}=0\,,$ (34) $\nabla\times\bf{E}+\dot{\bf{B}}=0\,,$ (35)
$\nabla.\left(\left(1-\frac{g_{DW}\phi^{2}}{f^{2}}\right)\bf{E}\right)=0\,,$
(36)
$\nabla\times\left(\left(1-\frac{g_{DW}\phi^{2}}{f^{2}}\right)\bf{B}\right)=\frac{\partial}{\partial
t}\left(\left(1-\frac{g_{DW}\phi^{2}}{f^{2}}\right)\bf{E}\right)\,.$ (37)

As in the previous section we neglect second derivatives and the square of
first derivatives. We can then obtain the wave equations for the electric and
magnetic fields (c.f. our earlier equations (21)-(22) for axion-like
interactions):

$\displaystyle\Box
E_{i}-\frac{g_{DW}}{f^{2}}\left(\partial_{t}\phi^{2}\partial_{t}E_{i}-\partial_{j}\phi^{2}\partial_{j}E_{i}\right)=0\,,$
(38) $\displaystyle\Box
B_{i}-\frac{g_{DW}}{f^{2}}\left(\partial_{t}\phi^{2}\partial_{t}B_{i}-\partial_{j}\phi^{2}\partial_{j}B_{i}\right)=0\,,$
(39)

where $j$ runs over spatial indices and there is an implied summation over
$j$.

Substituting in a plane-wave ansatz for the $E$ and $B$ fields,

$\mathbf{E}=\mathbf{E_{0}}\,e^{-i\omega t+i{\bf
k}.{\mathbf{r}}}\,,\quad\mathbf{B}=\mathbf{B_{0}}\,e^{-i\omega t+i{\bf
k}.{\mathbf{r}}}\,,$ (40)

gives the dispersion relation,

$\omega^{2}-k^{2}=i\frac{g_{DW}}{f^{2}}\left(\omega\partial_{t}\phi^{2}+\mathbf{k}.\nabla\phi^{2}\right).$
(41)

Solving for $\omega$ yields:

$\omega=\lvert
k\rvert+i\frac{g_{DW}}{2kf^{2}}\left(\omega\partial_{0}\phi^{2}+\mathbf{k}.\nabla\phi^{2}\right).$
(42)

Substituting this into the usual plane-wave ansatz we obtain

$\mathbf{E}\,=\,\mathbf{E_{0}}\,e^{\,\frac{g_{DW}}{2kf^{2}}\left(\omega\partial_{0}\phi^{2}+\mathbf{k}.\nabla\phi^{2}\right)\,t}e^{-i\lvert
k\rvert t+i\mathbb{k}\cdot\mathbb{r}}.$ (43)

It can clearly be seen that this is an ordinary plane-wave with an
exponentially growing amplitude. To aid in our treatment of cavities we employ
the same approach as in the previous sections and calculate the phase velocity
as:

$v_{P}=1+i\frac{g_{DW}}{2k^{2}f^{2}}\left(\omega\partial_{0}\phi^{2}+\mathbf{k}.\nabla\phi^{2}\right).$
(44)

Integrating as in the previous section gives

$2L=t-t_{0}+i\frac{g_{DW}}{2kf^{2}}\,\left(\phi_{\rm cl}^{2}(vt)-\phi_{\rm
cl}^{2}(vt_{0})\right),$ (45)

or

$t_{0}=t-2L+i\frac{g_{DW}}{2kf^{2}}\,\left(\phi_{\rm cl}^{2}(vt)-\phi_{\rm
cl}^{2}(vt_{0})\right).$ (46)

Note that here we cannot speak of a phase difference anymore but substituting
this solution into a plane-wave ansatz ($\mathbf{E}\propto e^{-i\omega t}$)
gives:

$\mathbf{E}\,=\,\mathbf{E_{0}}\,e^{\,\frac{g_{DW}}{2f^{2}}\phi_{\rm
cl}^{2}\rvert_{t-\tau_{s}}^{t}}$ (47)

or equivalently, the power in the detector,

$P\,=\,|E_{x}^{2}+E_{y}^{2}|\,\propto\,e^{\,\frac{g_{DW}}{f^{2}}\phi_{\rm
cl}^{2}\rvert_{t-\tau_{s}}^{t}}.$ (48)

In the presence of a domain wall, the interaction term in Eq. (32) gives an
exponentially growing power (although as we will discuss shortly, $g_{DW}$
must be small so the exponential enhancement will not be all that large). In
order to examine this, in Fig. 7 we plot $P/P_{0}-1$ where $P_{0}$ is the
observed power in the absence of domain walls (we subtract one in order to
better see the small deviations in the power). Also note that this extra
energy is supplied/absorbed (depending on the sign of $\Delta\phi^{2}$) by the
kinetic energy of the domain wall.

In an interferometer such as those considered in the previous sections, the
light in the two different arms will be affected equally (and the two
different polarisations no longer obey different dispersion relations). This
will produce a fluctuation in the power output at the detector but it will in
general not be possible to distinguish this signal from one such as in Section
2. However a coupling of this form will uniquely produce a power change in a
single beam of light since the change in power output is an intrinsic property
of the domain wall and not a consequence of a phase change. This would be one
possible method of determining the form of the coupling.

Before presenting results let us first examine some constraints on the value
of $g_{DW}$. Adding an interaction term of the form Eq. (32) leads to a change
in the QED coupling constant:

$\alpha\rightarrow\frac{\alpha}{1-g_{DW}\,\frac{\phi^{2}}{f^{2}}}.$ (49)

Since we can expect $\phi/f$ to be $\mathcal{O}(1)$, this can result in very
large changes to the value of the fine structure constant. Furthermore, since
a domain wall is an "event" rather than a continuously varying field, this
would lead to the value of the fine structure constant on Earth changing value
every time a domain wall passed. Clearly this would be catastrophic; to avoid
this issue we choose $g_{DW}$ so that the change in $\alpha$ falls within the
current uncertainties of its measurement. Since $\alpha$ is currently known to
a precision of approximately one part in $10^{10}$ we take $g_{DW}=10^{-10}$.
As in the previous section these results are independent of the direction of
domain wall travel.

Another way to avoid the change of the fine structure constant would be to
replace $\phi^{2}/f^{2}$ by some function which vanishes far from the domain
wall. Taking inspiration from Section 2 we choose:

$\mathcal{L}_{int}=\frac{g_{DW}}{4}\sin^{2}\left(\frac{N_{A}\phi}{f}\right)F_{\mu\nu}F^{\mu\nu}$
(50)

This would result in $\alpha$ changing only for a very short time as the
domain wall passed. It is possible that the current methods used to measure
$\alpha$ would be sensitive to these variations but we will interpret any null
result as a constraint on the domain wall event rate and not on $g_{DW}$.
Since the derivation of Eq. (48) made no assumptions about the form of $\phi$,
to obtain a new formula, one simply replaces all occurences of
$\phi^{2}/f^{2}$ by $\sin^{2}\left(\frac{N_{A}\phi}{f}\right)$. This gives

$P\propto e^{\,g_{DW}\sin^{2}\left(\frac{N_{A}\phi_{\rm
cl}}{f}\right)\big{\rvert}_{t-\tau_{s}}^{t}}.$ (51)

However we must also examine the assumption that second derivatives are
negligible. Referring to Eq. 29 we arrive at the requirement that

$4m\frac{N_{A}}{N_{\phi}}\frac{\operatorname{sech}\left(mx\right)}{\tan\left(4\frac{N_{A}}{N_{\phi}}\tan^{-1}\left(e^{mx}\right)\right)}\ll\omega$
(52)

To avoid a singularity, we must have
$4\frac{N_{A}}{N_{\phi}}\tan^{-1}\left(e^{mx}\right)<\pi$. As
$0<\tan^{-1}\left(e^{mx}\right)<\pi/2$ this is equivalent to requiring
$\frac{N_{A}}{N_{\phi}}\leq\frac{1}{2}$. Given also our requirement that
$\sin^{2}\left(\frac{N_{A}\phi_{\rm cl}}{f}\right)\rightarrow 0$ as
$x\rightarrow\infty$, the only value of $\frac{N_{A}}{N_{\phi}}$ which
satisifes both requirements is $\frac{N_{A}}{N_{\phi}}=\frac{1}{2}$.
Subsituting this value into Eq. 52 and using some trigonometric identities we
obtain the requirement that

$8m\tanh\left(mx\right)\ll\omega$ (53)

which is satisified for the range of $m$ and $\omega$ considered in this
paper. Results are presented in Fig. 7.

Figure 7: The value of $\frac{P}{P_{0}}-1$ plotted for $v=3\times
10^{-3},\,\,N_{\phi}=4,\,L=4\,\textrm{km},\,\mathcal{F}=450,\,g_{DW}=10^{-10}$
and $m=0.1$ neV in the case where the domain wall is coupled to the photon
kinetic term via (left) (32), (right) (50). On the right we additionally have
$N_{A}=2$. Here $P$ represents the observed power and $P_{0}$ is the power
which would be observed in the absence of a domain wall (for a single beam).
Note that the result is independent of the value of $f$.

It should be noted that a change in the fine structure constant would lead to
a change in the Bohr radius and thus the length of various objects as shown in
[22, 23], which would produce a signal independently of the effect calculated
above. Here we will briefly derive the effect on the power at the detector
when one accounts for both effects, i.e. simultaneously the power shift coming
from the intrinsic interaction with the domain wall and the phase shift caused
by the change in the Bohr radius.

The change in length is given by [22]:

$\delta L\approx-L\frac{\delta\alpha}{\alpha}.$ (54)

For the range of velocities and masses considered in this paper we can take
the change in length to be adiabatic. The length of the interferometer cavity
is insulated from changes in length by the pendulum suspension system used in
LIGO [22] and we need only worry about the change in size of the beam splitter
and of the mirrors which will cause a difference in optical path length. The
change in optical path length is given by [22]:

$\displaystyle\delta\left(L_{x}-L_{y}\right)\approx-\sqrt{2}nl\frac{\delta\alpha}{\alpha}+\Delta
w_{\textrm{x-mirror}}-\Delta w_{\textrm{y-mirror}}\approx$ (55)
$\displaystyle\sqrt{2}nlg_{DW}\frac{\phi^{2}}{f^{2}}\bigg{\rvert}_{\rm{beam-
splitter}}+wg_{DW}\left(\frac{\phi^{2}}{f^{2}}\bigg{\rvert}_{\textrm{y-mirror}}-\frac{\phi^{2}}{f^{2}}\bigg{\rvert}_{\textrm{x-mirror}}\right)$
(56)

where $n$ is the refractive index of the beam-splitter (changes in $n$ induced
by the domain wall are negligible), $l$ is the thickness of the beam-splitter
and $w$ is the width of the mirrors. This must be summed over the number of
round trips each photon takes, for each trip evaluating $\phi$ at the position
of the beam splitter/mirrors respectively. This effect turns out to be
completely dominant over the effect of the change in power. In Fig. 8 we plot
the phase shift caused by this effect.

Figure 8: The value of $\Delta\varphi$ plotted for
$\lambda=1064\,\textrm{nm},\,\alpha=\pi/2.2,\,v=3\times
10^{-3},\,n=1.45,\,\omega=1\,\textrm{eV},\,N_{\phi}=4,\,L=4\,\textrm{km},\,\mathcal{F}=450,\,g_{DW}=10^{-10},l=6\,\textrm{cm},\,w=20\,\textrm{cm},\,m=0.1$
neV in the case where the domain wall is coupled to the photon kinetic term
via (left) (32), (right) (50) where we also set $N_{A}=2$. Here we plot the
phase difference caused by the change in size of the beam splitter and mirrors
(the power fluctuation caused by the intrinsic interaction is negligible).
Note that we obtain a permanent phase shift in the left plot due to the
permanent change in the value of the fine structure constant after the passage
of the domain wall.

Finally for completeness we here derive the formula for the power observed at
the detector in the presence of both effects, although as already noted, the
effect of the term in Eq. 48 is negligible compared to the phase shift induced
by the change in Bohr radius. The power observed by the detector in the
presence of a phase difference between the two beams is

$P=\lvert A+Ae^{i\Delta\varphi}\rvert^{2}=2\lvert
A\rvert^{2}\left(1+\cos\left(\Delta\varphi\right)\right)$ (57)

where $A$ is the amplitude of each beam. If we now include the effect of the
intrinsic photon-domain wall interaction then the amplitude of each beam is
also enhanced by a factor $e^{\,\frac{g_{DW}}{2f^{2}}\phi_{\rm
cl}^{2}\rvert_{t-\tau_{s}}^{t}}$ so

$\displaystyle P=\lvert Ae^{\,\frac{g_{DW}}{2f^{2}}\phi_{\rm
cl}^{2}\rvert_{t-\tau_{s}}^{t}}+Ae^{i(\Delta\varphi_{I}+\Delta\varphi_{DW})+\frac{g_{DW}}{2f^{2}}\phi_{\rm
cl}^{2}\rvert_{t-\tau_{s}}^{t}}\rvert^{2}=$ (58) $\displaystyle 2\lvert
A\rvert^{2}e^{\,\frac{g_{DW}}{f^{2}}\phi_{\rm
cl}^{2}\rvert_{t-\tau_{s}}^{t}}\left(1+\cos\left(\Delta\varphi_{I}+\Delta\varphi_{DW}\right)\right)$
(59)

where $\Delta\varphi_{DW}$ is the phase shift caused by the change in length
of the beam splitter and mirrors and $\Delta\varphi_{I}$ is a phase shift in
the absence of domain walls which is often included in the experimental design
as, for various technical reasons, it makes detection easier (for details see
[14]).

The phase shift is given approximately by

$\Delta\varphi\approx\sqrt{2}nlg_{DW}\left(\frac{\mathcal{F}}{2\pi}\right)\left(\frac{2\pi}{\lambda}\right)$
(60)

as $\frac{\mathcal{F}}{2\pi}$ gives the mean number of round trips that light
makes in the interferometer. Note that here we have neglected the effect
coming from the change in width of the mirrors but as the effects are of
comparable magnitude this should not have a large effect on our estimate and
we have also taken
$\Delta\left(\frac{\phi^{2}}{f^{2}}\right)=\Delta\sin^{2}\left(\frac{N_{A}\phi}{f}\right)=\mathcal{O}\left(1\right)$.
As usual for LIGO to be sensitive to such a coupling we would require
$\Delta\varphi\gtrsim 10^{-10}$.

##  5 Conclusion

We briefly discuss here the prospect of detection at future, planned
interferometers such as LISA. We again consider the three couplings
separately. For the case of the coupling considered in Sec. 2, the longer
interferometer arms proposed in LISA will, in general, increase the magnitude
of the signal as the integrals in Eqs. 12,13 will grow larger, although as
with gravitational waves there is a limit to how much sensitivity can be
gained purely by increasing the size of the interferometer. Similar
conclusions apply to the coupling considered in Sec. 4.

For the case of the coupling considered in Sec. 3, it should be noted that
proposed experiments such as LISA would not be sensitive to this form of
coupling since they do not use polarised light. However more generally
considering the case of a larger interferometer, we point out at the end of
that section that there is a maximum phase shift which can be induced by the
domain wall. Therefore increasing the arm length (or equivalently increasing
the storage time), is really probing lower values of $m$.

The other main type of gravitational wave detector is a resonant bar detector.
Since the domain wall model in this paper couples only to the photon it would
not induce excitations of a bar as required by a resonant mass detector.
However the changes in the fine structure constant could produce a signal in
the case of the coupling considered in Sec. 4. The change in length of the bar
is, in principle, large enough to be detected, however the detection at a
resonant bar is complex and depends on the frequency of the signal. Hence
whether there is a realistic chance of detection is a more involved question
requiring more detailed calculation which we leave for future work. The
detection of domain walls using magnetometers was also considered in [13],
although in this paper the domain wall was coupled to fermions rather than
photons.

In conclusion, we have examined the effect of domain walls on the dispersion
relation of light and examined the effect of this on a signal in a Michelson
interferometer. In the case of an axion-like coupling, a standard Michelson
interferometer using unpolarised or linearly polarised light is insensitive to
this effect and we proposed to use instead detectors designed to search for
axions via a similar mechanism. The signals calculated here are within reach
of current gravitational wave detectors and we have examined the parameter
space for regions which current detectors woud be sensitive to.

## Acknowledgements

This work is supported by the STFC under consolidated grant ST/P001246/1 (VVK)
and an STFC studentship (DLM).

## References

  * [1] J. Jaeckel, V. V. Khoze, and M. Spannowsky, “Hearing the signal of dark sectors with gravitational wave detectors,” Phys. Rev. D 94 no. 10, (2016) 103519, [arXiv:1602.03901 [hep-ph]].
  * [2] P. Sikivie, “Of Axions, Domain Walls and the Early Universe,” Phys. Rev. Lett. 48 (1982) 1156–1159.
  * [3] LIGO Scientific, Virgo Collaboration, B. P. Abbott et al., “Observation of Gravitational Waves from a Binary Black Hole Merger,” Phys. Rev. Lett. 116 no. 6, (2016) 061102, [arXiv:1602.03837 [gr-qc]].
  * [4] G. Branco, P. Ferreira, L. Lavoura, M. Rebelo, M. Sher, and J. P. Silva, “Theory and phenomenology of two-Higgs-doublet models,” Phys. Rept. 516 (2012) 1–102, [arXiv:1106.0034 [hep-ph]].
  * [5] G. Belanger, B. Dumont, U. Ellwanger, J. Gunion, and S. Kraml, “Global fit to Higgs signal strengths and couplings and implications for extended Higgs sectors,” Phys. Rev. D 88 (2013) 075008, [arXiv:1306.2941 [hep-ph]].
  * [6] R. A. Battye, A. Pilaftsis, and D. G. Viatic, “Domain wall constraints on two-Higgs-doublet models with $Z_{2}$ symmetry,” Phys. Rev. D 102 no. 12, (2020) 123536, [arXiv:2010.09840 [hep-ph]].
  * [7] R. A. Battye, A. Pilaftsis, and D. G. Viatic, “Simulations of Domain Walls in Two Higgs Doublet Models,” JHEP 01 (2021) 105, [arXiv:2006.13273 [hep-ph]].
  * [8] G. Altarelli and F. Feruglio, “Discrete Flavor Symmetries and Models of Neutrino Mixing,” Rev. Mod. Phys. 82 (2010) 2701–2729, [arXiv:1002.0211 [hep-ph]].
  * [9] H. Baer, K.-Y. Choi, J. E. Kim, and L. Roszkowski, “Dark matter production in the early Universe: beyond the thermal WIMP paradigm,” Phys. Rept. 555 (2015) 1–60, [arXiv:1407.0017 [hep-ph]].
  * [10] K. Nakayama, F. Takahashi, and N. Yokozaki, “Gravitational waves from domain walls and their implications,” Phys. Lett. B 770 (2017) 500–506, [arXiv:1612.08327 [hep-ph]].
  * [11] J. Jaeckel, S. Schenk, and M. Spannowsky, “Probing Dark Matter Clumps, Strings and Domain Walls with Gravitational Wave Detectors,” [arXiv:2004.13724 [astro-ph.CO]].
  * [12] K. Saikawa, “A review of gravitational waves from cosmic domain walls,” Universe 3 no. 2, (2017) 40, [arXiv:1703.02576 [hep-ph]].
  * [13] M. Pospelov, S. Pustelny, M. Ledbetter, D. Jackson Kimball, W. Gawlik, and D. Budker, “Detecting Domain Walls of Axionlike Models Using Terrestrial Experiments,” Phys. Rev. Lett. 110 no. 2, (2013) 021803, [arXiv:1205.6260 [hep-ph]].
  * [14] M. Maggiore, Gravitational Waves. Vol. 1: Theory and Experiments.  Oxford Master Series in Physics. Oxford University Press, 2007.
  * [15] LIGO Scientific Collaboration, J. Aasi et al., “Advanced LIGO,” Class. Quant. Grav. 32 (2015) 074001, [arXiv:1411.4547 [gr-qc]].
  * [16] D. Harari and P. Sikivie, “Effects of a Nambu-Goldstone boson on the polarization of radio galaxies and the cosmic microwave background,” Phys. Lett. B 289 (1992) 67–72.
  * [17] W. DeRocco and A. Hook, “Axion interferometry,” Phys. Rev. D 98 no. 3, (2018) 035021, [arXiv:1802.07273 [hep-ph]].
  * [18] D. Blas, A. Caputo, M. M. Ivanov, and L. Sberna, “No chiral light bending by clumps of axion-like particles,” Phys. Dark Univ. 27 (2020) 100428, [arXiv:1910.06128 [hep-ph]].
  * [19] CAST Collaboration, K. Zioutas et al., “First results from the CERN Axion Solar Telescope (CAST),” Phys. Rev. Lett. 94 (2005) 121301, [arXiv:hep-ex/0411033].
  * [20] L. Maiani, R. Petronzio, and E. Zavattini, “Effects of Nearly Massless, Spin Zero Particles on Light Propagation in a Magnetic Field,” Phys. Lett. B 175 (1986) 359–363.
  * [21] Y. Liao, “Impact of spin-zero particle-photon interactions on light polarization in external magnetic fields,” Phys. Lett. B 650 (2007) 257–261, [arXiv:0704.1961 [hep-ph]].
  * [22] H. Grote and Y. V. Stadnik, “Novel signatures of dark matter in laser-interferometric gravitational-wave detectors,” Phys. Rev. Res. 1 no. 3, (2019) 033187, [arXiv:1906.06193 [astro-ph.IM]].
  * [23] Y. V. Stadnik and V. V. Flambaum, “Enhanced effects of variation of the fundamental constants in laser interferometers and application to dark matter detection,” Phys. Rev. A 93 no. 6, (2016) 063630, [arXiv:1511.00447 [physics.atom-ph]].