# Light scattering from an atomic gas under conditions of quantum degeneracy

V.M. Porozova  Department of Theoretical Physics, St.-Petersburg State
Polytechnic University, 195251, St.-Petersburg, Russia  L.V. Gerasimov
Faculty of Physics, M.V. Lomonosov Moscow State University, Leninskiye Gory
1-2, 119991, Moscow, Russia  M.D. Havey  Department of Physics, Old Dominion
University, Norfolk, VA 23529, United States  D.V. Kupriyanov
<EMAIL_ADDRESS>Department of Theoretical Physics, St.-Petersburg State
Polytechnic University, 195251, St.-Petersburg, Russia

###### Abstract

We consider a quantum theory of elastic light scattering from a macroscopic
atomic sample existing in the Bose-Einstein condensate (BEC) phase. The
dynamics of the optical excitation induced by an incident photon is influenced
by the presence of incoherent scattering channels. For a sample of sufficient
length the excitation transports as a polariton wave and the propagation
Green’s function obeys the scattering equation which we derive. The polariton
dynamics could be tracked in the outgoing channel of the scattered photon as
we show via numerical solution of the scattering equation for one-dimensional
geometry. The results are analyzed and compared with predictions of the
conventional macroscopic Maxwell theory for light scattering from a non-
degenerate atomic sample of the same density and size.

###### pacs:

42.50.Ct, 42.50.Nn, 42.50.Gy, 34.50.Rk

##  I Introduction

Light scattering from ultracold atomic systems existing under conditions of
quantum degeneracy is a challenging and intriguing issue for both quantum
optics and atomic physics. Together, investigation of these combined fields is
practically important for developing various quantum interface protocols
between light and matter subsystems. Although light scattering from either
degenerate Bose or Fermi gases is of strong interest, we consider in the
current context the degenerate Bose gas only, which is most typical for
alkali-metal systems. The superposed light and matter wave propagating as a
single quantum optical excitation through a Bose-Einstein condensate (BEC)
phase had been predicted in Politzer  even before BEC had been created in the
laboratory. Since the first successful experimental realizations of BEC in
alkali-metal systems reported in Cornell ; Ketterle , evident signatures of
cooperative dynamics in light scattering from the condensate have been
observed in a series of experiments. These include manifestation of
superradiant behavior of Rayleigh scattering in Pritchard ; Schneble ;
Hilliard , formation of superfluid vortexes induced by coherent optical
processes in Cornell99 ; Dalibard00 ; Phillips06  and spin vortexes in
Brachmann11 , and optical control of the BEC phase transition with Faraday
imaging technique in Sherman16 . The strong coherent coupling of light with a
sample led to the condensate fragmentation Pritchard ; Schneble ; Hilliard
and explanation of such a quite non-trivial optomechanical effect has been
attempted in Schneble  in terms of a Kapitza-Dirac diffraction phenomenon.

The above experiments have encouraged development of theoretical insights
towards deeper understanding and precise description of light scattering under
conditions of quantum degeneracy and from BEC in particular. The basic concept
of a master equation for the order parameter suggests a relevant approach
based on time dependent generalization of the non-linear Schrödinger (Gross-
Pitaevskii) equation Wright97 ; WillHoll99 ; Dowling05 . The coherent effects
of conversion of either linear or angular momentum from light to the
condensate are associated with a stimulated Raman process mediating the
dynamics of the order parameter Phillips06 ; Dowling05 . The superadiant
properties of the Rayleigh scattering, observed in a BEC, was explained by
making use of the effective Hamiltonian approach via the mechanism of
cooperative emission induced by a coherent classical pump in MooreMeystre99 ;
ZobayNikol06 ; Trifonov ; Avetisyan13 ; Avetisyan15 .

In the present report we are focusing on a microscopic quantum theory of a
single photon scattering towards ab-initio description of elastic light
scattering from a macroscopic atomic sample existing in the quantum degenerate
BEC phase. Following the second quantized formalism, Bogolubov theory
Bogoliubov47  and Gross-Pitaevskii model Gross ; Pitaevskii  we introduce a
set of coupled and closed diagram equations for the polariton propagator
contributing to the $T$-matrix and scattering amplitude. Our approach allows
us to follow important density corrections to the quasi-energy structure
caused by static interaction and radiation losses associated with the
incoherent scattering. We are aiming to test validity of the conventional
macroscopic Maxwell description for the quantum degenerate gas as well as to
follow possible deviations with light scattering from a non-degenerate atomic
sample of the same density and size.

This paper is organized as follows. In Section II we develop our general
theoretical framework of light scattering from a quantum degenerate atomic
gas. In Section III we derive the basic scattering equation via the Feynman
diagram method (briefly explained in Appendix A) and discuss general
properties of the Green’s function (polariton propagator) responsible for
transporting an optical excitation in a BEC sample. In Section IV we present
the results of our numerical simulations for light scattering in a one-
dimensional geometry; the calculational scheme is detailed in Appendix B. In
Section V we make some concluding remarks.

##  II The scattering problem under conditions of quantum degeneracy

The quantum-posed description of the photon scattering problem is based on the
formalism of the $T$-matrix, which is defined by

$\hat{T}(E)=\hat{V}+\hat{V}\frac{1}{E-\hat{H}}\hat{V}$ (1)

where $\hat{H}$ is the total system Hamiltonian consisting of the unperturbed
part $\hat{H}_{0}$ and an interaction term $\hat{V}$ such that
$\hat{H}=\hat{H}_{0}+\hat{V}$. The energy argument $E$ is an arbitrary complex
parameter in Eq. (1), which approaches a real physical value upon constructing
the scattering amplitude. The scattering process, evolving from an initial
state $|i\rangle$ to the final state $|f\rangle$, is expressed by the
following relation between the differential cross section and the scattering
amplitude, given by the relevant $T$-matrix element considered as a function
of the initial energy $E_{i}$:

${\,d\sigma}_{i\to f}=\frac{{\cal
V}^{2}{\omega^{\prime}}^{2}}{\hbar^{2}c^{4}{(2\pi)}^{2}}{|T_{g^{\prime}\mathbf{e}^{\prime}\mathbf{k}^{\prime};\,g\mathbf{e}\mathbf{k}}(E_{i}+i\,0)|}^{2}\,d\Omega$
(2)

Here the initial state $|i\rangle=|g;\mathbf{e},\mathbf{k}\rangle$ is
specified by the incoming photon s wave vector $\mathbf{k}$, frequency
$\omega\equiv\omega_{k}=ck$, and polarization vector $\mathbf{e}$, and the
atomic system populates a particular collective ground state $|g\rangle$. In
our case $|g\rangle=|\mathrm{BEC}\rangle^{N}$ initially performs a collective
state of $N$ atoms in the BEC phase. The final state
$|f\rangle=|g^{\prime};\mathbf{e}^{\prime},\mathbf{k}^{\prime}\rangle$ is
specified by a similar set of quantum numbers, excepting that
$|g^{\prime}\rangle$ can be now a disturbed condensate state for inelastic
channels, and the solid angle $\Omega$ is directed along the wave vector of
the outgoing photon $\mathbf{k}^{\prime}$. The presence of a quantization
volume $\cal V$ in this expression is caused by the second quantized structure
of the interaction operators. The optical theorem links the total cross
section with the diagonal $T$-matrix element

$\sigma_{\mathrm{tot}}=-\frac{2{\cal V}}{\hbar
c}\,\mathrm{Im}\,T_{ii}(E_{i}+i\,0)$ (3)

which gives a convenient tool for the cross-section evaluation via calculation
of only one $T$-matrix element for the elastic forward scattering.

In the second quantized representation the interaction term $\hat{V}$ in Eq.
(1), taken in the dipole long wavelength approximation ChTnDpRcGr ;
BerstLifshPitvsk ; Kupriyanov17 , is given by

$\hat{V}=-\sum\limits_{n}\int\,d^{3}r\left[d^{\mu}_{nm}\hat{E}_{\mu}({\mathbf{r}})\hat{\Psi}^{{\dagger}}_{n}(\textbf{r})\hat{\Psi}_{m}(\textbf{r})+h.c.\right]$
(4)

where $d^{\mu}_{nm}$ is the matrix element of the $\mu$-th vector component of
an atomic dipole moment, where $n$ and $m$ respectively specify the excited
and ground states of the atom. $\hat{E}_{\mu}({\mathbf{r}})$ is the $\mu$-th
vector component of the electric field operator and for sake of generality we
use co/contravariant notation for the vector and tensor indices. The operators
$\hat{\Psi}_{m}(\textbf{r})$ and $\hat{\Psi}_{n}^{\dagger}(\textbf{r})$ are
the second quantized annihilation and creation operators of an atom at
position $\mathbf{r}$ respectively in the ground and excited states. We will
further consider a BEC consisting of the simplest two-level atoms with a
${}^{1}S_{0}$ ground state and ${}^{1}P_{1}$ excited state such that quantum
numbers $n=0,\pm 1$ and $m=0$ respectively denote the single atom angular
momentum projection of the excited and the ground states.

In accordance with the general concept of quantum degeneracy for the system
ground state existing in the BEC phase at zero-temperature, see
BerstLifshPitvsk , we accept

$\hat{\Psi}_{0}(\textbf{r}){|\mathrm{BEC}\rangle}^{N}=\Xi(\textbf{r}){|\mathrm{BEC}\rangle}^{N-1}$
(5)

where $\Xi(\textbf{r})$ is the order parameter (often termed the
”wavefunction”) of the condensate. We consider the BEC as a macroscopic object
such that the order parameter is insensitive to any small variation of the
number of particles in the condensate. Then the scattering amplitude,
expressed by ”on-shell” $T$-matrix elements contributing to Eqs. (2) and (3)
for the scattering of an incident photon of frequency $\omega$ to the outgoing
photon of frequency $\omega^{\prime}$, is given by

$\displaystyle T_{fi}(E)=\frac{2\pi\hbar(\omega^{\prime}\omega)^{1/2}}{{\cal
V}}\iint\,d^{3}r^{\prime}\,d^{3}r\sum\limits_{n^{\prime},n}$
$\displaystyle\times$
$\displaystyle{(\mathbf{d}\cdot\mathbf{e}^{\prime})}^{*}_{n^{\prime}0}{(\mathbf{d}\cdot\mathbf{e})}_{0n}\,{\mathrm{e}}^{-i\mathbf{k}^{\prime}\mathbf{r}^{\prime}+i\mathbf{k}\mathbf{r}}\
\Xi^{*}(\mathbf{r}^{\prime})\,\Xi(\mathbf{r})$ $\displaystyle\times$
$\displaystyle\left(-\frac{i}{\hbar}\right)\int_{0}^{\infty}dt\;{\mathbf{e}}^{\frac{i}{\hbar}(E-E_{0}^{N-1}+i0)t}\;i\,G_{n^{\prime}n}(\mathbf{r}^{\prime},t;\mathbf{r},0)$

