\documentclass[fontsize=11pt,paper=a4,open=any,
twoside=no,toc=listof,toc=bibliography,headings=optiontohead,
captions=nooneline,captions=tableabove,english,DIV=12,numbers=noenddot,final,parskip=false,
headinclude=true,footinclude=false,BCOR=0mm]{scrartcl}
\pdfvariable suppressoptionalinfo 512\relax
\synctex=1

\author{Valentin Boettcher}
\usepackage{hirostyle}
\usepackage{hiromacros}
\addbibresource{references.bib}

\title{Input-Output Theory for Modulated Optical Fibre Resonators}
\date{2023}
\graphicspath{{graphics}}

\newcommand{\inputf}[0]{\ensuremath{\mathrm{in}}}
\newcommand{\outputf}[0]{\ensuremath{\mathrm{out}}}
\usetikzlibrary{math}
% \usetikzlibrary{external}
% \tikzexternalize[prefix=tikz/]
\usepackage{pgfplots}

\begin{document}
\maketitle
\tableofcontents

\section{Microscopic Derivation}
\label{sec:micr-deriv}
The setup we are describing consists of a general driven photonic
system \(A\) and a transmission line \(B\). The \(A\) system is
considered to have the Hamiltonian
\begin{equation}
  \label{eq:1}
  H_{A}=H_{0}+V(t) = ∑_{j,β;i,α} \pqty{H_{0}}_{i,α;j,β}a_{j,β}^†a_{i,α}=  ∑_{m} ω^{0}_{m} c_{m}^†c_{m} + V(t),
\end{equation}
where \(\comm{a_{i,α}}{a^†_{j,β}}=δ_{ij}δ_{αβ}\).  We assume that the
system \(A\) consists of several distinct resonators/cavities indexed
by the first index on the \(a^†\), who each have their own lengths
\(L_{A,i}\) and eigen-momenta \(k_{i,α} = 2πα/L_{A,i}\) with
\(α\in\ZZ\).

The eigenmodes of the system \(c_{m}\) are linear combinations of the
bare modes in the photonic system where we have
\begin{equation}
  \label{eq:43}
  c_{m} = ∑_{i,α} T^\ast_{i,α;m}a_{i,α},
\end{equation}
where \(T_{i,α;m}\) is the matrix whose rows are the normalized
eigenvectors of the matrix \(\pqty{H_{0}}_{i,α;j,β}\).


We designate the bare modes of the
EM field that are actually in contact with the transmission line as
the modes with subsystem index \(i=i_{0}\) which is suppressed for
clarity in all expressions concerning that subsystem. We find modes
\(a_{β}\) for the electric field in the subsystem in contact with the
transmission line
\begin{equation}
  \label{eq:4}
  E_{A}(x,t)= \iu \sqrt{\frac{\hbar}{2ε_{0}n_{A}^{2} L_{A,\perp}^{2} L_{A}}} ∑_{β}
  \sqrt{ω_{k_β}} \pqty{a_{β}(t)
    \eu^{\iu k_{β} x } - a_{β}^†(t)  \eu^{-\iu k_{β} x}},
\end{equation}
where \(L_{A,\perp}\) is a length scale that can be interpreted as the
diameter of the transmission line~\cite{Jacobs} and \(L_{A}\) is the
length of the cavity/resonator that hosts the electric field.  The
modes have wave numbers \(k_{β} = 2πβ/L_{A}\) for \(β \in \ZZ\) and
frequencies \(ω_{k_β} = c \abs{k_{β}}/n_{A}\), where \(n_{A}\) is the
refractive index inside the cavity. For simplicity we set \(\hbar
 = 1\) such that energy is measured in units of frequency.

\subsection{Transformations and Rotating Frames}
\label{sec:rotating-frames}

The bare \(a_{β}\) modes are linear combinations of the \(c_{m}\) and
can be related through
\begin{equation}
  \label{eq:5}
  a_{β} = ∑_{m} U_{βm} c_{m},
\end{equation}
where \(U_{βm} = T_{i_{0},β;m}\) is a not necessarily square matrix
that obeys the unitarity relation \(U U^† = \id\).  Transitioning into
a rotating frame with respect to \(H_{0}\) and employing the rotating
wave approximation removes all but the slowest-oscillating rotating
terms from the interaction
\begin{multline}
  \label{eq:12}
  h_{m}(t) = c_{m}(t)\eu^{\iu \pqty{ω^{0}_{m}-ε_{m}}t} \equiv c_{m}(t)\eu^{\iu \tilde{ω}^{0}_{m}t}
  \\\implies H_{A} \to \tilde{H}_{A}=
  ∑_{mn}\pqty{V_{mn}(t) + ε_{m}δ_{mn}}  \eu^{\iu t (ω^{0}_{n}-ω^{0}_{m})}\eu^{-\iu
      (ε_{n}-ε_{m})t}
  {c}_{m}^†{c}_{n} \approx ∑_{mn}\pqty{V^{0}_{mn}+ ε_{m}δ_{mn}} {h}_{m}^†{h}_{n},
\end{multline}
where
\(\abs{ε_{m}}\ll\abs{ω^{0}_{m}}\)\(\abs{ε_{m}-ε_{n}}\ll\abs{ω^{0}_{m}
  - ω^{0}_{n}}\) are the detunings of the drive with respect to the
energy levels of \(H_{0}\). \emph{This constitutes our target
  Hamiltonian which we can control through the modulation of
  \(V(t)\)}
\begin{equation}
  \label{eq:109}
  H^{T}_{mn}= {V^{0}_{mn} + δ_{mn}{ε_{m}-\iu η_{m}}}.
\end{equation}

Due to the coupling to the transmission line we will find that the
equation of motion for the \(h_{m}\) becomes non-unitary with a
damping term
\begin{equation}
  \label{eq:33}
  \iu \dot{h}_{m} = ∑_{n}\bqty{V^{0}_{mn} + δ_{mn}{ε_{m}-\iu η_{m}}}h_{n}.
\end{equation}

We can subsequently find a (non-unitary) transformation that diagonalizes
the RWA interaction
\begin{equation}
  \label{eq:30}
  ∑_{mn}\pqty{O^{-1}}_{γm}\bqty{V^{0}_{mn} + δ_{mn}\pqty{ε_{m}-\iu
      η_{m}}}O_{nγ'} = \pqty{ω_{γ}-\iu λ_{γ}} δ_{γ,γ'}.
\end{equation}
For \(η_{m}=0\) the columns of \(O\) are the normalized eigenvectors
of \(V_{mn}^{0}=\mel{m}{V^{0}}{n}\). So if \(\ket{ψ_{j}}\) is an
eigenvector of \(V\), then \(\braket{i}{ψ_{j}} = O_{ij}\)
\footnote{This is just a reminder for Valentin who can't seem to keep
  this in his head.}. For finite \(η_{m}\) we will find that the
eigenvalues will feature an imaginary part and \(O_{mγ}\) is no longer
unitary, except for the case where all \(η_{m}\) are the same. This
situation occurs if there are other dominating sources of loss such
that the coupling to the transmission line is not a factor or if we
couple to a sufficiently narrow range of modes so that the variation
of damping rates becomes negligible.

Transforming the \(h_{m}\) according to
\begin{equation}
  \label{eq:13}
  d_{γ} = ∑_{n}O^{-1}_{γn} h_{n} = ∑_{n}O^{-1}_{γn} \eu^{\iu
    \tilde{ω}^{0}_{n} t}{c}_{n}  \implies \iu \dot{d}_{γ} = ω_{γ}d_{γ}
\end{equation}
leaves us with a very simple equation of motion.


In summary, the bare modes of the resonators are denoted by
\(a_{j,α}\) where \(j\) refers to the resonator and \(α\) labels the
mode within that resonator. The eigenmodes \(c_{m}\) of the coupled
oscillators obeying \(H_{0}\) are related to the bare modes
\(α_{j,α}\) by \cref{eq:43}. The frame in which the equations of
motion for the modes reflects the target Hamiltonian
\(H^{T}\) are reached through a transformation into a
rotating frame in \cref{eq:12} leading to the \(h_{n}\) modes. The
resulting equations of motion for the \(h_{n}\) can be decoupled by
the transformation \(O_{mγ}\) giving the eigenmodes of the target
Hamiltonian \(d_{γ}\) including damping.

It is important to keep in mind that the actual observables are the
\(α_{β} = α_{i_{0},β}\) which couple to the transmission line, wheras
the modes in the rotating frame \(h_{n}\) correspond to the amplitudes
that evolve according to the target Hamiltonian.

Let us list the relation between the \(a\), \(c\), \(h\) and \(d\) operators
for later reference
\begin{equation}
  \label{eq:67}
  \begin{array}{@{}l| c c c c@{}}
    & α_{β} & c_{m} & h_{n} & d_{γ}\\
    \midrule
    a_{β} &   1   &  U_{βm}\equiv  T_{i_{0},β;m}  &  U_{βm}\eu^{-\iu \tilde{ω}^0_mt}  &  ∑_{m}U_{βm}O_{mγ}\eu^{-\iu \tilde{ω}^0_mt}    \\
    c_{m} &  U^{-1}_{mβ} = U^\ast_{βm}  &  1     &  δ_{mn}\eu^{-\iu \tilde{ω}^0_mt}  &  O_{mγ}\eu^{-\iu \tilde{ω}^0_nt}  \\
    h_{n} &  U^\ast_{βn}\eu^{\iu \tilde{ω}^0_nt}  &   δ_{nm}\eu^{\iu \tilde{ω}^0_mt}  &  1    &  O_{nγ} \\
    d_{γ} &  ∑_{m} U^\ast_{βm} O^{-1}_{γm} \eu^{\iu \tilde{ω}^0_mt}   &   O^{-1}_{γn} \eu^{\iu \tilde{ω}^0_mt}  &  O^{-1}_{γm}  &  1.  \\
  \end{array}
\end{equation}
The quantity \(x_{i}\) in each row is obtained from the quantity
\(y_{j}\) heading each row through the transformation \(A_{ij}(t)\) in
each cell by \(x_{i} = ∑_{j}A_{ij}(t)y_{j}\).

\subsection{Coupling to the Transmission Line}
\label{sec:coupl-transm-line}

The transmission line is considered to only
have one polarization direction and one dimension of
propagation, so that the vector potential is effectively scalar and we
have
\begin{equation}
  \label{eq:2}
  E_{B}(x, t) = \iu\sqrt{\frac{\hbar}{2ε_{0}n_{B}^{2}
      (2π)^{3}L_{B,\perp}^{2}}}  ∫{\sqrt{ω^{B}_{k}}} \pqty{b_{k}(t)
    \eu^{\iu k x } - b_{k}^†(t)  \eu^{-\iu k x}}\dd{k},
