Published in: Far East Journ. of Appl. Math., 116, N3, (2023), pp. 215-227.

#  Boundary values of analytic functions

Alexander G. Ramm
Department of Mathematics, Kansas State University,
Manhattan, KS 66506, USA
<EMAIL_ADDRESS>
http://www.math.ksu.edu/ $\sim$ ramm

###### Abstract

MSC: 30E25, 30B30 Key words: boundary values of analytic functions

Let $D$ be a connected bounded domain in ${\mathbb{R}}^{2}$, $S$ be its
boundary which is closed, connected and smooth. Let $\Phi(z)=\frac{1}{2\pi
i}\int_{S}\frac{f(s)ds}{s-z}$, $f\in L^{1}(S)$, $z=x+iy$. Boundary values of
$\Phi(z)$ on $S$ are studied. The function $\Phi(t)$, $t\in S$, is defined in
a new way. Necessary and sufficient conditions are given for $f\in L^{1}(S)$
to be boundary value of an analytic in $D$ function. The Sokhotsky-Plemelj
formulas are derived for $f\in L^{1}(S)$.

##  1 Introduction

Let $D$ be a connected bounded domain on the complex plane, $S$ be its
boundary, which is closed and $C^{1,a}$-smooth, $0<a\leq 1$. Consider an
analytic function in $D$, defined as

$\Phi(z)=c\int_{S}\frac{f(s)}{s-z}ds,\quad c:=\frac{1}{2\pi i}.$ (1)

We assume that $f\in L^{1}(S)$. This is the basic new assumption: in the
literature it was assumed that $f$ is H$\ddot{o}$lder-continuous or $f\in
L^{p}(S),\,\,p>1$, see [2], [7]. In [1] there is a result for $f\in L^{1}(S)$,
which we mention in Remark 3.

It is of great interest to have a proof of the Sokhotsky-Plemelj formulas for
$f\in L^{1}(S)$ and of the the relation (9) for $f\in L^{1}(S)$, see below.
The contents of Chapter 2 in [2] is based on the relation (9), which is proved
in [2] under the assumption that $f$ is H$\ddot{o}$lder-continuous, $f\in
H^{\mu}(S)$, $\mu\in(0,1]$. We prove that this relation holds almost
everywhere (a.e.) in the sense of the Lebesgue measure on $S$. The Sokhotsky-
Plemelj formulas for $f\in L^{1}(S)$ are important in applications to singular
integral operators and boundary value problems, [2], [8].

By $D^{+}$ we denote $D$, by $D^{-}=D^{\prime}$ we denote
${\mathbb{R}}^{2}\setminus\bar{D}$, the $\bar{D}$ is the closure of $D$. The
function $\Phi$ is analytic in $D$ and in $D^{\prime}$, $\Phi(\infty)=0$.
Boundary value of $\Phi$ on $S$ has been studied by many authors, [1]–[9]. In
[2] and [8] it is assumed that $f\in H^{\mu}(S)$. In [7] and in [6] it is
assumed that $f\in L^{p}(S)$, $p>1$. The case $p=1$ is not discussed in these
books. In [1] there are results related to $p=1$ related to the Cauchy
principal value definition of $\Phi(t)$. We give a new definition of the
singular integral $\Phi(t)$, see (12).

Define $\psi(z)$ by formula (2). Denote the limiting values of $\Phi(z)$ and
$\psi(z)$, when $z\to t\in S$ along the non-tangential to $S$ directions, by
$\Phi^{+}(t)$ and $\psi^{+}(t)$ when $z\in D,\,z\to t$, and by $\Phi^{-}(t)$
and $\psi^{-}(t)$ when $z\in D^{\prime},\,z\to t$. By $\Phi(t)$ and $\psi(t)$
we denote the values of $\Phi(z)$ and $\psi(z)$ when $z=t\in S$. Let us assume
that $z\to t$ along the unit normal $N_{t}$ to $S$, where the normal is
directed out of $D$.

We are interested in the following question:

How to characterize functions $f\in L^{1}(S)$ which are boundary values of
analytic functions in $D$?

We answer this question in Theorem 1 (see below).

Rewrite (1) as

$\Phi(z)=f(t)\nu(z)+\psi(z),\quad\psi(z):=c\int_{S}\frac{f(s)-f(t)}{s-z}ds,\quad
t\in S,\quad\nu(z):=c\int_{S}\frac{ds}{s-z},$ (2)

and

$\nu(z)=\begin{cases}1,\quad z\in D,\\\ \frac{1}{2},\quad z\in S,\\\ 0,\quad
z\in D^{\prime}.\end{cases}$ (3)

One has

$\Phi^{+}(t)=\lim_{z\to t,z\in D}\Phi(z)=f(t)+\psi^{+}(t)$ (4)

where $\psi^{+}(t)=\lim_{z\to t,z\in D}\psi(z)$ and

$\Phi^{-}(t)=\lim_{z\to t,z\in D^{\prime}}\Phi(z)=\psi^{-}(t).$ (5)

If $t\in S$, then one gets (see equation (3), the line $z\in S$) :

$\Phi(t)=\frac{f(t)}{2}+\psi(t):=\frac{f(t)}{2}+c\int_{S}\frac{f(s)-f(t)}{s-t}ds.$
(6)

The $\psi(t)$ is the value of $\psi(z)$ at $z=t$. The $\psi(t)$ and $\Phi(t)$
are understood as in Remark 3, see below.

If some equation holds almost everywhere with respect to the Lebesgue measure
on $S$, then we write that this equation holds $a.e.$

From formulas (3)–(6) one derives:

$\Phi^{+}(t)-\Phi^{-}(t)=f(t)+\psi^{+}(t)-\psi^{-}(t)\,\,a.e.,\quad\Phi^{+}(t)+\Phi^{-}(t)=f(t)+\psi^{+}(t)+\psi^{-}(t)\,\,a.e.$
(7)

In Lemma 1, see below, we prove that $\psi^{+}(t)=\psi^{-}(t)=\psi(t)\,\,a.e.$
Therefore, formula (7) can be rewritten as:

$\Phi^{+}(t)-\Phi^{-}(t)=f(t)\,\,a.e.,\quad\Phi^{+}(t)+\Phi^{-}(t)=f(t)+2\psi(t)\,\,a.e.$
(8)

We will need three lemmas. Lemma 1 is proved in Section 2.

Lemma 1. If $f\in L^{1}(S)$ and $S$ is $C^{1,a}-$smooth, $0<a\leq 1$, then

$\psi^{+}(t)=\psi^{-}(t)=\psi(t)\,\,\,a.e.$ (9)

Let $s,t\in S$

Lemma 2. One has

$\frac{1}{s-t\pm i0}=\frac{1}{s-t}\mp i\pi\delta(s-t),$ (10)

where $\delta(s-t)$ is the delta-function.

Proof. Formula (10) is understood as follows. Let $f\in L^{1}(S)$, $\phi\in
H^{\mu}(S)$, $N_{t}$ be a unit normal to $S$ directed out of $D$. Let

$\lim_{\epsilon\to+0}\int_{S}dt\phi(t)\int_{S}\frac{f(s)ds}{s-t-i\epsilon
N_{t}}:=\int_{S}dt\phi(t)\int_{S}\frac{f(s)ds}{s-t-i0}.$

Then

$\int_{S}dt\phi(t)\int_{S}\frac{f(s)ds}{s-t-i0}=\int_{S}dt\phi(t)\int_{S}\frac{f(s)ds}{s-t}+i\pi\int_{S}\phi(t)f(t)dt,\quad\forall\phi\in
H^{\mu}(S),$ (11)

and

$\int_{S}dt\phi(t)\int_{S}\frac{f(s)ds}{s-t}=\int_{S}dsf(s)\int_{S}\frac{\phi(t)dt}{s-t},\quad\forall\phi\in
H^{\mu}(S).$ (12)

Any $f\in L^{1}(S)$ is uniquely identified with the linear functional in
$L^{1}(S)$ defined as follows: $(f,\phi)=\int_{S}f(s)\phi(s)ds$,
$\forall\phi\in H^{\mu}(S)$, because the set $\\{\phi\\}|_{\forall\phi\in
H^{\mu(S)}}$ is dense in $L^{1}(S)$.

It is known, see [2], p. 52, that
$B\phi:=\frac{1}{i\pi}\int_{S}\frac{\phi(s)ds}{s-t}\in H^{\mu}(S)$, provided
that $\phi\in H^{\mu}(S)$ and $\mu\in(0,1)$. Moreover, the range of $B$ is
equal to $H^{\mu}(S)$ when $\phi$ runs through $H^{\mu}(S)$, see [2], p. 178.
This follows from the inversion formulas: if $\phi,\psi\in H^{\mu}(S)$ and
$B\phi=\psi$, then $\phi=B\psi$, see [8], p. 115, so $B^{2}=I$, the identity
operator. Therefore, the right side of (12) defines $f\in L^{1}(S)$ uniquely.
The integrals on the left and right sides of (12) are equal to the integral
$\int_{S\times S}\frac{f(s)\phi(t)dsdt}{s-t}$, which is absolutely convergent
on $S\times S$. We have proved formula (10) with the minus sign. Similarly one
proves this formula with the plus sign. Lemma 2 is proved. $\Box$

One can interpret formula (11) as follows: if $f\in L^{1}(S)$, then

$\int_{S}\frac{f(s)ds}{s-t-i0}=\int_{S}\frac{f(s)ds}{s-t}+i\pi f(t)\quad
a.e.,$ (13)

where the first integral on the right is understood as in formula (12), see
also Remark 3 below, and $f(t)$ is defined $a.e$: almost everywhere on $S$
with respect to the Lebesgue measure on $S$.

Usually the delta-function is defined as a linear continuous functional on the
space of continuous functions, $\int_{S}f(s)\delta(s-t)ds=f(t)$, $f\in C(S)$.
One may define delta-function $\delta(s-t)$ on $L^{1}(S)$ as a linear bounded
functional on $L^{1}(S)$, $\int_{S}f(s)\delta(s-t)ds=f(t)\quad a.e.$, $f\in
L^{(}S)$. In place of $f(t)$ one has to put $f(t)\,a.e.$ because $f\in
L^{1}(S)$ is defined a.e. with respect to the Lebesgue measure on $S$. By
Lusin’s theorem, see [10], p. 157, if $f\in L^{1}(S)$, then $f$ is continuous
on a subset $S_{\epsilon}$ of $S$, $|measS-measS_{\epsilon}|<\epsilon$, where
$\epsilon>0$ can be chosen arbitrarily small and $measS$ is the Lebesgue
measure of the set $S$. One can consider in this setting the delta-function
$\delta(s-t)$ as a kernel of the identity operator in $L^{1}(S)$.

In [4], p. 83, there is a formula $\frac{1}{x-i0}=\frac{1}{x}+i\pi\delta(x)$
understood in the sense of distributions. The formula in Lemma 2 is of the
similar type.

From Lemma 1 and formula (7) one derives the Sokhotsky-Plemelj formulas for
$f\in L^{1}(S)$, see Lemma 3. These formulas are derived in [2] and [8] under
the assumption that $f$ is H$\ddot{o}$lder-continuous. Under such an
assumption, these formulas hold everywhere, not almost everywhere.

Lemma 3. If $f\in L^{1}(S)$, then the Sokhotsky-Plemelj formulas hold:

$\Phi^{+}(t)=\Phi(t)+\frac{f(t)}{2}\,\,a.e.,\quad\Phi^{-}(t)=\Phi(t)-\frac{f(t)}{2}\,\,a.e.$
(14)

Proof of Lemma 3. Formulas (14) follow from formulas (4)–(9). Indeed, by
formulas (7) and (9) one concludes that $\Phi^{+}(t)-\Phi^{-}(t)=f(t)\,a.e.$
and $\Phi^{+}(t)+\Phi^{-}(t)=f(t)+2\psi(t)\,a.e.$ Since
$\Phi(t)=\Phi^{+}(t)-\frac{f(t)}{2}=\Phi^{-}(t)+\frac{f(t)}{2}$, formulas (14)
follow. Lemma 3 is proved. $\Box$

Let $z\in D$. The following question is of interest:

When is the boundary value $\Phi^{+}(t)$ of $\Phi(z)$ on $S$ equal to $f$
a.e.?

From formula (4) the answer follows immediately:

$\Phi^{+}(t)=f(t)$ if and only if (iff) $\psi^{+}(t)=0\,\,a.e.$

Equation (14) yields another necessary and sufficient condition:

$\Phi^{+}(t)=f(t)$ iff $\Phi(z)=0,\,\,z\in D^{\prime}$ and
$f(t)=\Phi(t)+\frac{f(t)}{2}$, or

$f(t)=\frac{1}{i\pi}\int_{S}\frac{f(s)}{s-t}ds\quad a.e.$ (15)

If equation (15) holds, then $\Phi^{-}(t)=0$ and, consequently,
$\Phi(z)=0\,\,\text{if}\,\,z\in D^{\prime}$. Equation (15) is a singular
integral equation for $f$. It is clear from this equation that the singular
integral operator in this equation maps every solution $f\in L^{1}(S)$ of this
equation into itself.

If one wants to formulate a necessary and sufficient condition for $f(s)\in
L^{1}(S)$ to be a boundary value of an analytic in $D^{\prime}$ function
$\Phi(z),\,\,\Phi(\infty)=0$, then an argument, similar to the above yields
the following conditions:

$f(t)=-\frac{1}{i\pi}\int_{S}\frac{f(s)}{s-t}ds\quad a.e.$ (16)

If equation (16) holds, then $\Phi^{+}(t)=0$ and, consequently,
$\Phi(z)=0\,\,\text{if}\,\,z\in D$.

Remark 1. Integral equation in formula (15) has infinitely many linearly
independent solutions: every boundary value of an analytic function $\Phi(z)$
in $D$ such that $\Phi(z)=0$ in $D^{\prime}$ solves this integral equation.

Remark 2. If $\Phi^{+}(t)=f(t)\,\,a.e.,f\in L^{1}(S)$, then $Bf\in L^{1}(S)$,
where $Bf:=\frac{1}{i\pi}\int_{S}\frac{f(s)}{s-t}ds\,\,a.e.$ If
$-\Phi^{-}(t)=f(t)\,\,a.e.,f\in L^{1}(S)$ and $\Phi(\infty)=0$, then $Bf\in
L^{1}(S)$. Since for some $f\in L^{1}(S)$ one does not have
$\Phi(z)=0,\,\,z\in D^{\prime}$, it follows that not every $f\in L^{1}(S)$ is
a boundary value of an analytic function in $D$.

Remark 3. Let us define $\Phi(t)$ as follows.