where $E_{0}^{N-1}$ is the initial energy of the condensate consisting of
$N-1$ particles. The internal dynamics of the scattering process is described
by a single optical excitation evolving in the condensate

$iG_{n^{\prime}n}(\mathbf{r}^{\prime},t^{\prime};\mathbf{r},t)=\langle\mathrm{BEC}|T\Psi_{n^{\prime}}(\mathbf{r}^{\prime};t^{\prime})\Psi_{n}^{{\dagger}}(\mathbf{r};t)|\mathrm{BEC}\rangle^{N-1}$
(7)

with projection onto the product of condensate and field vacuum states such
that entirely

$|\mathrm{BEC}\rangle^{N-1}\equiv|\mathrm{BEC}\rangle^{N-1}_{\mathrm{Atoms}}\times|0\rangle_{\mathrm{Field}}$
(8)

Eq. (7) defines the time ordered (causal) Green’s function (propagator)
associated with the polariton-type quasi-particle excitation superposed
between field and atom and propagating through the condensate consisting of
$N-1$ particles. The operators contributing to the polariton propagator are
the original atomic operators transformed in the Heisenberg representation and
dressed by the interaction process. In the matrix element of the $T$-matrix in
the form (LABEL:2.6) the outer operators $\hat{V}$ in its basic definition (1)
are disclosed in the rotating wave approximation (RWA). Such an assumption is
surely valid as far as we are interested in near resonant scattering when both
the frequencies $\omega$ and $\omega^{\prime}$ are close to the frequency of
the atomic transition $\omega_{0}$.

The Green’s function (7), rewritten in the interaction representation, can be
expanded in the perturbation theory series

$\displaystyle iG_{n^{\prime}n}(\mathbf{r}^{\prime},t^{\prime};\mathbf{r},t)$
(9)
$\displaystyle=\langle\mathrm{BEC}|\hat{S}^{-1}T\left[\Psi^{(0)}_{n^{\prime}}(\mathbf{r}^{\prime};t^{\prime})\Psi^{(0)\dagger}_{n}(\mathbf{r};t)\hat{S}\right]|\mathrm{BEC}\rangle^{N-1}$
$\displaystyle=\langle\mathrm{BEC}|T\left[\Psi^{(0)}_{n^{\prime}}(\mathbf{r}^{\prime};t^{\prime})\Psi^{(0)\dagger}_{n}(\mathbf{r};t)\hat{S}\right]|\mathrm{BEC}\rangle^{N-1}$

where in the interaction representation the $\Psi$-operators are superscripted
by zero-indices. We consider the condensate itself as a stable system, which
should not be modified by the interaction (4) without its advanced
perturbation by an incoming photon. This should be justified by the
requirement that the evolution operator

$\hat{S}=T\exp\left[-\frac{i}{\hbar}\int_{-\infty}^{\infty}\hat{V}^{(0)}(t)\,\mathrm{e}^{-0\cdot|t|}\right]$
(10)

does not change the BEC state such that
$\hat{S}|\mathrm{BEC}\rangle^{N-1}=|\mathrm{BEC}\rangle^{N-1}$. Although this
requirement seems as evidently accepted in assumptions of the RWA let us make
an important remark concerning its applicability.

The condensate, considered as a physical object, is not an ideal gas. The
small but physically important difference
$E_{0}^{N}-E_{0}^{N-1}=\varepsilon_{0}\equiv\mu_{c}+E_{0}$ gives a binding
energy for adding a particle into an atomic ensemble, which incorporates the
chemical potential $\mu_{c}$ and the internal ground state energy $E_{0}$ of a
single atom. The latter could be set as zero but in our derivation it is
convenient to leave $E_{0}$ as a physical parameter. For the quantum
degenerate gas, consisting of not extremely dense and weakly interacting atoms
and fairly described in a framework of the Gross-Pitaevskii model Gross ;
Pitaevskii , the following inequality is fulfilled

$\mu_{c}\lesssim\frac{\hbar^{2}k_{0}^{2}}{2{\mathrm{m_{A}}}}\ll\hbar\gamma$
(11)

where $k_{0}\equiv\lambdabar_{0}^{-1}$ is wave number for a resonant photon,
${\mathrm{m_{A}}}$ is the atomic mass and $\gamma$ is the natural spontaneous
decay rate for the upper state of the atom. In accordance with the model, see
BerstLifshPitvsk , the chemical potential for a homogeneous BEC is given by

$\mu_{c}=n_{0}\\!\int\\!d^{3}r\;U(r)>0$ (12)

where $U(r)$ is an interaction potential in the system of two atoms and
$n_{0}$ is the atomic density. The subtle point is that the interaction $U(r)$
incorporates both the short range repulsive part and the long range attractive
dipole-dipole polarization interactions. The latter is also known as the Van-
der-Waals interaction and the related asymptotic behavior of the potential
$U(r)$ is constructed in the second order of the same Hamiltonian (4) but with
keeping the terms beyond and alternative to the RWA concept.

The conflicting situation with double accounting of the interaction
Hamiltonian (4) can be resolved once we pay attention that the Van-der-Waals
interaction is meaningful on a distance of an atomic scale $r\sim O(1)a_{0}$,
where $a_{0}$ is the Bohr radius, but the optical coupling experiences the
distances $r\sim\lambdabar_{0}\gg a_{0}$ . That means that there is no
intersection in the diagram representation of $U(r)$ with those, which are
induced by the evolution operator (10), and which couples a pair of distant
atoms where one is always excited. In this case the evolution operator(9)
indeed does not affect the condensate state and the second line in Eq. (9) is
valid beyond the restrictions of the RWA approach as far as the internal
interaction in the atomic ensemble is weak and can be safely separated from
the optical excitation dynamics mediated by the scattering process.

Inequality (11) provides us the chemical potential as the smallest parameter
of the theory and is fulfilled up to the densities
$n_{0}\lambdabar_{0}^{3}\gtrsim 1$. This is a typical condition with
considering a condensate consisting of alkali-metal atoms. From the physical
point of view that means that we consider the BEC in conditions close to an
ideal gas and assume the matrix elements in (4), as well as the atomic energy
structure in the perturbation theory expansion, the same as for independent
atoms. Nevertheless we do not ignore the gas non-ideality and the interatomic
interaction $U(r)$ in the ground state as far as it is crucially important for
proper description of the general behavior of the order parameter
$\Xi=\Xi(\mathbf{r},t)$ under the framework of the Gross-Pitaevskii model with
including superfluidity as the main macroscopic quantum property of the
condensate. In its main approximations our consideration is applicable up to
the bound of $\mu_{c}\lesssim\hbar\gamma$

##  III Dynamics of the optical excitation in the condensate

###  III.1 Diagrammatic representation

The polariton propagator (9) can be expanded in the perturbation theory series
and the appearing terms can be regrouped with the Feynman diagram method. The
basic elements and definitions are listed in Appendix A. As far as the
considered interaction processes are primary developing in a near resonance
conditions we follow the RWA approach with keeping leading expansion terms.
Eventually the polariton propagator can be constructed as a dressed Green’s
function of an excited atom and obeys the following Dyson-type diagram
equation

(13)

where it is visualized as a doubly-straight line. This corresponds to the fact
that the original atomic propagator is assumed as ”dressed” here by all the
interaction processes. The inward and outward vertical arrows image the order
parameters and form the self-energy part responsible for coherent conversion
of the excitation between the field, in which free dynamics is expressed by an
undressed wavy line, and an atom subsequently recovered in the condensate
phase. However, as consistent with this diagram equation, the above coherent
process partly degrades because of interaction with the vacuum modes when the
excited atom emits a photon spontaneously and escapes coherent dynamics with
further drifting through the condensate as a spectator.

The latter process contributes in (13) by an incomplete polariton propagator,
which is imaged by a straight solid line in the diagrams and obeys the
following Dyson-type equation

(14)

which should be considered together with the equation for the dressed field
propagator

(15)

These two diagram equations are closed with respect to each other and
reproduce the self-consistent dynamics of an atomic dipole interacting with
its environment similar to the conditions in a disordered atomic gas. Indeed
any optical excitation created from the condensate has a chance to be
incoherently re-emitted into the vacuum modes and transfer the atom, emitting
such photon, out of the condensate phase. That is just described by the
incomplete polariton propagator (14) having a similar diagrammatic
representation as of atomic excitation in a disordered gas. Such an incoherent
scattering induces losses and leads to degradation of coherent dynamics
supported by the self-energy operator in Eq. (13). In a natural optical
association this process introduces the dielectric permittivity constructed
similarly as in a disordered atomic gas of the same density.

###  III.2 Incoherent losses and the dielectric permittivity of the condensate

The self-energy part in (15) (polarization operator) emphasizes the coherent
structure of the matter state considered in conditions of quantum degeneracy.
Nevertheless for an infinite, and locally homogeneous and isotropic medium,
which physically requires that the sample size as well as inhomogeneity scale
of the order parameter $\Xi(\mathbf{r})$ would be comparable or longer than
the radiation wavelength, the solution of Eq. (15) is expected to be similar
to the case of a disordered atomic gas of the same density. Indeed, both the
vertices in the self-energy part of Eq. (15) are linked by the propagator (14)
in which the respective resonant excitation degrades on a time scale of
natural decay when the excited atom can drift the distance much less than its
radiation wavelength. Thus both the vertices are taken in proximal spatial
points such that order parameter actually contributes to Eq. (14) as the local
atomic density $n_{0}(\mathbf{r})=|\Xi(\mathbf{r})|^{2}$. With this
simplification we can construct solution of equations (14) and (15) as for
infinite, homogeneous and isotropic medium in closed analytical form and
compare the result with similar performance of incoherent scattering process
developing in disordered atomic gas.

####  III.2.1 Analytical performance

For sake of convenience and for further derivation we switch the primed and
unprimed arguments and indices in notations of the Green’s functions, see as
an example Eq. (9), accordingly our definitions of Appendix A for the
undressed functions. In stationary and homogeneous conditions these functions
depend only on the difference between their spatial and time arguments. Then
we can make a Fourier transform for the ”dressed” Green’s functions,
constructed by the diagram method, and define

$\displaystyle{\cal D}^{(E)}_{\mu\mu^{\prime}}(\mathbf{k},\omega)$ (18)
$\displaystyle=\\!\int\\!d^{3}\\!R\int^{\infty}_{-\infty}\\!\\!d\tau\,\mathrm{e}^{i\omega\tau-i\mathbf{k}\cdot\mathbf{R}}\\!\left.{\cal
D}_{\mu\mu^{\prime}}^{(E)}(\mathbf{R},\tau)\right|_{\scriptsize{\begin{array}[]{c}\mathbf{R}=\mathbf{r}-\mathbf{r}^{\prime}\\\
\tau=t-t^{\prime}\end{array}}}$

for the photon propagator, fulfilling equation (15), and

$\displaystyle G^{(\gamma)}_{nn^{\prime}}(\mathbf{p},E)$ (21)
$\displaystyle=\\!\int\\!d^{3}\\!R\int^{\infty}_{-\infty}\\!\\!d\tau\,\mathrm{e}^{\frac{i}{\hbar}E\tau-\frac{i}{\hbar}\mathbf{p}\cdot\mathbf{R}}\\!\left.G_{nn^{\prime}}^{(\gamma)}(\mathbf{R},\tau)\right|_{\scriptsize{\begin{array}[]{c}\mathbf{R}=\mathbf{r}-\mathbf{r}^{\prime}\\\
\tau=t-t^{\prime}\end{array}}}$

for the incomplete polariton propagator, fulfilling equation (14). The
superscript $\gamma$ is added for associating such a propagator with
excitation dynamics mediated by spontaneous scattering processes. In
representation (18) we assume the ”dressed” positive frequency component of
the vacuum Green’s function (51) with $\omega>0$ and the equivalence between
the causal and retarded-type definitions for this case.

In the Fourier representation equation (15) can be straightforwardly resolved
with respect to the incomplete polariton propagator

$\displaystyle{\cal D}^{(E)}_{\mu\mu^{\prime}}(\mathbf{k},\omega)$
$\displaystyle=$
$\displaystyle-\frac{4\pi\hbar\omega^{2}}{\omega^{2}_{k}-\epsilon(k,\omega)\omega^{2}}\left[\delta_{\mu\mu^{\prime}}-c^{2}\,\frac{k_{\mu^{\prime}}k_{\mu}}{\epsilon(k,\omega)\omega^{2}}\right]$
(22) $\displaystyle\approx$
$\displaystyle-\frac{4\pi\hbar\omega^{2}}{\omega^{2}_{k}-\epsilon(\omega)\omega^{2}}\left[\delta_{\mu\mu^{\prime}}-c^{2}\,\frac{k_{\mu^{\prime}}k_{\mu}}{\epsilon(\omega)\omega^{2}}\right]$

where

$\epsilon(k,\omega)=1-\frac{4\pi}{\hbar}\,d^{2}_{0}\,n_{0}\,G^{(\gamma)}(\hbar
k,\hbar\omega+\varepsilon_{0})$ (23)

Here $d_{0}$ is the modulus of the transition dipole moment (the same for all
the transitions), $n_{0}=|\Xi|^{2}=\mathrm{const}$ is the density of atoms,
and for an isotropic medium with degenerate excited state
($E_{n}=\mathrm{const}_{n}$) we have

$G^{(\gamma)}_{nn^{\prime}}(\mathbf{p},E)\;=\;\delta_{nn^{\prime}}\,G^{(\gamma)}(p,E)$
(24)

With taking into account inequality (11) we expect negligible deviation in
(23) from the limit of immobile atoms and approximate
$\epsilon(k,\omega)\approx\epsilon(0,\omega)\equiv\epsilon(\omega)$, which
justifies the second line in Eq. (22). Equation (22) (as well as similar
tensor relations found later in the paper) is performed for Cartesian
components $\mu,\mu^{\prime}=x,y,z$, but for the case of spherical components
$\mu,\mu^{\prime}=0,\pm 1$ one has to change $\delta_{\mu\mu^{\prime}}\to
g_{\mu\mu^{\prime}}=(-)^{\mu}\delta_{\mu,-\mu^{\prime}}$.

The obtained result looks similar to that of a conventional medium beyond
quantum degeneracy. As we can see, with reference to BerstLifshPitvsk , such a
type of ”photon Green’s function in a medium” can be associated with a
fundamental solution of the macroscopic Maxwell equations where
$\epsilon(\omega)$ is the dielectric permittivity of the medium. However in
the case of quantum degeneracy both the excitations in the field and matter
subsystems, i.e. photon and excited atom, transport through the sample in a
superposed polariton mode, as suggested by the complete graph equation (13).
Although the association with a conventional medium is not intrinsically
consistent we shall call $\epsilon(\omega)$ as a dielectric permittivity of
the condensate with having in mind in such analogy that it is constructed with
involving only the contribution of over-condensate excitations created in the
incoherent scattering process.

Equation (14), decoded in the Fourier representation, contains the field
Green’s function (22) contributing to the self-energy part in the form of the
convolution integral with atomic propagator, see clarifying comment in
Appendix A. As far as recovering of the incoherent losses as well as
interaction with the quantized continuum are mostly important for near
resonant conditions, we can expect that in the integral evaluation, the
internal arguments are varied in sufficiently broad domains but located near
$\omega\sim\omega_{0}$ and $k\sim k_{0}=\omega_{0}/c$, where
$\omega_{0}=(E_{n}-E_{0})/\hbar$ is the atomic transition frequency.
Considering the field Green’s function as an analytical function of detuning
$\Delta=\omega-\omega_{0}$ in the complex half-plane where
$\mathrm{Im}[\Delta]>0$ the integral over $\omega$ (approximated as integral
over $\Delta$ in infinite limits) can be reliably reproduced by the residue at
the pole point $\omega_{E}=(E-E_{0})/\hbar$ (where
$\Delta\to\Delta_{E}=(E-E_{n})/\hbar$). In such an estimate we can safely
ignore the small pole displacement associated with the Doppler shift as a
negligible relativistic-type correction to the remaining integral evaluated
over $\mathbf{k}$-variable.

In these assumptions, equation (14) reads

$\left[E-\frac{p^{2}}{2\mathrm{m_{A}}}-E_{n}-\Sigma^{(\gamma)}(p,E)\right]\,G^{(\gamma)}(p,E)=\hbar$
(25)

and the self-energy part $\Sigma^{(\gamma)}(p,E)$ is expressed by the sum

$\Sigma^{(\gamma)}(p,E)=\Sigma^{(\mathrm{st})}(p,E)+\Sigma^{(\mathrm{rad})}(p,E)$
(26)

where the first term is given by

$\displaystyle\Sigma^{(\mathrm{st})}(p,E)$ $\displaystyle=$
$\displaystyle\frac{4\pi}{3}\int\frac{d^{3}k}{(2\pi)^{3}}\frac{d^{2}_{0}}{\epsilon\left(k,\omega_{E}\right)}$
(27) $\displaystyle\approx$
$\displaystyle\frac{4\pi}{3}\int\frac{d^{3}k}{(2\pi)^{3}}\frac{d^{2}_{0}}{\epsilon\left(\omega_{E}\right)}$

and can be associated with the interaction of the dipole with its own field in
the environment of the over-condensate medium, created in the incoherent
excitation process, see our comment above. The second term is given by

$\displaystyle\Sigma^{(\mathrm{rad})}(p,E)=-\frac{8\pi}{3}d^{2}_{0}\int\frac{d^{3}k}{(2\pi)^{3}}\frac{\omega_{E}^{2}}{c^{2}k^{2}-\epsilon\left(k,\omega_{E}\right)\omega_{E}^{2}}$
(28)
$\displaystyle\approx-\frac{8\pi}{3}d^{2}_{0}\int\frac{d^{3}k}{(2\pi)^{3}}\frac{\omega_{E}^{2}}{c^{2}k^{2}-\epsilon\left(\omega_{E}\right)\omega_{E}^{2}}$

and reveals radiation back action of the incoherent emission on the dipole’s
dynamics.

Equations (25)-(28) and (23) entirely construct one closed but quite
complicated self-consistent equation for the incomplete propagator
$G^{(\gamma)}(p,E)$, which has nonlinear and integral form. However the
equation can be essentially simplified with applying faithful approximation,
expressed by the second lines in Eqs. (27) and (28), which assumes that in Eq.
(25) the kinetic energy term for $p\sim\hbar k_{0}$ is small in comparison
with the self-energy part. As we have pointed out above, this is justified by
inequality (11). In this approximation the dielectric permittivity
$\epsilon(\omega)$ as well as the function $G^{(\gamma)}(p,E)$ (with
$p\sim\hbar k_{0}$) can be found in analytical form once we resolve the
problem with divergencies existing in both the contributions to the self-
energy part (26).

####  III.2.2 Renormalization of the self-energy divergences

Let us express contribution (27) in the following form

$\Sigma^{(\mathrm{st})}(p,E)\sim-\mathbf{d}\cdot\mathbf{E^{(\mathrm{vac})}}(\mathbf{0})-\mathbf{d}\left[\mathbf{E^{(\mathrm{med})}}(\mathbf{0})-\mathbf{E^{(\mathrm{vac})}}(\mathbf{0})\right]$
(29)

where we assumed that an atomic dipole $\mathbf{d}$ is located at the origin
of the coordinate frame and the diverging integral (27) was converted to the
dipole’s infinite electric field $\mathbf{E^{(\mathrm{med})}}(\mathbf{0})$ in
the medium with dielectric constant $\epsilon$. We also subtracted and added
the same quantity existing in vacuum with $\epsilon=1$. The vacuum term means
the dipole self-action i.e. an artificial object of the theory, which reveals
incorrectness of the dipole gauge on the distances comparable with atomic
scale. The infinite energy, associated with this term, should be incorporated
into the physical energy of the excited atom as internal energy of the point-
like dipole particle. Then the second term in Eq. (29) is a physical quantity
showing how the dipole self-action is modified in the environment of other
dipoles. One expects that the incoherent scattering is a locally cooperative
process and the selected dipole is indistinguishable from other proximal
dipoles responding the driving field of an exciting photon. Then in accordance
with the arguments performed in Refs. Kupriyanov17 ; OurPRA2009 ; Javanainen97
we can accept the standard Lorentz-Lorenz interpretation of the field and
energy shift, associated with static interaction of a collection of proximal
dipoles

$\displaystyle\mathbf{E^{(\mathrm{med})}}(\mathbf{0})-\mathbf{E^{(\mathrm{vac})}}(\mathbf{0})$
$\displaystyle\to$ $\displaystyle\frac{4\pi}{3}n_{0}\mathbf{d}$
$\displaystyle\Sigma^{(\mathrm{st})}(p,E)\to-\frac{4\pi}{3}n_{0}d_{0}^{2}$
(30)

where we substituted $\mathbf{d}^{2}\to d_{0}^{2}$.

The second contribution (28) can be interpreted as an interaction with the
quantized vacuum continuum manifestable via the radiation Lamb-shift and
spontaneous decay rate. With taking $\epsilon=1$ the integral transforms to

$\displaystyle\Sigma^{(\mathrm{vac})}(p,E)=-\left.\frac{8\pi}{3}d^{2}_{0}\int\frac{d^{3}k}{(2\pi)^{3}}\frac{\omega_{E}^{2}}{c^{2}k^{2}-\omega_{E}^{2}-i0}\right|_{E\sim
E_{n}}$ (31) $\displaystyle\Rightarrow\hbar\Delta_{L}-\frac{i\hbar}{2}\gamma$

where $\Delta_{L}\to\infty$ is the vacuum Lamb-shift further renormalized and
incorporated into the atomic energy $E_{n}$, ”dressed” by interaction with the
vacuum modes. The regularized integral given by the difference of (28) and
(31) becomes converging and reproducible by residues at its pole points.
Eventually we arrive at the following renormalization of the radiation
correction to the self-energy part

$\Sigma^{(\mathrm{rad})}(p,E)\Rightarrow-\frac{i\hbar}{2}\sqrt{\epsilon\left(\omega_{E}\right)}\gamma$
(32)

which contains both the radiation damping and energy shift modified by the
radiation coupling with the over-condensate environment.

####  III.2.3 The incomplete propagator in closed form

With substituting renormalized self-energy parts (30) and (32) into (26) and
(25) and in accordance with our definition of the dielectric permittivity
given by Eq. (23) with $k\to 0$ we obtain the following equation

$\epsilon\left(\omega-\frac{\mu_{c}}{\hbar}\right)=\frac{\omega-\omega_{0}-\displaystyle\frac{8\pi}{3\hbar}n_{0}d_{0}^{2}+\frac{i}{2}\sqrt{\epsilon(\omega)}\gamma}{\omega-\omega_{0}+\displaystyle\frac{4\pi}{3\hbar}n_{0}d_{0}^{2}+\frac{i}{2}\sqrt{\epsilon(\omega)}\gamma}$
(33)

In the left-hand side the frequency argument of the permittivity is displaced
by the chemical potential $\mu_{c}$. This emphasizes the fact that for a
single optical excitation from the condensate the extra action is needed,
which is a meaningful part of binding energy $\varepsilon_{0}=E_{0}+\mu_{c}$
given by the chemical potential. Although in our model this displacement is
rather small it recognizes a qualitatively important extension up to the case
of a strongly non-ideal gas. But with neglecting it, we obtain an equation for
the dielectric permittivity identical to an atomic ensemble consisting of cold
disordered and randomly distributed atomic dipoles, see Ref. OurPRA2009 .

Equation (33) can be analytically solved and its solution can be applicable
for the case of a inhomogeneous medium if the density
$n_{0}=n_{0}(\mathbf{r})$ and order parameter $\Xi(\mathbf{r})$ are varied on
a spatial scale comparable with the radiation wavelength or longer. Then
equation (25) suggests the following approximate form in the mixed space-
frequency representation

$\displaystyle\left[E+\frac{\hbar^{2}}{2\mathrm{m_{A}}}\triangle-
E_{n}+\frac{4\pi}{3}\,n_{0}(\mathbf{r})\,d_{0}^{2}\right.$ (34)
$\displaystyle\left.+\frac{i\hbar}{2}\sqrt{\epsilon\left(\mathbf{r},\omega_{E}\right)}\;\gamma\right]G^{(\gamma)}(\mathbf{r},\mathbf{r}^{\prime};E)=\hbar\,\delta(\mathbf{r}-\mathbf{r}^{\prime})$

where we parametrized the dielectric constant
$\epsilon=\epsilon(\mathbf{r},\omega)$ by its spatial dependence. Indeed, in
this equation $G^{(\gamma)}(\mathbf{r},\mathbf{r}^{\prime};E)$, considered as
a function of $\mathbf{r}-\mathbf{r}^{\prime}$, transports a single photon
excitation, created from the immobile condensate, from point
$\mathbf{r}^{\prime}$ to point $\mathbf{r}$, which degrades on a spatial scale
sufficiently less than $\lambdabar_{0}=k_{0}^{-1}$. Thus equation (34) accepts
only proximal spatial arguments
$\mathbf{r}\sim\mathbf{r}^{\prime}\sim(\mathbf{r}+\mathbf{r}^{\prime})/2$
where $n_{0}=n_{0}(\mathbf{r})$ is approximately constant.

We have constructed the incomplete polariton propagator (14) in the form,
which is similar to the complete excited state propagator of a single atom in
a disordered atomic gas of the same density. Such an analogy, emphasizing the
similarity in spontaneous scattering from both the systems, was expectable and
prefaced this part of our derivation. Nevertheless, as was pointed out above,
the analogy is not so straightforward and in the conditions beyond the Gross-
Pitaevskii model (i.e. for a non-ideal quantum gas with strong internal
coupling) it could appear important deviations in description of such
physically different systems.

###  III.3 The complete polariton propagator

With decoding the diagram equation (13) for the complete polariton propagator
we extend spontaneous dynamics, described by Eq. (34), with involving the
process of coherent conversion of the excitation between field and condensate

$\displaystyle\left[E+\frac{\hbar^{2}}{2\mathrm{m_{A}}}\triangle-
E_{n}+\frac{4\pi}{3}\,n_{0}(\mathbf{r})\,d_{0}^{2}\right.$ (35)
$\displaystyle\left.+\frac{i\hbar}{2}\sqrt{\epsilon\left(\mathbf{r},\omega_{E}\right)}\;\gamma\right]G_{nn^{\prime}}(\mathbf{r},\mathbf{r}^{\prime};E)$
$\displaystyle-\sum_{n^{\prime\prime}}\int\\!\\!d^{3}\\!r^{\prime\prime}\;\Sigma_{nn^{\prime\prime}}^{(\mathrm{c})}(\mathbf{r},\mathbf{r}^{\prime\prime};\\!E)\;G_{n^{\prime\prime}n^{\prime}}(\mathbf{r}^{\prime\prime}\\!,\mathbf{r}^{\prime};\\!E)$
$\displaystyle=\hbar\,\delta_{nn^{\prime}}\delta(\mathbf{r}\\!-\\!\mathbf{r}^{\prime})$

The kernel of the respective integral self-energy operator (with simplifying
argument superscripted from double prime to single prime) is given by

$\displaystyle\Sigma_{nn^{\prime}}^{(\mathrm{c})}(\mathbf{r},\mathbf{r}^{\prime};E)$
$\displaystyle=$
$\displaystyle\frac{1}{\hbar}\sum_{\mu\mu^{\prime}}\,\Xi(\mathbf{r})\,\Xi^{\ast}(\mathbf{r}^{\prime})\;d^{\mu}_{n0}d^{\mu^{\prime}}_{0n^{\prime}}$
(36) $\displaystyle\times$ $\displaystyle
D_{\mu\mu^{\prime}}^{(E)}\left(\mathbf{r}-\mathbf{r}^{\prime},\omega_{E}-\frac{\mu_{c}}{\hbar}\right)$

where the vacuum field Green’s function, expressed by the wavy line in the
diagram equation (13) and defined by Eqs. (45) and (51), contributes here in
the mixed space-frequency representation

$\displaystyle
D^{(E)}_{\mu\mu^{\prime}}(\mathbf{R};\omega)\\!\\!=\\!\\!-i\\!\int^{\infty}_{-\infty}\\!\\!d\tau\,\mathrm{e}^{i\omega\tau}\\!\left.\langle
TE^{(0)}_{\mu}(\mathbf{r},t)\,E^{(0)}_{\mu^{\prime}}(\mathbf{r}^{\prime},t^{\prime})\rangle\right|_{\scriptsize\begin{array}[]{c}\tau\\!=\\!t\\!-\\!t^{\prime}\\\
\mathbf{R}\\!=\\!\mathbf{r}\\!-\\!\mathbf{r}^{\prime}\end{array}}$
$\displaystyle=$
$\displaystyle-\hbar\frac{|\omega|^{3}}{c^{3}}\left\\{i\frac{2}{3}h^{(1)}_{0}\left(\frac{|\omega|}{c}R\right)\delta_{\mu\nu}\right.$
$\displaystyle+$
$\displaystyle\left.\left[\frac{X_{\mu}X_{\mu^{\prime}}}{R^{2}}-\frac{1}{3}\delta_{\mu\mu^{\prime}}\right]ih^{(1)}_{2}\left(\frac{|\omega|}{c}R\right)\right\\}$
(40)

Here the averaging is over the vacuum state and $h^{(1)}_{L}(\ldots)$ with
$L=0,2$ are the spherical Hankel functions of the first kind.

The derived equation (35) traces the dynamics of a single particle excitation
in the condensate with the assumption that the order parameter, density
distribution, dielectric permittivity, etc. have a smooth profile on a
mesoscopic scale, similar to the conventional macroscopic Maxwell theory. It
visualizes as a Schrödinger-type equation for an excited atom propagating in
space and modified by interacting with the environment. Here the kinetic
energy term is actually responsible for negligible drift of the excitation
during the decay time when the transferred momentum of the polariton is
limited by the value of $\hbar k_{0}$ in its order of magnitude. Nevertheless
the optical excitation itself can propagate through the sample with much
faster speed with approaching to speed of light, which can be demonstrated via
solution of equation (35) in the limit of infinite and homogeneous medium.

For an infinite, homogeneous and isotropic medium the solution of Eq. (35) can
be found in the reciprocal space as a linear combination of the transverse and
the longitudinal components with respect to the momentum argument. In our
further estimates in this subsection, with considering the internal binding
energy in the condensate as weak, we will ignore the chemical potential of the
condensate as a negligible quantity in comparison with the basic spectral
parameters such as spontaneous decay rate and recoil energy, see (11). Then
the Fourier components of the complete polariton propagator can be expanded as
follows

$\displaystyle G_{nn^{\prime}}(\mathbf{p},E)$ $\displaystyle=$ $\displaystyle
G_{\parallel}(p,E)\,\frac{p_{n}p_{n^{\prime}}}{p^{2}}$ (41) $\displaystyle+\
G_{\perp}(p,E)\\!\left[\delta_{nn^{\prime}}-\frac{p_{n}p_{n^{\prime}}}{p^{2}}\right]$

where, in accordance with the selection rules for the dipole moment operators
in Eq. (36), we associated the vector indices (in Cartesian frame) in the
quasi-particle momentum $\mathbf{p}$ with the quantum numbers of the atomic
excited state.

The longitudinal and transverse components of the polariton propagator are
respectively given by

$\displaystyle G_{\parallel}(p,E)\\!$ $\displaystyle=$
$\displaystyle\\!\hbar\left[E-E_{n}-\frac{p^{2}}{2\mathrm{m_{A}}}-\frac{8\pi}{3}n_{0}\,d_{0}^{2}\right.$
$\displaystyle\hskip 34.14322pt\left.+\
\frac{i\hbar}{2}\sqrt{\epsilon(\omega_{E})}\;\gamma\right]^{-1}$
$\displaystyle G_{\perp}(p,E)\\!$ $\displaystyle=$
$\displaystyle\\!\hbar\left[E-E_{n}-\frac{p^{2}}{2\mathrm{m_{A}}}+\frac{4\pi}{3}\,n_{0}\,d_{0}^{2}\right.$
(42)
$\displaystyle\left.+\frac{i\hbar}{2}\sqrt{\epsilon(\omega_{E})}\;\gamma-\;\frac{4\pi\,n_{0}\,d_{0}^{2}\,\omega_{E}^{2}}{\left(\omega_{E}^{2}-c^{2}p^{2}/\hbar^{2}\right)}\right]^{-1}$

With approaching the point of atomic resonance $E\to E_{n}$ the optical
excitation shows behavior associating with that of a non-condensed disordered
atomic gas. The collective dipole polarization is driven by the propagating
field and the environment of proximal dipoles induces the well known static
Lorentz-Lorentz red shift $-4\pi\,n_{0}\,d_{0}^{2}/3$ from the atomic
resonance as it contributes to the transverse part of the propagator. However
unlike a disordered gas there is an extra static frequency shift, induced by
the polarization interaction with the condensate background, which is given by
the last term in the right-hand side of Eq. (42). Indeed, considering the
quasi-particle as immobile with negligible momentum $p\ll\hbar\omega/c$ the
dependence on $E$ vanishes and this part of the interaction also becomes
static. In this limit the transverse component of the polariton propagator
coincides with its longitudinal part, such that the excitation process becomes
isotropic with positive static shift $8\pi\,n_{0}\,d_{0}^{2}/3$.

The spectral behavior of the polariton propagator in the form (41), (42)
consists of two branches. One is an atom-type excitation near atomic resonance
$E\sim E_{n}$, on which we have commented above. Another resonance exists in
the transverse part of the polariton propagator and is located near the energy
$E\sim E_{0}+c\,p$, which is a pole feature of the last term in the
denominator of the transverse component $G_{\perp}(p,E)$. This resonance
describes the optical excitation propagating through the sample with near
speed of light and creates the photon-type polariton branch. The detail
discussion of spectral behavior of the polariton modes in the infinite and
homogeneous medium is performed in Ref. Ezhova16 .

In general, with an inhomogeneous configuration with the order parameter of
arbitrary profile, equation (35) accepts only numerical solution. In the next
section we present such a solution in a one-dimensional geometry and compare
the results with predictions of conventional macroscopic Maxwell theory.

##  IV Results

Degenerate quantum gases have unique properties and are of particular interest
in reduced spatial dimensionality Salomon ; Levin . This motivates us to
consider initially our results for several instances of a one dimensional
model. Further, equation (35) is quite difficult for numerical solution in a
general three-dimensional configuration. Below we perform results of our
numerical simulations for a one-dimensional model expressed in terms of
transmission and reflection of light from a slab atomic sample, where atoms
can exist in either a quantum degenerate phase or as a disordered classical
gas. The considered geometries are shown in Fig. 1 for three tested
configurations - (a) a uniform slab of BEC with constant density, (b) an
inhomogeneous distribution parameterized by the order parameter with a cosine
profile, (c) interference of two matter waves for two BEC segments counter-
propagating through each other. In the last case, as we show, such an internal
motion of the overlapping condensate fragments can crucially modify the light
scattering process.

Figure 1: (Color online) Geometry of the considered one-dimensional scattering
process: (a) a uniform BEC slab of depth $L$ and with the order parameter
$\Xi=\sqrt{n_{0}}=\mathrm{const}_{z}$, where $n_{0}$ is the density of atoms;
(b) an inhomogeneous distribution parametrized by the order parameter
$\Xi(z)=\sqrt{n_{0}}\cos(\pi z/L)$; (c) interference of two matter waves
created by the BEC sample (b) split in two fragments, see Eq. (44) and
explanation in the text.

###  IV.1 Smooth profile of the order parameter

Any testable profile of the order parameter should be consistent with the
physical model of the condensate and, in the case of weak internal coupling,
performs as a possible solution of the Gross-Pitaevskii equation Pitaevskii ;
Gross . In the macroscopic limit any homogeneous spatial profile of the order
parameter can be suggested as an example of a Thomas-Fermi-type approximate
solution, for which the shape can be fitted by varying the trapping potential.
This approximation works for the condensate confined with an atomic trap where
the period of free oscillation is longer than $2\pi/\mu_{c}$ with $\mu_{c}$
estimated (in a homogeneous limit) by Eq. (12) and it is based on priority of
internal interaction. But even in the case of an ideal gas with $\mu_{c}\to 0$
the order parameter of a quite general profile can be accepted as well, but in
this case as the ground state eigenfunction of the stationary single particle
Schrödinger equation in the trap potential.

As a first example, let us consider the case of a homogeneous degenerate
quantum gas filling a slab of depth $L$ with the order parameter given by
$\Xi=\sqrt{n_{0}}=\mathrm{const}_{z}$, which is shown in Fig. 1(a). In a one-
dimensional geometry, with applying the Fourier transform, the scattering
equations (35) can be rewritten as an infinite set of the algebraic equations,
see Appendix B for derivation details. The obtained system of algebraic
equations can be numerically solved, which give us the spectra of transmission
$\mathrm{T}(\omega)$ and reflection $\mathrm{R}(\omega)$. The same quantities
can be independently constructed via solution of the macroscopic Maxwell
equations, see LaLfVIII , and they are given by

$\displaystyle\mathrm{T}(\omega)$ $\displaystyle=$
$\displaystyle\left|\frac{2\sqrt{\epsilon(\omega)}}{2\sqrt{\epsilon(\omega)}\cos\psi(\omega)-i(1+\epsilon(\omega))\sin\psi(\omega)}\right|^{2}$
$\displaystyle\mathrm{R}(\omega)$ $\displaystyle=$
$\displaystyle\left|\frac{\sin\left[\psi(\omega)\right]}{\sin\left[\psi(\omega)-i\ln\frac{1-\sqrt{\epsilon(\omega)}}{1+\sqrt{\epsilon(\omega)}}\right]}\right|^{2}$
(43)

where $\psi(\omega)=L\sqrt{\epsilon(\omega)}\omega/c$. With substituting here
the dielectric permittivity (33) (with canceled chemical potential) we arrive
at the result predicted for a macroscopic disordered gas, see OurPRA2009 .

In Fig. 2 we compare the spectra of light transmission through and reflection
from the condensate and disordered atomic gas of the same density
$n_{0}\lambdabar_{0}^{3}\sim 0.05$ and in the geometry of Fig. 1(a)The inset
shows the dielectric permittivity given by solution of Eq. (33). Since an
optical excitation from the condensate changes its energy the excitation
spectrum of non-ideal degenerate quantum bosonic gas is red-shifted from
atomic resonance by the value of the chemical potential. The shift is small
and seems negligible as far as the condition (11) is normally fulfilled for
any dipole-type transition and in alkali-metal systems in particular. Thus we
could safely ignore this shift with constructing the susceptibility for the
condensate as the solution of Eq. (33). Nevertheless, we leave it in our
reproduction of the spectral responses as far as such a red shift is a
physical effect and can be visible in the transmission and reflection spectra.
The red shift has been observed in the transmission spectrum of a BEC
consisting of helium atoms on a spectrally narrow dipole forbidden magnetic-
type transition NRV16 .

Figure 2: (Color online) The spectral dependences of transmission (upper
panel) and reflection (lower panel) calculated as a solution of the scattering
equations (35) vs. comparative solution of the Maxwell equations (43) in a
one-dimenional geometry for a homogeneous medium with a slab geometry shown in
Fig. 1(a). The graphs are plotted as a function of detuming
$\Delta=\omega-\tilde{\omega}_{0}$ from the displaced resonance frequency
$\tilde{\omega}_{0}=\omega_{0}-\mu_{c}/\hbar$, see text. The results are
performed for different sample depths $L$, scaled by the wavelength
$\lambda_{0}$ at the atomic resonance, and for the density
$n_{0}\lambdabar_{0}^{3}\sim 0.05$. The reflection spectra for different $L$
are unresolved in the graph with the plotted precision. The inset shows the
dielectric permittivity of the sample
$\epsilon(\omega)=\epsilon^{\prime}(\omega)+i\epsilon^{\prime\prime}(\omega)$
given by solution of Eq. (33) as a function of $\omega-\omega_{0}$. Both the
rounds of calculations give identical results, and to show this in the example
of $L=10\lambda_{0}$ we additionally indicate (by dotted curve) the prediction
of the macroscopic Maxwell theory.

Surprisingly, but this global offset of the spectral profile is only one
difference between the transmission and reflection spectra of degenerate and
non-degenerate atomic gases. To demonstrate this we plotted the graphs as a
function of detuning $\Delta=\omega-\tilde{\omega}_{0}$, where
$\tilde{\omega}_{0}=\omega_{0}-\mu_{c}/\hbar$, and where we additionally
displaced the spectra of a disordered gas on $\mu_{c}/\hbar$. We have obtained
excellent, i.e. point by point, coincidence of degenerate and non-degenerate
spectra despite the fact that they were calculated via solution of
exceptionally different equations. The small deviation for reflection near its
resonant point is a result of additional boundary contributions ignored in the
Fourier transformation of the Laplace operator to the algebraic form of Eq.
(68) and this incorrectness, as we have verified, softens in the macroscopic
limit $L/\lambda_{0}\to\infty$. The reflection itself is weak but not
negligible and results from the scattering from the sample edges and is
enhanced by interference effect. Such an excellent coincidence of two
independent rounds of calculations clearly indicates that for light scattering
from an ensemble of atoms, with uniform density distribution, the optical
response of the system is insensitive to either classical or quantum nature of
statistical averaging.

This can be confirmed by similar calculations performed for the order
parameter with a trigonometric profile $\Xi(z)=\sqrt{n_{0}}\cos(\pi z/L)$ ( in
geometry of Fig. 1(b)), and the results are shown in Fig. 3. For this case we
make additional simplifications with expanding
$\sqrt{\epsilon\left(z,\omega\right)}$ in a Taylor series near the vacuum
point $\epsilon=1$ and with keeping only the forwardly propagating wave in the
macroscopic Maxwell description of the problem. Again the calculations show
good (within the made approximations) agreement between both the approaches.
We used the same peak density $n_{0}\lambdabar_{0}^{3}\sim 0.05$ and the same
sample depths as in the plots of Fig. 2. In the case of smoothed sample bounds
with density profile $n_{0}(z)=n_{0}\cos^{2}(\pi z/L)$ the backward scattering
is expected as many orders of magnitude weaker process because of vanishing
boundary contributions. The latter can be seen via negligible response of the
reflected light as follows from the calculation data shown in the lower panel
of Fig. 3.

Figure 3: (Color online) Same as in Fig. 2 but for the density distribution
parametrized by the order parameter $\Xi(z)=\sqrt{n_{0}}\cos(\pi z/L)$ for a
slab geometry shown in Fig. 1(b). In the case of a smooth sample boundary the
backward scattering reveals a many orders of magnitude weaker process than in
the case of sharp boundaries.

###  IV.2 Interference of two counter-propagating BEC fragments

Finally, let us consider the experimental configuration when, as a result of
coherent interaction with light, a BEC sample is fractured into a number of
macroscopic pieces Schneble ; Hilliard . To simplify the complicated
experimental picture we model the process by the presence of only two
fragments counter-propagating with respect to each other in their center-of-
momentum reference frame. The considered configuration is shown in Fig. 1(c)
and reveals a strong density oscillation associated with interference of the
order parameters (matter wave-packets) from the BEC pieces in the area of
their overlapping. The existence of such a fringe structure of the density
distribution has been directly observed as an effect of interference of two
condensates in experiment Ketterle97 . The spatial phase matching condition,
determined by internal relative motion of the fragments, crucially affects the
scattering process. Indeed, the wave length of the oscillation is determined
by the speed of relative motion and after accumulation of essential linear
momentum from light can exceed a scale of the light wavelength. Then such a
density grating should lead to strong Bragg diffraction and, as we show by our
numerical simulations below, to significant enhancement of the backward
scattering.

The process can be described by the order parameter of the following spatial
profile

$\displaystyle\Xi(z)$ $\displaystyle=$
$\displaystyle\sqrt{2n_{0}}\cos\left(\frac{\pi z}{L}\right)\,\cos\left(\Delta
q\,z\right)$ (44) $\displaystyle=$
$\displaystyle\sqrt{\frac{n_{0}}{2}}\cos\left(\frac{\pi
z}{L}\right)\,\mathrm{e}^{i\Delta
q\,z}+\sqrt{\frac{n_{0}}{2}}\cos\left(\frac{\pi
z}{L}\right)\,\mathrm{e}^{-i\Delta q\,z}$ $\displaystyle\equiv$
$\displaystyle\Xi_{+}(z)+\Xi_{-}(z)$

which is constructed as an ideal overlap of two matter wave-packets associated
with the condensate fragments of identical shape and size counter-propagating
with respect to each other with the relative linear momentum $2\hbar\Delta q$
per atom. Let us make a clarifying comment concerning the validity and
consistency of the suggested profile as a solution of the time dependent
Gross-Pitaevskii equation.

Both of the partial contributions $\Xi_{+}(z)$ and $\Xi_{-}(z)$ are
representative solutions of the order parameter equation, for example, in the
Thomas-Fermi approximation. That can be justified via transforming dynamical
description of any of the wave-packets to that reference frame where the
particular fragment is motionless and then we arrive to the configuration
considered in the previous subsection. But the entire process of expansion and
fragmentation of the condensate, modeled by (44), can be imagined only after
the BEC is released from the trap and it results from both the external
disturbance and internal interaction processes. The superposed state (44) can
physically model the complicated dynamics of the condensate fragmentation once
we ignore the weak non-ideality of the atomic gas in comparison with the
kinetic energy associated with the relative motion of the fragments, see
inequality (11). This can be fulfilled for quite high relative speed with
$\Delta\,q\gg 1/L$ and $2\hbar^{2}\Delta q^{2}/\mathrm{m_{A}}>\mu_{c}$. Then
the factor ”$\cos\left(\Delta q\,z\right)$” is a strongly oscillating function
of $z$, which implies its averaging $\langle\cos^{2}\Delta q\,z\rangle\to 1/2$
in the normalization of the order parameter by a total number of particles.
Then expansion (44) corresponds to beginning of the splitting process of the
released matter wave $\Xi(z)=\sqrt{n_{0}}\cos(\pi z/L)$, as shown in Fig.
1(c), in two separated wave-packets $\Xi_{+}(z)$ and $\Xi_{-}(z)$ propagating
in opposite directions.

Figure 4: (Color online) Same as in Fig. 2 but for the density distribution
parametrized by the order parameter $\Xi(z)=\sqrt{2n_{0}}\cos\left(\pi
z/L\right)\,\cos\left(\Delta q\,z\right)$ with $2\Delta q=k_{0}$ for the
geometry shown in Fig. 1(c). Both the forward and backward scattering have a
clear signature of the coherent enhancement due to the effect of the Bragg
diffraction. In the upper panel the dotted curve indicates the reference
transmission spectrum for $L=10\lambda_{0}$ with smoothed profile of the order
parameter and corresponds to the configuration of a disordered atomic gas.

In Fig. 4 we show the spectra of transmission and reflection for the order
parameter with the spatial profile given by Eq. (44). It is expected that for
a classical disordered gas any internal motion of its macroscopic fragments
with a rather slow relative speed would not modify the scattering process at
all. As an example, such an expansion with a relative speed given by the
recoil limit $\sim\hbar k_{0}/\mathrm{m_{A}}$ would induce only a negligible
Doppler shift between the spectral outputs from both the fragments. But in the
case of BEC such an internal motion dramatically modifies the scattering
process. As pointed out above, the spatial modulation of the order parameter
initiates a mechanism of the Bragg diffraction and scattering on the spatially
oscillating density. As a consequence, this leads to strong enhancement of the
backward scattering and it is manifestable in an abrupt structure of the
transmission spectrum as well. The strongest scattering is observed for the
modulation wave number $\Delta q=k_{0}$ when the condensate expands with the
relative speed $v_{0}=2\hbar k_{0}/\mathrm{m_{A}}$. As follows from the
dependencies of Fig. 4 this effect experiences as a broader spectral domain as
the sample spatial scale is longer.

In Fig. 5 we reproduce the dependence of the reflection coefficient as a
function of $2\pi/\Delta q$ for different sample depths $L$. As can be seen
from these graphs, the reflection always has a local maxima at the points
$\Delta q=2\pi/\lambda_{0},\,2\pi/2\lambda_{0}\ldots$. This is optimal
condition for manifestation of the Bragg diffraction, which creates the
oppositely propagating polariton wave via scattering of the impinging wave on
periodic structure. As a consequence of the Bragg-type scattering an
additional amount of linear momentum transfers to the condensate and enforces
its fragmentation. So the Bragg diffraction also results in a certain
optomechanical action on the system and accordingly leads to kinematic
entanglement of the spatially structured BEC, see Hilliard .

In our calculation model we can describe such an effect of optomechanical
interface primary for the backward and forward scattering channels.
Nevertheless, in experiment Schneble  the fragmentation was observed for the
scattering directions orthogonal to the incident light along the major axis of
an ellipsoid-shaped condensate sample. The observed effect had been associated
in Schneble  with the Kapitza-Dirac phenomenon of the matter wave scattering
on the spatial structure created by an electromagnetic wave. In this sense, we
can point out that in the case of excitation of a BEC sample by an external
light pulse, consisting of many photons, the entire dynamics apparently
results from several physical processes, which includes internal interactions,
disturbance of the matter-wave (order parameter) by external driving field and
formation of the polariton structure by the optical excitation. Then the Bragg
scattering reveals a coherent mechanism for rearranging of photon-type
polariton waves (see Section III.3) propagating in different directions. The
coherently scattered photons emerge the sample with indicating prior
propagation directions of these waves.

Figure 5: (Color online) Reflection coefficient for the order parameter of
Fig. 4, at the point of atomic resonance, plotted as a function of
$2\pi/\Delta q$ for different sample depths $L$.

##  V Conclusion

In this paper we have developed a formalism of the microscopic quantum
scattering theory directed towards ab-initio description of the elementary
process of a single photon scattering from a quantum degenerate atomic gas.
The gas exists in the BEC phase parametrized by the order parameter introduced
in the framework of the Gross-Pitaevskii model. The main mathematical object
of our calculational approach is the single particle Green’s function
(propagator) tracking the propagation of a specific polariton wave through the
condensate. The polariton is created as a quantum superposed state between the
photon and condensate.

The polariton propagation is disturbed by the process of incoherent scattering
and its entire dynamics is described by the closed scattering equation for the
complete polariton propagator as we derived. The crucial difference with the
light propagation through a disordered and non-degenerate atomic gas is that
in the considered case the atomic medium represents a coherent matter wave
strongly rejecting its classical interpretation. The conventional vision of
the macroscopic Maxwell description of the electromagnetic wave in a bulk
medium seems insufficient and can be even incorrect in some situations.

To clarify the above point we have solved the derived scattering equations in
a one-dimensional geometry and compared the result with predictions of the
conventional macroscopic Maxwell theory for the disordered atomic gas of the
same density and size as the BEC sample. For steady state conditions and
uniform distribution of the order parameter we obtained identical results for
the transmission and reflection spectra for both the approaches. Nevertheless
we observe a significant difference once the BEC is fractured into a number of
the interfering matter wave fragments. In the latter case the scattering
process evolves towards conditions of Bragg diffraction, which strongly
affects the process and can coherently redirect the propagating polariton wave
in the backward or other directions associated with the condensate
fragmentation.

## Acknowledgements

This work was supported by the Russian Foundation for Basic Research under
Grant 15-02-01060. We also acknowledge financial support by the National
Science Foundation under Grant No. NSF-PHY-1606743.

##  Appendix A Overview of the diagram approach

Below we introduce basic elements of the diagram equations, which are
constructed and discussed in the main text. We follow standard definitions and
rules of the microscopic version of the Feynman diagram method, as described
in Ref. BerstLifshPitvsk , but revise it for a non-relativistic dipole-type
coupling of light and atoms, see Ref. Kupriyanov17 . The expansion of the
evolution operator (10) in the Green’s function (9) generates the sequence of
expectation values of the various operator products, which after a set of
transpositions and with the aid of the Wick theorem can be regrouped to the
results visualized by diagram images. The diagrams consist of the objects
listed below.

The undisturbed causal-type electric field Green’s function is defined via
transposition of the field operators in any pair product from chronologically
$T$-ordered to normally $N$-ordered form

$\displaystyle
iD_{\mu\mu^{\prime}}^{(E)}(\mathbf{r},t;\mathbf{r}^{\prime},t^{\prime})$ (45)
$\displaystyle=T\\!\left[E^{(0)}_{\mu}\\!(\mathbf{r};t)E^{(0)}_{\mu^{\prime}}\\!(\mathbf{r}^{\prime};t^{\prime})\right]\\!-\\!N\\!\left[E^{(0)}_{\mu}\\!(\mathbf{r};t)E^{(0)}_{\mu^{\prime}}\\!(\mathbf{r}^{\prime};t^{\prime})\right]$

It can be linked with a fundamental object of quantum electrodynamics namely
with the causal-type photon propagator

$D_{\mu\mu^{\prime}}^{(E)}(\mathbf{r},t;\mathbf{r}^{\prime},t^{\prime})=\left.\frac{1}{c^{2}}\frac{\partial^{2}}{\partial
t\partial
t^{\prime}}D_{\mu\mu^{\prime}}^{(c)}(\mathbf{r},t;\mathbf{r}^{\prime},t^{\prime})\right|_{\scriptsize{\begin{array}[]{c}\mathbf{r}\neq\mathbf{r}^{\prime}\\\
\mathrm{or}\\\ t\neq t^{\prime}\end{array}}}$ (46)

where we follow gradient invariance of the theory and fix the propagator by a
vanishing scalar potential such that $\mu,\,\mu^{\prime}=x,y,z$. With
simplifying notation for each argument $\mu,\mathbf{r},t\to x$ and
$\mu^{\prime},\mathbf{r}^{\prime},t^{\prime}\to x^{\prime}$ the electric field
Green’s function is imaged by a wavy line

$iD^{(E)}(x,x^{\prime})\Leftrightarrow\begin{array}[c]{c}\scalebox{1.0}{\includegraphics*{eqA_3.eps}}\end{array}$
(47)

where the ending indices are often omitted in graph equations. This function
depends only on the difference of its spatial and time arguments and its
Fourier image is given by

$\displaystyle D^{(E)}_{\mu\mu^{\prime}}(\mathbf{k},\omega)$ (51)
$\displaystyle=\\!\int\\!d^{3}\\!R\int^{\infty}_{-\infty}\\!\\!d\tau\,\mathrm{e}^{i\omega\tau-i\mathbf{k}\cdot\mathbf{R}}\\!\left.D_{\mu\mu^{\prime}}^{(E)}(\mathbf{R},\tau)\right|_{\scriptsize{\begin{array}[]{c}\mathbf{R}=\mathbf{r}-\mathbf{r}^{\prime}\\\
\tau=t-t^{\prime}\end{array}}}$
$\displaystyle=\\!\frac{4\pi\hbar\omega^{2}}{\omega^{2}-\omega^{2}_{k}+i0}\left[\delta_{\mu\mu^{\prime}}-c^{2}\,\frac{k_{\mu^{\prime}}k_{\mu}}{\omega^{2}}\right]$

where $\omega_{k}=ck$.

The electric field Green’s function is expressed via solution of the
microscopic Maxwell equations with a point-like dipole source and for
$\omega>0$ coincides with positive-frequency component of the retarded-type
fundamental solution of these equations
$D^{(R)}_{\mu^{\prime}\mu}(\mathbf{k},\omega)$

$\left.D^{(E)}_{\mu\mu^{\prime}}(\mathbf{k},\omega)\right|_{\omega>0}=\frac{\omega^{2}}{c^{2}}\left.D^{(R)}_{\mu\mu^{\prime}}(\mathbf{k},\omega)\right|_{\omega>0}$
(52)

The positive frequency domain is only important in the RWA approach and in
this approximation it is convenient to add an arrow in the diagram (47) for
indicating creation and annihilation events of a virtual photon at the edging
points of the wavy line.

The undisturbed atomic Greens’s function is defined via transposition of the
atomic operators in any pair product from chronologically ordered to normally
ordered form. For operators of the excited state this reads

$\displaystyle
iG_{nn^{\prime}}^{(0)}(\mathbf{r},t;\mathbf{r}^{\prime},t^{\prime})$ (53)
$\displaystyle=T\\!\left[\Psi^{(0)}_{n}(\mathbf{r};t)\Psi^{(0)\dagger}_{n^{\prime}}(\mathbf{r}^{\prime};t^{\prime})\right]\\!-\\!\Psi^{(0)\dagger}_{n^{\prime}}(\mathbf{r}^{\prime};t^{\prime})\Psi^{(0)}_{n}(\mathbf{r};t)$

and similarly with replacement $n,n^{\prime}\to m\\!=\\!m^{\prime}\\!=0$ for
operators of the ground state. With simplifying notation for each argument
$n,\mathbf{r},t\to x$ and $n^{\prime},\mathbf{r}^{\prime},t^{\prime}\to
x^{\prime}$ the atomic Green’s function is imaged by an arrowed straight line

$iG^{(0)}(x,x^{\prime})\Leftrightarrow\ \
\begin{array}[c]{c}\scalebox{1.0}{\includegraphics*{eqA_7.eps}}\end{array}$
(54)

where the ending indices are often omitted in graph equations. This function
also depends only on the difference between its spatial and time arguments and
its Fourier image is given by

$\displaystyle G^{(0)}_{nn^{\prime}}(\mathbf{p},E)$ (57)
$\displaystyle=\\!\int\\!d^{3}\\!R\int^{\infty}_{-\infty}\\!\\!d\tau\,\mathrm{e}^{\frac{i}{\hbar}E\tau-\frac{i}{\hbar}\mathbf{p}\cdot\mathbf{R}}\\!\left.G_{nn^{\prime}}^{(0)}(\mathbf{R},\tau)\right|_{\scriptsize{\begin{array}[]{c}\mathbf{R}=\mathbf{r}-\mathbf{r}^{\prime}\\\
\tau=t-t^{\prime}\end{array}}}$ $\displaystyle=$
$\displaystyle\delta_{nn^{\prime}}\,\frac{\hbar}{E-p^{2}/{2\mathrm{m_{A}}}-E_{n}+i0}$
(58)

where $\mathrm{m}_{A}$ is the atomic mass and the internal atomic state is
assumed to be degenerate such that $E_{n}=\mathrm{const}_{n}$.

The atomic Green’s function is expressed by the fundamental solution (atomic
propagator) of the Schrödinger equation for a free atom which describes
propagation of an atomic wave initially localized in a certain spatial point.
As follows from (53) this function vanishes if $t<t^{\prime}$ such that the
causal-type atomic propagator is identical to the retarded-type propagator.

There are different diagram vertices indicating optical interactions of
different types. If a virtual photon interacts with an atom, which is also
presented as a virtual object in a diagram, then in the RWA we associate the
process with the following two vertexes

$\displaystyle\frac{i}{\hbar}\,d^{\mu}_{nm}\ $ $\displaystyle\Leftrightarrow$
$\displaystyle\frac{i}{\hbar}\,d^{\mu}_{mn}\ $ $\displaystyle\Leftrightarrow$
(59)

If a similar process is developing with condensate particles we associate it
with the vertexes of another type

$\displaystyle\frac{i}{\hbar}\,d^{\mu}_{nm}\,\Xi(\mathbf{r})\,\mathrm{e}^{-\frac{i}{\hbar}\epsilon_{0}t}\
$ $\displaystyle\Leftrightarrow$
$\displaystyle\frac{i}{\hbar}\,d^{\mu}_{mn}\,\Xi^{\ast}(\mathbf{r})\,\mathrm{e}^{+\frac{i}{\hbar}\epsilon_{0}t}\
$ $\displaystyle\Leftrightarrow$ (60)

which describe either excitation of an atom from the condensate phase (upper
diagram) or its recovering in the condensate phase (lower diagram). The
detailed specification of vertexes is usually unnecessary and often omitted if
it does not confuse interpretation of the diagram.

In the original representation each vertex corresponds to the integral over
respective spatial and time variables and each contributing line is decoded in
accordance with (47) and (54). In the stationary and homogeneous conditions
after Fourier transform, the external lines are decoded in accordance with
(51) and (58) but internal lines, when they shape a loop, contribute as
convolution-type integrals over reciprocal variables such as energy
(frequency) and momentum (wave vector), with conserving total energy and
momentum transporting by the diagram. For more details we readdress reader to
Refs.BerstLifshPitvsk ; Kupriyanov17 .

##  Appendix B One-dimensional scattering

Consider the scattering problem for a slab geometry of an atomic medium,
homogeneous and infinite in the plane transverse to the wave vector of the
incident photon. In this case the $T$-matrix element, given by Eq. (LABEL:2.6)
and selected for either forward or backward elastic scattering channels, is
given by

$\displaystyle T_{i^{\prime}i}(E)=\frac{2\pi\omega}{{\cal
L}}\iint\,dz^{\prime}\,dz\sum\limits_{n^{\prime},n}{(\mathbf{d}\cdot\mathbf{e})}^{*}_{n^{\prime}0}{(\mathbf{d}\cdot\mathbf{e})}_{0n}$
(61)
$\displaystyle\times{\mathbf{e}}^{-ik^{\prime}z^{\prime}+ikz}\;\Xi^{*}(z^{\prime})\,\Xi(z)\,G_{n^{\prime}n}(z^{\prime},z;E-E_{0}^{N-1})$

where the output frequency and polarization are unchanged such that
$\omega^{\prime}=\omega$ and $\mathbf{e}^{\prime}=\mathbf{e}$, and we
redefined $f=i^{\prime}$ with emphasizing the physical equivalence of initial
and final states in one dimensional scattering process. All the integrands are
considered as functions of longitudinal coordinates $z,z^{\prime}$ and the
polariton propagator is proportional to a $\delta$-function of transverse
coordinates $x,y$ and $x^{\prime},y^{\prime}$, see Eq. (35). The integral
evaluated in the transverse plane over variables $dxdy$ and
$dx^{\prime}dy^{\prime}$ cancels out the area scale ${\cal L}_{x}{\cal L}_{y}$
in the normalization volume ${\cal V}={\cal L}_{x}{\cal L}_{y}{\cal L}_{z}$
and we denoted ${\cal L}_{z}={\cal L}$.

Let us express the $S$-matrix components via $T$-matrix

$\displaystyle S_{i^{\prime}i}$ $\displaystyle=$
$\displaystyle\delta_{i^{\prime}i}-i\,\frac{{\cal L}}{\hbar
c}\,T_{i^{\prime}i}(E_{i}+i0)$ (62)

In a one-dimensional geometry for non-degenerate ground state of the
degenerate quantum gas the light scattering can be described by coefficients
of transmission $\mathrm{T}(\omega)$, reflection $\mathrm{R}(\omega)$ and
losses $\mathrm{L}(\omega)$, which are subsequently given by

$\displaystyle\mathrm{T}(\omega)$ $\displaystyle=$
$\displaystyle\left.{\left|S_{i^{\prime}i}\right|^{2}}\right|_{k^{\prime}=k>0}$
$\displaystyle\mathrm{R}(\omega)$ $\displaystyle=$
$\displaystyle\left.{\left|S_{i^{\prime}i}\right|^{2}}\right|_{k^{\prime}=-k<0}$
$\displaystyle\mathrm{L}(\omega)$ $\displaystyle=$ $\displaystyle
1-\mathrm{T}(\omega)-\mathrm{R}(\omega)$ (63)

and can be found via solution of the simplified equations (35)-(40) as we show
below.

Consider the example of the slab with the order parameter
$\Xi(z)=\sqrt{n_{0}}=\mathrm{const}_{z}$ inside the medium. In this case the
integral equation (35) can be transformed to the set of algebraic equations
via spatial Fourier transform with periodic boundary conditions on the sample
bounds. The azimuthal symmetry justifies the diagonal structure of the
polariton propagator

$G_{nn^{\prime}}(z,z^{\prime};E)=\delta_{nn^{\prime}}\,G(z,z^{\prime};E)$ (64)

Then, with the assumption that the origin of the coordinate frame is located
in the middle point and $z\in(-L/2,L/2)$, where $L$ is the sample length, it
can be expanded as

$\displaystyle G_{ss^{\prime}}(E)$ $\displaystyle=$
$\displaystyle\frac{1}{L}\,\iint\limits_{-L/2}\limits^{L/2}dzdz^{\prime}\mathrm{e}^{-ik_{s}z+ik_{s^{\prime}}z^{\prime}}G(z,z^{\prime};E)$
$\displaystyle G(z,z^{\prime};E)$ $\displaystyle=$
$\displaystyle\frac{1}{L}\,\sum_{s,s^{\prime}}\mathrm{e}^{ik_{s}z-ik_{s^{\prime}}z^{\prime}}G_{ss^{\prime}}(E)$
(65)

where $k_{s}=2\pi s/L$ and $k_{s^{\prime}}=2\pi s^{\prime}/L$ with
$s,s^{\prime}=0,\pm 1,\pm 2,\ldots$. The Green’s function (64) contributes to
the transmission amplitude (62) at a specific energy argument
$E-E_{0}^{N-1}=E_{i}-E_{0}^{N-1}=\hbar\omega+E_{0}^{N}-E_{0}^{N-1}=\hbar\omega+\varepsilon_{0}$
and we denote

$\left.{G_{ss^{\prime}}(E)}\right|_{E=\hbar\omega+\varepsilon_{0}}\equiv
G_{ss^{\prime}}(\omega)$ (66)

and consider the Fourier components as functions of the frequency of the
incident photon. Then the $S$-matrix elements (62) can be expressed as follows

$\displaystyle S_{i^{\prime}i}=\delta_{i^{\prime}i}-\frac{8\pi i\omega}{L\hbar
c}\,n_{0}d_{0}^{2}\,\sum_{s^{\prime},s}$ (67)
$\displaystyle\times\frac{\sin\left(k^{\prime}-k_{s^{\prime}}\right)\frac{L}{2}}{k^{\prime}-k_{s^{\prime}}}\;\frac{\sin\left(k-k_{s}\right)\frac{L}{2}}{k-k_{s}}\;G_{s^{\prime}s}(\omega)$

where $k=\omega/c$ and $k^{\prime}=\pm\omega/c$.

With substituting (64) and applying transforms (65) to Eq. (35), considered in
a one-dimensional configuration, we arrive at the following system of
algebraic equations

$\displaystyle\left[\omega-\tilde{\omega}_{0}+\frac{\hbar\,k_{s}^{2}}{2\mathrm{m_{A}}}+\frac{4\pi}{3\hbar}\,n_{0}\,d_{0}^{2}+\frac{i}{2}\sqrt{\epsilon\left(\omega\right)}\;\gamma\right]G_{ss^{\prime}}(\omega)$
(68)
$\displaystyle-\sum_{s^{\prime\prime}}\Sigma_{ss^{\prime\prime}}^{(c)}(\omega)G_{s^{\prime\prime}s^{\prime}}(\omega)=\delta_{ss^{\prime}}$

where
$\tilde{\omega}_{0}=(E_{n}-\varepsilon_{0})/\hbar=(E_{n}-E_{0}-\mu_{c})/\hbar\equiv\omega_{0}-\mu_{c}/\hbar$
with same $E_{n}$ for all the upper state Zeeman sublevels. We approximated
$\epsilon(\omega+\mu_{c}/\hbar)\approx\epsilon(\omega)$, see Eq. (33) and
related comment.

The matrix of the self-energy part is given by

$\displaystyle\Sigma_{ss}^{(c)}(\omega)=\frac{4\pi}{\hbar}n_{0}d_{0}^{2}\,\frac{\omega^{2}}{\omega^{2}-c^{2}\,k_{s}^{2}}$
(69) $\displaystyle-\frac{4\pi
i}{\hbar}n_{0}d_{0}^{2}\,\frac{\omega}{cL}\frac{\displaystyle\frac{\omega^{2}}{c^{2}}+k_{s}^{2}}{\left(\displaystyle\frac{\omega^{2}}{c^{2}}-k_{s}^{2}\right)^{2}}\left[1-\exp\left(i\frac{\omega}{c}L\right)\right]$

for $s^{\prime\prime}=s$ and

$\displaystyle\Sigma_{ss^{\prime\prime}}^{(c)}(\omega)\\!\\!$ $\displaystyle=$
$\displaystyle\\!\\!-(-)^{s-s^{\prime\prime}}\frac{4\pi
i}{\hbar}n_{0}d_{0}^{2}\,\frac{\omega}{cL}\frac{\displaystyle\frac{\omega^{2}}{c^{2}}+k_{s}k_{s^{\prime\prime}}}{\left(\displaystyle\frac{\omega^{2}}{c^{2}}-k_{s}^{2}\right)\\!\\!\left(\displaystyle\frac{\omega^{2}}{c^{2}}-k_{s^{\prime\prime}}^{2}\right)}$
(70) $\displaystyle\times$
$\displaystyle\left[1-\exp\left(i\frac{\omega}{c}L\right)\right]$

for $s^{\prime\prime}\neq s$. For a sample of infinite length $L\to\infty$
Eqs. (66), (68)-(70) reproduce the transverse component of the polariton
propagator in an infinite and uniform medium, see Eq. (42), and in this case
the scattering process manifests itself mainly via the incoherent channels.

For the sample of finite length the system (68) consists of an infinite number
of equations. Nevertheless it can be numerically solved with cutoff by a
limited number of the involved equations. With increasing of this number the
iterative process becomes internally converging and approaching the exact
solution. The performed calculation scheme can be straightforwardly
generalized if the order parameter is non-uniform and described by
trigonometric functions such as $\Xi(z)\sim\cos(\pi z/L)$,
$\Xi(z)\sim\mathrm{e}^{i\kappa_{1}z}\cos(\pi
z/L)+\mathrm{e}^{i\kappa_{2}z}\cos(\pi z/L)$ etc., which we have considered in
our numerical simulations.

## References

## References

  * (1) H.D. Politzer Phys. Rev. A 43, 6414(R) (1991).
  * (2) M.H. Anderson, J.R. Ensher, M.R. Matthews, C.E. Wieman, and E.A. Cornell, Science 269, 198 (1995).
  * (3) K.B. Davis, M.-O. Mewes, M.R. Andrews, N.J. van Druten, D.S. Durfee, D.M. Kurn, and W. Ketterle, Phys. Rev. Lett. (PRL) 75 3969 (1995).
  * (4) S. Inouye, A.P. Chikkatur, D.M. Stamper-Kurn, J. Stenger, D.E. Pritchard, W. Ketterle, Science 285, 571 (1999).
  * (5) D. Schneble, Y. Torii, M. Boyd, E.W. Streed, D.E. Pritchard, W. Ketterle, Science 300, 475 (2003).
  * (6) A. Hilliard, F. Kaminski, R. le Targat, C. Olausson, E.S. Polzik, and J.H. Müller, Phys. Rev. A 78, 051403(R) (2008); J. H. M ller, D. Witthaut, R. le Targat, J. J. Arlt, E. S. Polzik, A. J. Hilliard, Journal of Modern Optics, Volume 63, Issue 18, 1886 (2016).
  * (7) M.R. Matthews, B.P. Anderson, P.C. Haljan, D.S. Hall, C.E. Wieman, and E.A. Cornell, Phys. Rev. Lett. 83, 2498 (1999).
  * (8) F. Chevy, K.W. Madison, and J. Dalibard, Phys. Rev. Lett. 85, 2223 (2000).
  * (9) M.F. Andersen, C. Ryu, Pierre Cladé, Vasant Natarajan, A. Vaziri, K. Helmerson, and W.D. Phillips Phys. Rev. Lett. 97, 170406 (2006).
  * (10) J.F.S. Brachmann, W.S. Bakr, J. Gillen, A. Peng, M. Greiner, Optics Express Vol. 19, No. 14 , 12984 (2011).
  * (11) M.G. Bason, R. Heck, M. Napolitano, O. El asson, R. Müller, A. Thorsen, W.-Z. Zhang, J. Arlt, J.F. Sherson, e-print arXiv:1607.02934 (2016).
  * (12) K.-P. Marzlin, W. Zhang, and E.M. Wright Phys. Rev. Lett. 79, 4728 (1997).
  * (13) J.E. Williams and M.J. Holland, Nature 401 568 (1999).
  * (14) K.T. Kapale and J.P. Dowling, Phys. Rev. Lett. 95, 173601 (2005).
  * (15) M.G. Moore and P. Meystre, Phys. Rev. Lett. 83, 5202 (1999).
  * (16) O. Zobay, G.M. Nikolopoulos, Phys. Rev. A.73, 013620 (2006).
  * (17) E.D. Trifonov, Laser Physics, 12, No. 1, pp. 211 216 (2002).
  * (18) Ya. A. Avetisyan, E.D. Trifonov, Phys. Rev. A 88, 025601 (2013).
  * (19) Ya.A. Avetisyan and E.D. Trifonov, Phys.-Usp. 58 286 (2015).
  * (20) N.N. Bogoliubov, J. Phys. (USSR), 11:23 (1947).
  * (21) L.P. Pitaevskii, J. Exptl. Theoret. Phys. (USSR), 40, 646-651 (1961).
  * (22) E.P. Gross, Nuovo Cimento 20 454 (1961); J. Math. Phys. 4 195 (1963).
  * (23) C. Cohen-Tannoudji, J. Dupont-Roc, G. Grynberg, Atom-Photon Interactions. Basic Processes and Applications (John Wiley $\&$ Sons, Inc., 1992).
  * (24) V.B. Beresteskii, E.M. Lifshits, L.P. Pitaevskii, Course of Theoretical Physics: Quantum Electrodynamics (Oxford: Pergamon Press, 1981); E.M. Lifshits, L.P. Pitaevskii, Course of Theoretical Physics: Statistical Physics. Part II (Oxford: Pergamon Press, 1980).
  * (25) D.V. Kupriyanov, I.M. Sokolov, M.D. Havey, Physics Reports 671, 1 (2017).
  * (26) I.M. Sokolov, M.D. Kupriyanova, D.V. Kupriyanov, M.D. Havey Phys. Rev. A 79 053405 (2009).
  * (27) J. Ruostekovski and J. Javanainen, Phys. Rev. A 56, 2056 (1997).
  * (28) V.M. Ezhova, L.V. Gerasimov and D.V. Kupriyanov Journal of Physics: Conference Series 769 012045 (2016).
  * (29) C. Salomon, G. Shlyapnikov, L.F. Cugliandolo, Many-Body Physics with Ultracold Gases, (Oxford: University Press, UK, 2013).
  * (30) Katheryn Levin, Alexander L. Fetter, and Dan M. Stamper-Kurn, Ultracold Bosonic and Fermionic Gases, (Elsevier: Oxford, UK, 2012).
  * (31) L.D. Landau, E.M. Lifshits Course of Theoretical Physics: Electrodynamics of continuous media (Oxford: Pergamon Press, 1981).
  * (32) R.P.M.J.W. Notermans, R.J. Rengelink, and W. Vassen, Phys. Rev. Lett. 117, 213001 (2016).
  * (33) M.R. Andrews, C.G. Townsend, H.-J. Miesner, D.S. Durfee, D.M. Kurn, W. Ketterle, Science, 275, 637 (1997).