\end{equation}
with \(\comm{b_{k}}{b_{q}^†}=δ(k-q)\), \(ω^{B}_{k} = c \abs{k}/n_{B}\)
with \(n_{B}\) being the refractive index of the fibre and
\(L_{B,\perp}\) being the perpendicular length scale as discussed
above. Note that the \(b_{k}\) here have dimensions of \(\sqrt{[L]}\)
as opposed to \(\sqrt{[t]}\), as is the usual convention in
input-output theory. If a stochastic theory is desired, the latter
convention is preferrable and can be obtained through substituting
\(k\to \pm ω/c n_{B}\) and rescaling
\(b_{k}\to b_{k}/ \sqrt{c n_{B}^{-1}}\).

An interaction between the transmission line and the
system \(A\) roughly inspired by coupled mode theory is
\begin{equation}
  \label{eq:3}
  H_{I} = g_{0} ∫ E_{A,+}(x,t)E_{B,-}(x,t) f(x) \dd{x} + \hc,
\end{equation}
where the subscripts \(\pm\) denote positive or negative frequency
portions of the fields and \(f(x)\) is a dimensionless weighting
function with compact support \([-Δx/2, Δx/2]\) whose maximum is
unity.  Coupling only the positive/negative parts simplifies the
calculations and is consistent with the later application of the
rotating wave approximation. A possible phase shift between the fields
has been absorbed into the definition of the creation and annihilation
operators.

Expanding the fields in \cref{eq:3} we obtain
\begin{equation}
  \label{eq:6}
  H_{I} = {g_{0}} \frac{\hbar Δx}{2 ε_{0}n_{A}n_{B} (2π)^{3}
    L_{A,\perp}L_{B,\perp}\sqrt{L_{A}}}  ∑_{β}∫
  \sqrt{ω^{B}_{k}ω_{k_{β}}}\,\tilde{f}(k-k_{β})\,  b^†_{k}
  a_{β} \dd{k} + \hc
\end{equation}
The Fourier transform of the weighting function
\begin{equation}
  \label{eq:7}
  \tilde{f}(k) = \frac{1}{Δx} ∫f(x)\eu^{-\iu k x} \dd{x}
\end{equation}
controls how ``far'' the interaction reaches in \(k\)-space. In the
extreme case \(Δx\to 0\) every \(b_{k}\) couples to every \(a_{β}\),
whereas for \(Δx\to ∞\) only modes with matching wave-numbers
couple. As the \(b_{k}\) will contain both the coherent drive with a
laser and the output field amplitudes it is desirable to have this
coupling to be as local in \(k\)-space as possible for targeted
control and precise readout. In the limit of weak coupling between
transmission line and system, which we will assume in a short while,
the rotating wave approximation will ensure that our result won't
depend significantly on the choice of \(f\).

The coupling constant \(g_{0}\) in \cref{eq:6} has the dimensions of
\([L]^{2}\times [ε_{0}]\). We define a new coupling constant that has
units of energy as
\begin{equation}
  \label{eq:8}
  g_{0} = g\frac{n_{A}n_{B}ε_{0} L_{A,\perp}L_{B,\perp} 2(2π)^{3}}{\hbar ω_{0}},
\end{equation}
where \(ω_{0}\) is a typical frequency\footnote{For example, the
  frequency of the drive laser.}.
Using this, \cref{eq:6} becomes
\begin{equation}
  \label{eq:9}
  \begin{aligned}
    H_{I} &= \frac{gΔx}{
     \sqrt{L_{A}}}  ∑_{β}∫
    G_{β}(k)  b^†_{k}
    a_{β} \dd{k} + \hc
    &G_{β}(k) &= \frac{\sqrt{ω^{B}_{k}ω_{k_β}}}{ω_{0}} \tilde{f}(k-k_{β}).
  \end{aligned}
\end{equation}
We note that for
a \(ω_{k_{β}}= ω_{0} + δω\) with \(δω \ll ω_{0}\) the coupling factor
\(G_{β}(k)\) only depends on the difference \(k-k_{β}\). By defining
\begin{equation}
  \label{eq:11}
  \mathcal{O(k)} = \frac{Δx}{\sqrt{L_{A}}} ∑_{β} G_β(k)a_{β} =
  \frac{Δx}{\sqrt{L_{A}}} ∑_{β,m} G_β(k)U_{β,m}c_{m}
\end{equation}
the interaction takes on the more familiar form
\begin{equation}
  \label{eq:14}
  H_{I} = {g}  ∫
   b^†_{k} \mathcal{O}(k)
  \dd{k} + \hc
\end{equation}

Changing variables from \(k\) to\footnote{This is a bit
  unconventional.} \(ω^{B}_{k}=k c / n_{B}\) in
\cref{eq:9} we obtain
\begin{equation}
  \label{eq:17}
  H_{I} = \frac{gΔx}{\sqrt{L_{A}}}  ∑_{β}∫_{-∞}^{∞}
   G'_{β}(ω)f^†_{ω}
  {a_{β}} \dd{ω} + \hc,