For any $\phi\in H^{\mu}(S)$, $\mu\in(0,1)$, we define $\Phi(t)$ by the
relation:

$\int_{S}\Phi(t)\phi(t)dt=\int_{S}f(s)\eta(s)ds,\quad\forall\phi\in
H^{\mu}(S),$ (17)

where $\eta(s):=c\int_{S}\frac{\phi(t)}{s-t}dt:=B\phi$.

If $\phi\in H^{\mu}(S)$, $\mu\in(0,1)$, then the function $\eta(s)\in
H^{\mu}(S)$, provided that $S$ is smooth. This is proved in [2], p. 53.

The singular integral equation

$\frac{1}{i\pi}\int_{S}\frac{\phi(s)ds}{s-t}=h(t),\quad\phi\in
H^{\mu}(S),\,\,\mu\in(0,1),$ (18)

can be solved uniquely and explicitly for any $h\in H^{\mu}(S)$:

$\frac{1}{i\pi}\int_{S}\frac{h(s)ds}{s-t}=\phi(t),$ (19)

and its solution belongs to $H^{\mu}(S)$. This is proved in [2], p. 66.

Consequently, any $h\in H^{\mu}(S)$, $\mu\in(0,1)$, can be represented by the
integral (18) with $\phi\in H^{\mu}(S)$. This is in contrast to the integral
equation (18) with $\phi\in L^{1}(S)$: not for every such $\phi$ the $h$
belongs to $L^{1}(S)$, see Remark 5 below.

In [1] it is proved that $\lim_{\epsilon\to 0}\Phi_{\epsilon}(t)$ exists for
a.e. $t\in S$ if $f\in L^{1}(S)$, where
$\Phi_{\epsilon}(t):=c\int_{|t-s|>\epsilon}\frac{f(s)}{s-t}ds$.

There exists a sequence $f_{n}\in
H^{\mu}(S),\,\,lim_{n\to\infty}\|f_{n}-f\|=0$, where the norm is
$L^{1}(S)$-norm. Clearly
$\int_{S}\phi(t)\Big{(}\int_{S}\frac{f_{n}(s)ds}{s-t}\Big{)}dt=\int_{S}f_{n}(s)\Big{(}\int_{S}\frac{\phi(t)dt}{s-t}\Big{)}ds,\quad\forall\phi\in
H^{\mu}(S).$ The limit of the right side, as $n\to 0$, exists, so the limit of
the left side exists and is equal to that of the right side. This leads to the
formula (20) below. An alternative derivation of this formula goes as follows.
Since $\phi(t)f(s)\in L^{1}(S\times S)$, the operator $\int_{S\times
S}\phi(t)f(s)|s-t|^{-1}dsdt$ is a bounded operator in $L^{1}(S\times S)$, see
[5], p.159. Consequently, the change of the order of integration in the double
integral:

$\int_{S}\phi(t)\Big{(}\int_{S}\frac{f(s)ds}{s-t}\Big{)}dt=\int_{S}f(s)\Big{(}\int_{S}\frac{\phi(t)dt}{s-t}\Big{)}ds,\quad\forall\phi\in
H^{\mu}(S),$ (20)

is legitimate.

Let us summarize the results:

Theorem 1. A necessary and sufficient condition for $f\in L^{1}(S)$ to be the
boundary value of an analytic in $D$ function $\Phi(z)$ is equation (15).

A necessary and sufficient condition for $f\in L^{1}(S)$ to be the boundary
value of an analytic in $D^{\prime}$ function $\Phi(z),\,\Phi(\infty)=0$, is
equation (16).

##  2 Proof of Lemma 1.

Write formula (2) for $\psi$ as

$\psi^{\pm}(t)=c\int_{S}\frac{f(s)-f(t)}{s-(t\pm i0)}ds.$ (21)

Using Lemma 2, one obtains:

$\psi^{\pm}(t)=c\int_{S}\frac{f(s)-f(t)}{s-t}ds\pm
c\int_{S}[f(s)-f(t)]\delta(s-t)ds=\psi(t)\,\,\,a.e.,$ (22)

because

$\int_{S}[f(s)-f(t)]\delta(s-t)ds=0\,\,\,a.e.$ (23)

Let us explain (23). The union of the discontinuity points of a summable
function $f$ and the points where the function $f$ is not defined has measure
zero. The $L^{1}(S)$ norm of a summable function is not changed if this
function is changed on the set of measure zero. This yields the conclusion
that equation (23) holds.

Usually the distribution $\delta(s-t)$ is defined as a linear bounded
functional on the space of continuous functions with the $\sup$-norm. On the
other hand, one can define $\delta(s-t)$ as the kernel of the identity
operator $I$ in the space $L^{p}(S),\,p\geq 1$. In this case, the relation
$\int_{S}f(s)\delta(s-t)ds=f(t)$ holds $a.e.$ The expression $f(t)$ for $f\in
L^{p}(S)$, $p\geq 1$, is not defined at every point $t$, but is defined almost
everywhere. The set of points at which $f$ is discontinuous, or is not defined
(including the points at which $f=\infty$) is of measure zero. Outside of this
set, that is, a.e., one has

$If=\int_{S}\delta(s-t)f(s)ds=f(t)\quad a.e.,$ (24)

where $I$ is the identity operator with the kernel $\delta(s-t)$. Formula (24)
gives another way to see that formula (23) is valid. Since in the general
theory of distributions the product of two distributions is not defined
(including the case when one of the distributions is a function, belonging to
$L^{p}(S)$, $p\geq 1$), formula (24) is useful for $f\in L^{p}(S)$, $p\geq 1$.

Finally, recall that for a continuous function $f$ on $S$ equation (23) is
obvious. If $f\in L^{1}(S)$ then the Luzin’s theorem, see, for example, [10],
p.157, says that for any $f\in L^{1}(S)$, where $S$ is bounded, and for an
arbitrary small $\epsilon>0$, there exists a closed subset $Q$ of $S$,
$measQ>measS-\epsilon$, on which $f$ is continuous. Since $\epsilon>0$ is
arbitrary small, equation (23) holds.

Lemma 1 is proved. $\Box$

Remark 4. Let $S=S_{0}$ be the unit circle on the complex plane and

$B_{0}f:=\frac{1}{i\pi}\int_{S_{0}}\frac{f(s)}{s-t}ds.$ (25)

For an arbitrary $f\in L^{1}(S_{0})$ the $B_{0}f$ may be more singular than
$f$.

Example 1. Let

$f:=\sum_{n=1}^{\infty}\frac{e^{in\phi}-e^{-in\phi}}{n}=2i\sum_{n=1}^{\infty}\frac{\sin
n\phi}{n}=2i\frac{\pi-\phi}{2},\quad 0<\phi<2\pi,$ (26)

see [3], p.52. The function (26) is a bounded function on the interval
$0\leq\phi\leq 2\pi$.

One can check, that

$B_{0}f:=\sum_{n=1}^{\infty}\frac{B_{0}e^{in\phi}-B_{0}e^{-in\phi}}{n}=\sum_{n=1}^{\infty}\frac{e^{in\phi}+e^{-in\phi}}{n}=2\sum_{n=1}^{\infty}\frac{\cos
n\phi}{n}=-\frac{1}{2}\ln(2-2\cos\phi),\quad 0<\phi<2\pi.$ (27)

see [3], p.52. The validity of the change of the order of application of
$B_{0}$ and the summation with respect to $n$ follows from the convergence of
the series $\sum_{n=1}^{\infty}\frac{\cos n\phi}{n}$.

We have used the formulas, see [7], p. 47,

$B_{0}e^{in\phi}=e^{in\phi},\quad n>0;\quad
B_{0}e^{-in\phi}=-e^{-in\phi},\quad n>0.$ (28)

Thus, $B_{0}f$ has a singularity at $\phi=0$ on the interval $0\leq\phi\leq
2\pi$, while $f$ is a bounded function.

Remark 5. One can check that if $f\in L^{1}({\mathbb{R}})$ and
$Sf:=\frac{1}{2\pi i}\int_{{\mathbb{R}}}\frac{f(s)ds}{s-t}$, then for some $f$
the $Sf$ is more singular than $L^{1}({\mathbb{R}})$. Indeed, consider $Sf$ as
a convolution in the sense of distributions.

Take the Fourier transform of $Sf$ and use the formula
$F(f*\frac{1}{s})=Ff(i\pi{sign}\xi)$, where $Ff=\int_{{\mathbb{R}}}e^{i\xi
s}f(s)ds$, $*$ denotes the convolution:
$f*h:=\int_{{\mathbb{R}}}f(s)h(t-s)ds$, $Ff*h=F(f)F(h)$. One has
$F\frac{1}{s}=i\pi{sign}\xi$, see [4]. Since the set $Ff$ is the set of
continuous functions of $\xi$ when $f$ runs through $L^{1}({\mathbb{R}})$, the
set of functions $\\{Ff\times{sign}\xi\\}$ is not the set of continuous
functions of $\xi$. Therefore, $Sf$ does not belong to $L^{1}({\mathbb{R}})$
for some $f\in L^{1}({\mathbb{R}})$.

##  3 Conclusion

Let $D$ be a connected bounded domain in ${\mathbb{R}}^{2}$, $S$ be its
boundary, which is closed, connected and smooth. Let
$\Phi=\Phi(z),\,\,z=x+iy,$ be an analytic function in $D$.

A necessary and sufficient condition for a function $f=f(s)\in L^{1}(S)$ to be
a boundary value of an analytic in $D$ function $\Phi(z)$ is equation (15).

A necessary and sufficient condition for a function $f=f(s)\in L^{1}(S)$ to be
a boundary value of an analytic in $D^{\prime}$ function $\Phi(z)$,
$\Phi(\infty)=0$, is equation (16).

## References

  * [1] A. Calderon, Cauchy integral on Lipschitz curves and related operators, Proc. Natl. Acad. Sci. USA, 74, N4, (1977), 1324-1327.
  * [2] F. Gahov, Boundary value problems, Nauka, Moscow, 1977. (in Russian)
  * [3] I. Gradshtein, I. Ryzhik, Tables of integrals, series and products, Gos. Izdat. Fiz.-Math. Lit., Moscow, 1962. (In Russian)
  * [4] I. Gel’fand, G. Shilov, Generalized functions, Vol. 1, Gos. Izdat. Fiz.-Math. Lit., Moscow, 1959. (In Russian)
  * [5] D. Gilbarg, N. Trudinger, Elliptic partial differential equations of second order, Springer-Verlag, Berlin, 1983.
  * [6] B. Khvedelidze, Linear discontinuous boundary value problems of function theory, singular integral equations and some applications, Trudy Tbilisskogo math. instituta Akad. Nauk Grusinskoi SSR, 23, (1956), 3-158.
  * [7] S. Mikhlin, S. Pr$\ddot{o}$ssdorf, Singular integral operators, Springer-Verlag, New York, 1986.
  * [8] N. Muskhelishvili, Singular integral equations, Nauka, Moscow, 1968. (In Russian)
  * [9] I. Privalov, Boundary values of univalent analytic functions, Gostekhizdat, Moscow, 1950. (In Russian)
  * [10] G. Shilov, Mathematical analysis, Fizmatgiz, Moscow, 1960


# Can we distinguish between black holes and wormholes by their Einstein-ring
systems?

Naoki Tsukamoto,111Electronic<EMAIL_ADDRESS>Tomohiro Harada
and Kohji Yajima  Department of Physics, Rikkyo University, Tokyo 171-8501,
Japan

###### Abstract

For the last decade, gravitational lensing in the strong gravitational field
has been studied eagerly. It is well known that, for the lensing by a black
hole, an infinite number of Einstein rings are formed by the light rays which
wind around the black hole nearly on the photon sphere, which are called
relativistic Einstein rings. This is also the case for the lensing by a
wormhole. In this paper, we study the Einstein ring and relativistic Einstein
rings for the Schwarzschild black hole and the Ellis wormhole, the latter of
which is an example of traversable wormholes of the Morris-Thorne class. Given
the configuration of the gravitational lensing and the radii of the Einstein
ring and relativistic Einstein rings, we can distinguish between a black hole
and a wormhole in principle. We conclude that we can detect the relativistic
Einstein rings by wormholes which have the radii of the throat $a\simeq 0.5$pc
at a galactic center with the distance $10$Mpc and which have $a\simeq 10$AU
in our galaxy using by the most powerful modern instruments which have the
resolution of $10^{-2}$arcsecond such as a $10$-meter optical-infrared
telescope. The black holes which make the Einstein rings of the same size as
the ones by the wormholes are galactic supermassive black holes and the
relativistic Einstein rings by the black holes are too small to measure with
the current technology. We may test the hypotheses of astrophysical wormholes
by using the Einstein ring and relativistic Einstein rings in the future.

###### pacs:

04.20.-q, 04.70.-s

††preprint:  RUP -11-5

##  I I. Introduction

Gravitational lensing is a very useful tool for astrophysics and cosmology. At
first the gravitational lensing mainly was investigated on a theoretical basis
in the weak gravitational field. Using the gravitational lensing, we determine
the cosmological constant, the distribution of dark matter and the Hubble
constant, the existence of extrasolar planets and so on (see Schneider et al.
Gravitational_lenses  and Perlick Perlick_2004_Living_Rev ; Perlick_2010  for
the detail of the gravitational lens, and references therein).

For the last decade, gravitational lensing in the strong gravitational field
has been studied eagerly (see Virbhadra and Keeton Virbhadra_Keeton_2008 ,
Virbhadra Virbhadra_2009 , Bozza Bozza_2010 , Bozza and Mancini
Bozza_Mancini_2012  and references therein). Frittelli et al.
Frittelli_Kling_Newman_2000 , Virbhadra and Ellis Virbhadra_Ellis_2000 ;
Virbhadra_Ellis_2002  and Bozza et al. Bozza_Capozziello_Iovane_Scarpetta_2001
studied the gravitational lensing in the strong field with the Schwarzschild
spacetime and found the infinite Einstein rings which are too close to each
other to separately resolve. In this paper, we call these rings relativistic
Einstein rings. The gravitational lensing in the strong field on the
spherically symmetric static spacetime was investigated by Bozza Bozza_2002 ,
Hasse and Perlick Hasse_Perlick_2001  and Perlick Perlick_2004_Phys_Rev_D .
They showed that the relativistic Einstein rings are formed not only in the
Schwarzschild spacetime but also in the other spherically symmetric static
spacetime.

General relativity permits nontrivial topology of the spacetime such as
wormhole spacetimes (see Visser Lorentzuan_Wormholes  for the details of
wormholes). Some hypotheses of astrophysical wormholes have been investigated
Harko_Kovacs_Lobo_2009 ; Abdujabbarov_Ahmedov_2009 ; Pozanenko_Shatskiy_2010 .
For example, Kardashev et al. suggest that some active galactic nuclei and
other compact astrophysical objects may be explained as wormholes
Kardashev_Novikov_Shatskiy_2007 . We may test these hypotheses by using the
gravitational lensing in the future.

Kim and Cho Kim_Cho_1994  and Cramer et al.
Cramer_Forward_Morris_Visser_Benford_Landis_1995  pioneered gravitational
lensing effects by wormholes. Since then, the gravitational lensing effects by
various wormholes have been investigated Rahaman_Kalam_Chakraborty_2007 ;
Safonova_Torres_2002 ; Safonova_Torres_Romero_2002 ;
Safonova_Torres_Romero_2001 ; Eiroa_Romero_Torres_2001 ;
Nandi_Zhang_Zakharov_2006 ; Cramer_Forward_Morris_Visser_Benford_Landis_1995 .

The Ellis spacetime which was investigated by Ellis Ellis_1973  is an example
of traversable wormholes of the Morris-Thorne class Morris_Thorne_1988 ;
Morris_Thorne_Yurtsever_1988 . The deflection angle of light in the Ellis
wormhole geometry was studied by Chetouani and Clement Chetouani_Clement_1984
and recently Nakajima and Asada Nakajima_Asada_2012 . The gravitational
lensing on the Ellis geometry was studied by Dey and Sen Dey_Sen_2008 , Abe
Abe_2010  and Toki et al. Toki_Kitamura_Asada_Abe_2011  in the weak
gravitational field and Perlick Perlick_2004_Phys_Rev_D , Nandi et al.
Nandi_Zhang_Zakharov_2006  and Tejeiro and Larranaga Tejeiro_Larranaga_2005
in the strong gravitational field.

Perlick Perlick_2004_Phys_Rev_D , Nandi et al. Nandi_Zhang_Zakharov_2006  and
Tejeiro and Larranaga Tejeiro_Larranaga_2005  pointed out that the qualitative
features of the gravitational lensing in the Ellis spacetime are very similar
to the ones in the Schwarzschild spacetime for their photon spheres and their
asymptotic flatness.

In this paper, we will consider the Einstein ring and relativistic Einstein
rings in the Ellis spacetime and the Schwarzschild spacetime, both of which
are static and spherically symmetric ones. We ask whether we can distinguish
the Einstein-ring systems on the Schwarzschild spacetime and on the Ellis
spacetime. To answer this question, we focus on the relations between the
Einstein ring and the relativistic Einstein rings.

This paper is organized as follows. In Sec. II we will review the deflection
angle on the Ellis wormhole spacetime. In Sec. III we give the radii of the
Einstein ring and relativistic Einstein rings in the Ellis spacetime. In Sec.
IV we will compare the Einstein ring and the relativistic Einstein rings in
the Ellis spacetime to the ones in the Schwarzschild spacetime. In Sec. V we
summarize and discuss our result. In this paper we use the units in which
$c=1$.

##  II II. Ellis wormhole spacetime and deflection angle

In this section, we review the deflection angle on the Ellis wormhole
spacetime Chetouani_Clement_1984 ; Nakajima_Asada_2012 . The line element in
the Ellis wormhole solution is written in the following form:

$ds^{2}=-dt^{2}+dr^{2}+(r^{2}+a^{2})d\Omega^{2},$ (1)

where $d\Omega^{2}=d\theta^{2}+\sin^{2}\theta d\phi^{2}$ and $a$ is a positive
constant. Introducing $\rho^{2}=r^{2}+a^{2}$, we can rewrite this into

$ds^{2}=-dt^{2}+\left(1-\frac{a^{2}}{\rho^{2}}\right)^{-1}d\rho^{2}+\rho^{2}d\Omega^{2},$
(2)

where $\rho=\pm a$ corresponds to the wormhole throat. The spacetime has the
Killing vectors $t^{\mu}\partial_{\mu}=\partial_{t}$ and
$\phi^{\mu}\partial_{\mu}=\partial_{\phi}$ for stationarity and axial
symmetry.

We can concentrate ourselves on the equatorial plane because of spherical
symmetry. Using the conservation of the energy $E\equiv-
g_{\mu\nu}k^{\mu}t^{\nu}$ and angular momentum $L\equiv
g_{\mu\nu}k^{\mu}\phi^{\nu}$ and $k^{\mu}k_{\mu}=0$, where $k^{\mu}$ is the
photon wave number, the photon trajectory is then given by

$\frac{1}{\rho^{4}}\left(\frac{d\rho}{d\phi}\right)^{2}=\frac{1}{b^{2}}\left(1-\frac{a^{2}}{\rho^{2}}\right)\left(1-\frac{b^{2}}{\rho^{2}}\right),$
(3)

where $b\equiv L/E$ is the impact parameter of the photon.

We can see that the photon is scattered if $|b|>a$, while it reaches the
throat if $|b|<a$. Since we are interested in the scattering problem, we
assume $|b|>a$. Using $u=1/\rho$, we find

$\left(\frac{du}{d\phi}\right)^{2}=\frac{1}{b^{2}}(1-a^{2}u^{2})(1-b^{2}u^{2}).$
(4)

Putting

$G(u)=a^{2}(a^{-2}-u^{2})(b^{-2}-u^{2}),$ (5)

the azimuthal angle $\phi$ can be given as a function of $u$ by

$\phi=\pm\int^{b^{-1}}_{u}\frac{du}{\sqrt{G(u)}}.$ (6)

Here we have set $\phi(b^{-1})=0$. The deflection angle $\alpha$ is then
calculated to give

$\alpha=2\int^{b^{-1}}_{0}\frac{du}{\sqrt{G(u)}}-\pi.$ (7)

In the present case, we find

$\int^{b^{-1}}_{0}\frac{du}{\sqrt{G(u)}}=\int^{\pi/2}_{0}\frac{d\theta}{\sqrt{1-\left(\frac{a}{b}\right)^{2}\sin^{2}\theta}}=K\left(\frac{a}{b}\right),$
(8)

where we have transformed $u=b^{-1}\sin\theta$ and $K(k)$ denotes the complete
elliptic integral of the first kind (for example, see
Handbook_of_Elliptic_Integrals_for_Engineers_and_Scientists ), which is
defined as

$K(k)=\int^{\pi/2}_{0}\frac{d\theta}{\sqrt{1-k^{2}\sin^{2}\theta}}.$ (9)

Hence, the deflection angle is given by

$\alpha=2K\left(\frac{a}{b}\right)-\pi.$ (10)

Since $K(k)$ admits a power series

$K(k)=\frac{\pi}{2}\sum_{n=0}^{\infty}\left[\frac{(2n-1)!!}{(2n)!!}\right]^{2}k^{2n},$
(11)

where $n!!$ denotes the double factorial of $n$ and $(-1)!!=1$, we get the
deflection angle

$\alpha=\pi\sum_{n=1}^{\infty}\left[\frac{(2n-1)!!}{(2n)!!}\right]^{2}\left(\frac{a}{b}\right)^{2n}.$
(12)

Thus, the deflection angle is approximately given in the weak-field regime
$|b|\gg a$ by

$\alpha\simeq\frac{\pi}{4}\left(\frac{a}{b}\right)^{2}.$ (13)

In general, the deflection angle is always greater than its weak-field
approximation and is diverging as $|b|\to a$.

##  III III. Einstein ring and relativistic Einstein rings of Ellis wormhole

In this section, we examine the diameter angles of the Einstein ring and the
relativistic Einstein rings on the Ellis spacetime. Now we will consider the
case that both the observer and the source object are far from the lensing
object, or $D_{l}\gg b$ and $D_{ls}\gg b$, where $D_{l}$ and $D_{ls}$ are the
separations between the observer and lens and between the lens and source,
respectively. The configuration of the gravitational lensing is given in Fig.
1.

Figure 1: The configuration of the gravitational lensing. The light rays
emitted by the source $S$ are deflected by the lens $L$ (a wormhole or a black
hole) and reach the observer $O$ with the angle of the lensed image $\theta$,
instead of the real angle $\phi$. $b$ and $\bar{\alpha}$ are the impact
parameter and the effective deflection angle, respectively. $D_{l}$ and
$D_{ls}$ are the separations between the observer and the lens and between the
lens and the source, respectively.

Then, the lens equation is given by

$D_{ls}\bar{\alpha}=D_{s}(\theta-\phi),$ (14)

where $\bar{\alpha}=(\alpha~{}\mbox{mod}~{}2\pi)$ is the effective deflection
angle, $\theta$ and $\phi$ are the angles of the lensed image and the real
image from the observer, respectively, and $D_{s}=D_{l}+D_{ls}$ is the
separation between the observer and source. Note that we have assumed
$|\bar{\alpha}|\ll 1$, $|\theta|\ll 1$ and $|\phi|\ll 1$. The deflection angle
can be expressed $\alpha=\bar{\alpha}+2\pi n$, where $n$ is a non-negative
integer, denoting the winding number of the light ray.

The ring image corresponds to the image angle $\theta$ for vanishing real
angle $\phi=0$. By the symmetry, the image is necessarily a ring with the
diameter angle $\theta$.

Since $b=D_{l}\theta$, we find that the ring image is given by

$\theta_{n}=\frac{a}{D_{l}}\frac{1}{k_{n}},$ (15)

where $k_{n}\in(0,1)$ is a unique root of the transcendental equation

$\displaystyle 2K(k)-\frac{\eta}{k}$ $\displaystyle=$
$\displaystyle\left(2n+1\right)\pi,$ (16) $\displaystyle\eta$ $\displaystyle=$
$\displaystyle\frac{D_{s}}{D_{l}D_{ls}}a.$ (17)

We should note that $2K(k)-\eta/k$ is monotonically increasing with respect to
$k$ and changes from $-\infty$ to $\infty$ as $k$ increases from $0$ to $1$.
The uniqueness of the root follows from the monotonicity. Moreover, we can
conclude that $k_{n}$ monotonically increases and approaches 1 as $n\to\infty$
and hence the image angle $\theta_{n}$ monotonically decreases and approaches
$a/D_{l}$.

In the weak-field regime $|b|\gg a$, the winding number $n$ should be $n=0$.
Using the deflection angle (13), we can solve the equation (16) approximately
and get the diameter angle of the Einstein ring

$\displaystyle\theta_{0}$
$\displaystyle\simeq\left(\frac{\pi}{4}\frac{D_{ls}}{D_{s}D_{l}^{2}}a^{2}\right)^{\frac{1}{3}}$
(18) $\displaystyle\simeq
2.0\,\textrm{arcsecond}\left(\frac{D_{ls}}{10\textrm{Mpc}}\right)^{\frac{1}{3}}\left(\frac{20\textrm{Mpc}}{D_{s}}\right)^{\frac{1}{3}}\left(\frac{10\textrm{Mpc}}{D_{l}}\right)^{\frac{2}{3}}\left(\frac{a}{0.5\textrm{pc}}\right)^{\frac{2}{3}}.\qquad$

This approximation is good for $D_{l}\gg a$ and $D_{ls}\gg a$. The relative
error is $\sim 10^{-2}$ for $a=0.5$pc and $D_{l}=D_{ls}=10$ Mpc.

In the especially strong-field regime, where the winding number $n$ becomes
$n\geq 1$, we can easily check that $a\simeq b$ or $k_{n}\simeq 1$ satisfies
the transcendental equation (16) in numerical calculations. Physically this
means that the light rays which wind around the wormhole nearly on the photon
sphere make the relativistic Einstein rings Bozza_2002 ;
Perlick_2004_Phys_Rev_D . Then the diameter angles of the relativistic
Einstein rings are approximately given by

$\displaystyle\theta_{n\geq 1}$ $\displaystyle\simeq$
$\displaystyle\frac{a}{D_{l}}$ (19) $\displaystyle\simeq$ $\displaystyle
1.0\times
10^{-2}\,\textrm{arcsecond}\left(\frac{10\textrm{Mpc}}{D_{l}}\right)\left(\frac{a}{0.5\textrm{pc}}\right).\qquad$

Regardless of the values of $D_{ls}$, $D_{l}$ and $a$, the relative error of
the above approximation to the direct numerical solution of the outermost
relativistic Einstein ring ($n=1$) is $\sim 10^{-3}$ and those of the other
relativistic Einstein rings ($n\geq 2$) are smaller than $10^{-5}$. This
implies that it is difficult to resolve each relativistic Einstein ring
separately.

Thus, we conclude that there is one Einstein-ring image and countably infinite
relativistic Einstein-ring images, the latter of which accumulate to form the
apparently single ring image of the throat with the diameter $a/D_{l}$. This
conclusion does not depend on the value of $\eta$.

If we are given the distance $D_{s}$ to the source from the observer, the
distance $D_{l}$ to the lens from the observer and the radius $\theta_{0}$ of
the Einstein ring, we can determine the radius of the throat $a$ from Eq.
(18). Then, we can use $\theta_{n\geq 1}$ (19) to test the assumption that the
lens object is a wormhole.

From Eqs. (18) and (19) we obtain the relation between $\theta_{0}$ and
$\theta_{n\geq 1}$ by

$\displaystyle\theta_{n\geq
1}\simeq\left(\frac{4}{\pi}\frac{D_{s}}{D_{ls}}\right)^{\frac{1}{2}}\theta_{0}^{\frac{3}{2}}.$
(20)

This relation generally holds in astrophysical situations, as long as $a\ll
D_{l}$ and $a\ll D_{ls}$ are satisfied. If the lens is identified with a
wormhole, we can even estimate the radius of the wormhole throat in terms of
$D_{s}$, $\theta_{0}$ and $\theta_{n}$ through

$a\simeq\theta_{n\geq 1}\left(1-\frac{4}{\pi}\theta_{n\geq
1}^{-2}\theta_{0}^{3}\right)D_{s}.$ (21)

This follows from Eqs. (19) and (20).

##  IV IV. Comparison between wormholes and black holes

In this section we compare the Einstein ring and the relativistic Einstein
rings in the Ellis spacetime to the ones in the Schwarzschild spacetime and
show that we can distinguish between black holes and wormholes.

Now, we will briefly review the deflection angle, Einstein ring and
relativistic Einstein rings for the Schwarzschild spacetime
Bisnovatyi_Tsupko_2008 ; Bozza_Capozziello_Iovane_Scarpetta_2001 ; Muller_2008
; Virbhadra_Ellis_2000 ; Virbhadra_Ellis_2002 ; Hasse_Perlick_2001  and
present the relation between $\theta_{0}$ and $\theta_{n\geq 1}$.

In the weak-field regime $b\gg r_{g}$, where $r_{g}=2GM$ is the Schwarzschild
radius of the black hole of mass $M$, the deflection angle is approximately
given by

$\alpha\simeq\frac{2r_{g}}{b}.$ (22)

The diameter angle of the Einstein ring is given by

$\displaystyle\theta_{0}$
$\displaystyle\simeq\sqrt{\frac{2D_{ls}}{D_{l}D_{s}}r_{g}}$ (23)
$\displaystyle\simeq
2.0\,\textrm{arcsecond}\left(\frac{D_{ls}}{10\textrm{Mpc}}\right)^{\frac{1}{2}}\left(\frac{M}{10^{10}M_{\odot}}\right)^{\frac{1}{2}}\left(\frac{10\textrm{Mpc}}{D_{l}}\right)^{\frac{1}{2}}\left(\frac{20\textrm{Mpc}}{D_{s}}\right)^{\frac{1}{2}}.\qquad$

We can determine $r_{g}$ in the same way as the radius of the throat $a$.

In the especially strong-field regime, where the winding number $n$ becomes
$n\geq 1$, the impact parameter $b$ that satisfies the lens equation should be
nearly the critical impact parameter $b\simeq(3\sqrt{3}/2)r_{g}$ (see
Bozza_Capozziello_Iovane_Scarpetta_2001 ; Bozza_2002 ; Virbhadra_Ellis_2000 ).
Then the diameter angles of the inseparable relativistic Einstein rings
$\theta_{n\geq 1}$ are given by

$\displaystyle\theta_{n\geq 1}$ $\displaystyle\simeq$
$\displaystyle\frac{3\sqrt{3}}{2}\frac{r_{g}}{D_{l}}$ (24)
$\displaystyle\simeq$ $\displaystyle 5.1\times
10^{-5}\,\textrm{arcsecond}\left(\frac{M}{10^{10}M_{\odot}}\right)\left(\frac{10\textrm{Mpc}}{D_{l}}\right).\qquad$

It is useful to remember that the leading term of the deflection angle in the
weak-field regime is the second order of the small amount $a/b$ on the Ellis
geometry (13), while it is the first order of the small amount $r_{g}/b$ on
the Schwarzschild geometry (22). So the relation between $\theta_{0}$ and
$\theta_{n\geq 1}$ for the Schwarzschild spacetime

$\displaystyle\theta_{n\geq
1}\simeq\frac{3\sqrt{3}}{4}\frac{D_{s}}{D_{ls}}\theta_{0}^{2}$ (25)

is different from that on the Ellis spacetime (20). Figure. 2 shows the angle
of the relativistic Einstein ring $\theta_{n\geq 1}$ versus the angle of the
Einstein ring $\theta_{0}$ for $D_{l}=D_{ls}=10$Mpc. Thus, we can distinguish
between black holes and wormholes in principle if we are given $D_{s}$,
$D_{l}$, $\theta_{0}$ and $\theta_{n\geq 1}$.

Figure 2: The angle of the relativistic Einstein ring $\theta_{n\geq 1}$
versus the angle of the Einstein ring $\theta_{0}$ for $D_{l}=D_{ls}=10$Mpc.
The broken (green) and solid (red) lines plot the cases where the lens objects
are a wormhole and a black hole, respectively.

If the lens is identified with a black hole, we can estimate the Schwarzschild
radius and,hence, the black hole mass by

$r_{g}\simeq\frac{2}{3\sqrt{3}}\theta_{n\geq
1}\left(1-\frac{3\sqrt{3}}{4}\theta_{0}^{2}\theta_{n\geq 1}^{-1}\right)D_{s}$
(26)

in terms of $\theta_{0}$, $\theta_{n\geq 1}$ and $D_{s}$. This follows from
Eqs. (24) and (25).

##  V V. Discussion and conclusion

It is well known that the qualitative features of the gravitational lensing on
the Ellis spacetime are very similar to the ones on the Schwarzschild
spacetime for their photon spheres and their asymptotic flatness
Perlick_2004_Phys_Rev_D ; Nandi_Zhang_Zakharov_2006 ; Tejeiro_Larranaga_2005 .
However, we realize that their quantitative features are very different due to
their different weak-field behaviors.

We consider the experimental situation where we know the separation $D_{s}$
between the observer and the source and the separation $D_{l}$ between the
observer and the lens. We assume that we do not know whether the lens object
is a black bole or a wormhole and do not know its parameter, i.e., the mass
$M$ or the radius $a$ of the throat in advance.

We need at least two observable quantities to determine whether the lens
object is a black hole or wormhole since the lens system has one parameter in
this situation. First, we observe an Einstein ring and determine the parameter
for both possibilities. Second, we observe relativistic Einstein rings and
tell the wormhole from the black hole. If the predicted relativistic ring
angles by the black hole and by the wormhole were of similar size, we could
not discern the difference. However, Eqs. (20), (25) and Fig. 2 show that we
do not confuse them.

We conclude that we can detect the relativistic Einstein rings by wormholes
which have $a\simeq 0.5$pc at a galactic center with the distance
$D_{l}=D_{ls}=10$Mpc and which have $a\simeq 10$AU in our galaxy with the
distance $D_{l}=D_{ls}=10$kpc using the most powerful modern instruments which
have the resolution of $10^{-2}$arcsecond such as a $10$-meter optical-
infrared telescope. Note that the corresponding black holes which have the
Einstein rings of the same size are galactic supermassive black holes with
$10^{10}M_{\odot}$ and $10^{7}M_{\odot}$, respectively, and that the
relativistic Einstein rings by these black holes are too small to measure with
the current technology.

In fact, our results imply that we can distinguish between slowly rotating
Kerr-Newmann black holes and the Ellis wormholes with their Einstein-ring
systems. This is because the leading term of the deflection angle for the
lensing by the Kerr-Newmann black holes in the weak-field regime is equal to
the one for the lensing by the Schwarzschild black holes, while the black hole
charge and small spin only slightly change the radii of the relativistic
Einstein rings Bozza_2002 ; Eiroa_Romero_Torres_2002 ; Sereno_2004 ;
Bozza_2003 ; Sereno_Luca_2006 . Moreover, this also suggests that it is much
more challenging to determine the charge and/or small spin of black holes than
to distinguish between black holes and the Ellis wormholes.

We assumed that the observer, the lensing object and the source object are
directly-aligned, though such a configuration is fairly rare. In general the
strong gravitational lensing effect is observed as broken-ring images which
are called relativistic images Virbhadra_Ellis_2000 . Therefore, more
realistic problem is to size the relativistic images. Our result suggests that
we can distinguish black holes and wormholes by using the relativistic images.
To observe the relativistic images is one of the challenging works with many
difficulties. Bozza et al. pointed out that relativistic images are always
very faint with respect to the weak field images
Bozza_Capozziello_Iovane_Scarpetta_2001 . The Very Large Telescope
Interferometer (VLTI) has high resolution Delplancke_Gorski_Richichi_2001 ;
Eckart_Bertram_et_al_2003  but it will not work because of this demagnifying
effect.

We also assumed point-like sources, although astrophysical sources have their
own size. If the source object is a galaxy, it may conceal the relativistic
Einstein rings, especially, in the case that the lens object is a black hole.
Testing some hypotheses of astrophysical wormholes by using the relativistic
Einstein rings and the Einstein ring is left as future work.

Tejeiro and Larranaga Tejeiro_Larranaga_2005  investigated the gravitational
lensing effect of the wormhole solution obtained by connecting the Ellis
solution as an interior region and the Schwarzschild solution as an exterior
region Lemos_Lobo_Quinet_2003 . They concluded that we cannot distinguish the
Schwarzschild black hole and the wormhole unless the discontinuity of the
magnification curve at the boundary is observed. This does not contradict our
results because their wormhole solution behaves as the Schwarzschild solution
in the weak-field regime.

Acknowledgements:

The authors would like to thank U. Miyamoto, H. Asada and F. Abe for valuable
comments and discussion. The work of N. T. and K. Y. was supported in part by
Rikkyo University Special Fund for Research. T. H. would also thank the
Astronomy Unit, Queen Mary, University of London for its hospitality. T. H.
was supported by the Grant-in-Aid for Young Scientists (B) (No. 21740190) and
the Grant-in-Aid for Challenging Exploratory Research (No. 23654082) for
Scientific Research Fund of the Ministry of Education, Culture, Sports,
Science and Technology, Japan. N. T. thanks the Yukawa Institute for
Theoretical Physics at Kyoto University, where this work progressed during the
YITP workshop YITP-W-12-08 on ”Summer School on Astronomy & Astrophysics
2012”.

## References

  * (1) P. Schneider, J. Ehlers and E. E. Falco, Gravitational Lenses (Springer-Verlag, Berlin, 1992).
  * (2) V. Perlick, Living Rev. Relativity 7, 9 (2004), http://relativity.livingreviews.org/Articles/lrr-2004-9.
  * (3) V. Perlick, arXiv:1010.3416.
  * (4) K. S. Virbhadra and C. R. Keeton, Phys. Rev. D 77, 124014 (2008).
  * (5) K. S. Virbhadra, Phys. Rev. D 79, 083004 (2009).
  * (6) V. Bozza, Gen. Relativ. Gravit. 42, 2269 (2010).
  * (7) V. Bozza and L. Mancini, Astrophys. J. 753, 56 (2012).
  * (8) S. Frittelli, T. P. Kling and E. T. Newman, Phys. Rev. D 61, 064021 (2000).
  * (9) K. S. Virbhadra and G. F. R. Ellis, Phys. Rev. D 62, 084003 (2000).
  * (10) K. S. Virbhadra and G. F. R. Ellis, Phys. Rev. D 65, 103004 (2002).
  * (11) V. Bozza, S. Capozziello, G. Iovane and G. Scarpetta, Gen. Relativ. Gravit. 33, 1535 (2001).
  * (12) V. Bozza, Phys. Rev. D 66, 103001 (2002).
  * (13) W. Hasse and V. Perlick Gen. Relativ. Gravit. 34, 415 (2001).
  * (14) V. Perlick, Phys. Rev. D 69, 064017 (2004).
  * (15) M. Visser, Lorentzian Wormholes: From Einstein to Hawking (American Institute of Physics, Woodbury, NY, 1995).
  * (16) T. Harko, Z. Kovacs and F. S. N. Lobo, Phys. Rev. D 79, 064001 (2009).
  * (17) A. A. Abdujabbarov and B. J. Ahmedov, Astrophys. Space Sci. 321, 225, (2009).
  * (18) A. Pozanenko and A. Shatskiy, arXiv:1007.3620.
  * (19) N. S. Kardashev, I. D. Novikov and A. A. Shatskiy, Int. J. Mod. Phys. D 16, 909 (2007).
  * (20) S. W. Kim and Y. M. Cho, in Evolution of the Universe and its Observational Quest (Universal Academy Press, Tokyo, 1994), p. 353.
  * (21) J. G. Cramer, R. L. Forward, M. S. Morris, M. Visser, G. Benford and G. A. Landis, Phys.Rev. D 51, 3117 (1995).
  * (22) K. K. Nandi, Y. Z. Zhang and A. V. Zakharov, Phys. Rev. D 74, 024020 (2006).
  * (23) M. Safonova, D. F. Torres and G. E. Romero, Mod. Phys. Lett. A 16, 153 (2001).
  * (24) E. Eiroa, G. E. Romero and D. F. Torres, Mod. Phys. Lett. A 16, 973 (2001).
  * (25) M. Safonova, D. F. Torres and G. E. Romero, Phys. Rev. D 65, 023001 (2001).
  * (26) M. Safonova and D. F. Torres, Mod. Phys. Lett. A 17, 1685 (2002).
  * (27) F. Rahaman, M. Kalam and S. Chakraborty, Chin. J. Phys. 45, 518 (2007).
  * (28) H. G. Ellis, J. Math. Phys. 14, 104 (1973).
  * (29) M. S. Morris and K. S. Thorne, Am. J. Phys. 56, 395 (1988).
  * (30) M. S. Morris, K. S. Thorne and U. Yurtsever, Phys. Rev. Lett. 61, 1446 (1988).
  * (31) L. Chetouani and G. Clément, Gen. Relativ. Gravit. 16, 111 (1984).
  * (32) K. Nakajima and H. Asada, Phys. Rev. D 85, 107501 (2012).
  * (33) T. K. Dey and S. Sen, Mod. Phys. Lett. A, 23, 953 (2008).
  * (34) F. Abe, Astrophys. J. 725, 787 (2010).
  * (35) Y. Toki, T. Kitamura, H. Asada and F. Abe, Astrophys. J. 740, 121 (2011).
  * (36) J. M. Tejeiro and E. A. Larranaga, Rom. J. Phys. 57, 736 (2012).
  * (37) P. F. Byrd and M. D. Friedman, Handbook of Elliptic Integrals for Engineers and Scientists (Springer-Verlag, New York, 1971).
  * (38) T. Müller, Phys. Rev. D 77, 124042 (2008).
  * (39) G. S. Bisnovatyi-Kogan and O. Yu. Tsupko, Astrophys. 51, 99 (2008).
  * (40) E. F. Eiroa, G. E. Romero and D. F. Torres, Phys. Rev. D 66, 024010 (2002).
  * (41) M. Sereno, Phys. Rev. D 69, 023002 (2004).
  * (42) V. Bozza, Phys. Rev. D 67, 103006 (2003).
  * (43) M. Sereno and F. De Luca, Phys. Rev. D 74, 123009 (2006).
  * (44) F. Delplancke, K. Gorski and A. Richichi, Astron. Astrophys. 375, 701 (2001).
  * (45) A. Eckart, T. Bertram, N. Mouawad, T. Viehmann, C. Straubmeier and J. Zuther, Ap&SS, 286, 269 (2003).
  * (46) J. P. S. Lemos, F. S. N. Lobo and S. Q. de Oliveira, Phys. Rev. D 68, 064004 (2003).


Multimode-cavity picture of non-Markovian waveguide QED

Dario Cilluffo 0000-0001-6862-0511
Institut für Theoretische Physik and IQST, Albert-Einstein-Allee 11, Universität Ulm, 89069 Ulm, Germany

Luca Ferialdi0000-0003-4016-3800
Universit$\grave{a}$  degli Studi di Palermo, Dipartimento di Fisica e Chimica – Emilio Segr$\grave{e}$, via Archirafi 36, I-90123 Palermo, Italy

G. Massimo Palma0000-0001-7009-4573
Universit$\grave{a}$  degli Studi di Palermo, Dipartimento di Fisica e Chimica – Emilio Segr$\grave{e}$, via Archirafi 36, I-90123 Palermo, Italy
NEST, Istituto Nanoscienze-CNR, Piazza S. Silvestro 12, 56127 Pisa, Italy

Giuseppe Calajò0000-0002-5749-2224
Istituto Nazionale di Fisica Nucleare (INFN), Sezione di Padova, I-35131 Padova, Italy

Francesco Ciccarello0000-0002-6061-1255
Universit$\grave{a}$  degli Studi di Palermo, Dipartimento di Fisica e Chimica – Emilio Segr$\grave{e}$, via Archirafi 36, I-90123 Palermo, Italy
NEST, Istituto Nanoscienze-CNR, Piazza S. Silvestro 12, 56127 Pisa, Italy

We introduce a picture to describe and intrepret waveguide-QED problems in the non-Markovian regime of long photonic retardation times resulting in delayed coherent feedback. The framework is based on an intuitive spatial decomposition of the waveguide into blocks. Among these, the block directly coupled to the atoms embodies an effective lossy multimode cavity leaking into the rest of the waveguide, in turn embodying an effective white-noise bath. The dynamics can be approximated by retaining only a finite number of cavity modes that yet eventually grows with the time delay. This description captures the atomic as well as the field's dynamics, even with many excitations, in both emission and scattering processes. As an application, we show that the recently identified non-Markovian steady states can be understood by retaining very few or even only one cavity modes.

Most quantum optics phenomena investigated so far occur under Markovian conditions (lack of memory effects). One major reason behind this is that atom-photon interaction is usually weak, while light travels very fast. This allows to neglect photonic time delays (retardation times), which is a tremendous simplification of the dynamics underpinning standard tools such as the Lindblad master equation and the input-output formalism <cit.>. Recent years have yet seen a growing attention to the non-Markovian regime of non-negligible photonic time delays, especially in the emerging area of waveguide Quantum ElectroDynamics (QED) which generally investigates the coherent interaction between quantum emitters and the one-dimensional (1D) field of a waveguide <cit.>. Such non-Markovian regime can today be accessed in some waveguide-QED experiments, through superconducting qubits coupled to surface acoustic waves <cit.> or slow-light modes (near a band edge) <cit.> and even cold atoms coupled to fiber-ring resonators <cit.>. While complicating the dynamics considerably, time delays can be leveraged for a variety of unprecedented phenomena and applications, such as persistent quantum beats <cit.>, stabilization of Rabi oscillations <cit.>, peculiar inelastic two-photon scattering <cit.>, generation of photonic cluster states <cit.>, excitation of dressed bound states in the continuum <cit.>, enhanced Dicke superradiance <cit.>, anomalous population trapping <cit.>, stationary oscillations of giant atoms <cit.>, enhanced energy-time entanglement <cit.>, improved single-photon sources <cit.>, genuinely non-Markovian steady states <cit.>.

Basic setup and idea. (a) A two-level atom is coupled to a semi-infinite waveguide, whose left end acts like a perfect mirror. We conveniently decompose the waveguide into a pair of blocks $A$ and $B$ mutually coupled with rate $\gamma$. Block $A$ is an effective multimode and intrinsically open cavity (of length $L\gtrsim x_0$) in contact with block $B$ (the latter one working as a semi-infinite waveguide itself). (b) Normal frequencies $\{\Omega_\nu\}$ of block $A$. The atom (energy $\omega_0$) resonates with mode of frequency $\Omega_0$ and is detuned from off-resonant modes $\Omega_{\pm 1}$, $\Omega_{\pm 2}$,...  The frequency spacing of $A$ modes $\Delta \Omega$ scales as $\sim 1/\tau$ with $\tau=v/x_0$ the time delay, hence the stronger the non-Markovian effects, the larger is the number of $A$ modes to account for. Notice that each frequency $\Omega_\nu$ has width $\gamma\sim 1/\tau$ comparable with $\Delta \Omega$, reflecting the open nature of block $A$.

From a merely computational viewpoint, it is well-established that such class of non-Markovian dynamics can be efficiently tackled via numerical methods <cit.> (moreover diagrammatic approaches were developed <cit.>).
Notwithstanding, mostly due the daunting complication of delayed (coherent) feedback, the underlying physics remains generally involved and non-trivial to interpret. In particular, to our knowledge, no simple picture to understand such dynamics was so far identified even in the weakly non-Markovian limit (the lowest non-trivial order in the characteristic time delay of the problem) where one could expect a relatively simple effective model capturing the essential physics to exist.

With the above motivations, this work introduces a physical picture for describing and interpreting non-Markovian waveguide-QED problems that is built upon a real-space decomposition of the waveguide into blocks. The essential idea is inspired by cavity QED, where non-negligible time delays occur if the atom couples significantly to many cavity modes (see <cit.>): the longer the delay (say the time taken by a photon to reach a cavity mirror), the more cavity modes need to be considered. Of course, there is no actual cavity in a waveguide. Yet, nothing prevents one from viewing the latter as a set of communicating blocks [see blocks $A$ and $B$ in <ref>(a)] and treating the finite-size block directly coupled to the atom [block $A$ in <ref>(a)] as an open cavity that leaks into the rest of the waveguide (this in turn described as a white-noise bath). If the size of the fictitious cavity [block $A$ in <ref>(a)] is chosen to be comparable with the characteristic length of the problem, then the frequency spacing between block-$A$ normal modes scales as the inverse of the characteristic time delay $\tau$.
Despite the relative arbitrariness of the block-$A$ boundaries, it turns out that this picture allows to define an effective Hamiltonian that fully captures both (driven) emission and photon scattering, including many-excitation dynamics where atomic non-linearities are important.

System and basic parameters.—To present our theory, we will consider the case study where a two-level quantum emitter (henceforth called “atom") is coupled to a semi-infinite waveguide [see <ref>(a)]. In spite of its apparent simplicity, this system hosts rich physics and is complex enough to show most salient effects of waveguide QED in the regime of long time delays (see <cit.>), which includes (via a suitable mapping) emission of a giant atom <cit.> and even some paradigmatic sub- and super-radiance phenomena <cit.>. The left end of the waveguide [see <ref>(a)] works as a perfect mirror placed at distance $x_0$ from the atom. The waveguide sustains a 1D field with linear dispersion $\omega=v k$ (with $v$ the photon group velocity and $k$ the wave vector). The atom's ground and excited states $|g\rangle$ and $|e\rangle$, respectively, are separated in frequency by $\omega_0=vk_0$; hence $k_0$ is the wave vector (modulus) of a photon resonant with the atom. Under weak coupling (enabling the rotating-wave approximation), the Hamiltonian reads ($\hbar=1$) <cit.>

\begin{align}\label{Hmodel}
H&=\omega_0\, \hat\sigma_+\hat\sigma_- \!-iv\!\!\int_{0}^{\infty}\!\! \!{\rm d} x\!\left[ \hat a^{\dagger}_R(x)\partial_x \hat a_R(x)\!-\! \hat a^{\dagger}_L(x)\partial_x \hat a_L(x)\right] \nonumber\\
&\quad+ g\! \int_{0}^{\infty} \!\!\!\!{\rm d} x \left[ \hat\sigma_+\left( \hat a_{ L}(x)+ \hat a_{ R}(x)\right)\!+\!{\rm H.c.}\right]\delta(x{-}x_0)\,,\!\!\!
\end{align}

with $\partial_x=\frac{d}{d x}$, $ \hat\sigma_-=\hat\sigma_+^\dag=|g\rangle\langle e|$,
$ \hat a_{R(L)}(x)$ the bosonic field operator annihilating a right-going (left-going) photon at position $x$ and $g$ the atom-photon coupling strength.
This model is not analytically solvable in general, except for single-excitation dynamics (such as spontaneous emission or single-photon scattering <cit.>).
The essential physical parameters are:

\begin{align}
\Gamma=\frac{2g^2}{v}\,,\,\,\,\,\tau= \frac{2 x_0}{v}\,,\,\,\,\,\phi=2k_0x_0\,\,.\label{param0}
\end{align}

Here, $\Gamma$ is the standard decay rate that the atom would have without the feedback effect of the mirror (as if the waveguide were infinite, instead of ending at $x=0$). Importantly,
$\tau$ (time delay) is the time taken by a photon resonant with the atom to travel (twice) the atom-mirror distance and  $\phi$ the corresponding accumulated phase. The strength of non-Markovian effects is measured by $\Gamma \tau$, quantifying how long is the time delay $\tau$ compared to the lifetime $1/\Gamma$. Thus, the Markovian regime occurs for $\Gamma \tau\ll1$, while for growing $\Gamma \tau$ non-Markovian effects get stronger and stronger.

Effective model.—We view the waveguide as two joint “blocks" [see <ref>(a)]: block $A$ (for $x\in [0,L[$ with $L> x_0$) and block $B$ (for $x\ge L$). Importantly, $A$ is the block directly coupled to the atom and has finite length. In contrast, block $B$ (uncoupled from the atom) has infinite length and can be seen itself as a semi-infinite waveguide, yet having its left edge at $x=L$ [see <ref>(a)]. Based on such block-decomposition of the waveguide, by calling $\hat\alpha_\nu$ and $\hat\beta_\omega$ respectively the normal-mode (bosonic) ladder operators of blocks $A$ and $B$, one can replace (<ref>) with the effective Hamiltonian <cit.>

\begin{align}
H_{\rm eff}=&\,\,\omega_0\, \hat\sigma_+\hat\sigma_-{+}\sum_\nu \Omega_{\nu}\, {\hat\alpha}^{\dagger}_{\nu}{\hat\alpha}_{\nu}{+}\int\!{\rm d}\omega\,\omega\, \hat\beta^\dag_\omega \hat\beta_\omega+ \sqrt{\frac{\gamma}{2\pi}}\,\sum_\nu\int\! {\rm d} \omega \,(  \hat\alpha_\nu^\dag  \hat\beta_\omega+{\rm H.c.})+\sum_\nu g_\nu \,( \hat\alpha_\nu^\dag \hat\sigma_-+{\rm H.c.})\label{Heff}\,,
\end{align}

\begin{align}
\Omega_\nu&=\omega_0+v \frac{\nu\pi}{L}\,,\,\,\gamma=\frac{2v }{L}\,,\label{param1}\\
g_\nu&=g\,(-1)^{\nu} \sqrt{\frac{2}{L}}\sin{\left(\frac{\nu\pi}{L} x_0+\frac{\phi}{2}\right)}\label{param2}\,
\end{align}

with $\nu$ running over all integers, while $\omega$ takes values throughout the real axis.

The second (third) term of Hamiltonian (<ref>) describes the free Hamiltonian of block $A$ (block $B$), where in particular $\Omega_\nu$ [(<ref>)] are the normal frequencies of block $A$ [see <ref>(b)]. The fourth term of (<ref>) couples block $A$ and block $B$ with a characteristic rate $\gamma$ given in (<ref>) [notice that $1/\gamma$ is the characteristic time taken by a photon to leak out of block $A$]. Finally, the last term describes the coupling between the atom and each block-$A$ normal mode with corresponding coupling strength $g_\nu$ [see <ref>(b)]. The expression (<ref>) of $g_\nu$ reflects the sinusoidal spatial shape of the $A$'s normal modes (just like a standard cavity-QED system). As a hallmark of the present framework, $L$ (length of block $A$) is a free parameter of the model except for two conditions: (i) it must be strictly greater than $x_0$ to ensure that block $A$ contains the atom, but in practice is required to be still comparable with $x_0$ (more on this later on); (ii) $L$ must be a multiple integer of $\lambda_0/2$ with $\lambda_0=2\pi/k_0$ the atomic wavelength. Condition $(ii)$ makes sure that there is a block-$A$ mode resonant with the atom: this mode (henceforth called “resonant mode") is labeled by $\nu=0$ [indeed (<ref>) yields $\Omega_{\nu=0}=\omega_0$].

Hamiltonian $H_{\rm eff}$ formally describes an effective cavity-QED system in that the atom is coupled to a multimode lossy cavity (block $A$) leaking into a white-noise photonic bath (embodied by block $B$). Notice that block $A$ is an intrinsically low-finesse cavity: indeed [(<ref>)]] $\Delta \Omega=\Omega_{\nu+1}-\Omega_{\nu}\sim \gamma$, the frequency spacing $\Delta\Omega$ between $A$ modes [<ref>(b)] is comparable with the loss rate $\gamma$ of block $A$. Physically, this stems from the inherently open nature of the fictitious cavity $A$ which fully lacks the mirror at $x=L$.

To end up with (<ref>), we resort to the standard discretization of a waveguide combined with weak coupling (the latter one allowing to linearize the dispersion law) <cit.>. Discretizing the system this way enables a clean definition of the two blocks, leading to a natural identification of their associated Hamiltonian and normal modes, whose continuous limit is eventually worked out.

Dependence on time delay.—
The essence of the present picture is that the dynamics of the joint system can be effectively described by replacing Hamiltonian (<ref>) with (<ref>), where the latter can be approximated by retaining only a finite number of $A$ modes which yet eventually grows with the time delay $\tau$. To see the last key property, recall that we require $L$ to be of the order of $x_0$ ($L=2 x_0$). Thus the spacing of $A$ modes $\Delta\Omega=\pi v/L$ is of the order of $ \tau^{-1}$ [see (<ref>), (<ref>) and (<ref>)(b)], the detuning from the atom of the $A$ modes with $|\nu|\ge 1$, henceforth called “off-resonant modes", scales as the inverse of time delay (recall that block $A$ is defined so as to ensure $\omega_{0}=\Omega_0$).
Accordingly, in the Markovian regime of vanishing $\tau$ all these off-resonant modes $|\nu|\ge 1$ are far-detuned from the atom so that only the resonant mode $\nu=0$ needs to be accounted for [<ref>(b)].

As the time delay grows up so that non-Markovian effects get increasingly important, more and more off-resonant modes must be retained in general. In other words, for given $\tau$, one neeeds to retain all the modes  $\nu=0,\,\pm 1,...,\pm N_A$ with $N_A$ eventually growing with $\tau$. Thus the multimode nature of block $A$ (instead of a more canonical one-mode cavity) reflects occurrence of retardation effects: while this is a well-known fact in standard cavity QED <cit.>, the present framework provides ground to take advantage of this property also in the study of waveguide-QED systems.

Testing the framework.—To check the effectiveness of the waveguide decomposition into blocks, in <ref> we set $L=2 x_0$, $\Gamma \tau\ge2$ (relatively long delay) and the representative phase $\phi=\pi/2$ [except for panel (b) where $\phi=2\pi$] for two paradigmatic dynamics: spontaneous emission [panels (a)-(b)] and scattering of a coherent-state wave packet [(c)-(d)].

In each case, the dynamics predicted by the effective model in (<ref>) by retaining only the $A$ modes $\nu=0, \pm1,...,\pm N_A$ is compared with the exact solution of (<ref>) obtained through either analytical methods when available [as in the single-excitation dynamics of panels (a)-(b)] or numerical simulations based on Matrix Product States <cit.> [panels (c)-(d)].

(a) and (b): Atomic excitation, $\langle e| \rho|e\rangle=\rho_{ee}$ with $\rho$ the atom's density matrix, when the atom is initially in state $\ket{e}$ and the field in the vacuum state: exact analytical solution <cit.> (black line) versus the approximated using (<ref>) by retaining only modes $\nu=0, \pm1,...,\pm N_A$ for $\phi=\pi/2$ (a) and $\phi=2\pi$ (b). We set $\Gamma \tau=2$. (c)-(d) Field's output intensity  (c) and two-photon correlation function  (d)  after the scattering of a left-incoming Gaussian coherent wave packet (see <cit.>). We set the wave packet's width to $W=2.5\Gamma$ and the average number of photons to $n_{\rm ph}=0.5$ with phase $\phi=\pi/2$ and delay $\Gamma \tau=4$. In (a)-(d), we set $L=2 x_0$.

As $N_A$ is made larger, the mismatch between exact and approximated dynamics gets smaller and smaller until becoming negligible, including the special value of phase $\phi=2\pi$ in <ref>(b) where it is known that the atom does not fully decay <cit.>. Notice that, unlike <ref>(a) and (b), the dynamics in panels (c)-(d)  involves  many excitations, providing evidence that the picture is effective even when the atom's intrinsic nonlinearity has substantial effects. This is especially striking in panel <ref>(d), reporting the two-photon correlation function of scattered light, which clearly shows a multi-photon peak [this adds to the delayed single-photon peak of the output intensity in panel (c)]. Analogous conclusions hold for different settings of the parameters.

We note that the convergence rate is generally affected by $\phi$ and $L$, a major reason being the sinusoidal dipendence on these parameters of the coupling strength $g_\nu$ [(<ref>)]. For analogous reasons, in general the convergence is not strictly monotonic, it can happen that $N_A+1$ modes perform as $N_A$ or even worse <cit.>. Notwithstanding, convergence eventually occurs because $|g_\nu|\le g \sqrt{2/L}$ for any $\nu$ while the detuning $|\Omega_\nu-\omega_0|$ grows linearly with $\nu$ [(<ref>)-(<ref>)]. Although somewhat implicit in the above, it is worth stressing here that the advantage of the framework is not to provide a fast computational numerical tool (where efficient techniques exist already) but rather an intuitive physical picture connecting non-Markovian waveguide QED with cavity QED. This is illustrated next with some important instances.

Markovian limit and Purcell effect.—
For $\Gamma\tau\ll 1$ we are in the Markovian regime: time delay is negligible, but the feedback provided by the mirror affects atomic emission resulting in a decay rate modulated by $\phi$ as $\Gamma'(\phi)=2\Gamma \sin^2\frac{\phi}{2}$ <cit.>, which was experimentally confirmed  <cit.>.
Thus emission can be either enhanced or suppressed; in particular $\Gamma'=2\Gamma$ for $\phi=(2m+1) \pi$ while $\Gamma'=0$ for $\phi=2m \pi$ ($m$ is an integer). Now, using our framework, we see that, due to $\Delta\Omega\sim 1/\tau$ [(<ref>)], for negligible $\tau$ the off-resonant modes of block $A$ are very far-detuned from the atom and thus can be neglected. Only the resonant mode $\nu=0$ thus needs to be accounted for. Its bandwidth is yet very large since we also have $\gamma\sim \Delta\Omega \sim 1/\tau$ [(<ref>)]. The system therefore reduces to an atom coupled to a standard one-mode cavity, but in the bad cavity limit. The corresponding atom-mode coupling strength is  $g_0=g\sqrt{\frac{2}{L}}\sin{\frac{\phi}{2}}$ [(<ref>)], which indeed vanishes for $\phi=2m\pi$ meaning that in this case the atom sits right on a cavity field's node and is thus unable to emit. Standard cavity-QED theory (bad-cavity limit) then predicts the decay rate $4g_0^2/\gamma$, which indeed exactly matches $\Gamma'(\phi)$. This shows that the action of the mirror, despite no actual cavity is present, can still be seen as a manifestation of the standard Purcell effect in a cavity.

Non-Markovian steady states. Each point of the plane has coordinates $(|\rho_{eg}|, \rho_{ee})$, with $\rho_{eg}$ ($\rho_{ee}$) the coherences (excited-state population) of a possible steady state $\hat\rho$ of a driven atom emitting into the semi-infinite waveguide. Emission into a Markovian bath [see (<ref>)] can only yield steady states within the shaded area <cit.>. The blue dots are exact steady states computed through MPS simulations for Rabi frequencies (in units of $\Gamma$) $\Omega_D=0.5, 1, 1.5,  2,  2.5,  3,  3.5,  4$ from bottom to top (as shown by the curved arrow) and for $\Gamma \tau=0.25$ and $\phi=\pi$. The other points (squares, rhombi and triangles) are approximated steady states calculated through the effective Hamiltonian (<ref>) by retaining $N_A$ modes of block $A$ for $L=x_0$.

Non-Markovian steady states.—Recently, Ask and Johansson showed that a driven atom in front of a mirror for non-negligible $\Gamma\tau$ can reach steady states unattainable in a standard Markovian bath <cit.>. Specifically, let

\begin{equation}
\dot{\hat\rho} = -i\left[\tfrac{1 }{2}\Omega_D\hat\s_x , \hat\rho\right] \kappa \,\mathcal{D}[\hat\sigma_-]\hat\rho  + \kappa_\phi \,\mathcal{D}[\hat{\sigma}_+ \hat{\sigma}_-] \label{ME}
\end{equation}

with $\mathcal{D}[A]\hat\rho=A\hat\rho A^\dag-\tfrac{1}{2}\{A^\dag A,\hat\rho\}$, be the standard Markovian master equation (ME) of an atom subject to a classical drive of Rabi frequency $\Omega_D$, pure dephasing with rate $\kappa_\phi$ and decaying into a Markovian bath with rate $\kappa$. An atom's state $\rho$ is fully specified by the (excited-state) population $\rho_{ee}$ and coherence $\rho_{eg}$. After a transient, the  emitter reaches a steady state obtained by imposing $\dot{\rho} =0$ in (<ref>). It turns out that, irrespective of $\Omega$ and $\kappa$, the steady state must lie within the shaded region bounded by the elliptical line shown in <ref>. This bound in particular holds for a driven atom in a semi-infinite waveguide in the Markovian regime $\Gamma\tau\ll1$ (see previous section). Yet, it can be violated in the non-Markovian regime: see the blue dots in <ref> computed via exact MPS numerical simulations <cit.> corresponding to the steady states occurring for $\Gamma \tau=0.25$, $\phi=\pi$ and (from bottom to top, see curved arrow) growing values of $\Omega_D$.

The steady states in <ref> can be well-approximated using our block-decomposition framework, as shown by <ref> where the agreement with the exact solution (blue dots) grows with  $N_A$. Remarkably, retaining even a single cavity mode ($N_A=0$; orange squares) provides an excellent quantitative approximation at low values of $\Omega_L/\Gamma$ and, as long as population $\rho_{ee}$ is concerned, even at larger ones; in any case, on a qualitative ground, it appears to capture most of the relevant physics. In this case, the emitter's steady state can be obtained through partial trace from the bipartite Markovian ME governing the evolution of $\hat\varrho$ (joint state of the atom and mode $\nu=0$)

\begin{align}
\dot{\hat\varrho} = -i[\tfrac{1}{2}\Omega_D\hat\s_x + g_0 (\hat\alpha^\dag_0 \hat\s_- + {\rm H.c.}), \hat\varrho]
+\gamma \mathcal{D}[\hat\alpha_0]\hat\varrho \,,
\end{align}

with $\gamma$ and $g_0$ given by (<ref>) and (<ref>). We point out that occurrence of population inversion at long times, $\rho_{ee}>1/2$, is a sufficient condition to reach non-Markovian steady states beyond the elliptical bound in <ref>. Occurrence of population inversion for a driven two-level system coupled to a lossy cavity mode is a well-established quantum optics effect  <cit.>, which highlights a further interesting connection between non-Markovian waveguide QED and cavity QED.

Conclusions.— We presented a picture to describe and interpret waveguide-QED dynamics with delayed feedback. After decomposing the waveguide into blocks, the block coupled to the atom is viewed as an open cavity leaking into the rest of the waveguide. The longer the time delay, the more modes of such open cavity generally need to be accounted for. The picture captures both the atom and field dynamics, even when many excitations are present.

The idea that atoms in waveguides could be understood in terms of effective cavities appeared several times (see <cit.>), sometimes relying on the well-known mirror-like behaviour of an atom <cit.>. In contrast, the cavity central to our framework relies on a fully transparent fictitious mirror. Yet, it is instrumental to the establishment of a sharp link between emission phenomena featuring delayed feedback and cavity-QED physics with the bonus that the picture can capture the waveguide-field dynamics as well.
Our theory can be seen to define a (so called) Markovian embedding (or dilation) in the following sense: one replaces an open system (the atom in our case) immersed in a non-Markovian bath with an enlarged open system (atom plus $A$ modes here) that is instead immersed in a Markovian bath. This is arguably the most common strategy to attack non-Markovian problems and is typically accomplished by adding auxiliary lossy modes to the open system <cit.>. Two remarkable features yet stand out in the present approach: (i) the open-cavity modes have a clear physical meaning and can be straightforwardly visualized as degrees of freedom taken out of the bath (the waveguide); (ii) besides the open system, the framework can describe as well the bath dynamics.
We anticipate that this work could offer a novel alternative approach to understanding and interpreting waveguide-QED phenomena within the relatively unexplored non-Markovian regime.

D.C. acknowledges support from the BMBF project PhoQuant (grant no. 13N16110) and the state of Baden-Württemberg through bwHPC and the German Research Foundation (DFG) through grant no. INST
40/575-1 FUGG (JUSTUS 2 cluster). G.M.P. acknowledges support from MUR under PRIN Project no. 2022FEXLYB
“Quantum Reservoir Computing" (QuReCo).
G.C. acknowledges that results incorporated in this standard have received funding from the  T-NiSQ consortium agreement financed by QUANTERA 2021. Numerical simulations were performed using QuTip <cit.> and mpnum <cit.>. F.C. acknowledges support from  European Union – Next Generation EU through Project Eurostart 2022 “Topological atom-photon interactions for quantum technologies" (MUR D.M. 737/2021) and through Project PRIN 2022-PNRR no. P202253RLY “Harnessing topological phases for quantum technologies".

Supplemental Material for
“Multimode-cavity picture of non-Markovian waveguide QED”

D.  Cilluffo,$^{1}$, L.  Ferialdi,$^{2}$, G. M. Palma,$^{2,3}$, G. Calajò,$^{4}$ and F. Ciccarello$^{2,3}$

${}^1$Institut für Theoretische Physik and IQST, Albert-Einstein-Allee 11, Universität Ulm, 89069 Ulm, Germany

${}^2$Universit$\grave{a}$  degli Studi di Palermo, Dipartimento di Fisica e Chimica – Emilio Segr$\grave{e}$, via Archirafi 36, I-90123 Palermo, Italy

${}^3$NEST, Istituto Nanoscienze-CNR, Piazza S. Silvestro 12, 56127 Pisa, Italy

${}^4$ Istituto Nazionale di Fisica Nucleare (INFN), Sezione di Padova, I-35131 Padova, Italy

(Dated: August 27, 2024)

§ DERIVATION OF THE EFFECTIVE HAMILTONIAN

Here we derive the effective Hamiltonian (<ref>) by first discretizing the waveguide, then decomposing it into blocks and finally taking the continuous limit. We start with a review of the discretization of a continuos waveguide.

§.§ Discretized waveguide: review

Consider a homogeneous coupled-cavity array with cavities located at positions $n$, whose free Hamiltonian has the usual tight-binding form

\begin{equation}
\hat 	H_F=\omega_c\sum_n 	\hat  a_n^\dag  	\hat a_n-J \sum_n  (	\hat  a_{n+1}^\dag  	\hat  a_n+{\rm H.c.})\,\label{HB-wav}\, ,
\end{equation}

with $\omega_c$ the frequency of each cavity and $J$ the cavity-cavity hopping rate.

Under periodic boundary conditions ($	\hat  a_1\equiv 	\hat  a_{N+1}$) the waveguide's normal modes  are the usual plane waves with dispersion law

\begin{equation}
\omega_k=\omega_{c}-2J\cos k\label{wk2}\, ,
\end{equation}

which gives rise to an energy band of width $4J$, where $k_m=2\pi m/ (N+1)$ with $m=1,...,N$. The corresponding group velocity is $v_k=d \omega/dk=2J  \sin k $, hence (far from the band edges)  the dispersion law can be linearized around a specific wavevector $k_0$ as $\omega_k\simeq (\omega_{c}-2J\cos k_0)+v (k-k_0) $ with $v=v_{k_0}$.

In the presence of an atom coupled to the waveguide at position $n_0$, the total Hamiltonian reads

\begin{equation}
\hat 	H=\omega_{0}\hat \sigma_{+}\hat \sigma_-+\hat H_F+\hat V_{AF} \,.\label{Htot}
\end{equation}

with the atom-field interaction Hamiltonian given by

\begin{equation}
\hat V_{AF}=g\, (\hat \sigma_{+} \hat a_{n_0}+{\rm H.c.})\label{vaf}\,\,.
\end{equation}

Note that the coupling-strength $g$ here is defined differently from the one in Eq. (<ref>) of the main text, and indeed they have different dimensions.
If the atom is tuned inside the photonic band ($\omega_{c}-2 J<\omega_{0}< \omega_{c}+2 J$), only photonic modes of frequency $\omega_k$ close to $\omega_0$ need to be considered, hence one can linearize the dispersion law around $\pm k_0$ with $k_0$ defined by $\omega_0=\omega_{k_0}$.

It can be shown that the real-space representation of the linearized field's Hamiltonian so obtained is analogous to the second term of (<ref>).

If we now consider open boundary conditions, such that $	\hat  a_0=	\hat  a_{N+1}=0$, the spectrum (<ref>) is unaffected but the normal modes now have sinusoidal shape

\begin{equation}\label{alpha1}
\hat{\alpha}_m=\sqrt{\frac{2}{N+1}}\sum_{n =1}^{N}\,\!\sin{\left(k_m n \right)} \,\hat{a}_n \,\,\,\,{\rm with}\,\,\,\,k_m=\frac{\pi m}{N+1}\,\,\,{\rm for} \,\,m=1,...,N\,.
\end{equation}

In the limit $N\rightarrow \infty$ we obtain a semi-infinite waveguide (or equivalently a waveguide subject to a hard-wall boundary condition at $n=0$).

§.§ Decomposition into blocks

Decomposition of the discretized waveguide into blocks $A$ and $B$ with $b_1\equiv a_{N_A+1}, ...,b_{N_B}\equiv a_{N}$. Note that $N_A\ge n_0$, meaning that block $A$ is directly coupled to the atom.

We start from Hamiltonian (<ref>) under open boundary conditions and split the discretized waveguide in two blocks (see <ref>): block $A$ for $1\le n\leq N_A$ with $N_A\ge n_0$ and block B for $N_A+1\leq n\leq N$ (so that block $A$ is the one of the two directly coupled to the atom). Corresponding to this block decomposition, we rearrange the free field's Hamiltonian as

\begin{align}
\hat H_F=\underbrace{\omega_c \sum_{n=1}^{N_A} \hat a^\dag_n\hat a_n-J \sum_{n =1}^{N_A-1} \left(\hat a_n  \hat a_{n +1}^\dag+{\rm H.c.}\right)}_{\hat H_A}+\underbrace{\omega_c \sum_{n=1}^{N_B} \hat b^\dag_n\hat b_n-J \sum_{n =1}^{N_B-1} \left(\hat b_n  \hat b_{n +1}^\dag+{\rm H.c.}\right)}_{\hat H_B}\,\,\underbrace{-J \left(\hat a_{N_A} \hat b_{1}^\dag+{\rm H.c.}\right)}_{\hat V_{AB}}\label{Hdec}
\end{align}

where $N_B=N-N_A$ is the size of block $B$. Notice that for convenience we have re-named the bosonic site modes of block $B$ as $\hat b_n{=}\hat a_{n-N_A}$. In (<ref>), $\hat H_A$ ($\hat H_B$) is the free Hamiltonian of block $A$ (block $B$). It is important to stress that the two blocks are mutually coupled with interaction Hamiltonian $\hat V_{AB}$ and coupling strength just equal to the photon hopping rate $J$: their coupling guarantees a photon coming from block $A$ to reach $B$ (or the other way around) without suffering any back-reflection.

The sine-shaped normal modes of block $A$ read [compare with (<ref>)]

\begin{equation}
\hat{\alpha}_m =\sqrt{\frac{2}{N_A+1}}\sum_{n =1}^{N_A}\,\!\sin{\left(k^A_m n \right)} \,\hat{a}_n \,\,\,\,\,\,{\rm with}\,\,\,\,k^A_m=\frac{\pi m}{N_A+1}\,\,\,{\rm for} \,\,m=1,...,N_A\,.
\end{equation}

In terms of these, the free Hamiltonian of block $A$ takes the diagonal form

\begin{equation} \label{H0f-diagonal}
\hat{H}_A\!=\!\sum_{m=1}^{N_A} \omega^A_{m}\, \hat{\alpha}^{\dagger}_{m}\hat{\alpha}_{m}\,\,,
\end{equation}

where $\omega^A_m=\omega_c-2J \cos{k^A_m}$.
Real-space block-A operators can be expressed in terms of these normal modes as

\begin{equation}\label{nba}
\hat a_n \!=\sqrt{\frac{2}{N_A+1}}\,\sum_{m=1}^{N_A}\,\!\sin{(k^A_m n )} \,\hat{\alpha}_m\,.
\end{equation}

This allows us to express even the atom-field interaction Hamiltonian (<ref>) in terms of the $A$'s normal modes as

\begin{align}
\hat V_{AF}&= \sum_{m=1}^{N_A} g_m\left( \hat\sigma_+\hat \alpha_m+{\rm H.c.}\right)\,,
\end{align}

where $g_m$ measures how strongly is the atom coupled to the $m$-th normal mode of block $A$

\begin{align}\label{gm}
g_m=g\, \sqrt{\frac{2}{N_A+1}}\sin{\left(k^A_m{n_0} \right)}\,.
\end{align}

Likewise, using the same decomposition (<ref>) for $n=N_A$, the coupling Hamiltonian between blocks A and B [(<ref>)] is arranged as

\begin{equation}
\hat V_{AB}=-J\sqrt{\frac{2}{N_A+1}}\sum_{m=1}^{N_A}\sin{\left(k^A_m N_A\right)} \,\hat{\alpha}_{m}\,\hat b_1^\dag=\sum_{m=1}^{N_A} \xi_m\,\left(\hat \alpha_m \hat b_1^\dag+{\rm H.c.}\right)\,,\label{vab1}
\end{equation}

where we defined

\begin{align}
\xi_m=J\sqrt{\frac{2}{N_A+1}}\,(-1)^{m}\sin\left(\frac{m\pi }{N_A+1}\right)\,,\label{xim1}
\end{align}

and used the identity

\begin{align}
\sin{\left(k^A_m \,N_A\right)}= -(-1)^m\sin\left(k^A_m\right)\,.
\end{align}

Analogously to block $A$, we can define normal modes also for block $B$ as [(<ref>)]

\begin{equation}
\hat{\beta}_m\!=\sqrt{\frac{2}{N_B+1}}\,\,\sum_{n =1}^{N_B}\,\!\sin{\left(k^B_m n \right)} \,\hat{b}_n\,\,\,\,{\rm with}\,\,\,\,k^B_m=\tfrac{\pi m}{N_B+1}\,\,\,{\rm for} \,\,m=1,...,N_B
\end{equation}

such that [(<ref>)] $\hat{H}_B\!=\!\sum_{m=1}^{N_A} \omega^B_{m}\, \hat{\beta}^{\dagger}_{m}\hat{\beta}_{m}$ with $\omega^B_m=\omega_c-2J \cos{k^B_m}$. The inverse transformation reads

\begin{equation}
\hat b_n \!=\sqrt{\frac{2}{N_B+1}}\,\sum_{m=1}^{N_B}\,\!\sin{\left(k^B_m n \right)} \,\hat{\beta}_m\,,
\end{equation}

and allows now to arrange the $A$-$B$ interaction Hamiltonian (<ref>) as

\begin{equation}
\hat V_{AB}=\sum_{m=1}^{N_A}\sum_{m'=1}^{N_B} \xi_m\chi_{m'}\,\left(\hat \alpha_m \hat \beta_{m'}^\dag+{\rm H.c.}\right)\,.
\end{equation}

\begin{align}
\chi_{m}=\sqrt{\frac{2}{N_B+1}}\sin{\left( \frac{m\pi}{N_B+1}  \right)}
\end{align}

To summarize, in terms of block-$A$ and block-$B$ normal modes the total Hamiltonian reads

\begin{align}
\hat H= \omega_{0}\hat \sigma_{+}\hat \sigma_-+\underbrace{\sum_{m=1}^{N_A} \omega^A_m\, \hat{\alpha}^{\dagger}_{m}\hat{\alpha}_{m}+\sum_{m=1}^{N_B} \omega^B_m\, \hat{\beta}^{\dagger}_{m}\hat{\beta}_{m}+\sum_{m=1}^{N_A}\sum_{m'=1}^{N_B} \xi_m\chi_{m'}\,\left(\hat \alpha_m \hat \beta_{m'}^\dag+{\rm H.c.}\right)}_{=\,\hat H_F}+\underbrace{\sum_{m=1}^{N_A} g_m\left( \hat\sigma_+\hat \alpha_m+{\rm H.c.}\right)}_{=\,\hat V_{AF}}\label{Ht}
\end{align}

where we recall that

\begin{align}\label{kkp}
\omega^A_m=\omega_c-2J \cos{k^A_m}\,,\quad\omega^B_m=\omega_c-2J \cos{k^B_m}\,.
\end{align}

§.§ Continuous limit and linearization

Since we are considering a semi-infinite waveguide, the length of block B must diverge as $N_B\rightarrow \infty$. Accordingly, block-$B$ normal ladder operators $\{{\hat\beta}_m\}$ become a continuum of singular bosonic modes $\{\hat\beta(k)\}$ with $0\le k <\pi$ fulfilling $[\hat\beta(k),\hat\beta^\dag(k')]=\delta(k-k')$. Specifically, $\hat \beta(k)$ is obtained as the continuous limit of the rescaled ladder operators $\hat\beta_m/\sqrt{\Delta k}$ with $\Delta k=2\pi/N_B$. This way, in Hamiltonian (<ref>) we can make the replacements

\begin{align}
\sum_{m=1}^{N_B} \omega^B_{m}\, \hat{\beta}^{\dagger}_{m}\hat{\beta}_{m}\rightarrow\int_0^\pi {\rm d}k\,\omega^B(k)\hat \beta^\dag(k)\hat \beta(k)\,,\,\,\,	\sum_{m'=1}^{N_B}\,\chi_{m'}\hat\beta_{m'} \rightarrow\,\int_0^\pi{\rm d} k\,\sqrt{\frac{2}{\pi}}\sin k \,\,\hat\beta(k)\,\,.\label{sinkk}
\end{align}

Here, in each sum we multiplied and divided the summand by $\Delta k$, expressed it in terms of $\hat\beta_m/\sqrt{\Delta k}$ and finally carried out the continuous limit thus turning sums into integrals over the first Brillouin zone.

We assume now that the atom is tuned on resonance with a specific normal mode of block $A$ whose wavevector labeled by $m=m_0$,
\begin{equation}
\omega_0=\omega_c-2J \cos{k^A_{m_0}}\,. %=\omega_{k_0}
\end{equation}

Also, to ensure the weak-coupling regime, we assume that $g\ll J$ and that $\omega_0\equiv \omega_{m_0}$ is sufficiently far from the band edges $\omega_{c}\pm 2J$ (where singularities occur). Accordingly, we can effectively approximate the dispersion law of block $A$ to the first order around $m=m_0$ as

\begin{align}
\omega^A_m&\simeq\omega_0+v\left(k_{m}-k_{m_0}\right)=\omega_0 +v \frac{\pi}{N_A+1}\nu\,\,,\label{w0k0}
\end{align}

where we used (<ref>) and replaced the effective group velocity $v= 2 J \sin k_{m_0}$. In the last identity, we introduced the integer number $\nu=m-m_0$ (taking both negative and positive values). For analogous reasons, since the coupling strength between the generic block-$A$ mode and block B is much smaller than the waveguide bandwidth in the thermodynamic limit, $\xi_m \ll J$ [(<ref>)], we linearize the dispersion relation of block B as $\omega^B(k)\simeq\omega_0+vk$, where $k$ now runs between $-\infty$ and $+\infty$ and approximate $\xi_m$ [(<ref>)] to the lowest order around $k_{m_0}$ obtaining

\begin{align}
\xi_m \simeq J\sqrt{\frac{2}{N_A+1}}\,(-1)^{m}\sin{\left(k_{m_0}\right)}=\frac{v}{2}\sqrt{\frac{2}{N_A+1}}\,(-1)^{\nu}\,.
\end{align}

§.§ Final continuous Hamiltonian

Putting everything together we get the total Hamiltonian

\begin{align}
\hat H_F&=\sum_{\nu=-\infty}^\infty\Omega_{\nu}\, \hat{\alpha}^{\dagger}_{\nu}\hat{\alpha}_{\nu}+\int_0^\pi {\rm d}k\,\omega(k)\hat \beta^\dag(k)\hat \beta(k)+\sum_{\nu=-\infty}^\infty \frac{v}{2}\sqrt{\frac{2}{N_A+1}}\,(-1)^{\nu}\,\left(\hat \alpha_\nu \int_0^\pi{\rm d} k\,\sqrt{\frac{2}{\pi}}\,\,\hat\beta^\dag(k)+{\rm H.c.}\right)\nonumber\\
\hat V_{AF}&=\sum_{\nu=-\infty}^\infty g_\nu  \hat\sigma_{+}\hat\alpha_\nu+{\rm H.c.}\,,
\end{align}

where $g_\nu $ is just (<ref>) expressed in terms of $\nu=m-m_0$ and $\phi=2 k_{m_0} n_0$.

\begin{align}
g_\nu&=g\, \sqrt{\frac{2}{N_A+1}}\sin{\left(\nu\frac{\pi}{N_A+1} n_0+ \frac{\phi}{2}\right)}\,\,,
\end{align}

with $\Omega_\nu=\omega_0+ vk_\nu$ and $\omega(k)=\omega_0+ vk$.
Each mode $\nu$ of the cavity is coupled to block-B with strength $\sim v/\sqrt{N_A}$, which for $N_A$ large enough will be far smaller than $v$. This and the fact that the coupling of each mode to block-B modes is flat (i.e., frequency-independent) yields that the integrals over $k$ can be extended to the entire real axis. By passing in addition to the frequency domain (using $\omega=v k$ and $\hat \beta(k)=\hat \beta (\omega)/\sqrt{v}$), we get the free-field Hamiltonian in the form

\begin{align}
\hat H_F&=\sum_{\nu} \Omega_{\nu}\, \hat{\alpha}^{\dagger}_{\nu}\hat{\alpha}_{\nu}+\int_{-\infty}^{\infty} {\rm d}\omega\,\omega\hat \beta^\dag(\omega)\hat \beta(\omega)+\sum_{\nu} \sqrt{\frac{v}{N_A+1}}\,\left(\hat \alpha_\nu \int_{-\infty}^\infty{\rm d} \omega\,\sqrt{\frac{1}{\pi}} \,\,\hat\beta^\dag(\omega)+{\rm H.c.}\right)\,,
\end{align}

where we also redefined the cavity modes as $\hat \alpha_{\nu}\rightarrow (-1)^\nu\hat \alpha_{\nu}$, in a way that factor $(-1)^{\nu}$ is now incorporated in the definition of the emitter-cavity.
The latest step is to take the continuous limit of the tight binding model, which essentially leads to the replacement $N_A+1\rightarrow L_A$, which leads to (<ref>)-(<ref>) in the main text.

§ SPONTANEOUS EMISSION IN THE MARKOVIAN REGIME

It is well-known that Hamiltonian (<ref>) implies that the atom's excitation amplitude $\epsilon(t)$ obeys the exact delay differential equation <cit.>

\begin{equation}
\dot \varepsilon(t)=-\frac{\Gamma}{2} \,\epsilon(t)+\frac{\Gamma}{2}e^{i\phi} \,\epsilon(t-\tau)\, \Theta (t-\tau)\,.\label{delay-de}
\end{equation}

For very short time delay (Markovian regime), we can replace $t-\tau$ $\simeq t$, so that (<ref>) reduces to

\begin{equation}
\dot \epsilon=-\frac{\Gamma}{2} \,\epsilon+\frac{\Gamma}{2}e^{i\phi} \,\epsilon=i \frac{\Gamma}{2}\sin \phi\, \epsilon- \frac{\Gamma}{2}\left(1-\cos \phi\right) \epsilon\,.\label{delay-de2}
\end{equation}

Hence, the excited-state population will decay as $|\epsilon|^2=e^{-\Gamma' t}$ with

\begin{equation}
\Gamma'=\Gamma \left(1-\cos \phi\right)=2 \Gamma \sin^2 \frac{\phi}{2}\label{Gap}\,.
\end{equation}

§.§ One block-$A$ mode

Consider the effective Hamiltonian (<ref>) and approximate it by retaining only the block-$A$ resonant mode $\nu=0$, which is justified in the limit of very short time delay. Also, we take $L=x_0$, hence [(<ref>)]

\begin{equation}
\end{equation}

Let the atom initially in state $\ket{e}$ with mode $\nu=0$ and all the modes of block $B$ initially in the vacuum state. Then the joint state of the atom and mode $\nu=0$ has the form

\begin{equation}
\ket{\Psi(t)}=\epsilon (t)\ket{e,0}+a_0 (t)\ket{g,1}
\end{equation}

with $\ket{e,0}$ ($\ket{g,1}$) the state where the atom is in the excited (ground) state while mode $\nu=0$ has zero (one) photon, where $\epsilon (t)$ and $a_0 (t)$ fulfill the differential system (we set in a rotating frame such that $\omega_0=0$)

\begin{align}
\dot\epsilon=-i g_0  \alpha_0\,\,\,\,\,	\dot\alpha_0=- \tfrac{\gamma}{2} \alpha_0-i g_0 \epsilon\nonumber
\end{align}

subject to the intial condition $\epsilon(0)=1$, $\alpha_0(0)=0$.

In the Laplace domain (variable $t$ replaced by $s$), the system reads

\begin{align}
s\tilde\epsilon-1=-i g_0  \tilde\alpha_0\,\,\,\,\, s\tilde\alpha_0=- \tfrac{\gamma}{2} \tilde\alpha_0-i g_0 \tilde\epsilon\label{epsa}\,.
\end{align}

In the limit of very short delay, we get $\gamma\gg g_0$ [(<ref>)-(<ref>)]; hence we can replace $\alpha(t)$ with its stationary value for given $\epsilon(t)$ which is equivalent to setting $s=0$ in the second identity of (<ref>). This yields $\tilde\alpha_0=-i \frac{2g_0}{\gamma}\tilde\epsilon$. Replacing in the equation for $\tilde{\epsilon}$ we get

\begin{equation}
\tilde\epsilon=\frac{1}{s+\Gamma'}\,,\label{eps2}
\end{equation}

where $\Gamma'=\frac{2 g_0^2}{\gamma}$ matches (<ref>).

§ TENSOR NETWORK SIMULATIONS

The exact results throughout the main text are obtained by leveraging the MPS formalism specifically adapted for photonic circuits featuring time delays, as illustrated in <cit.>.
We consider the reduced quantum state of the emitter $A$ and the non-Markovian bath $F$ (field in the region within the atom and the mirror) as initially uncorrelated, i.e. $\rho_{AF}(t=0) = \rho_A (0) \otimes \rho_F (0)$. For a single emitter coupled to a semi-infinite waveguide, a simplification of the problem's geometry arises by transforming the configuration into an equivalent one with a chiral infinite waveguide and the emitter coupled to the waveguide at two points separated by a distance of $2 x_0$ <cit.>.
In the interaction picture with respect to the bath and the emitter free Hamiltonians,  (<ref>) is can be recast as <cit.>:

\begin{equation}
\hat H (t) = \hat{H}_{\rm dr} + g  (\,\hat \sigma^{+}_A ( {\hat b}_{t} + e^{i 2 \omega_0 \tau} \, {\hat b}_{t-2\tau} )+ {\rm H.c.} )\,\label{Vt1}\,
\end{equation}

where we are including a classical driving Hamiltonian on the atom $\hat{H}_{\rm dr} = \Omega (\sigma^{+}_A + \sigma^{-}_A)$ and ${\hat b}_{t}$ are the time-domain ladder operator of the chiral bath.
The overall state $\rho_{AF}$ is evolved according to a stroboscopic map with discrete time steps $\Delta t$ chosen to be small compared to the relevant frequencies of the system, consistent with the approach used in quantum collision models <cit.>.
Hence we define the discrete-time propagator

\begin{align}
\label{s_eq:prop}
\hat{U}_n = \exp{-i \hat{H}_n \Delta t}\,,
\end{align}

\begin{align}
\hat{H}_n = \hat{H}_{\rm dr} + \frac{g}{\sqrt{\Delta t}}  (\,\hat \sigma^{+}_A ( {\hat b}_{n} + e^{i \phi} \, {\hat b}_{n-\ell} )+ {\rm H.c.} ) \,,
\end{align}

where we introduced the discrete bosonic noise operators $\hat{b}_n = \frac{1}{\sqrt{\Delta t}} \int_{t_{n-1}}^{t_{n}} \hat{b}(s) ds$ and $\ell=2 \tau/{\Delta t}$.
The discretization of the interaction reflects in the representation of the environment as a chain of $\ell$ quantum harmonic oscillators and the joint state of atom and environment as

\begin{align}
\rho_{AF} = \sum_{\overline{i},\overline{i}'} c_{\overline{i},\overline{i}'} |\overline{i} \rangle \langle \overline{i}' |
\end{align}

with the basis state $|\overline{i}\rangle = | i_A, i_{1}, \ldots, i_{{\ell}} \rangle$, where the numbers identify the $\ell$ oscillators.
This representation can be reformulated using singular value decomposition between each possible bipartition of the chain, yielding a Matrix Product Operator (MPO)

\begin{align}
\rho_{AF} = \sum_{\overline{i},\overline{i}'} \sum_{\overline{\kappa}} A^{i_A,i'_A}_{\kappa_1}
A^{i_{1},i'_{1}}_{\kappa_1,\kappa_2}\ldots A^{i_{\ell},i'_{\ell}}_{\kappa_{\ell-1}} |\overline{i} \rangle \langle \overline{i}' |\,.
\label{s_eq:mpo}
\end{align}

Here, the indices $\overline{i}$ and $\overline{i}'$ iterate over the computational basis of each subsystem (physical indices), and the contracted indices $\overline{\kappa}$ (virtual or bond indices) run from $0$ to a maximum value $D_{\rm max}$ called bond dimension, capturing correlations between the sites <cit.>.
Truncating the bond dimension up to a certain threshold corresponds to discarding the smallest singular values in the decomposition mentioned above. This approximation proves highly effective, especially when the subsystems exhibit weak correlations, significantly reducing the computational resources required to manage the state. Any operator can be expressed as a Matrix Product Operator (MPO) and operates on the state by contracting the relevant physical indices.
After the action of a non-local operator, the bond dimension between two sites involved in the evolution typically increases. Therefore, a compression step, i.e., truncation of the bond dimension, is always performed to effectively manage the dimension of the tensor network.
In particular, the propagator (<ref>) features terms acting on $A$ and on two oscillators (the $1$st and the $\ell$th) at the same time. Given the high cost of compression for long-range interactions, we mitigate this by simplifying the propagator to the application of nearest neighbor unitary operations. This is achieved through a suitable swap scheme that maintains the physics unaltered <cit.>. The complete process of evolution-update is illustrated in <ref> a-d, utilizing the Penrose notation for tensors.

Update of the MPO representing the joint state of the atom and the non-Markovian bath according to (<ref>). The tensor in red represents $A^{i_A,i'_A}$, indicating the atom's state. The blue tensors correspond to the chain's harmonic oscillators, while the green ones represent the first and $\ell$th oscillators of the chain, which interact with the atom. Each physical leg denotes the pair of vector spaces corresponding to each tensor. Consequently, the operators acting on the chain should be understood as the corresponding maps, or equivalently, the operators acting on the local purification form matrix product state (PMPS) of the system's density matrix <cit.>. a) After the $n-1$th step the joint atom-bath system and the oscillator $1$ on the left are uncorrelated. b) The $\ell$th oscillator is put on the right of the atom through a sequence of swap operations (not depicted). c) the tensor network corresponding to the propagator (<ref>) acts on the systems $1, A$ and $\ell$. d) After the interaction, the $\ell$th oscillator is traced, $A$ and 1 swap their positions, and a new uncorrelated oscillator is added to the left. A compression follows after each transformation of the chain.
e) Initial state for the scattering dynamics. A chain of oscillators encoding the state of the incoming pulse is placed on the left of the joint atom-bath system.

This scheme can be directly applied to spontaneous-emission dynamics without modification.
For scattering dynamics, an additional step is necessary to represent the incident pulse because the incoming field constitutes a correlated state across multiple oscillators.
Let $m$ be the number of such oscillators. Thus we define the discrete n-particle wavepacket operator as

\begin{align}
\hat{\psi}_n^\dag = (\sum_{i=0}^m \sqrt{\Delta t} \,\xi_i \, \hat{a}_i^\dag)^n \,,
\label{s_eq:n_ph_op}
\end{align}

where $\xi_i$ is a discrete sample of the Gaussian amplitude

\begin{align}\label{eq:gaussian}
\xi(t)_{t_0} = \left(\frac{W^2}{2 \pi}\right)^{1/4} \exp \left\{\frac{1}{4} W^2 (t-t_0)^2 \right\}\,,
\end{align}

with $W$ and $t_0$ the frequency bandwidth and the center of the pulse respectively.
The coherent-state pulse with average photon number $|\alpha|^2$ reads <cit.>

\begin{align}
\ket{\psi}_{\rm pulse} = e^{-|\alpha|^2/2}\sum_{n=0}^{\infty} \frac{\alpha^n}{n!} \hat{\psi}_n^\dag \ket{0}\,.
\end{align}

This state can be generated as a Matrix Product Operator (MPO) by applying the tensor network corresponding to (<ref>) to a chain of $m$ oscillators initially in the vacuum state. For low powers, the first terms of the sum above are sufficient to represent the pulse.
The atomic population, $\varepsilon$,  the output field intensity$I_{\rm out}$, and the two-photon correlation function,
$G^{(2)}$,  are obtained through measurements on the reduced state of the atom, $\hat\rho_A (t) = {\rm Tr}_{\rm F}\{ \hat\rho_{\rm AF}(t)\}$, and of the $\ell$th oscillator, $\hat\rho_\ell (t) = {\rm Tr}_{\rm AF \backslash \ell}\{ \hat\rho_{\rm AF}(t)\}$ respectively:

\begin{align}
\varepsilon (t) &= {\rm Tr}\{ \hat\rho_A (t) \hat\s_+ \hat\s_- \} \,,\\
I_{\rm out} (t) &= {\rm Tr}\{ \hat\rho_\ell (t) \hat{a}^{\dagger}_{\ell}(t) \hat{a}_{\ell}(t) \} \,,\\
G^{(2)} (t) & = {\rm Tr}\{ \hat\rho_\ell (t) \hat{a}^{\dagger}_{\ell}(t) \hat{a}^{\dagger}_{\ell}(t) \hat{a}_{\ell}(t) \hat{a}_{\ell}(t) \} \,.
\end{align}

§ INPUT OUTPUT FORMALISM FOR THE EFFECTIVE MODEL

To simulate the scattering of a photonic wavepacket with the emitter within our effective model we treat the block B waveguide modes, $\hat\beta_\omega$, as a Markovian bath. In this way is possible to eliminate these degrees of freedom obtaining the following master equation for an open multimode cavity QED system

\begin{equation}\label{eq:Me_scattering}
\dot{\hat\rho} = -i[H_{\rm CQED}+H_D, \hat\rho]
+\sqrt{\gamma} \mathcal{D}[\hat A]\hat\rho \,,
\end{equation}

\begin{equation}
H_{\rm CQED}=\omega_0\, \hat\sigma_+\hat\sigma_-{+}\sum_\nu \Omega_{\nu}\, {\hat\alpha}^{\dagger}_{\nu}{\hat\alpha}_{\nu}+\sum_\nu g_\nu \,( \hat\alpha_\nu^\dag \hat\sigma_-+{\rm H.c.})
\end{equation}

is the system Hamiltonian with $\Omega_{\nu}$ and $g_\nu$ being defined in the main text. In Eq. (<ref>) we introduced the usual Linbladian dissipator
$\mathcal{D}[\hat A]\rho=\hat A\rho\hat A^\dag-\{\hat A^\dag \hat A,\rho\}/2$
applied to the collective mode operator $\hat A=\sum_\nu\hat\alpha_{\nu}$ with $\gamma=2 v/L$ being the decay rate of block A into the block B bath. In Eq. (<ref>) we included   a coherent  driving term of the modes $\alpha_\nu$:

\begin{equation}
H_D=\sqrt{\gamma}\sum_{\nu}[E_{\rm in}(t)e^{i\omega_{\rm in} t}\hat\alpha^\dagger_\nu +\rm H.c]
\end{equation}

where $\omega_{\rm in}$  is  the frequency   of the driving  input field.
The shape of the input pulse is determined by the field amplitude $E_{\rm in}(t)$, normalized with respect to the number of photons $n_{\rm ph}$, $\int dt |E_{\rm in}(t)|^2=n_{\rm ph}$. For the scattering process discussed in the main text, we employ the same Gaussian pulse shape as specified in Eq. (<ref>),  $E_{\rm in}(t):=\xi(t)_{t_0}$.
The open-system dynamics of the multimode  cavity QED system, described by Eq. (<ref>), is simulated using a quantum trajectories approach, averaging over $n_t=4000$ trajectories <cit.>.
Once solved the system dynamics
the output field can be reconstructed using the following input-output equation <cit.>:

\begin{equation}
\hat E_{\rm out}(t) =\hat E_{\rm in}(t)+i\sqrt{\gamma}\sum_\nu\hat\alpha_{\nu}(t).
\end{equation}

The intensity of the output field and the two-photon correlation function can then be computed in terms of the modes of the effective cavity and read $I_{\rm out}(t)=\langle \hat\alpha^{\dagger}_{\nu}(t)\hat\alpha_{\nu}(t)\rangle$ and $G^{(2)}(t)=\langle \hat\alpha^{\dagger}_{\nu}(t)\hat\alpha^{\dagger}_{\nu}(t)\hat\alpha_{\nu}(t)\hat\alpha_{\nu}(t)\rangle$, respectively.