\end{equation}
where \(f_{ω}=\sqrt{\frac{n_{B}}{c}}b_{\frac{ω n_{B}}{c}}\) with
\(\comm{f_{ω}}{f_{ω'}^†}=δ(ω-ω')\) and
\(G'_{β}(ω)=G_{β}\pqty{\frac{ω n_{B}}{c}}\).


\subsection{Rotating Wave and First Markov Approximation}
\label{sec:rotating-wave-first}
Following the route taken in \cite{Jacobs}, the next step would be to
transition into a rotating frame so that
\(\tilde{H}_{A}=\tilde{H}_{B}=0\) and apply the rotating wave
approximation. Here, the rotating terms that would occur have the
frequencies of the form \(ω^{0}_m + ω_{γ}\) which are not guaranteed
to be spaced sufficiently far apart for the RWA to
apply\footnote{Consider, for example the SSH model where the
  \(k\)-space density can be arbitrarily high depending on the length
  of the chain.}. We therefore work in the frame of the
\(h_{m}\) and \(\tilde{f}_{ω} = f_{ω}\eu^{\iu \abs{ω}t}\) to
obtain
\begin{equation}
  \label{eq:10}
  \tilde{H}_{I}= \frac{gΔx}{
    \sqrt{L_{A}}}  ∑_{β,m}∫
  G'_{β}(ω)  \eu^{-\iu
    (\tilde{ω}^{0}_{m}-\abs{ω}) t}
  U_{β,m} \tilde{f}_k^†h_{m} \dd{ω} + \hc.
\end{equation}


\begin{figure}[H]
  \centering
  {\fontsize{8pt}{1em}
  \input{graphics/rwa_illustr.pdf_tex}}
\caption{\label{fig:rwa_illustr} In the rotating wave approximation
  the bare frequencies of the resonator only couple to the
  transmission line in frequency sub-intervals
  \([ω_{m}-λ_{m}, ω_{m}+λ_{m}]\).  A second effect that comes into
  play is the geometrically induced coupling amplitude
  \(\tilde{G'}_{m}(ω)\), which is visualized around \(ω_{m}\) under
  the assumption that all \(G'_{β}\) that enter into
  \(∑_{β}U_{βm}G'_{β}\) share a similar profile which is a valid
  assumption in the applications considered \cref{sec:appl-non-mark}.}
\end{figure}
For \(g \ll \tilde{ω}_{m}^{0}\) each \(h_{m}\) in \cref{eq:10} only
interacts with non-overlapping sub-intervals
\([\tilde{ω}^{0}_{m}-λ_{m}, \tilde{ω}^{0}_{m}+λ_{m}]\) of the transmission frequency axis
(rotating wave approximation) with \(g\ll λ_{m} \ll \tilde{ω}_{m}^{0}\). This
situation is illustrated in \cref{fig:rwa_illustr}. Also, the coupling
amplitude \(G_{β}(ω)\) is local in frequency space and can assist the
RWA depending on the choice of parameters and how close the
\(\tilde{ω}^{0}_{m}\) are to the \(ω_{k_{β}}\), i.e. how local in
frequency space \(U_{βm}\) is. We obtain
\begin{equation}
  \label{eq:16}
  \tilde{H}_{I}\approx \frac{gΔx}{
    \sqrt{L_{A}}}  ∑_{β,m}∫_{\tilde{ω}^{0}_{m}-λ_{m}}^{\tilde{ω}^{0}_{m}+λ_{m}}
  \eu^{-\iu
    (\tilde{ω}^{0}_{m}-\abs{ω}) t}
  U_{β,m} \pqty{G'_{β}(ω)  \tilde{f}_{ω}^† + G'_{β}(-ω)
    \tilde{f}_{-ω}^†}h_{m} \dd{ω} + \hc
\end{equation}
For any finite \(Δx\) and
\(\tilde{ω}_{0}^{m},ω_{k_{β}}\gg \frac{2πc}{Δx n_{A}}\) we can assume
\begin{equation}
  \label{eq:44}
  G'_{β}\pqty{-\sgn(β)  ω}\approx 0
\end{equation}
in \cref{eq:16}.

As each \({h}_{m}\) is now interacting with non-overlapping
transmission-line field modes, we can introduce a separate field for
each \({h}_{m}\) that commutes with all other fields and extend
the integration bounds to infinity again\footnote{This is called the
  ``First Markov Approximation'' in \refcite{Gardiner1985}.}.
Care has to be taken to maintain consistency with \cref{eq:44},
\begin{equation}
  \label{eq:16}
  \tilde{H}_{I}= \frac{gΔx}{
    \sqrt{L_{A}}}  ∑_{β,m}∫_{0}^{∞}
  \eu^{-\iu
    (\tilde{ω}^{0}_{m}-\abs{ω}) t}
  U_{β,m} G'_{β}(\sgn({β})ω)  \tilde{f}^{m,†}_{\sgn({β})ω}{h}_{m} \dd{ω} + \hc
\end{equation}
which becomes\footnote{A lot of discussion for a simple result :).}
\begin{equation}
  \label{eq:18}
  H_{I}= ∑_{m}∫_{-∞}^{∞}
  \tilde{G}_{m}(k) {b}^{m,†}_{k}{c}_{m} \dd{k}
\end{equation}
upon transitioning out of the rotating frame with
\begin{equation}
  \label{eq:15}
  \tilde{G}_{m}(k) =
\frac{gΔx}{
  \sqrt{L_{A}}}  ∑_{β\gtrless 0}U_{β,m} G_{β}(k)δ_{\sgn(β),\sgn(k)}
\end{equation}

The equation of motion
for the transmission line modes become
\begin{gather}
  \iu\dot{b}^{m}_{k} = ω^{B}_{k} b_{k}^{m,†} +
  \tilde{G}_{m}(k) \eu^{-\iu \tilde{ω}^{0}_{m}t}h_{m}\\
    \label{eq:19}
  \implies b^{m}_{k}(t) = b^{m}_{k}(0) \eu^{-\iu ω_{k}^{B}t} -\iu
  \tilde{G}_{m}(k) ∫_{0}^{t}\eu^{-\iu
    ω_{k}^{B}(t-s)} \eu^{-\iu ω_{k}^{B}s}h_{m}(s)\dd{s}.
\end{gather}
The equation of motion for \(h_{m}\) is
\begin{equation}
  \label{eq:21}
  \iu\dot{h}_{m} = ∑_{n}\bqty{V^{0}_{mn} + ε_{m}δ_{nm}}{h}_n +
  \underbrace{\eu^{\iu \tilde{ω}_{m}^{0}t}∫_{-∞}^{∞}\tilde{G}^\ast_{m}(k)
    b_{k}^{m}(t)\dd{k}}_{\equiv I}.
\end{equation}
Further inspection of the rightmost term in \cref{eq:21} yields
\begin{equation}
  \label{eq:22}
  \begin{aligned}
    I &= \eu^{\iu \tilde{ω}_{m}^{0}t} ∫_{-∞}^{∞}\tilde{G}^\ast_{m}(k)
        b_{k}^{m}(t)\dd{k} \\
      &= ∫_{-∞}^{∞}\tilde{G}^\ast_{m}(k)
        b_{k}^{m}(0)\eu^{-\iu ({ω}^{B}_{k} - \tilde{ω}^{0}_{m})t}\dd{k} -\iu  ∫_{0}^{t}∫_{-∞}^{∞}\abs{\tilde{G}_{m}(k)}^{2}
        {h}_{m}(s)\eu^{-\iu {ω}^{B}_{k}(t-s)} \eu^{\iu
        \tilde{ω}^{0}_{m}(t-s)}\dd{k}\dd{s}\\
       &=II + III.
  \end{aligned}
\end{equation}
As the RWA limits the interaction of each \(h_{m}\) to an narrow
frequency/momentum band, we assume that \(\tilde{G}_{m}\) is
approximately constant close to the resonance frequencies \([\tilde{ω}^{0}_{m}-λ_{m}, \tilde{ω}^{0}_{m}+λ_{m}]\) 
\begin{equation}
  \label{eq:23}
  \begin{aligned}
    \tilde{G}_{m}(k) &\approx
    δ_{m}\tilde{G}_{m}\pqty{\sgn(k) ω_{m}^{0}\frac{n_{B}}{c n_{A}}} =
    δ_{m}\frac{gΔx}{\sqrt{L_{A}}}∑_{β}U_{βm}G_{β}\pqty{\sgn(k) ω_{m}^{0}\frac{n_{B}}{c
    n_{A}}} δ_{\sgn(β),\sgn(k)} \\
    &\equiv ∑_{β}g^{0}_{β}\sqrt{ω^{B}_{k}} U_{βm}δ_{\sgn(β),\sgn(k)} \equiv g_{m, \sgn(k)}\sqrt{ω^{B}_{k}} 
  \end{aligned}
\end{equation}
in the interval (see \cref{eq:16}) where \(δ_{m}\) is a possible
scaling factor to better approximate \(\tilde{G}_{m}(k)\) as a
constant in \cref{eq:16}.

Using this in \(I\), we obtain
\begin{equation}
  \label{eq:24}
  \begin{aligned}
    II &= {\eu^{\iu ω_{m}^{0}t}} \bqty{g_{m,+}^\ast ∫_{0}^{∞}\sqrt{ω^{B}_{k}}b^{m}_{k}(0)\eu^{-\iu
    ω^{B}_{k}t}\dd{k} + g_{m,-}^\ast∫_{-∞}^{0}\sqrt{ω^{B}_{k}}b^{m}_{k}(0)\eu^{-\iu
         ω^{B}_{k}t}\dd{k}}\\
    &\equiv {\eu^{\iu ω_{m}^{0}t}}\pqty{
      g_{m,+}^\ast b_{\inputf,+}^{m}(t) + g_{m,-}^\ast b_{\inputf,-}^{m}(t)},
  \end{aligned}
\end{equation}
where \(b_{\inputf,+(-)}^{m}(t)\) is identified as the
right(left)-moving input field and is proportional to the annihilation
part of the electric field. The second part of \cref{eq:22} becomes
\begin{equation}
  \label{eq:25}
  III= -\iu ∫_{0}^{t}\eu^{\iu \tilde{ω}^{0}_{m}(t-s)}{h}_{m}(s)
  \bqty{ \abs{g_{m,+}}^{2} ∫_{0}^{∞}ω^{B}_{k}\eu^{-\iu ω^{B}_{k}(t-s)} \dd{k} + \abs{g_{m,-}}^{2}  ∫_{-∞}^{0}ω^{B}_{k}\eu^{-\iu ω^{B}_{k}(t-s)} \dd{k}}\dd{s}.
\end{equation}
Now we use the identity
\begin{equation}
  \label{eq:26}
  ∫_{0}^{∞}ω^{B}_{k}\eu^{-\iu ω^{B}_{k}(t-s)} \dd{k} = -\iu ∂_{s}\frac{n_{B}}{c}
  \bqty{\mathcal{P}\frac{-i}{t-s} + π δ(t-s)},
\end{equation}
but neglect the principal value, as it leads only to rapidly
oscillating terms that are inconsistent with the RWA, to obtain
\begin{equation}
  \label{eq:27}
  \begin{aligned}
    III&= -2\iu η_{m}∫_{0}^{t}\eu^{\iu
    \tilde{ω}^{0}_m(t-s)}\bqty{\pqty{\tilde{ω}^{0}_{m}+\iu ∂_{s}}\tilde{c}_{m}(s)}
         δ(t-s)\dd{s}\\
       &= -\iu η_{m} h_{m}(t) + \frac{η_{m}}{\tilde{ω}^{0}_{m}}
  \dot{h}_{m}(t) \overset{η_{m}\ll \tilde{ω}_{m}^{0},\, V\ll H_{0}}{\approx}-\iu η_{m} h_{m}(t),
  \end{aligned}
\end{equation}
where the factor \(1/2\) in the last equality stems from the fact that
we only use half of the delta function and
\begin{equation}
  \label{eq:45}
  η_{m}\equiv π\frac{\tilde{ω}^{0}_{m} n_{B}}{c}\bqty{\abs{g_{m,-}}^{2}+\abs{g_{m,+}}^{2}}.
\end{equation}
We have also neglected the correction to the eigen-energy
\(η_{m}/\tilde{ω}^{0}_{m}\) for the same reason we neglected the
principal value.

Note that \cref{eq:45} is an incoherent sum of the couplings to the
right moving and left moving fields in the transmission line.
Altogether we arrive at
\begin{equation}
  \label{eq:28}
  \iu\dot{h}_{m} = {∑_{n}\underbrace{\bqty{V^{0}_{mn} +
        \pqty{ε_{m}-\iu η_{m}}δ_{nm}}}_{\equiv H^{T}_{mn}}{h}_n +
    {\eu^{\iu \tilde{ω}_{m}^{0}t}}
    ∑_{σ=\pm}g_{m,σ}^\ast b_{\inputf,σ}^{m}(t)} .
\end{equation}
The usual situation is that \(b^{m}_{\inputf, -} = 0\) and we can
restrict ourselves to the coupling to the right-moving input field.

\subsection{Input-Output Relation and further Simplifications}
\label{sec:input-outp-relat}
Integrating \cref{eq:19} over all \(k\) yields
\begin{equation}
  \label{eq:29}
  \begin{aligned}
    {b_{\outputf}^{m}(x,t)} &\equiv
  ∫ \sqrt{ω^{B}_{k}} b_{k}^{m}(t) \eu^{\iu k
                                                     t}\dd{k}\\
    &= b_{\inputf}^{m}(x, t) -\iu
      g_{m,\sgn(x)}\frac{\tilde{ω}^{0}_{m}π n_{B}}{c}
      {h}_{m}(τ(x,t))\eu^{-i \tilde{ω}^{0}_{m}τ(x,t)}Θ(τ(x,t)),
  \end{aligned}
\end{equation}
which is the input-output relation with the retarded time
\begin{equation}
  \label{eq:20}
  τ(x,t)=t - \frac{\abs{x}n_{B}}{c}.
\end{equation}
The coupling constant accounts for the direction of propagation and
the time argument is properly retarded. We defined
\begin{equation}
  \label{eq:48}
  b_{\inputf}^{m}(x,t) = ∫\sqrt{ω^{B}_{k}} b_{k}^{m}(0)\eu^{\iu \pqty{kx -
    ω_{k}^{B}t}}\dd{k}
\end{equation}
used that
\begin{equation}
  \label{eq:42}
  ∫_{0}^{∞}ω^{B}_{k}\eu^{-\iu ω^{B}_{k}(t-s)}\eu^{\pm\iu k x} \dd{k} =
  -\iu ∂_{s}\frac{n_{B}}{c}
  \bqty{\mathcal{P}\frac{-i}{t-s \pm \frac{x n_{B}}{c}} + π δ\pqty{t-s\mp
    \frac{x n_{B}}{c}}},
\end{equation}
along with similar arguments to the above.

The case of \(x=0\) is recovered by defining
\begin{equation}
  \label{eq:47}
  \lim_{x\to0} g_{m,\sgn(x)=0} = \frac{1}{2} \pqty{g_{m,+} + g_{m,-}},
\end{equation}
which amounts to taking half of each delta function in
\cref{eq:42}. It shall be noted, that it is physical to assume
\(x>0\), as we necessarily measure outside the fibre-coupler between
transmission line and resonator. By neglecting the \(k\)-dependence of
the coupling in \cref{eq:23} through invocation of the RWA we have
effectively ignored length \(Δx\), but to maintain consistency with
\cref{eq:44} we should assume it to be finite.
We can also neglect the retardation if \(x / v_{g}\) is
much smaller than a typical timescale we're interested in.


To integrate \cref{eq:28}, we
first diagonalize \(V^{0}_{mn} + δ_{mn}\pqty{ε_{m}-i η_{m}}\)
\begin{equation}
  \label{eq:65}
  V^{0}_{mn} + δ_{mn}\pqty{ε_{m}-i η_{m}} \to ∑_{γγ'}
  δ_{γγ'}(ω_{γ}-\iu λ_{γ})
\end{equation}
using the transformation \(O_{nγ}\) (see \cref{sec:rotating-frames})
and find
\begin{equation}
  \label{eq:32}
  \dot{d}_{γ} = ∑_{m}O^{-1}_{γm}\dot{h}_{m} =
  -\iu\bqty{\pqty{ω_{γ} - \iu λ_{γ}}d_{γ} +
    ∑_{σ=\pm}∑_{m}\pqty{O^{-1}}_{γm}{g_{m,σ}^\ast } \eu^{\iu \tilde{ω}_{m}^{0}t}
    b_{\inputf,σ}^{m}(t)}.
\end{equation}

We now introduce some additional simplifications beginning with
equating all input fields \(b_{\inputf}^{m}\). This is allowed, as we
will transition to the classical picture later, where the commutation
relations do not matter. We also assume that we're working in a region
in \(m\) space, where the \(\sqrt{ω_{0}}g_{β}^{0}\approx \sqrt{κ}\) and
\(\tilde{ω}^{0}_{m}\approx{ω_{0}}\), where \(ω_{0}\) is a typical
frequency in the input field, can be assumed to be approximately
constant. With these considerations in mind we can simplify
\cref{eq:45,eq:32} to
\begin{equation}
  \label{eq:64}
  η_{m}=\abs{κ}\frac{πn_{B}}{c}∑_{σ=\pm,β,β'}U_{βm}U^\ast_{β'm}δ_{\sgn(β),σ} δ_{\sgn(β'),σ}
\end{equation}
and
\begin{gather}
  \label{eq:34}
  \dot{d}_{γ} =
  -\iu\bqty{\pqty{ω_{γ}-\iu λ_{γ}}d_{γ} + \sqrt{κ^\ast} ∑_{σ=\pm}
    U^{σ}_{γ}(t) {b_{\inputf}(t)}}\\
  U^{σ}_{γ}(t) = ∑_{m,β} δ_{\sgn({β}),σ}U^\ast_{βm}\pqty{O^{-1}}_{γm}\eu^{\iu \tilde{ω}_{m}^{0}t}.
\end{gather}

These simplifications still capture the essence of the physics, as
demonstrated in the current long-range SSH experiment.

We can now proceed to integrate \cref{eq:34} to obtain
\begin{equation}
  \label{eq:36}
  d_{γ}(t)= d_{γ}(0) \eu^{-\pqty{\iu ω_{γ} + {λ}_{γ}}t} -
  \frac{i}{\sqrt{κ}} Σ_{σ=\pm} ∫_{0}^{t}χ_{γ}(t-s) U^{σ}_{γ}(s)
  {b_{\inputf,σ}(t)}\dd{s}
\end{equation}
with
\begin{equation}
  \label{eq:37}
  χ_{γ}(t) = \abs{κ} \eu^{-\pqty{\iu ω_{γ} + λ_{γ}}t}.
\end{equation}

When constructing the total output field, we have to remember how the
separate fields \(b_{\outputf,m}\) came about. We assumed that each
\(c_{m}\) only interacted with a finite range of modes (see
\cref{eq:16}) in the transmission line and then just extended the
resulting sub-fields back to full independent fields for
simplicity. Now, we have to perform the reverse process, which amounts
to summing together all system (resonator) contributions in
\cref{eq:34} as these only excite the sub-fields and we can safely
glue them back together. To be consistent, we have to sum together the
finite ranges of the input fields which amounts to having \emph{one}
whole copy of the input field.
This leads us to
\begin{equation}
  \label{eq:38}
  {b_{\outputf}(x,t)} \equiv
  b_{\inputf}(x, t) -i  θ(τ(x,t)) \frac{\sqrt{κ} πn_{B}}{c}
  ∑_{γ}\bqty{U^{\sgn(x)}_{γ}\pqty{τ(x,t)}}^\ast d_{γ}(τ(x,t))
\end{equation}

Transitioning to expectation values and using \(\ev{d_{γ}(0)}=0\) we
find
\begin{equation}
  \label{eq:39}
  \ev{{b_{\outputf}(x,t)}} =
  \ev{b_{\inputf}(x,t)} - ∑_{σ=\pm}∫_{0}^{τ(x,t)}χ_{\sgn(x),σ}(τ(x,t),s) \ev{b_{\inputf,σ}(s)} \dd{s}
\end{equation}
with the time non-local susceptibility for the left and right moving
input fields
\begin{equation}
  \label{eq:40}
  χ_{δ,σ}(t,s) = \frac{π n_{B}}{c}Θ(t) ∑_{γ}\pqty{U^{δ}_{γ}(t)}^\astχ_{γ}(t-s)U^{σ}_{γ}(s).
\end{equation}

For an input field with no left-moving components and a measurement
position \(x>0\) we have
\begin{equation}
  \label{eq:31}
  \ev{{b_{\outputf}(x>0,t)}} =
  \ev{b_{\inputf}(x,t)} -∫_{0}^{τ(x,t)}χ_{++}(τ(x,t),s) \ev{b_{\inputf}(s)} \dd{s}.
\end{equation}
with \(b_{\inputf}(s) = b_{\inputf,+}(s) + b_{\inputf,-}(s) =
b_{\inputf,+}(s)\). 

\subsection{Langevin-Equations for Lossy Oscillators}
\label{sec:lang-equat-lossy}

In the above we have assumed that \(H_{0}\) is hermitian. This,
however, ceases to be the case when we assume some a-priori
phenomenological decay in the bare components of the system and we
cannot write \(H_{0}=H^{H}_{0} - \iu η \id\) with \(H^{H}_{0}\)
hermitian. To retain consistency, the decay rates have to be
introduced on the level of the equations of motion of the mode
operators \(a_{i,α}\) after deriving them from the hermitian
Hamiltonian. The equations of motion can then still be decoupled by
diagonalizing the non-hermitian that includes the phenomenological
decay.

We find\footnote{Assuming that the non-hermiticity is small enough for
the matrix to remain diagonalizable.}
\begin{equation}
  \label{eq:53}
  ∑_{iα;jβ}\pqty{T^{-1}}_{m;i,α} \pqty{H_{0}}_{i,α;j,β}T_{j,β;n} =
  \pqty{ω_{m}^{0}-\iu η_{m}^{0}}δ_{nm},
\end{equation}
where \(T\) is the matrix whose rows are the eigenvectors of
\(H_{0}\). Note that \(T\) is not unitary anymore.  For notational
convenience we will write \(T^{-1}_{m;i,α}\) instead of
\(\pqty{T^{-1}}_{m;i,α}\) and use explicit fractions if we want to
express the multiplicative inverse. The mode operators transform as
\begin{equation}
  \label{eq:60}
  c_{m} = ∑_{i,α} T^{-1}_{m;i,α}a_{i,α},
\end{equation}
which are \emph{not} to be identified with bosons anymore, as the
non-unitarity of \(T\) breaks the bosonic commutation
relations. Again, we express the modes that are in contact with the
transmission line as \(a_{α}=a_{i_{0},α}\) and find
\begin{equation}
  \label{eq:69}
  α_{α} = ∑_{α} T_{i_{0},α;m}c_{m} \equiv ∑_{α}U_{αm} c_{m}.
\end{equation}
For convenience we define
\begin{equation}
  \label{eq:70}
  U^{-1}_{mα}\equiv T^{-1}_{m;i_{0}α}.
\end{equation}

The modulation term \(V\) transforms as
\begin{equation}
  \label{eq:74}
  V_{mn}=∑_{iα;jβ}\pqty{T^{-1}}_{m;i,α} V_{i,α;j,β}T_{j,β;n},
\end{equation}
and is no longer hermitian.

We start by writing down the equations of motion for the original
modes, assuming \(H_{0}\) to be hermitian, introduce the non-hermitian
terms and express everything in terms of the \(c_{m}\) using
\(T\). Subsequently, we change into a rotating frame
\begin{equation}
  \label{eq:66}
  h_{m} = c_{m}\eu^{\iu \tilde{ω}^{0}_{m}t},
\end{equation}
rotating away only the unitary evolution. Applying the rotating wave
and first Markov approximations works out precisely as in
\cref{sec:rotating-wave-first}.


To account for non-unitarity we have to make the following
replacements along the way
\begin{align}
  \label{eq:68}
  \tilde{G}_{m}(k) &\rightarrow \tilde{G}_{m}(k)= \frac{gΔx}{\sqrt{L_{A}}} ∑_{β} U_{βm}
                     G_{β}(k) δ_{\sgn(β),\sgn(k)}\\
  \tilde{G}^\ast_{m}(k) &\rightarrow\tilde{G}^{-1}_{m}(k) = \frac{g^\astΔx}{\sqrt{L_{A}}} ∑_{β} U^{-1}_{mβ}
                          G^\ast_{β}(k) δ_{\sgn(β),\sgn(k)}\\
  g_{m,σ}&\rightarrow g_{m,σ}=∑_{β}g^{0}_{β} U_{βm}δ_{\sgn(β),σ}\\
  g^\ast_{m,σ}&\rightarrow \bar{g}_{m,σ}=∑_{β}\pqty{g^{0}_{β}}^\ast U^{-1}_{mβ}δ_{\sgn(β),σ}\\
\end{align}
which gives us
\begin{align}
  \label{eq:72}
  η_{m}=\frac{\tilde{ω}_{m}^{0}π n_{B}}{c} ∑_{σ} g_{mσ}\bar{g}_{mσ},
\end{align}
which might have an imaginary part.

This leaves us with
\begin{equation}
  \label{eq:73}
  \iu\dot{h}_{m} = {∑_{n}\Bqty{V^{0}_{mn} + \bqty{ε_{m}-\iu \pqty{η_{m} +
    η_{m}^{0}}}δ_{nm}}{h}_n +
    {\eu^{\iu \tilde{ω}_{m}^{0}t}}
    ∑_{σ=\pm}\bar{g}_{m,σ} b_{\inputf,σ}^{m}(t)},
\end{equation}
where the \(η_{m}\) might shift the energies \(ε_{m}\)
slightly. 

Diagonalizing
\begin{equation}
  \label{eq:77}
  ∑_{mn}O^{-1}_{γ'm}\bqty{∑_{n}\Bqty{V^{0}_{mn}+ \bqty{ε_{m}-\iu
        \pqty{η_{m}^{0}+η_{m}}}δ_{nm}}h_{m}}O_{nγ} = \pqty{ω_{γ}-\iu λ_{γ}}δ_{γ,γ'}
\end{equation}
and defining
\begin{equation}
  \label{eq:78}
  d_{γ} = ∑_{n}O^{-1}_{γn}h_{n} \implies h_{n}=∑_{γ}O_{nγ}d_{γ}
\end{equation}
will give us the equivalent of \cref{eq:32}
\begin{equation}
  \label{eq:80}
  \dot{d}_{γ}=-\iu\pqty{ω_{γ}+\sqrt{κ^\ast}∑_{σ=\pm}U^{σ}_{γ}{b_{\inputf,σ}}}d_{γ}
  - λ_{γ}d_{γ}
\end{equation}
where we have set \(\sqrt{ω_{0}}g_{β}^{0}=\sqrt{κ}\) and defined
\begin{equation}
  \label{eq:82}
  U^{σ}_{γ} =
  ∑_{mβ}\eu^{\iu\tilde{ω}^0_mt}O^{-1}_{γm}U^{-1}_{mβ}δ_{\sgn(β),σ}.
\end{equation}
This also simplifies \cref{eq:64} to
\begin{equation}
  \label{eq:88}
  η_{m}=\abs{κ}\frac{πn_{B}}{c}∑_{σ=\pm,β,β'}U_{βm}U^{-1}_{mβ}δ_{\sgn(β),σ} δ_{\sgn(β'),σ}.
\end{equation}

Further defining
\begin{align}
  \label{eq:83}
  \bar{U}^{σ}_{γ}&=∑_{mβ}\eu^{-\iu\tilde{ω}^0_mt}O_{mγ}U_{βm}δ_{\sgn(β),σ}\qq{and}&
  χ_{γ}&=\abs{κ}\eu^{-\pqty{\iu ω_{γ}+λ_{γ}}t},
\end{align}
we obtain
\begin{equation}
  \label{eq:86}
  \ev{{b_{\outputf}(x,t)}} =
  \ev{b_{\inputf}(x,t)} - ∑_{σ=\pm}∫_{0}^{τ(x,t)}χ_{\sgn(x),σ}(τ(x,t),s) \ev{b_{\inputf,σ}(s)} \dd{s}
\end{equation}
with the time non-local susceptibility for the left and right moving
input fields
\begin{equation}
  \label{eq:87}
  χ_{δ,σ}(t,s) = \frac{π n_{B}}{c}Θ(t) ∑_{γ}\bar{U}^{δ}_{γ}(t)χ_{γ}(t-s)U^{σ}_{γ}(s).
\end{equation}

These equations are essentially the same as \cref{eq:39,eq:40},
accounting for the non-unitary transformations and the apriori decay
rates when diagonalizing the equations of motion for the \(h_{m}\).

For completeness, we give the equivalent of \cref{eq:67} for the non-unitary case
\begin{equation}
  \label{eq:35}
  \begin{array}{@{}l| c c c c@{}}
    & α_{β} & c_{m} & h_{n} & d_{γ}\\
    \midrule
    a_{β} &   1   &  U_{βm}\equiv  T_{i_{0},β;m}  &  U_{βm}\eu^{-\iu \tilde{ω}^0_mt}  &  ∑_{m}U_{βm}O_{mγ}\eu^{-\iu \tilde{ω}^0_mt}    \\
    c_{m} &  U^{-1}_{mβ}\equiv  T^{-1}_{m;i_{0},β} &  1     &  δ_{mn}\eu^{-\iu \tilde{ω}^0_mt}  &  O_{mγ}\eu^{-\iu \tilde{ω}^0_nt}  \\
    h_{n} &  U^{-1}_{nβ}\eu^{\iu \tilde{ω}^0_nt}  &   δ_{nm}\eu^{\iu \tilde{ω}^0_mt}  &  1    &  O_{nγ} \\
    d_{γ} &  ∑_{m} U^{-1}_{mβ} O^{-1}_{γm} \eu^{\iu \tilde{ω}^0_mt}   &   O^{-1}_{γn} \eu^{\iu \tilde{ω}^0_mt}  &  O^{-1}_{γm}  &  1  \\
  \end{array}.
\end{equation}

\section{Application to the Non-Markovian Quantum Walk}
\label{sec:appl-non-mark}
The experimental setup for implementing the non-Markovian quantum walk
discussed in~\cite{Ricottone2020,Kitagawa2010} is illustrated in
\cref{fig:schematic}. The abstract system introduced in
\cref{sec:micr-deriv} is replaced by a small \(S\) and a large \(B\)
fibre loop with lengths \(L_{B}\gg L_{S}\). The resonant modes of the
loops do have the free spectral ranges
\begin{equation}
  \label{eq:41}
  \begin{aligned}
    Ω_{B} &= \frac{2πc}{n_{B}} & Ω_{S} &= \frac{2πc}{L_{S}},
  \end{aligned}
\end{equation}
where \(n\) is the refractive index of the respective fibres.

Attaching the transmission line to the smaller loop has the advantage
that we only excite and detect what we will later identify as the
\(A\) level of the non-Markovian quantum walk in momentum space.

We assume that the loops share some common eigenfrequency
\begin{equation}
  \label{eq:46}
  ω_{s} = n_{0}^{S}Ω_{S} = n_{0}^{B}Ω_{B}\iff \frac{n_{0}^{B}}{L_{B}}
  = \frac{n_{0}^{S}}{L_{S}}
\end{equation}
with \(n,m\gg 1\) and \(ω_{s}\) close to the frequency of the
laser. We choose to index the eigenfrequencies of the loops relative
to \(ω_{s}\) so that
\begin{equation}
  \label{eq:49}
  \begin{aligned}
    ω_{n}^{X} &= ω_{s} + n Ω_{X} & k_{n}^{X} = \frac{2π}{L_{B}}
                                   n_{0}^{X} + \frac{2π}{L_{B}} n = k_{0} + \frac{2π}{L_{B}} n,
  \end{aligned}
\end{equation}
where \(X=S,B\).

\begin{figure}[H]
  \centering
  \includegraphics[width=.9\textwidth]{walker_setup}
  \caption{\label{fig:schematic}Schematic of the experimental setup
    and the different frames that the input-output theory deals
    with. Two fibre loops \(B\) (Big) and \(S\) (small) are being
    coupled with strength \(δ\).  The golden parts correspond to the
    bare modes (\(ω_{m}^{B},ω_{s}\)) being coupled to form the
    eigenmodes with frequencies \(ω_{m}^{0}\) of the composite \(B+S\)
    system. We then transition into a rotating frame (Blue) with
    respect to \(H_{0}\) and the \(ω^{0}_{m}\) and consider a
    modulation through the EOM to give us an effectively
    time-independent target Hamiltonian \(H^{T}_{mn}\). Diagonalizing
    this Hamiltonian gives us the energies \(ω_{λ}\) that we can
    detect when measuring the transmission (Green). Damping terms have
    been left out of the picture for simplicity.}
\end{figure}
These loops are coupled to each other with amplitude \(δ\eu^{\iu ϕ}\) which can
possibly contain a phase due to the choice of the coordinate
origin. Let us now assume that
\begin{equation}
  \label{eq:50}
  \frac{Ω_{S}}{Ω_{B}} = 2N+1
\end{equation}
with \(N\in \NN\). We denote by \(a\) the annihilation operator for the
mode with \(ω=ω_{s}\) in the small loop and suppress all other modes
in the small loop, as we will not populate them. The annihilation
operators \(f_{n}\) destroy modes with frequencies \(ω^{B}_{n}\) in
the big loop, where we limit \(n\) to the range \([-N, N]\) for the
same reasons as above.
This leads to the Hamiltonian
\begin{equation}
  \label{eq:51}
  H_{A}= ω_{s}a^†a + ∑_{n=-N}^{N} ω_{n}^{B}f_{n}^†f_{n} + δ \pqty{\eu^{\iu
    ϕ}f_{0}^†a + \eu^{-\iu ϕ}a^†f_{0}} + ∑_{mn}J_{mn}(t) f_{m}^†f_{n},
\end{equation}
where \(J_{mn}(t)\) is the coupling mediated by the EOM in the big
loop. The phase \(ϕ=k_{0}L_{S}/2\) of the coupling stems from the
location of the fibre coupler along the loop but can potentially also
contain other contributions.

To be bring the Hamiltonian into the form \cref{eq:1}, we define
\begin{equation}
  \label{eq:52}
  \begin{aligned}
    c_{\pm} &= \frac{1}{\sqrt{2}}\pqty{f_{0}\eu^{-\iu ϕ} \pm a} &
                                                                  c_{n\neq
                                                                  0}&=f_{n}\\
    ω_{\pm}^{0}&=ω_{s}\pm δ & ω_{n\neq 0}&= ω_{s} + Ω_{B} n\\
    V_{\pm,\pm}&=\frac{J_{00}}{2} & V_{\pm, n\neq0} &=
                                                      \frac{J_{0n}}{\sqrt{2}}\eu^{\iu
                                                      ϕ}\\
    V_{n\neq \pm, m\neq \pm} &= J_{nm}
  \end{aligned}
\end{equation}
and obtain
\begin{equation}
  \label{eq:54}
  H_{A} = ∑_{n} ω_{n}^{0}c_{n}^†c_{n} + ∑_{nm} V_{nm}(t) c_{n}^†c_{m},
\end{equation}
where the index \(n\) can take on the values in
\(\pqty{[-N,N]\setminus \{0\}} \cap \{+, -\}\). The spectra of the
coupled and uncoupled systems are visualized in \cref{fig:spectra}.


Upon transitioning into the rotating frame and applying the RWA we
arrive at a target Hamiltonian
\begin{equation}
  \label{eq:111}
  H^{T}_{nm}=V^{0}_{n,m} + \pqty{ε_{m} -\iu \pqty{η^{0}_{m} + η_{m}}
    δ_{nm}},
\end{equation}
where we can control the \(V^{0}_{n,m}\) through the drive amplitudes
and the \(ε_{m}\) through the drive
detunings.

The choice of \(δ\) will be discussed in
\cref{sec:choice-hybridization} and the origin of the \(η^{0}_{m}\),
as well as the effects of asymmetric loss in the two loops are
investigated in \cref{sec:effects-asymm-damp}.
\Cref{sec:rotat-wave-inter} is concerned with relating the
\(V^{0}_{nm}\) with the drive amplitudes, frequencies and phases,
whereas \cref{sec:steadyst-transm,sec:steady-state-transm} will deal
with the transmission through transmission lines attached to either
the small or the big loop.

\tikzmath{
  \Nmodes = 5;
  integer \NmodesBetween;
  \NmodesBetween = \Nmodes - 1;
  \delta = .2;
}
\begin{figure}[p]
  \centering
\begin{tikzpicture}[y=-1.5cm]
  \draw[->] (-1.5, -1) -- node[left] {\(ω\)} ++(0,2);

  \foreach \y in {-\Nmodes,...,\Nmodes}
  \draw  (0, \y) node[left] {\(ω^{B}_{\y}\)} -- ++(1, 0);

  \foreach \y/\name in {-\Nmodes/-1,0,\Nmodes/1}
  \draw (1.1, \y) -- ++(1, 0) node[right] {\(ω^{S}_{\name}\)};
\end{tikzpicture}
\begin{tikzpicture}[y=-1.5cm]
  \foreach \y in {-\NmodesBetween,...,-1}
  \draw  (0, \y) node[left] {\(ω^{0}_{\y}\)} -- ++(1, 0);

  \foreach \y in {1,...,\NmodesBetween}
  \draw  (0, \y) node[left] {\(ω^{0}_{\y}\)} -- ++(.5, 0) node(mc\y){}
  -- ++(.5, 0) node(mr\y){};

  \foreach \y in {-\Nmodes,\Nmodes} {
    \foreach \sub/\name in {-\delta/-, \delta/+}
    \draw  (0, \sub+\y) --  ++(.5, 0)  -- ++(.5, 0);
  }

  \foreach \y in {0} {
    \foreach \sub/\name in {-\delta/-, \delta/+}
    \draw (0, \sub+\y) node[left] {\(ω^{0}_{\name}\)} --  ++(.5, 0) node(pmc\name){} -- ++(.5, 0) node[right] (pmr\name){};
  }

  \foreach \y in {-\Nmodes,0,\Nmodes}
  \draw[dashed,color=gray]  (0, \y) -- ++(1, 0);

  \draw[<->] (pmc+) -- node[right] {\(2δ\)} (pmc-);
  \draw[<->] (pmc+) -- node[right] {\(Ω_{B} - δ\)} (mc1);
  \draw[<->] (mc1) -- node[right] {\(Ω_{B}\)} (mc2);
\end{tikzpicture}
\caption{\label{fig:spectra}The spectra of the uncoupled loops (left) \(B,S\) and the resulting
  spectrum after coupling the loops (right) for \(N=2\) and \(δ=Ω_{B}/5\).}
\end{figure}

\subsection{The Choice of the Hybridization Amplitude}
\label{sec:choice-hybridization}

It remains to be discussed which choice of \(δ\) is suitable. By
modulating \(V(t)\) and applying the RWA we can couple certain levels
in the spectrum of the system. To be able to couple one level to many
and implement the non-Markovian quantum walk, we have to select out
unique level-spacings. At the same time, we want to maximize the
frequency of the residual rotating terms.

A suitable mode for this one-to-many coupling is the \(c_{+}\) mode
with frequency \(ω_{+}^{0}=ω_{s}+δ\).  We identify the \(A\) site of
the non-Markovian Quantum Walk in \(k\)-space with \(c_{+}\) and the
\(j\)th bath level with \(c_{j}\) (\(j\in [1, N]\)).

Let us denote the frequency
difference of the \(m\)th and \(n\)th mode by
\(Δ_{nm}=ω^{0}_{n}-ω^{0}_{m}\).
For \(n>m\) we have
\begin{align}
  \label{eq:56}
  Δ_{+,-}&=2δ\\
  \label{eq:57}Δ_{n>0,-} &= n Ω_{B} + δ\\
  \label{eq:58}Δ_{n>0,+} &= n Ω_{B} - δ\\
  \label{eq:59}Δ_{n>0,m>0} &= (n-m) Ω_{B}.
\end{align}
To couple exactly one other mode with \(n>0\) to the \(+\) mode, the
RWA requires that \(\abs{Δ_{n,+}}\neq \abs{Δ_{kl}}\) for \(k\neq n\), \(l\neq +\).
This requirement yields the following restrictions on the value of
\(δ\)
\begin{align}
  \label{eq:61}
  δ&\neq \frac{n}{2}Ω_{B} & \text{\cref{eq:58,eq:57}}\\
  δ&\neq \frac{n}{3}Ω_{B} & \text{\cref{eq:58,eq:56}}\\
  δ&\neq {n}Ω_{B} & \text{\cref{eq:58,eq:59}\footnote{Duh!}}.
\end{align}
To maximize the residual rotating terms, the minimum of the
\cref{eq:56,eq:57,eq:58,eq:59} has to be maximized for \(δ \in [0, Ω_{B}]\)
\begin{equation}
  \label{eq:62}
  Δ_{\max}\equiv \max_{δ}Δ_{\min}(δ) = \max_{δ}\min\Bqty{2δ, \abs{Ω_{B}-δ}, \abs{Ω_{B}-3δ},
    \abs{2Ω_{B}-3δ}, \abs{Ω_{B}-2δ}}.
\end{equation}
We find that \(Δ_{\max}=2Ω_{B}/5\) for
\begin{equation}
  \label{eq:63}
  δ_{\mathrm{opt}}=Ω_{B}/5,
\end{equation}
as can be ascertained from \cref{fig:delta_choice}.
\begin{figure}[H]
  \centering
  \begin{tikzpicture}
    \begin{axis}[
      scale only axis=true,
      width=.8\columnwidth,
      height=.2\columnwidth,
      xmin = 0, xmax = 1,
      ymin = 0, ymax = .5,
      axis lines* = left,
      xtick = {0}, ytick = \empty,
      clip = false,
      xtick={},ytick={},
      minor tick num=5,
      grid=both,
      grid style={line width=.1pt, draw=gray!10},
      major grid style={line width=.2pt,draw=gray!50},
      axis lines=middle,
      ylabel = {\(Δ_{\min}/Ω_{B}\)},
      xlabel = {\(δ/Ω_{B}\)},
      x label style={at={(axis description cs:0.5,-0.2)},anchor=north},
      y label style={at={(axis description cs:-0.06,.5)},rotate=90,anchor=south},
      ]
      \addplot[domain = 0:1, restrict y to domain = 0:1, samples =
      1000]{min(2*x, 1-x, abs(1-3*x), abs(2-3*x),
        abs(1-2*x))};
      \addplot[color = black, mark = *, only marks, mark size = 3pt]
      coordinates {(.2, .4)};
      \addplot[color = black, dashed, thick] coordinates {(.2, 0) (.2,
        .4) (0, .4)};

      \addplot[color = gray, mark = *, only marks, mark size = 3pt]
      coordinates {(.25, .25)};
      \addplot[color = gray, dashed, thick] coordinates {(.25, 0) (.25,
         .25) (0, .25)};
    \end{axis}
  \end{tikzpicture}
  \caption{\label{fig:delta_choice} The minimal rotating frequencies
    \cref{eq:62} for the range of possible \(δ\). The black marker
    highlights \(δ=Ω_{B}/5\) and the grey marker marks \(δ=Ω_{B}/4\).}
\end{figure}




\subsection{Effects of Asymmetric Damping}
\label{sec:effects-asymm-damp}

In the above, we have not accounted for the apriori damping in the
fibre loops. As they are of vastly differing lengths, it is sensible
to expect different damping rates for each.

To account for this, we modify \cref{eq:51} to
\begin{equation}
  \label{eq:89}
  \begin{aligned}
    H_{A}= ω_{s}a_{S,0}^†a_{S,0} &+ ∑_{n=-N}^{N} ω_{n}^{B}a_{B,n}^†a_{B,n} + δ \pqty{\eu^{\iu
                                   ϕ}a_{B,0}^†a_{S,0} + \eu^{-\iu ϕ}a_{S,0}^†a_{B,0}} + ∑_{mn}J_{mn}(t) a_{B,m}^†a_{B,n}
    \\&- \iu η_{S}a_{S,0}^†a_{S,0}  - \iu \pqty{η_{S}+2Δδ} ∑_{n}a_{B,n}^†a_{B,n}
  \end{aligned}
\end{equation}
with the damping rates \(η_{S}\) and the damping asymmetry
\(\abs{Δ}< 1\).


The diagonalizing transformation is then
{\renewcommand\arraystretch{1.5}
\begin{equation}
  \label{eq:90}
  \begin{aligned}
    T_{i,0;\mp}  &\triangleq \frac{1}{\sqrt{2}}
  \begin{pmatrix}
    -\eu^{-\iu ψ} & \eu^{\iu ψ} \\
    \eu^{\iu ϕ} & \eu^{\iu ϕ}
  \end{pmatrix}
  & T^{-1}_{\mp;i,0} &\triangleq
                    \frac{1}{\sqrt{2}\cos(ψ)}
                    \begin{pmatrix}
                      -1 & \eu^{-\iu (ϕ-ψ)}\\
                      1 & \eu^{-\iu (ϕ+ψ)}
                    \end{pmatrix}
  \end{aligned}
\end{equation}}
with \(i=B,S\) and
\begin{equation}
  \label{eq:93}
  \eu^{\iu ψ} \equiv \iu Δ + \sqrt{1-Δ^{2}}  \implies ψ = \sin^{-1}(Δ).
\end{equation}
For all other matrix elements \(T\) equals the
identity
\begin{equation}
  \label{eq:95}
  T_{i,α;m} = T^{-1}_{m;i,α} = δ_{i,B} δ_{αm}
\end{equation}
for \(m\neq\pm\).


With this, we have for \(m=\pm,-N,-N+1,\ldots,-1,1,2,\ldots,N\)
\begin{equation}
  \label{eq:94}
  \begin{aligned}
    ω^{0}_{m\neq \pm} &= ω_{n}^{B} & ω^{0}_{\pm} &= ω_{s} \pm δ \cos(ψ)\\
    η^{0}_{m\neq\pm} &= η_{B} = η_S+2δΔ & η^{0}_{\pm} &=  η_{S} + δΔ =
                                                      \frac{η_{S}+ η_{B}}{2}.
  \end{aligned}
\end{equation}
Note that the eigenenergies of the \(\pm\) modes are slightly shifted. 

Next, we compute the transformation of the interaction
\(T^{-1}VT\) (see \cref{eq:74}) and find
\begin{equation}
  \label{eq:96}
  \begin{aligned}
    V_{\pm,n\not\in \{+,-\}} &= \frac{J_{0,n}}{\sqrt{2}\cos(ψ)}\eu^{-\iu(ϕ\pm ψ)}  &
    V_{n\not\in \{+,-\},\pm}&= \frac{J_{n,0}}{\sqrt{2}}\eu^{\iu ϕ} & V_{n\neq \pm,
                                                         m\neq \pm} &=
                                                         J_{nm} \\
    V_{+,-}&=\frac{J_{00}}{2\cos(ψ)}\eu^{-\iu ψ} & V_{-,+} &=\frac{J_{00}}{2\cos(ψ)}\eu^{\iu ψ}.
  \end{aligned}
\end{equation}
Evidently, the matrix \(T^{-1}VT\) is not hermitian except in the
limit \(Δ\to 0\).

In all of the above one can set \(\cos(ψ)=1\) to only account for
\(Δ\) in leading order.


\subsection{Rotating-Wave Interaction}
\label{sec:rotat-wave-inter}
The modulation term in \cref{eq:89} can be written as
\begin{equation}
  \label{eq:103}
  V_{nm}(t) = \hat{V}_{nm} f(t),
\end{equation}
where \(f(t)\) is proportional to the voltage applied to the
EOM. Note that \(V_{nm}\) is already expressed in the basis of the
\(c_{m}\) (see \cref{eq:96}).

Let us now assume that the voltage modulation takes the form of
\begin{equation}
  \label{eq:104}
  f(t)=∑_{j}\frac{\hat{f}_{j}}{2} \sin\qty[\pqty{\hat{ω}_{j}+\hat{δ}_{j}}t
  + \varphi_{j}] = -\iu ∑_{j}\hat{f}_{j}\bqty{\eu^{\iu
      \pqty{\hat{ω}_{j}+\hat{δ}_{j}}t}\eu^{\iu \varphi} - \eu^{-\iu
      \pqty{\hat{ω}_{j}+\hat{δ}_{j}}t}\eu^{-\iu \varphi_{j}}}.
\end{equation}

Transforming into the rotating frame of the \({h}_{m}\) (see
\cref{eq:66}), we have
\begin{equation}
  \label{eq:105}
  \tilde{V}_{mn}=-\iu ∑_{j}\hat{V}_{mn}\hat{f}_{i}\bqty{\eu^{i \pqty{\hat{ω_{j}} +\hat{δ}_{j} -
        \pqty{ω^{0}_{n}-ω^{0}_{m} + ε_{m}-ε_{n}}}t} \eu^{\iu \varphi_{j}}
    - \eu^{-i \pqty{\hat{ω_{j}} +
        \hat{δ}_{j}+\pqty{ω^{0}_{n}-ω^{0}_{m} + ε_{m}-ε_{n}}}t} \eu^{-\iu \varphi_{j}}}.
\end{equation}

As discussed in \cref{sec:choice-hybridization}, we want to couple the
\(+\) and \(m\neq \pm\). Comparing with \cref{eq:12}, we set
\begin{align}
  \label{eq:107}
  \hat{ω_{j}}&=ω^{0}_{j} - ω^{0}_{+} = - \pqty{ω^{0}_{-j} - ω^{0}_{-}}& \hat{δ}_{j} = ε_{+}-ε_{j}
\end{align}
where \(0<j\leq N\) and obtain from \cref{eq:105} by neglecting all
counter-rotating terms
\begin{equation}
  \label{eq:108}
  \begin{aligned}
    \tilde{V}_{+,j} &= -\iu \hat{f}_{j} \hat{V}_{+,j} \eu^{i\varphi_{j}}
                  \eu^{\iu \pqty{ε_{+}-ε_{j}}t}
    & \tilde{V}_{j,+} &= \iu \hat{f}_{j} \hat{V}_{j,+}
                    \eu^{-\iu\varphi_{j}} \eu^{\iu
                    \pqty{ε_{j}-ε_{+}}t}\\
    \tilde{V}_{-j,-} &= -\iu \hat{f}_{j} \hat{V}_{-j,-} \eu^{i\varphi_{j}}
                  \eu^{\iu \pqty{ε_{+}-ε_{j}}t}
    & \tilde{V}_{-,-j} &= \iu \hat{f}_{j} \hat{V}_{-,-j}
                    \eu^{-\iu\varphi_{j}} \eu^{\iu
                    \pqty{ε_{j}-ε_{+}}t}.
  \end{aligned}
\end{equation}
Note, that \cref{eq:108} implies \(ε_{-} = -ε_{+}\) and
\(ε_{-j} = -ε_{j}\).

Comparing \cref{eq:108} with \cref{eq:96} and identifying
\(J_{mn}=f(t)\hat{J}_{mn}\) we get
\begin{align}
  \label{eq:113}
  \begin{aligned}
    V^{0}_{+,j} &= -\iu \hat{f}_{j}
                  \frac{\hat{J}_{0,j}}{\sqrt{2}\cos(ψ)}
                  \eu^{-i\pqty{ϕ+ ψ-\varphi_{j}}}
    & V^{0}_{j,+} &= \iu \hat{f}_{j} \frac{\hat{J}_{0,j}}{\sqrt{2}}
                    \eu^{i\pqty{ϕ-\varphi_{j}}}\\
    V^{0}_{-j,-} &= -\iu \hat{f}_{j} \frac{\hat{J}_{-j,0}}{\sqrt{2}}
                    \eu^{i\pqty{ϕ+\varphi_{j}}}
    & V^{0}_{-,-j} &= \iu \hat{f}_{j} \frac{\hat{J}_{0,-j}}{\sqrt{2}\cos(ψ)}
                  \eu^{-i\pqty{ϕ- ψ+\varphi_{j}}}.
  \end{aligned}
\end{align}

\subsection{Steady-state Transmission on the Small Loop}
\label{sec:steadyst-transm}

We can now proceed to calculate the steady state transmission for the
transmission line attached to the small loop (see
\cref{fig:schematic}).

As our input field has only right-moving components and our detector
is situated to the right of the fibre-coupler, we can set \(δ,σ=+\) in
\cref{eq:87} and will suppress those indices in the following.

Using \cref{sec:effects-asymm-damp}, we identify
\begin{equation}
  \label{eq:97}
  \begin{aligned}
    T_{i_{0},β;m} &\to T_{S,0;\pm} &\implies U_{β,m} &\to U_{0,m} =
                                                       T_{S,0;m} =
                                                       \frac{1}{\sqrt{2}}
                                                       \eu^{\iu m
                                                       ψ}\pqty{δ_{m,+}-δ_{m,-}}\\
    && U^{-1}_{m,β}&\to U^{-1}_{m,0} =
                     T^{-1}_{m;S,0}=\frac{1}{\sqrt{2}\cos(ψ)}\pqty{δ_{m,+}-δ_{m,-}}
  \end{aligned}
\end{equation}
which yields using \cref{eq:88}
\begin{equation}
  \label{eq:98}
  η_{m} = \frac{\abs{κ} πn_{T}}{2c\cos(ψ)} \pqty{δ_{m,+}
    \eu^{\iu ψ} + δ_{m,-}\eu^{-\iu ψ}},
\end{equation}
where \(n_{T}\) is the refractive index of the transmission line
(replace \(n_{B}\to n_{T}\) in \cref{sec:micr-deriv}).

To obtain the \(ω_{γ}\), \(λ_{γ}\) and \(O_{m,γ}\) we have to
diagonalize the target Hamiltonian
\begin{equation}
  \label{eq:99}
  H^{T}_{nm}=V^{0}_{n,m} + \pqty{ε_{m} -\iu \pqty{η^{0}_{m} + η_{m}} δ_{nm}}.
\end{equation}
Both apriori loss \(η_{m}^{0}\) and the loss induced by the coupling
to the transmission line \(η_{m}\), as well as the detunings \(ε_{m}\)
enter the final interaction hamiltonian.

As
\begin{equation}
  \label{eq:100}
  \Im η_{\pm}= \pm \frac{\abs{κ} πn_{T}}{2c\cos(ψ)}Δ,
\end{equation}
we obtain a correction to the on-site energies \(ε_{\pm}=\pm ε_{A}\),
whereas
\begin{equation}
  \label{eq:101}
  \Re η_{\pm} = \frac{\abs{κ} πn_{T}}{2c}
\end{equation}
is independent of \(Δ\) which is consistent with the symmetry of
\(η^{0}_{\pm}\). Note however that \cref{eq:100} is of very small
magnitude as it is the product of two quantities that are small
compared to \(Ω_{B}\) and \(1\) respectively.


In light of \cref{eq:113} we may choose\footnote{The right hand sides
  denote quantities from the previous set of notes.}
\begin{gather}
  \label{eq:114}
    ε_{+} = ε_{A}- \Im η_{+}= ε_{A} - \frac{\abs{κ}
      πn_{T}}{2c\cos(ψ)}Δ\\
    \begin{aligned}
    ε_{j} &= ω_{j},  & \varphi_{j} &= 0, & \hat{f}_{i}&= \frac{\sqrt{2}}{\hat{J}_{0,i}}η_{j}
    \end{aligned}
\end{gather}
so that \(V^{0}_{n,m}\) most closely resembles the target
Hamiltonian. The drive phases \(\varphi_{j}\) can also be set to
\(\varphi_{j}=ϕ + \frac{π}{2}\), to remove the phase in the
interaction if it is known. The phase in the interaction does not
influence the observable \(ρ_{A}\). However it does influence the
interference with a reference light beam. Also, the magnitude of
\cref{eq:100} is likely negligible. If it is required, we can
determine it by choosing \(N=1\) and measuring the eigenenergies with
the result obtained below.

To calculate the susceptibility (see \cref{eq:87}), we evaluate
\begin{align}
  \label{eq:116}
  U_{γ} &= ∑_{σ=\pm} \eu^{\iu \tilde{ω}_{σ} t}O^{-1}_{γσ}
  T^{-1}_{σ;S0} =  ∑_{σ=\pm} \frac{σ
          }{\sqrt{2}\cos(ψ)} \eu^{\iu \pqty{ω_{σ}^{0}-ε_{σ}} t}O^{-1}_{γσ}\\
  \bar{U}_{γ}&= ∑_{σ=\pm} \eu^{-\iu \pqty{ω_{σ}^{0}-ε_{σ}}
               t}O_{σγ}T_{S0;σ} = ∑_{σ=\pm} \frac{σ \eu^{\iu σ ψ}}{\sqrt{2}}\eu^{-\iu \pqty{ω_{σ}^{0}-ε_{σ}}t}O_{σγ}
\end{align}
with
\begin{equation}
  \label{eq:117}
  \tilde{ω}^0_{\pm} \equiv ω_{σ}^{0}-ε_{σ} = ω_{S} \pm δ\cos(ψ) - ε_{A}
  + \frac{\abs{κ}
    πn_{T}}{2c\cos(ψ)}Δ.
\end{equation}

Finally we arrive at
\begin{equation}
  \label{eq:118}
  \begin{aligned}
    χ(t,s) &= χ_{0}Θ(s)
    ∑_{γ,σ,σ'}\eu^{{-\iu\bqty{\tilde{ω}^0_{σ}t -
    \tilde{ω}^0_{σ'}s + ω_{γ}(t-s)} -
             λ_{γ}(t-s)}}O_{σ,{γ}}O^{-1}_{γ,σ'}\frac{σ σ' \eu^{\iu σ
             ψ}}{\cos(ψ)}\\
    &= θ(s) χ_{1}(t-s) + χ_{2}(t,s).
  \end{aligned}
\end{equation}
where
\begin{equation}
  \label{eq:119}
  χ_{0}=\frac{\abs{κ}π n_{B}}{2c}.
\end{equation}
The stationary and non-stationary susceptibilities work out to be
\begin{align}
  \label{eq:121}
  χ_{1}(t) &=  χ_{0}
    ∑_{γ,σ}\eu^{{-\bqty{\iu\pqty{\tilde{ω}^0_{σ} + ω_{γ}} +
             λ_{γ}}t}}O_{σ,{γ}}O^{-1}_{γ,σ}\frac{\eu^{\iu σ
             ψ}}{\cos(ψ)}\\
  χ_{2}(t,s) &= -χ_{0} ∑_{γσ} \eu^{-\bqty{\iu\pqty{\tilde{ω}^0_{σ} +
               ω_{γ}} + λ_{γ}}t}  \eu^{\bqty{\iu\pqty{\tilde{ω}^0_{\bar{σ}} +
               ω_{γ}} + λ_{γ}}s} O_{σ,{γ}}O^{-1}_{γ,\bar{σ}},
\end{align}
where \(\bar{σ}=-σ\). As \(V^{0}_{mn}\) decomposes into two blocks,
where modes with \(n=-,-1,-2,\dots,-N\) don't couple to modes with
\(n=+,1,2,\dots,N\) so that we obtain two non-overlapping sets of
eigenstates, we have
\begin{equation}
  \label{eq:122}
  O_{σ,γ}=0 \wedge O^{-1}_{\bar{σ},γ} =0\; \forall σ=\pm\implies   χ_{2}(t,s) = 0
\end{equation}
and the non-stationary contribution to the susceptibility
vanishes. Persistent oscillations in the output intensity therefore
would likely signal a breakdown of the RWA.

We can now proceed to calculate the response of the system to a
coherent input beam with frequency \(ω\) in the limit of
\(t\gg λ_{γ}\)
\begin{equation}
  \label{eq:123}
  ∫_{0}^{t}\eu^{-\iu ω s} χ(t-s)\dd{s} = \eu^{\iu ω t} χ_{0} ∑_{σγ}\frac{O_{σ,{γ}}O^{-1}_{γ,σ}}{\iu\pqty{\tilde{ω}^0_{σ} + ω_{γ}-ω} +
    λ_{γ}}\frac{\eu^{\iu σψ}}{\cos(ψ)}\equiv \eu^{-\iu ω t}F_{S}(ω).
\end{equation}

Both the magnetic and the electric field are proportional to
\(\Im b_{\outputf}\) so that the absolute value of the pointing
vector, e.g. the intensity, averaged over the oscillation period of
the input light becomes
\begin{equation}
  \label{eq:128}
  \bar{I}  = I_{0}\pqty{1-2\Re{F_{S}(ω)} + \abs{F_{S}(ω)}^{2}} \approx
  I_{0}\pqty{1-2\Re{F_{S}(ω)}} ,
\end{equation}
where we have used \(b_{\inputf}=b_{0}\eu^{-\iu ωt}\).

For \(Δ=0\) we have
\begin{equation}
  \label{eq:130}
  \Re{F_{S}(ω)} = χ_{0} ∑_{σγ}\frac{λ_{γ}\abs{O_{σ,γ}}^{2}}{\pqty{\tilde{ω}^0_{σ} + ω_{γ}-ω}^{2} +
    λ_{γ}^{2}},
\end{equation}
whereas \(Δ\neq 0\) will very slightly shift the peak locations and
influence the peak heights. We also see, that we only have a good
signal on the states that have some overlap with the small loop. The
locations of the peaks (dips in the transmission) are the eigenenrgies
\(ω_{γ}\) of the target Hamiltonian \(H^{T}\) relative to the
eigenenrgies of the unmodulated system and the drive detunings
\(\tilde{ω}^0_{σ}=ω^{0}_{σ}-ε_{σ}\). Wheras the \(ω^{0}_{σ}\) define
the ``zero'' of energy, the shifts by \(ε_{σ}\) can be interpreted as
an effect similar to the AC stark shift that arises due to the drive
detuning. Another effect of the drive, namely persistent oscillations
of the output intensity, is not observed here, as we don't couple the
\(σ=\pm\) states with the drive.  Finally, the peak heights and widths are
controlled by the loss rates \(λ_{γ}\).

\subsection{Steady-state Transmission on the Big Loop}
\label{sec:steady-state-transm}

To probe the structure of the states in bath, we have to probe the big
loop.

Analogous to \cref{sec:steadyst-transm} we can also obtain the
transmission for a transmission line (inclusive laser) attached to the
big loop. Note, that the phase \(ϕ\) is now \(ϕ=k_{0}L_{B}/2\),
whereas all the other model parameters retain their meaning.  The main
difference to the calculations for the small loop is the number of
modes that interact with the transmission line, leading to
non-vanishing stationary (low-frequency) oscillations in the output
intensity.

Just as in \cref{sec:steadyst-transm}, we begin by identifying the
relation between the bare modes in the big loop and the eigenmodes of
the unmodulated system.
\begin{equation}
  \label{eq:55}
  \begin{aligned}
    T_{i_{0},β;m} &\to T_{B,β;\pm} &\implies U_{β,m}  &= T_{S,0;m} = δ_{β0}
                                                       \frac{\eu^{\iu ϕ}}{\sqrt{2}}
                                                       ∑_{σ=\pm} δ_{m,σ}
    + ∑_{j\neq 0}δ_{βj}δ_{mj}\\
    && U^{-1}_{m,β} &=
                      T^{-1}_{m;S,0}=δ_{β,0}\frac{\eu^{-\iu ϕ}}{\sqrt{2}\cos(ψ)}
                      ∑_{σ=\pm} δ_{m,σ}\eu^{-\iu σ ψ} +∑_{j\neq 0}δ_{βj}δ_{mj}
  \end{aligned}
\end{equation}
where \(j\in [-N,N] \setminus 0\) and \(β\in [-N,N]\).

The decay rate introduced by the coupling to the transmission line
works out to be
\begin{equation}
  \label{eq:71}
  \begin{aligned}
    η_{\pm} = \frac{\abs{κ} πn_{T}}{2c\cos(ψ)} \eu^{\mp \iu ψ}
  \end{aligned}
\end{equation}
where the sign in the exponent is \emph{inverted} compared to
\cref{eq:98}. Therefore the effective on-site energy for the \(\pm\)
states will be shifted in the inverse direction. Note however, that it
is not guaranteed that \(κ\), \(η_{S}\), and \(Δ\) will be the same as
in \cref{sec:steadyst-transm}.

Using the fact that either \(O_{σ,γ}=0\) or \(O^{-1}_{\bar{σ},γ} =0\)
for any value of \(σ\) we find for the transmission
\begin{equation}
  \label{eq:79}
  T_{B}(ω,t) = θ(s) χ_{0}' ∑_{σ=\pm}\bqty{T^{B}_{σ,σ}(ω) + ∑_{n\neq\pm}\pqty{T^{B}_{σ,n}(ω,t) + T^{B}_{n,σ}(ω,t)} + ∑_{n,m\neq\pm} T^{B}_{σ,m,n}(ω,t)},
\end{equation}
with
\begin{equation}
  \label{eq:91}
  \begin{aligned}
      χ_{0}' &= \abs{κ'}\frac{π n_{T}}{c} & \tilde{ω}_{m} = ω^{0}_{m} - ε_{m}
  \end{aligned}
\end{equation}
and
\begin{align}
  \label{eq:85}
  T^{B}_{σ,σ}(ω) &= ∑_{γ}\frac{O_{σ,{γ}}O^{-1}_{γ,σ}\eu^{\iu σψ}}{2\cos(ψ)}\frac{1}{\iu\pqty{\tilde{ω}_{σ} + ω_{γ}-ω} +
                   λ_{γ}}\\
  \label{eq:92}
  T^{B}_{σ,n}(ω,t) &=δ_{σ,\sgn(n)}\frac{O_{σ,{γ}}O^{-1}_{γ,n}\eu^{\iu ϕ}}{\sqrt{2}}
                                        \frac{\eu^{-\iu(\tilde{ω}_{σ}-\tilde{ω}_{n})t}}{\iu\pqty{\tilde{ω}_{n}+ ω_{γ}-ω} +λ_{γ}}\\
  \label{eq:102}
  T^{B}_{n,σ}(ω,t) &= δ_{σ,\sgn(n)}\frac{O_{n,{γ}}O^{-1}_{γ,σ}\eu^{-\iu (ϕ+ψ)}}{\sqrt{2}}
                                        \frac{\eu^{-\iu(\tilde{ω}_{n}-\tilde{ω}_{σ})t}}{\iu\pqty{\tilde{ω}_{σ}+ ω_{γ}-ω} + λ_{γ}}\\
  \label{eq:106}
  T^{B}_{σ,m,n}(ω,t) &= δ_{σ,\sgn(n)} δ_{σ,\sgn(m)}O_{nγ}O^{-1}_{γm} \frac{\eu^{\iu(\tilde{ω}_{m}-\tilde{ω}_{n})t}}{\iu\pqty{\tilde{ω}_{m}+ ω_{γ}-ω} + λ_{γ}}.
\end{align}

The stationary transmission peaks around \(\tilde{ω}_{\pm}\) and has
subpeaks shifted by \(ω_{γ}\) just as in \cref{sec:steadyst-transm},
where the peak height is roughly proportional to the overlap of the
\(γ\) and \(\pm\) states. In the regime we're interested in, there
will only be one state with substantial overlap with the \(\pm\)
states (a.k.a. the \(A\) site). In the same frequency region,
\cref{eq:102} will also have peaks. Those peaks however will be
suppressed, as their height is proportional to the overlap of the
\(\pm\) states and the \(n\neq \pm\), i.e. the \(A\)-site, and bath
sates with the eigenstate \(γ\). The frequency of the steady state
oscillations of \cref{eq:102} allows to tune the relative energies of
the \(A\) sites and the bath site. The same signal may be retrieved
more cleanly from \cref{eq:92}, where the peaks are situated around
\(\tilde{ω}_{n}\). The transmission component \cref{eq:106} will only
be significant if \(m=n\), as the \(γ\) states that don't overlap with
the \(\pm\) states are almost identical to the \(n\neq\pm\),
i.e. bath, states. In this case the transmission does not exhibit
oscillations making the signal from \cref{eq:92} even clearer.
Comparing \cref{eq:92} and \cref{eq:102}, we can extract the damping
asymmetry \(Δ = \sin(ψ)\).

Time-averaging \cref{eq:79} leaves us with the stationary transmission
\begin{equation}
  \label{eq:110}
  T_{B}(ω) = θ(s) χ_{0}' ∑_{σ=\pm}\bqty{T^{B}_{σ,σ}(ω) + ∑_{n}T^{B}_{σ,n,n}(ω)}.
\end{equation}

\newpage
\printbibliography{}
\end{document}


%%% Local Variables:
%%% mode: latex
%%% TeX-master: t
%%% TeX-output-dir: "output"
%%% TeX-engine: luatex
%%% jinx-languages: "en_CA"
%%% End: