#  Light-matter interaction at rough surfaces:
a morphological perspective on laser-induced periodic surface structures

Vladimir Yu. Fedorov<EMAIL_ADDRESS>Laboratoire Hubert Curien,
Université Jean Monnet, Saint-Etienne, France  Jean-Philippe Colombier
<EMAIL_ADDRESS>Laboratoire Hubert Curien,
Université Jean Monnet, Saint-Etienne, France

###### Abstract

We use ab-initio electromagnetic simulations to investigate light absorption
by rough surfaces in the context of the formation of laser-induced periodic
surface structures. Our approach involves modeling a realistic rough surface
using a statistical description of its continuous height distribution via a
corresponding correlation function. We study the influence of incident light
polarization and various statistical properties of surface roughness, such as
root-mean-square height and correlation length, on the distribution of
absorbed laser energy. By analyzing light absorption in different layers of
the surface selvedge, we elucidate how different features of surface
morphology influence the shape of the resulting periodic surface structures.
We show that circularly polarized laser pulses are highly sensitive to initial
or progressively developing asymmetries in surface roughness.

##  I Introduction

When an intense ultrashort laser pulse hits the surface of a solid, it
triggers a chain of events, each unfolding on a specific time scale. In the
initial stage occurring on a femtosecond time scale, within the pulse
duration, a strong laser field excites electrons in the material; light
scatters by surface roughness and, in case of metallic media, generates
surface plasmons. Then, on a picosecond time scale, after the passage of the
laser pulse, the electron-phonon thermalization stage begins. At this stage,
hot excited electrons, colliding with atoms of the material, transfer their
energy to the lattice, thereby heating the material. Next, on a time scale
from a few picoseconds to nanoseconds, hydrodynamic processes come into play.
During this period, material melting, cavitation and explosive boiling occur,
accompanied by convection and hydrodynamic instabilities that ultimately shape
the surface as the material solidifies. Finally, approximately a millisecond
later (corresponding to a typical kilohertz laser repetition rate) the
succeeding pulse within the pulse train arrives and triggers the same chain of
events again. During this iterative process, each subsequent laser pulse
interacts with the surface that has been modified by the previous pulses. This
interaction forms a feedback loop, which, after multiple cycles, leads to the
emergence of periodic ripples known as laser-induced periodic surface
structures (LIPSS) Bonse2017 ; Graef2020 ; Bonse2020 ; Rudenko2023 .

The formation of LIPSS is a universal phenomenon: they can be produced on
almost all types of materials (dielectrics, semiconductors, metals). By fine-
tuning parameters of LIPSS, we can change surface properties and thereby
functionalize the irradiated material. For example, since LIPSS act as a
diffraction grating, by adjusting their period we can select which wavelengths
of light will be reflected and thus change the color of materials Vorobyev2008
; Dusser2010 ; Li2015 . In addition, with LIPSS we can alter the wettability
of the surface and create hydrophobic or hydrophilic coatings Zorba2008 ;
Kietzig2009 ; Kunz2018 . Also, using LIPSS we can control the friction
coefficient of materials Yu1999 ; Bonse2018 ; Kunz2020  or impart
antibacterial properties to their surfaces Jalil2020 ; Sotelo2023 ; Outon2024
.

Depending on the period, LIPSS can be divided in two main classes: LIPSS whose
period is close to the laser wavelength are referred to as low spatial
frequency LIPSS (LSFLs) and those whose period is significantly less than the
laser wavelength are referred to as high spatial frequency LIPSS (HSFLs)
Bonse2017 ; Graef2020 ; Rudenko2023 . Among the two, HSFLs are of particular
interest because they enable surface nanostructuring with a resolution well
below the optical diffraction limit Stoian2020 . Interestingly, to date HSFLs
have only been achieved using ultrashort femtosecond or picosecond laser
pulses.

The orientation of LIPSS depends on material properties and can be either
parallel or perpendicular to the laser polarization. The determination of
which orientation prevails occurs at the initial stage of laser-matter
interaction where anisotropy of light scattering by inevitable surface defects
leads to symmetry breaking. The distribution of scattered light intensity, and
consequently, the spatial modulation of local energy (fluence) imprinted on
the material surface through absorption, is dictated by the surface roughness.
The inhomogeneous distribution of absorbed light energy acts as a seed,
setting on the amplification process through the inter-pulse feedback loop,
thereby shaping the final morphology of LIPSS. Thus, an accurate description
of surface morphology and its influence on scattering and absorption of
incoming light is of paramount importance for understanding the physics
underlying the LIPSS generation.

Although there are a number of experiments on the influence of surface quality
on LIPSS formation, a thorough and comprehensive study in this area is clearly
lacking. It was shown that the surface roughness directly affects the period
of LIPSS and thus can alter certain surface functionality Fuentes2019 ;
Sotelo2023 . The formation of consistent LIPSS can be significantly enhanced
by polishing the surface in the direction perpendicular to the expected
orientation of LIPSS Ardron2014 . In contrast, deposition of thin films
through pulse vapor deposition provides a flexible method for preventing
aligned scratch and adjusting the initial nanostructures of the topography,
which in turn influences the feedback process significantly Prudent2021 .

In theoretical studies, the importance of surface roughness in explaining the
origin of LIPSS has become clear starting from the pioneering works of Sipe et
al. Sipe1983  who linked the appearance of LIPSS with the interference of the
incident laser radiation with surface electromagnetic waves generated by
scattering at the rough surface. Later, the crucial role of light interaction
with surface inhomogeneities in the formation of LIPSS has been confirmed
using ab-initio electromagnetic simulations, free from the limitations of the
Sipe’s theory Skolski2012 ; Skolski2014 ; Zhang2015 ; Zhang2020 . However, in
all of these simulations the surface roughness was considered simply as a
discrete binary mask — a collection of identical randomly distributed sub-
wavelength scatterers on top of a flat surface. Meanwhile, any realistic rough
surface possesses a continuous distribution of heights, which requires
description in terms of statistical properties like the height probability
distribution and the height correlation function. In this paper we apply such
statistical description of surface roughness in order to simulate the
scattering and subsequent absorption of ultrashort laser pulses by rough
surfaces. We study how various statistical properties of the surface roughness
affect the distribution of the absorbed light energy. In particular, we
demonstrate that different layers of surfaces with a continuous height
distribution scatter light in different directions. This observation allows us
to trace which surface layers are responsible for the formation of certain
features of the LIPSS shape. Finally, we explore how light polarization
affects the distribution of absorbed energy on rough surfaces. We show that
circularly polarized laser pulses allow to detect initial or gradually
developing asymmetries in surface roughness.

##  II Model of the surface roughness

A rough surface can be described by the roughness function $R(x,y)$, which
defines the random deviations of the surface height relative to some reference
plane Ogilvy1992 . With its help, the statistical properties of the surface
roughness can be expressed using the correlation function $C(X,Y)$ of the form

$\displaystyle C(X,Y)=\frac{\langle R(x,y)R(x+X,y+Y)\rangle}{\sigma^{2}},$ (1)

where $\sigma=\sqrt{\langle R^{2}\rangle}$ is the root-mean-square (rms)
height and $\langle\dots\rangle$ denotes the spatial averaging. The
correlation function $C(X,Y)$ describes the spatial coherence between
different surface points separated by a distance $d=\sqrt{X^{2}+Y^{2}}$.
Typical cases are the Gaussian correlation function

$\displaystyle C(X,Y)=\sigma^{2}\exp\left(-\frac{X^{2}+Y^{2}}{\xi^{2}}\right)$
(2)

and the exponential correlation function

$\displaystyle
C(X,Y)=\sigma^{2}\exp\left(-\frac{\sqrt{X^{2}+Y^{2}}}{\xi}\right)$ (3)

with $\xi$ being the correlation length.

In turn, the power spectral density $W(K_{x},K_{y})$ characterizes the average
distribution of the spatial frequency components of randomly fluctuated
surface profile. It is related to the correlation function $C(X,Y)$ by the
Fourier transform as follows:

$\displaystyle W(K_{x},K_{y})=$
$\displaystyle\qquad\frac{1}{(2\pi)^{2}}\iint_{-\infty}^{\infty}C(X,Y)e^{-i(K_{x}X+K_{y}Y)}dXdY,$
(4)

where $K_{x}$ and $K_{y}$ are the spatial frequencies.

To generate a rough surface numerically, we introduce a spectral function
$F(K_{x},K_{y})$ which we obtain from the power spectral density
$W(K_{x},K_{y})$ by modifying its amplitude and phase by a sequence of
normally distributed random numbers Thorsos1988 ; Mack2013 . This new function
allows us to calculate the roughness $R(x,y)$ using the following inverse
Fourier transform:

$\displaystyle
R(x,y)=\iint_{-\infty}^{\infty}F(K_{x},K_{y})e^{i(K_{x}x+K_{y}y)}dK_{x}dK_{y}.$
(5)

Function $F(K_{x},K_{y})$ is constructed as follows Thorsos1988 ; Mack2013 .
For all positive frequencies $K_{x}>0$, $K_{y}>0$ we define $F(K_{x},K_{y})$
as

$\displaystyle
F(K_{x},K_{y})=\sqrt{L_{x}L_{y}W(K_{x},K_{y})}~{}\frac{\eta_{1}+i\eta_{2}}{\sqrt{2}}$
(6)

and for zero and Nyquist frequencies ($K_{x}=0$, $K_{y}=0$ or
$K_{x}=\max(K_{x})$, $K_{y}=\max(K_{y})$) as

$\displaystyle F(K_{x},K_{y})=\sqrt{L_{x}L_{y}W(K_{x},K_{y})}~{}\eta_{3}.$ (7)

Here $L_{x}$ and $L_{y}$ are the lengths of $x$ and $y$ numerical domains
while each of $\eta_{1}$, $\eta_{2}$ and $\eta_{3}$ denotes a sequence of
normally distributed random numbers in the range $[0,1]$ with zero mean and
unity standard deviation. Since the surface roughness $R(x,y)$ must be a real
valued function, $F(K_{x},K_{y})$ have to satisfy the following additional
conditions:

$\displaystyle F(-K_{x},-K_{y})$ $\displaystyle=F^{*}(K_{x},K_{y}),$ (8)
$\displaystyle F(-K_{x},K_{y})$ $\displaystyle=F^{*}(K_{x},-K_{y}).$ (9)

In other words, $F(K_{x},K_{y})$ must be ”conjugate symmetric” about the
origin, i.e. the reflection of $F(K_{x},K_{y})$ at any point relative the
origin must be its complex conjugate.

Thus, in order to obtain the surface roughness function $R(x,y)$ with the
desired statistical properties we first select a correlation function
$C(X,Y)$. Then, using the Fourier transform, we calculate the power spectral
density $W(K_{x},K_{y})$. After, we construct the function $F(K_{x},K_{y})$
and, finally, applying the inverse Fourier transform to $F(K_{x},K_{y})$, we
calculate $R(x,y)$.

Figure 1 shows a typical surface roughness $R(x,y)$ generated by the above
approach for the Gaussian, Eq. (2), and exponential, Eq. (3), correlation
functions. For the both cases we use the same seed for the random numbers, the
rms height $\sigma=50$ nm and consider two different correlation lengths,
$\xi=100$ and 200 nm. In Figs. 1(a,b,c) we see that for the Gaussian
correlation function the correlation length $\xi$ determines the average scale
of surface inhomogeneities. In turn, Figs. 1(d,e,f) show that for the
exponential correlation function a faster decay of the spatial coherence
results in a surface with very rapidly varying small-scale perturbations
superimposed on the Gaussian roughness. In this case the correlation length
$\xi$ determines the amount of these perturbations and the sharpness of their
peaks.

Figure 1:  Surface roughness $R(x,y)$ for the Gaussian (a,b,c) and exponential
(d,e,f) correlation functions with the rms height $\sigma=50$ nm and two
different correlation lengths $\xi=100$ (b,e) and 200 nm (c,f). The line plots
(a,d) show the corresponding cross-sections at $y=0$.

##  III Numerical model

To study the interaction of ultrashort laser pulses with rough surfaces we
numerically solve the following system of Maxwell equations using the finite-
difference time-domain method (FDTD) Taflove2005 :

$\displaystyle\vec{\nabla}\times\vec{E}=-\frac{\partial\vec{B}}{\partial
t},\quad\vec{\nabla}\times\vec{H}=\frac{\partial\vec{D}}{\partial t},$ (10)

where $\vec{E}(x,y,z,t)$ and $\vec{H}(x,y,z,t)$ are the electric and magnetic
field vectors, $\vec{B}=\mu_{0}\vec{H}$, with $\mu_{0}$ being the vacuum
permeability, and
$\widetilde{D}=\varepsilon_{0}\varepsilon(\omega)\widetilde{E}$, where
$\widetilde{\hphantom{E}}$ means the temporal spectrum, while
$\varepsilon_{0}$ is the vacuum permittivity and $\varepsilon(\omega)$ is the
frequency-dependent permittivity of the medium.

In our simulations we use the calculation grid with the sizes
$L_{x}=L_{y}=12.2$ $\mu$m and $L_{z}=0.9$ $\mu$m in the $x$, $y$ and $z$
directions, respectively, with the corresponding step sizes $\Delta x=\Delta
y=10$ nm and $\Delta z=5$ nm. To avoid nonphysical reflections, at each end of
the grid we place convolutional perfectly matched layers of 0.1 $\mu$m
thickness.

As the incoming laser radiation we use a plane wave propagating along $z$
direction, whose temporal envelope $A(t)$ is given by a one period of
$\sin^{2}$ function: $A(t)=\sin^{2}(\pi t/4\tau_{0})$, where $\tau_{0}$ is the
laser pulse duration. In the case of a laser pulse polarized in the $x$
direction, we excite the $x$ component of the electric field:
$E_{x}(t)=A(t)\cos(\omega_{0}t)$, where $\omega_{0}=2\pi c_{0}/\lambda_{0}$ is
the central frequency with $\lambda_{0}$ and $c_{0}$ being the central
wavelength and the speed of light in vacuum. To model a circularly polarized
laser pulse we also add the $y$ component of the field with the relative phase
shift of $\pi/2$: $E_{y}(t)=A(t)\sin(\omega_{0}t)$. To preserve the pulse
energy, for circularly polarized pulses we divide $E_{x}$ and $E_{y}$ by
$\sqrt{2}$. In our simulations we use ultrashort laser pulses with the
duration $\tau_{0}=50$ fs and the central wavelength $\lambda_{0}=1.03$
$\mu$m.

At the distance of 0.6 $\mu$m below the radiation plane (as measured relative
to the reference plane of the roughness function $R(x,y)$), we place a semi-
infinite stainless steel medium with a rough surface. To model the dispersive
response of stainless steel we apply the auxiliary differential equation
method Taflove2005  and assume the Drude permittivity
$\varepsilon(\omega)=1-\omega_{p}^{2}/(\omega^{2}+i\omega\gamma)$ with the
plasma frequency $\omega_{p}=19.2\times 10^{15}$ 1/s and the damping rate
$\gamma=9.15\times 10^{15}$ 1/s Rudenko2019 . Under these parameters the
complex refractive index $n=n^{\prime}+in^{\prime\prime}$ of stainless steel
at central frequency $\omega_{0}$ has the real and imaginary parts equal to
$n^{\prime}=3.02$ and $n^{\prime\prime}=3.51$. The corresponding skin depth is
$l=1/(2n^{\prime\prime}\omega_{0}/c_{0})=23.34$ nm.

In our studies we are mainly interested in the distribution of the laser pulse
energy absorbed on the surface: the surface areas that have absorbed
sufficiently large amount of energy will later be extruded after the
thermalization and hydrodynamic stages of the laser-matter interaction.
Therefore any asymmetry in the distribution of the absorbed energy serves as a
seed for LIPSS growth. To estimate the amount of energy absorbed on the
surface let us first consider the power $P$ delivered by the laser pulse per
unit time per unit volume of the material: $P=\vec{E}\cdot\vec{J}$, where
$\vec{E}$ is the electric field of the laser pulse and $\vec{J}$ is the
current of the excited electrons Griffiths1999 . For the current $\vec{J}$ we
can write the Omh’s law in the form $\vec{J}=\sigma\vec{E}$, where the
conductivity $\sigma$, taken at the central frequency $\omega_{0}$, can be
expressed through the imaginary part of the permittivity
$\varepsilon=\varepsilon^{\prime}+i\varepsilon^{\prime\prime}$ as
$\sigma=\varepsilon_{0}\omega_{0}\varepsilon^{\prime\prime}$. Finally, we can
write $\varepsilon^{\prime\prime}=2n^{\prime}n^{\prime\prime}$, where
$n^{\prime}$ and $n^{\prime\prime}$ are the real and imaginary parts of the
refractive index $n=n^{\prime}+in^{\prime\prime}=\sqrt{\varepsilon}$.
Combining all these expressions together, we obtain
$P=2\varepsilon_{0}n^{\prime}n^{\prime\prime}\omega_{0}E^{2}$. Then, to obtain
the spatial distribution of the energy $Q$ delivered by the laser pulse per
unit volume we have to integrate $P$ over time:
$Q(x,y,z)=2\varepsilon_{0}n^{\prime}n^{\prime\prime}\omega_{0}\int_{-\infty}^{\infty}E(x,y,z,t)^{2}dt$.
Thus, the electric field $\vec{E}$, obtained during the FDTD simulations,
allows us to calculate the desired distribution of the absorbed energy $Q$ in
the entire material below its surface.

##  IV Results and discussion

First, let us study the distribution of the absorbed energy $Q(x,y,z)$ in
different layers of the rough surface exposed to an $x$-polarized laser pulse.
To model the rough surface we take the Gaussian correlation function with the
rms height $\sigma=50$ nm and the correlation length $\xi=100$ nm. Figure 2(a)
shows the $y=0$ cross-section of the corresponding roughness function
$R(x,y)$. In turn, Figs. 2(b,c,d) show the $xy$ slices of $Q$ in three
different surface layers $z=+70$, 0, -70 nm (see the dashed lines in Fig.
2(a)). Note that in all absorbed energy maps presented in the paper, in order
to better examine small details, we plot only the central part of the
calculation grid with the size 4 by 4 $\mu$m, while the full grid size is 12
by 12 $\mu$m. The white areas in Figs. 2(b,c,d) correspond to the regions
outside of the material, where $Q=0$. We observe that, as a result of near
field enhancement, the maxima of the absorbed energy $Q$ are concentrated in
the vicinity of the sharp edges of the surface.

Figures 2(f,g,h) show the spatial spectra of the $Q$ slices in the spatial-
frequency coordinates $k_{x}$ and $k_{y}$ normalized by the wave number
$k_{0}=\omega_{0}/c_{0}$. These spatial spectra give us an idea of typical
sizes and orientation of the absorbed energy spots. The spectral components of
$Q$ at spatial frequencies close to or less than $k_{0}$ correspond to the
large spots whose characteristic sizes are comparable to or exceeding the
laser wavelength $\lambda_{0}$. They are responsible for the formation of
LSFLs. In contrast, the spectral components of $Q$ at spatial frequencies
significantly exceeding $k_{0}$ correspond to the small-scale sub-wavelength
spots responsible for generation of HSFLs. For example, in Figs. 2(f,g,h) we
see that the highest frequency components of the spectra lie in the region of
$5k_{0}$. Therefore, we can conclude that in all three considered layers the
distribution of the absorbed energy $Q$ has morphological features with a
minimum size approximately five times smaller than the laser wavelength
$\lambda_{0}$.

Figure 2:  (a) The $y=0$ cross-section of the rough surface. The dashed lines
mark $z=\pm 70$ nm levels. The surface has the Gaussian correlation function
with the rms height $\sigma=50$ nm and the correlation length $\xi=100$ nm.
(b,c,d) The distribution of absorbed energy $Q$ in three surface layers
$z=+70$, 0, -70 nm. (e) The distribution of absorbed energy $Q_{\text{i}}$
integrated over all surface layers along $z$; the double arrow shows the
direction of laser polarization. (f–i) The spatial power spectra of the
corresponding $Q$ distributions. The arrows in (i) mark the characteristic
spectral patterns known as ”type-r” and ”type-s” features.

Figures 2(f,g,h) also show that the spectral distributions of the absorbed
energy $Q$ in different layers of the surface are not isotropic. The spectrum
corresponding to $z=+70$ nm layer is elongated in the $k_{y}$ direction (see
Fig. 2(f)), while the spectrum corresponding to $z=-70$ nm layer is elongated
in the perpendicular $k_{x}$ direction (see Fig. 2(h)); the spectrum for $z=0$
layer represents the transition between the two previous ones. Therefore,
taking into account that smaller features in the spatial domain correspond to
higher spectral frequencies, we can conclude that the spots of absorbed energy
$Q$ in different surface layers have, in average, different orientation. At
the top edge of the surface the spots of the absorbed energy are mainly
stretched in the $x$ direction, parallel to the laser polarization, while at
the bottom edge they line up in the perpendicular $y$ direction. In order to
explain this result we refer to our recent study where we showed that a single
nanobump or nanocavity, irradiated by linearly polarized laser light, both
have an asymmetric scattering pattern similar to that of a dipole Rudenko2019
. In particularly, we showed that nanobumps and nanocavities scatter light
predominantly in mutually perpendicular directions: nanobumps along the laser
polarization and nanocavities in the perpendicular one (see Fig. 3 in
Rudenko2019 ). In turn, as shown in Fig. 1(a), any rough surface can be
considered as a collection of peaks (bumps) and valleys (cavities). Therefore,
at the upper edge of the surface, consisting of peaks, light is mainly
scattered parallel to the direction of polarization, while at the lower edge
of the surface, populated by valleys, the main direction of light scattering
is perpendicular to the polarization direction. Thus, the upper and lower
edges of the rough surface tend to form LIPSS with mutually perpendicular
orientations.

The final shape of LIPSS results from the accumulation of absorbed energy in
all surface layers. Therefore, the ultimate quantity of relevance for
predicting light-surface coupling is the absorbed energy $Q_{\text{i}}$
integrated over all $z$ layers:
$Q_{\text{i}}(x,y)=\int_{-\infty}^{\infty}Q(x,y,z)dz$. Figure 1(e) shows the
distribution of this integrated absorbed energy $Q_{\text{i}}$ where we see
that the regions of high losses are mainly elongated along $x$ direction.
Therefore, we can expect that in this particular case the resulting LIPSS will
be oriented parallel to the laser polarization. In turn, in Fig. 1(i) we see
that the spatial spectrum of $Q_{\text{i}}$ has a well-recognizable shape and
contains spectral features known as ”type-r” and ”type-s” Bonse2020 . The
above analysis of absorbed energy spectra across various surface layers leads
to the conclusion that type-r features primarily originate from the upper edge
of the surface, where the scattering of light is govern by surface peaks.
Conversely, type-s features predominantly arise from the lower layers of the
surface, where light scatters by the surface valleys.

To understand how different types of correlation functions affect light
absorption, let us consider two rough surfaces with the Gaussian, Eq. (2), and
exponential, Eq. (3), correlation functions. Figure 3 shows the difference
between the corresponding distributions of the integrated absorbed energy
$Q_{\text{i}}$ and their spatial spectra. We can see that compared to the case
of the Gaussian correlation function, the distribution of $Q_{\text{i}}$ at
the surface with the exponential correlation function has much more fine
details. This difference in the absorbed energy distributions is reflected in
their spatial spectra, where we see that although both spectra are
qualitatively very similar, the one corresponding to the exponential
correlation function has much higher spectral components in the $k_{y}$
direction. This observation is expected since, as demonstrated earlier, the
higher frequencies along $k_{y}$ arise from light scattered at the peaks in
the top edge of the rough surface. As depicted in Fig. 1, for the surface with
the exponential correlation function, these peaks are more concentrated and
smaller in scale.

Figure 3:  The distributions of integrated absorbed energy $Q_{\text{i}}$
(a,b) and their spectra (c,d) for rough surfaces described by the Gaussian
(a,c) and exponential (b,d) correlation functions with the rms height
$\sigma=50$ nm and correlation length $\xi=100$ nm. The double arrow shows the
direction of laser polarization.

Next, we consider a surface with the Gaussian correlation function and examine
how changes in its rms height $\sigma$ and correlation length $\xi$ affect the
absorbed energy $Q_{\text{i}}$. Figure 4 shows the distributions of
$Q_{\text{i}}$ and their spectra for different values of $\sigma$ and $\xi$.
In Fig. 4(b) we can see that the increase of the correlation length $\xi$
leads to the compression of the spectra in the $k_{y}$ direction. This
compression corresponds to an increase in the scale of the absorbed energy
spots in the $y$ direction (see Fig. 4(a)). Interestingly, at $\xi=200$ nm,
the spectra become equally elongated in the $k_{x}$ and $k_{y}$ directions
which means that the average sizes of the $Q_{\text{i}}$ spots in the $x$ and
$y$ directions become the same (see the right columns in Fig. 4(a) and (b)).
Thus, with the increase of the correlation length, we observe a transition
from the distribution of absorbed energy oriented along the laser polarization
to a speckle-like pattern.

Meanwhile, Fig. 4 shows that the increase of the rms height $\sigma$ leads to
an increase in the total absorption: the spatial distributions of
$Q_{\text{i}}$ and their spectra become more intense at higher $\sigma$
values. Additionally, in Fig. 4(b) we see that with the increase of $\sigma$
the spectral features of type-s (see Fig. 2(i)) gradually spread out in the
$k_{y}$ direction. In the case of the smallest $\xi=50$ nm and the largest
$\sigma=100$ nm the type-s features completely disappear. With such $\sigma$
and $\xi$, the surface consists of very narrow (small correlation length) and
very deep (large rms height) inhomogeneities. As a result, little light
penetrates to the bottom of such elongated pits where it can be scattered by
the surface valleys in the right direction.

Figure 4:  The distributions of integrated absorbed energy $Q_{\text{i}}$ (a)
and their spectra (b) for rough surfaces described by the Gaussian correlation
function with different values of the rms height $\sigma$ and the correlation
length $\xi$. The double arrow shows the direction of laser polarization.

Now, let us consider a surface with the Gaussian correlation function of the
form

$\displaystyle
C(X,Y)=\sigma^{2}\exp\left(-\frac{X^{2}}{\xi_{x}^{2}}-\frac{Y^{2}}{\xi_{y}^{2}}\right)$
(11)

with the rms height $\sigma=50$ nm and two different correlation lengths
$\xi_{x}$ and $\xi_{y}$ along the $x$ and $y$ directions. Here we still use
$x$-polarized incident laser pulses. The correlation function given by Eq.
(11) describes a surface with initially broken symmetry where the scale of the
inhomogeneities in one direction is larger then in another. In other words,
such surface initially have a ripple-like structure in the $x$ or $y$
direction. It may correspond to a surface subjected to mechanical polishing
along one specific direction or it may represent a surface that has already
been irradiated by a small number of preceding laser pulses. Figure 5 shows
the distributions of integrated absorbed energy $Q_{\text{i}}$ together with
their spectra for different values of $\xi_{x}$ and $\xi_{y}$.

Let us first take a look at the case with the largest $\xi_{x}=200$ nm and
$\xi_{y}$ decreasing from 200 to 50 nm (the last columns in Fig. 5). This case
describes a transition from a surface with isotropic roughness to one where
surface inhomogeneities, being compressed in the $y$ direction, turn out to be
elongated in the $x$ direction, parallel to the polarization. The spectra in
Fig. 5(b) show that in this case the decrease in $\xi_{y}$ leads to higher
spectral components in the $k_{y}$ direction, which in the real domain
corresponds to the compression of the absorbed energy spots along $y$.

Now, let us consider the scenario where $\xi_{y}=200$ nm is the largest and
$\xi_{x}$ decreases from 200 to 50 nm (the first rows in Fig. 5). Here we see
a transition from the isotropic roughness to the one where the surface
inhomogeneities, being compressed in the $x$ direction, become elongated along
$y$, in the direction perpendicular to the laser polarization. As in the
previous case, the spectra in Fig. 5(b) show that decreasing $xi_{x}$ leads to
the appearance of higher frequency components, but this time in the $y$
direction. The corresponding distribution of the absorbed energy becomes
elongated along $y$.

The rest of Fig. 5 confirms the previously observed trends: as surface
irregularities extend further along one of the $x$ or $y$ directions, the
absorbed energy spots become increasingly elongated in the corresponding
direction. Thus, any asymmetry in surface roughness, whether it originated as
an inherent feature of the surface or arose after laser exposure, will be
magnified by subsequent laser pulses. This mechanism constitutes the
foundation of the inter-pulse feedback loop, which ultimately leads to the
formation of LIPSS from even the slightest initial symmetry breaking.
Moreover, the random nature of the surface roughness ensures that no two LIPSS
patterns are perfectly alike. This highlights the individuality of each
imprint left by exposure to light, evoking to the distinctiveness of a
fingerprint, which offers anti-counterfeiting abilities.

Figure 5:  The distributions of integrated absorbed energy $Q_{\text{i}}$ (a)
and their spatial spectra (b) for rough surfaces described by an asymmetric
Gaussian correlation function with the fixed rms height $\sigma=50$ nm and
different correlation lengths $\xi_{x}$ and $\xi_{y}$ along $x$ and $y$
directions. The double arrow shows the direction of laser polarization.

Finally, we study how the polarization of the incoming light affects the
distribution of the integrated absorbed energy $Q_{\text{i}}$. Figure 6 shows
the distributions of $Q_{\text{i}}$ and their spectra for $x$-polarized,
$y$-polarized and circularly polarized laser pulses. In Figs. 6(a,b) and (d,e)
we see that in the case of $x$\- and $y$-polarized laser pulses the shapes of
the $Q_{\text{i}}$-distributions and their spectra are the same but rotated by
90 degrees relative to each other. In other words, as expected, the
orientation of LIPSS follows the direction of the laser polarization. However,
Figs. 6(c) and (f) show that in the case of circular polarization the
distribution of $Q_{\text{i}}$ does not show any preferred orientation and the
corresponding spatial spectrum is symmetrical about the origin. This can be
explained by the fact that for circularly polarized light, when the electric
field vector rotates, it continuously passes through all possible
orientations. The absence of a preferred orientation in the spectral domain
means that the distribution of the absorbed energy in the real space is
isotropic. In Fig. 6(c) we see that the spots of the absorbed energy line up
in a circle-like pattern. Therefore, we can expect that the resulting LIPSS
will exhibit the same radial symmetry.

Figure 6:  The distributions of integrated absorbed energy $Q_{\text{i}}$
(a,b,c) and their spectra (d,e,f) for $x$-polarized (a,d), $y$-polarized (b,e)
and circularly polarized (c,f) laser pulses. The arrows show the direction of
each laser polarization.

However, what happens if a circularly polarized laser pulse interacts with a
surface that initially has some asymmetry in its roughness? In this case, the
breaking of symmetry should occur not due to a dedicated direction of laser
polarization (circular polarization does not have a preferred direction), but
due to the properties of the surface itself. In order to answer this question
let us examine a rough surface with an asymmetric Gaussian correlation
function defined by Eq. (11). We fix the rms height $\sigma=50$ nm and the $y$
correlation length $\xi_{y}=100$ nm, and consider the $x$ correlation length
values $\xi_{x}=100$, 200, and 400 nm. While $\xi_{x}=100$ nm defines a
surface with an isotropic roughness, $\xi_{x}=200$ and 400 nm values
correspond to surfaces whose inhomogeneities are stretched in the $x$
direction by two and four times, respectively. Figure 7 shows the resulting
distributions of the integrated absorbed energy $Q_{\text{i}}$ and their
spectra. We can see that the increase of the $x$ correlation length leads to
the elongation of the absorbed energy spots in the $x$ direction. Meanwhile,
the corresponding spectra compress in the $k_{x}$ direction, remaining the
same in the $k_{y}$ direction. This means that the increase of the correlation
length in one direction can change the size of the absorbed energy spots only
in that same direction. Moreover, we see that the compression of spectra in
the $k_{x}$ direction is directly proportional to the increase in $\xi_{x}$:
twice larger $\xi_{x}$ corresponds to twice narrower spectrum. Therefore, by
using circularly polarized laser pulses, we can effectively test the asymmetry
of surface roughness.

Figure 7:  The distributions of integrated absorbed energy $Q_{\text{i}}$
(a,b,c) and their spectra (d,e,f) for a circularly polarized laser pulse
interacting with rough surfaces described by an asymmetric Gaussian
correlation function with the fixed rms height $\sigma=50$ nm and $y$
correlation length $\xi_{y}=100$ nm but with different $x$ correlation
lengths: $\xi_{x}=100$ nm (a,d), $\xi_{x}=200$ nm (b,e) and $xi_{x}=400$ nm
(c,f). The arrow show the direction of laser polarization.

Interestingly, the compression of the absorbed energy spectra occurs only due
to high-frequency components, responsible for HSFLs formation. In turn, the
low-frequency components responsible for LSFLs (the central ring structures in
Figs. 7(d–f) at frequencies close to $k_{0}$) do not change with the
increasing asymmetry of the surface. Taking into account that HSFLs need
ultrashort laser pulses to be formed, we can conclude that sensing of the
surface roughness asymmetries using circular polarization would be possible
only with picosecond or femtosecond laser pulses.

##  V Conclusion

In conclusion, we use ab-initio electromagnetic simulations to study the
interaction of ultrashort laser pulses with rough surfaces. We investigate the
spatial distribution of the absorbed laser energy in order to better
understand the underlying physics of LIPSS formation. In contrast to existing
studies, to model realistic rough surfaces we apply a statistical description
of distribution of their heights. Such approach allows us to introspect the
distribution of the absorbed laser energy across the surface selvedge. As a
result, we discovered that the high-frequency (type-s) and low-frequency
(type-r) components of LIPSS morphology emerge in distinct surface layers:
type-s features predominantly appear in the upper edge, while type-r features
are mainly present in the lower edge of the surface. We investigate how
different statistical properties of surface roughness affect the distribution
of absorbed laser energy. Our findings demonstrate that as the correlation
length increases, the pattern of absorbed energy spots transitions from the
one aligned with the laser polarization to a speckle-like distribution.
Moreover, we demonstrate that an increase in the rms height of surface
roughness leads to the suppression of the type-r LIPSS features up to their
complete disappearance. We show that the lack of the dedicated direction of
the electric field vector in circularly polarized pulses leads to the
formation of radially symmetric distribution of the absorbed energy. The
resulting LIPSS can be expected to inherit this symmetry and form a pattern
without any preferential orientation. The symmetry breaking induced by the
roughness under circularly polarized laser pulses can be exploited to detect
and measure asymmetries in surface roughness. In the case of circular
polarization, changes in the scale of surface inhomogeneities lead to a
directly proportional deformation of the absorbed energy spectrum. Our
findings allow us to deepen the understanding of the physics behind the
formation of LIPSS and offer valuable guidance for precision laser-induced
surface fabrication, with potential applications ranging from tailored
photonics devices to advanced materials design.

###### Acknowledgements.

This work was supported by the LABEX MANUTECH-SISE (ANR-10-LABX- 0075) of
Université de Lyon, within the Plan France 2030 operated by the French
National Research Agency (ANR).

## References

  * (1) J. Bonse, S. Hohm, S. V. Kirner, A. Rosenfeld, and J. Kruger, “Laser-induced periodic surface structures— a scientific evergreen,” IEEE Journal of Selected Topics in Quantum Electronics, vol. 23, no. 3, 2017.
  * (2) S. Gräf, “Formation of laser-induced periodic surface structures on different materials: fundamentals, properties and applications,” Advanced Optical Technologies, vol. 9, no. 1–2, pp. 11–39, 2020.
  * (3) J. Bonse and S. Gräf, “Maxwell meets Marangoni — a review of theories on laser‐induced periodic surface structures,” Laser & Photonics Reviews, vol. 14, no. 10, p. 2000215, 2020.
  * (4) A. Rudenko and J.-P. Colombier, How Light Drives Material Periodic Patterns Down to the Nanoscale, pp. 209–255.  Springer International Publishing, 2023.
  * (5) A. Y. Vorobyev and C. Guo, “Colorizing metals with femtosecond laser pulses,” Applied Physics Letters, vol. 92, no. 4, 2008.
  * (6) B. Dusser, Z. Sagan, H. Soder, N. Faure, J.-P. Colombier, M. Jourlin, and E. Audouard, “Controlled nanostructrures formation by ultra fast laser pulses for color marking,” Optics Express, vol. 18, no. 3, p. 2913, 2010\.
  * (7) G. Li, J. Li, Y. Hu, C. Zhang, X. Li, J. Chu, and W. Huang, “Femtosecond laser color marking stainless steel surface with different wavelengths,” Applied Physics A, vol. 118, no. 4, pp. 1189–1196, 2015.
  * (8) V. Zorba, E. Stratakis, M. Barberoglou, E. Spanakis, P. Tzanetakis, and C. Fotakis, “Tailoring the wetting response of silicon surfaces via fs laser structuring,” Applied Physics A, vol. 93, no. 4, 2008.
  * (9) A.-M. Kietzig, S. G. Hatzikiriakos, and P. Englezos, “Patterned superhydrophobic metallic surfaces,” Langmuir, vol. 25, no. 8, pp. 4821–4827, 2009.
  * (10) C. Kunz, F. Müller, and S. Gräf, “Multifunctional hierarchical surface structures by femtosecond laser processing,” Materials, vol. 11, no. 5, p. 789, 2018.
  * (11) J. Yu and Y. Lu, “Laser-induced ripple structures on Ni–P substrates,” Applied Surface Science, vol. 148, no. 3–4, pp. 248–252, 1999.
  * (12) J. Bonse, S. Kirner, M. Griepentrog, D. Spaltmann, and J. Krüger, “Femtosecond laser texturing of surfaces for tribological applications,” Materials, vol. 11, no. 5, p. 801, 2018.
  * (13) C. Kunz, J. Bonse, D. Spaltmann, C. Neumann, A. Turchanin, J. F. Bartolomé, F. A. Müller, and S. Gräf, “Tribological performance of metal-reinforced ceramic composites selectively structured with femtosecond laser-induced periodic surface structures,” Applied Surface Science, vol. 499, p. 143917, 2020.
  * (14) S. A. Jalil, M. Akram, J. A. Bhat, J. J. Hayes, S. C. Singh, M. ElKabbash, and C. Guo, “Creating superhydrophobic and antibacterial surfaces on gold by femtosecond laser pulses,” Applied Surface Science, vol. 506, p. 144952, 2020.
  * (15) L. Sotelo, T. Fontanot, S. Vig, P. Herre, P. Yousefi, M. H. Fernandes, G. Sarau, G. Leuchs, and S. Christiansen, “Influence of initial surface roughness on LIPSS formation and its consecutive impact on cell/bacteria attachment for TiAl6V4 surfaces,” Advanced Materials Technologies, vol. 8, no. 12, 2023.
  * (16) J. Outón, M. Carbú, M. Domínguez, M. Ramírez-del Solar, G. Alba, M. Vlahou, E. Stratakis, V. Matres, and E. Blanco, “Size matters: how periodicity and depth of LIPSS influences E. coli adhesion on ferritic stainless steel,” Applied Surface Science, vol. 663, p. 160225, 2024.
  * (17) R. Stoian and J.-P. Colombier, “Advances in ultrafast laser structuring of materials at the nanoscale,” Nanophotonics, vol. 9, no. 16, pp. 4665–4688, 2020.
  * (18) Y. Fuentes-Edfuf, J. A. Sánchez-Gil, C. Florian, V. Giannini, J. Solis, and J. Siegel, “Surface plasmon polaritons on rough metal surfaces: Role in the formation of laser-induced periodic surface structures,” Acs Omega, vol. 4, no. 4, pp. 6939–6946, 2019.
  * (19) M. Ardron, N. Weston, and D. Hand, “A practical technique for the generation of highly uniform LIPSS,” Applied surface science, vol. 313, pp. 123–131, 2014.
  * (20) M. Prudent, F. Bourquard, A. Borroto, J.-F. Pierson, F. Garrelie, and J.-P. Colombier, “Initial morphology and feedback effects on laser-induced periodic nanostructuring of thin-film metallic glasses,” Nanomaterials, vol. 11, no. 5, p. 1076, 2021.
  * (21) J. E. Sipe, J. F. Young, J. S. Preston, and H. M. van Driel, “Laser-induced periodic surface structure. I. Theory,” Physical Review B, vol. 27, no. 2, pp. 1141–1154, 1983.
  * (22) J. Z. P. Skolski, G. R. B. E. Römer, J. V. Obona, V. Ocelik, A. J. Huis in ’t Veld, and J. T. M. De Hosson, “Laser-induced periodic surface structures: Fingerprints of light localization,” Physical Review B, vol. 85, p. 075320, 2012.
  * (23) J. Z. P. Skolski, G. R. B. E. Römer, J. Vincenc Obona, and A. J. Huis in ’t Veld, “Modeling laser-induced periodic surface structures: Finite-difference time-domain feedback simulations,” Journal of Applied Physics, vol. 115, p. 103102, 2014.
  * (24) H. Zhang, J.-P. Colombier, C. Li, N. Faure, G. Cheng, and R. Stoian, “Coherence in ultrafast laser-induced periodic surface structures,” Physical Review B, vol. 92, p. 174109, 2015.
  * (25) H. Zhang, J.-P. Colombier, and S. Witte, “Laser-induced periodic surface structures: Arbitrary angles of incidence and polarization states,” Physical Review B, vol. 101, p. 245430, 2020.
  * (26) J. A. Ogilvy, Theory of wave scattering from random rough surfaces.  Bristol, England: Institute of Physics Publishing, reprinted ed., 1992\.
  * (27) E. I. Thorsos, “The validity of the Kirchhoff approximation for rough surface scattering using a Gaussian roughness spectrum,” The Journal of the Acoustical Society of America, vol. 83, no. 1, pp. 78–92, 1988.
  * (28) C. A. Mack, “Generating random rough edges, surfaces, and volumes,” Applied Optics, vol. 52, no. 7, p. 1472, 2013.
  * (29) A. Taflove and S. C. Hagness, Computational electrodynamics.  Norwood, MA: Artech House, 3rd ed., 2005.
  * (30) A. Rudenko, C. Mauclair, F. Garrelie, R. Stoian, and J. Colombier, “Light absorption by surface nanoholes and nanobumps,” Applied Surface Science, vol. 470, pp. 228–233, 2019.
  * (31) D. J. Griffiths, Introduction to electrodynamics, ch. 8.  New Jersey: Prentice Hall, 1999.


# Entanglement in Resonance Fluorescence

Juan Camilo López Carreño<EMAIL_ADDRESS>Institute of Theoretical
Physics, University of Warsaw, ul. Pasteura 5, 02-093, Warsaw, Poland
Santiago Bermúdez Feijoo  Departamento de Física, Universidad Nacional de
Colombia, Ciudad Universitaria, K. 45 No. 26–85, Bogotá D.C., Colombia
Magdalena Stobińska  Faculty of Mathematics, Informatics and Mechanics,
University of Warsaw, ul. Banacha 2, 02-097 Warsaw, Poland

(August 28, 2024)

###### Abstract

Particle entanglement is a fundamental resource upon which are based many
quantum technologies. However, the up-to-now best sources of entangled photons
rely on parametric down-conversion processes, which are optimal only at
certain frequencies, which rarely match the energies of condensed-matter
systems that can benefit from entanglement. In this Article, we show a way to
circumvent this issue, and we introduce a new source of entangled photons
based on resonance fluorescence delivering photon pairs as a superposition of
vacuum and the Bell state $\left|\Phi^{-}\right\rangle$. Our proposal relies
on the emission from the satellite peaks of a two-level system driven by a
strong off-resonant laser, whose intensity controls the frequencies of the
entangled photons. Furthermore, the degree of entanglement can be optimized
for every pair of frequencies, thus demonstrating a clear advantage over
existing technologies. Finally, we illustrate the power of our novel source of
entangled single-photon pairs by exciting a system of polaritons and showing
that they are left in a maximally entangled steady state.

Resonance fluorescence, the interaction between an artificial atom and
coherent light, has been the subject of fundamental research from the early
stages of quantum optics 1, 2. In particular, the observation of photon
antibunching from this interaction 3 paved the way for investigations
regarding the quantum character of light, and the phenomena that has enabled,
for instance, the pursue of single-photon transistors 4. In turn, single
photon emission is a key element for many quantum information technologies 5,
6, allowing the possibility to design protocols for, e.g., quantum
teleportation 7 and quantum cryptography 8. Usually, the source of light is a
laser (with a well defined energy) which effectively matches a single energy
transition of the artificial atom, and therefore, in practice, one deals with
the excitation of a so-called “two-level system” (2LS). Thus, the laser can
only induce a single excitation in the atom at the time, and which has led
resonance fluorescence to be regarded as an ultrabright source of quantum
light 9, 10, 5, with high single-photon purity 11, 12. However, recent
investigations that analysed the luminescence with spectral resolution 13,
have found that the emission form a 2LS actually consists of multiple highly-
correlated photons 14, 15, 16, 17, 18, 19, 20.

Figure 1: (Color online). The Mollow triplet driven out of resonance. (a)
Emission spectrum when the laser is resonant (dashed) and detuned (solid) from
the 2LS. The emission lines are identified with the four possible transitions
(two of them are degenerate) between the dressed states (shown in the inset).
(b) Scheme of our proposed source of entangled photon pairs emitted from the
sidebands of the Mollow triplet. (c) Energy transitions enabling the emission
from the 2LS as a function of the driving intensity. In the detuned case
(solid, dark lines), and as opposed to the resonant one (dashed, light lines),
the satellite peaks become the dominant feature of the spectrum, and the
dynamics is given by transitions that change the quantum state of the 2LS,
i.e., by transitions of the type
$\left|\pm\right\rangle\rightarrow\left|\mp\right\rangle$ (shown in red and
blue). For the figure we used $\gamma_{\sigma}$ as the unit,
$\Omega/\gamma_{\sigma}=4$ (marked in panel (c) as an horizontal gray line)
and $(\omega_{\sigma}-\omega_{\mathrm{L}})/\gamma_{\sigma}=25/2$.

Such a multi-photon structure is particularly revealed when the intensity of
the driving laser is strong and the 2LS enters into the so-called _Mollow
regime_ 21, in which the emission spectrum of the 2LS consists of a triplet,
as illustrated in dashed lines in Fig. 1(a). Notably, although the photons
emitted from the 2LS are perfectly antibunched, selecting particular frequency
regions of the emission allows to unveil a richer landscape of photon
correlations 14, 15, 16—ranging from antibunching to superbunching statistics,
passing through thermal and uncorrelated light. In fact, the statistical
variability of the photons emitted by resonance fluorescence allows to design
exotic sources of light 22, excite other optical targets 23, 24, 25, and
perform the so-called _Mollow spectroscopy_ 23, whereby, e.g., the internal
structure of complex and highly-dissipative quantum systems can be probed with
a minimal amount of photons, namely, two.

In this Article, we approach a different type of quantum correlations, namely
entanglement, which has been observed in Resonance Fluorescence when either a
collection of atoms is considered 26 or its biexciton structure is taken into
account 27 and the 2LS is coupled to a microcavity 28, 29, 9, 30, 31, 32.
Here, we demonstrate that time-frequency entanglement can be extracted from
resonance fluorescence, specifically, from the emission of a Mollow triplet,
without the need of coupling it to other objects or taking into account any
additionally internal structures of the 2LS. Instead, one only needs to
include the observation of the emission into the description. In fact, we show
that when the 2LS is driven out of resonance, the emission from the sidebands
of the triplet behave as a heralded source of entangled photon pairs, as
sketched in Fig. 1(b).

An important advantage that our proposal has over existing technologies, e.g.,
sources based on parametric downconversion processes is that the entangled
photons pairs emitted by our source follow an antibunched (and not an
uncorrelated) statistics. Thus, while in sources based on SPDC a higher
intensity of driving produces higher-order processes, our source produces more
pairs of entangled photons. Finally, noting that the Mollow regime has been
successfully realized in a variety of systems, including quantum dots 33, 34,
35, 36, 37, 38, molecules 39, cold atom ensembles 40, confined single atoms
41, photonic chips 42, and superconducting qubits 43, 44, 45 (and 2LSs can
also be constructed in other platforms including superconducting circuits 46,
47, 48, 49, 50, 51, 52 and photonic structures 53, 54, 55), our source is able
to operate on a wide gamut of frequencies and to interface with, e.g.,
condensed-matter systems that can benefit from entangled excitation 56, thus
making resonance fluorescence a compelling alternative to the existing sources
of entangled photons.

The rest of this Article is organized as follows: We first demonstrate that
energy-time entangled between photons emitted from the sidebands of the 2LS is
unveiled simply by including the observation of the light into the description
of our system. Next, we use a quantum Monte Carlo experiment to demonstrate
that, as a consequence of the detuned excitation, the entangled photon pairs
are heralded. Finally, we show that our source is able to drive complex
condensed-matter systems (e.g., exciton-polaritons 57) into a maximally
entangled steady state, despite them being immersed within a highly-
dissipative environment.

Figure 2: (Color online). Heralded single photons observed through a
frequency-resolved quantum Monte Carlo experiment. (a, b) Crossed-correlations
between photons detected with frequencies $\omega_{+}$ and $\omega_{-}$,
showing the agreement with theoretical prediction (solid blue lines) and the
quantum Monte Carlo experiment (black bars). When the driving laser is
resonant to the 2LS, the correlation function is completely symmetric (cf.
panel a). When the driving is taken out of resonance, the shape of the
correlation function resembles a $\lambda$, indicating that the emission of a
photon with frequency $\omega_{+}$ heralds the emission of a photon with
frequency $\omega_{-}$ (cf. panel b). (c) Ratio
$r(n,\tau)=p(n,\tau,\Delta)/p(n,\tau,0)$ of the probability to detect one
(blue) and two or more (red) photons after detecting the heralding photon,
when the driving is made out of resonance and in resonance. Taking the laser
out of resonance enhances the single-photon heralding probability within a
time-window $\tau^{\ast}$ of almost two lifetimes of the 2LS, as indicated by
the dashed vertical line. Additionally, the probability to herald two or more
photons is suppressed when the driving is done out of resonance. For the
figures we used $\gamma_{\sigma}$ as the unit, and the parameters that
optimize entangled emission from the sidebands of the triplet (cf. Section II
of the Supplemental Material); namely $\Omega/\gamma_{\sigma}=1$,
$\Gamma/\gamma_{\sigma}=1$, and $\Delta=1.85\Omega$.

## Results

Figure 3: (Color online). Characterisation of our source of entangled pairs of
single photons. (a) Maximum logarithmic negativity that can be reached as a
function of the linewidth of the detectors and the intensity of the driving
laser, which indicates that the optimal condition for extracting entanglement
is to have narrow detectors observing a Mollow triplet with a very large
splitting (namely, the bottom right corner of the panel). (b) The optimum
detuning at which the 2LS must be driven depends both on the intensity of the
driving and the linewidth of the detectors or, equivalently, the optical
targets receiving the entangled photons. (c) Violation of the CSI showing the
regions where the negativity obtained in panel (a) is not an artifact of the
competition between the detectors for the photons of the sidebands. (d) The
ratio between the linewidth of the detectors and the splitting between the
sidebands. It shows that the best configurations are those for which this
ratio is as small as possible, and that it has to be less than 1 to violate
the CSI. (e) Emission rate from the detectors, thus completing the mapping of
the quality and the brightness of the source of entangled single photons based
on the Mollow triplet. The red, diagonal, dotted line on each panel indicates
the boundary above which the emission does not violate the CSI.

Measurement of the photons
A key aspect of quantum mechanics is that measurements affect the quantum
state of the system under observation. Thus, a correct description of the
emission from the 2LS should also incorporate in the dynamics the effect of
the observation by a physical detector, namely one able to observe photons of
a given frequency $\tilde{\omega}$ with a finite linewidth $\tilde{\Gamma}$.
Such a description can be done with the theory of frequency-resolved
correlations 13, whereby detectors are considered as quantum mechanical
objects that receive the excitation from the source of light without returning
any feedback. In practice, this can be achieved by coupling the emitter and
the detectors with a vanishing strength 13, or by using the formalism of
_cascaded systems_ 58, 59. While these two approaches yield exactly the same
normalized correlations (i.e., $n^{\textrm{th}}$-order correlation functions)
60, the latter also provides the correct un-normalized correlations (e.g., the
mean population or the emission rate of the system), and therefore it is the
implementation that we will use throughout this Article. Using the master
equation outlined in the Methods section, we gain access to the correlations
between photons from the Mollow triplet emitted at frequencies $\omega_{1}$
and $\omega_{2}$. In particular, we are interested in the second-order
correlation between photons emitted at the frequencies of the satellite peaks
(with energies $\omega_{\pm}$, c.f. the scheme in Fig. 1, associated to
bosonic annihilation operator $a_{1}$ and $a_{2}$) separated by a time $\tau$,
i.e., we compute
$g_{12}^{(2)}(\tau)=\langle{a_{1}^{\dagger}}{a_{2}^{\dagger}}(\tau)a_{2}(\tau)a_{1}\rangle/(\langle{a_{1}^{\dagger}}a_{1}\rangle\langle{a_{2}^{\dagger}}a_{2}\rangle)$
for the case where $\omega_{1}=\omega_{+}$ and $\omega_{2}=\omega_{-}$. Note
that while the inclusion of detection into the description of the emission
leads to a loss of antibunching of the signal 61, it also unveils the multi-
photon structure behind the dynamics of the 2LS 16. When the latter is driven
resonantly by a laser, its $g_{12}^{(2)}(\tau)$ is completely symmetric, as
shown in Fig. 2(a). Such a shape is an indication that the emission from the
high-energy sideband (with frequency $\omega_{+}$) can occur either before or
after the emission from the low-energy peak (with frequency $\omega_{-}$);
i.e., there is not a causal relation between consecutive emissions. Such a
symmetry is broken when the driving laser (with frequency
$\omega_{\mathrm{L}}$) is taken out of resonance from the 2LS (with natural
frequency $\omega_{\sigma}$). Figure 2(b) shows the case for
$\Delta\equiv(\omega_{\sigma}-\omega_{\mathrm{L}})=1.85\Omega$, where it is
clear that the emission from the high-energy sideband occurs _after_ the
emission of a photon from the low-energy one (the reverse situation occurs if
the detuning becomes negative). While the type of correlation shown in Fig.
2(b) is indicative of a heralded emission, it does not necessarily imply that
the emission consists of single photons. To show that this is in fact the
case, we performed a frequency-resolved quantum Monte Carlo simulation 60 of
the emission from the two sidebands of the Mollow triplet to compare the cases
when the 2LS is driven in and out of resonance. Thus, from the simulations we
are able to obtain the probability $p(n,\tau,\Delta)$ to detect $n$ photons of
energy $\omega_{+}$ within a time-window $\tau$ after a photon of energy
$\omega_{-}$ has been measured, provided that the detuning between the 2LS and
the laser is $\Delta$. Figure 2(c) shows the ratios
$r(n,\tau)=p(n,\tau,\Delta)/p(n,\tau,0)$, which illustrate that the detuning
enhances the probability to detect one photon by about 20% (blue line) while
simultaneously decreases the probability to detect two or more photons (red
line). Note that the ratio $r(1,\tau)$ becomes less than one for time windows
larger than $\tau^{\ast}\approx 2/\gamma_{\sigma}$, shown as a vertical dashed
line in Fig. 2(c). Thus, the heralded photons are more likely to be emitted
within the time window $\tau\leq\tau^{\ast}$ when the 2LS is driven out of
resonance. Together, Fig. 2(b) and (c) are the evidence that demonstrates that
a 2LS driven out of resonance by a laser is a source of heralded single
photons, whose frequencies correspond to the energies of the sidebands of the
Mollow triplet.

Source of entangled photon pairs
Now that we have established that the emission from the satellite peaks of the
detuned Mollow triplet is composed of highly correlated pairs of single
photons, and that their emission can be observed in a heralded fashion, we
investigate another type of quantum correlation: entanglement. It has been
theoretically predicted 17 and experimentally observed 15 that pairs of
photons emitted from the Mollow spectrum at various frequencies violate the
Cauchy-Schwarz inequality (CSI), which can only happen in systems displaying
entanglement 62. Thus, in our Article we quantify entanglement through the so-
called _logarithmic negativity_ $\mathcal{N}(\rho)$ 63, 64, 65, which is an
entanglement monotone 66 that quantifies the degree to which the partial
transposition of the quantum state violates the criterion of positivity 67.

Independently of the detuning between the 2LS and the laser that takes it into
the Mollow regime, there are always pairs of frequencies $\tilde{\omega}_{1}$
and $\tilde{\omega}_{2}$ for which the CSI is violated 17. However, although
the latter is an indication of entanglement, the logarithmic negativity only
becomes nonzero when the laser is taken out of resonance from the 2LS. Figure
3(a) shows the maximum $\mathcal{N}(\rho)$ that can be extracted from photons
emitted at the sidebands, depending on both the intensity of the driving
$\Omega$ and the linewidth of the detector $\Gamma$. Here, each point is
obtained for the optimum detuning $\Delta_{\mathrm{opt}}$ between the laser
and the 2LS, which we display in Fig. 3(b). For the largest part of the
figure, the detuning that optimizes the entanglement is around $\Delta\sim
2\Omega$. However, in the upper left corner of the panel we find a region for
which the optimal condition is found near resonance. However, looking at the
map of the CSI violation, shown in in Fig. 3(c), we find that such a region is
compatible with a classical state, as the ratio $R$ falls below one (cf. Sec.
II of the Supplemental Material). For visual aid, we have added the $R=1$
contour as a red dashed line in all the panels. Note that the behaviour of the
CSI, together with the artificial disruption in the logarithmic negativity, is
a consequence of the competition of the detectors for photons emitted from the
two sidebands, as quantified by the ratio $\Gamma/\omega_{+}$, which we show
in Fig. 3(d): the CSI violation starts precisely when this ratio becomes
strictly less than one, and the photons emitted from the lateral peaks become
distinguishable in frequency. Finally, in Figure 3(e) we show the emission
rate of the detectors, i.e., $I=\Gamma\langle{a^{\dagger}}a\rangle$, as a
function of the intensity of the laser and the linewidth of the detectors (the
detuning between the 2LS and the laser is taken as in panel (b)), thus making
evident the interplay between the quantity and the quality of the signal: the
highest degree of entanglement is found when the peaks are very well separated
from each other (i.e, $\Omega/\gamma_{\sigma}\gg 1$) and the linewidth of the
observer is narrow (with $\Gamma/\gamma_{\sigma}$ as small as possible), which
in turn comes with the price that such a narrow linewidth decreases the
emission rate of the source. In fact, in the configurations realized in the
bottom right corner of the panels of Fig. 3, the quantum state of the pairs of
entangled single photons is described—with a 97$\%$ fidelity—as the
superposition of the vacuum and the state
$\tilde{\rho}\equiv\left|\Phi^{-}\right\rangle\left\langle\Phi^{-}\right|$.
The latter is given by the Bell state

$\left|\Phi^{-}\right\rangle=\frac{1}{\sqrt{2}}\left(\left|0,0\right\rangle-\left|1,1\right\rangle\right)\,,$
(1)

with a purity of $91.6\%$, and its contribution to the full quantum state
ranges from 0 to 0.6$\%$. Such a small contribution indicates that, although
the photons are maximally entangled, one needs to wait for them to be emitted,
as discussed above. Furthermore, we find that entanglement is spoiled as the
linewidth of the detectors becomes large as compared to the emission lines of
the triplet. This is because wide detectors effectively erase the spectral
information of the photons, and the emission from the sidebands becomes
indistinguishable. However, one can overcome such an issue by driving the 2LS
deeper into the Mollow regime, i.e., increasing the intensity of the laser,
and thus taking the satellite peaks further away from each other. Notably, the
pair of parameters $(\Omega\,,\Delta)$ that optimize the logarithmic
negativity between the photons from the sidebands of the triplet do not
optimize the violation of the CSI nor maximize their second-order correlation
function (cf. Section II of the Supplemental Material).

Entangling polaritons
A direct application of the results presented in the previous section is the
excitation of one of the most ubiquitous systems in condensed matter physics;
namely a pair of coupled bosonic fields. While the latter can represent a
large variety of quantum systems, in the following we will associate them to
the so-called exciton-polaritons (henceforth, simply _polaritons_); which are
pseudo-particles arising from the strong coupling between a photon and an
exciton, either within a semiconductor microcavity 57 or on an organic sample
68, 69. The dynamics of polaritons driven by resonance fluorescence is given
by the master equation outlined in the Methods section, and we assume that the
polaritons are in the strong coupling regime, in which the energy states
become dressed and the light emitted is observed at the frequencies of the
lower- and upper-polariton branches (cf. the full derivation in Sec. III of
the Supplemental Material).

Figure 4: (Color online). Maximally entangled polaritons. (a) Concurrence
between the lower- and the upper-polariton branches, reaching up to 91%. (b)
Quantum state tomography of the quantum state of polaritons. (c) Same as panel
(b) after a post-selection process removing the vacuum. The remaining quantum
state has a 0.997 fidelity with a superposition of vacuum and the Bell state
$\left|\\!\left|\Psi^{-}\right\rangle\\!\right\rangle$. The figure has been
done using the parameters that optimize the entanglement of the photons
emitted by the source, namely $\Omega/\gamma_{\sigma}=4.9$,
$\Delta/\gamma_{\sigma}=8.92$ and letting the photon and the exciton be in
resonance with the higher- and lower-energy sideband of the triplet,
respectively. For panels (b) and (c) we have also used the parameters that
yield the maximum concurrence between the polaritons, i.e.,
$\Gamma_{a}/\gamma_{\sigma}=10$, $g/\gamma_{\sigma}\approx 300$ and
$\Gamma_{b}\ll\Gamma_{a}$.

Figure 4(a) shows the concurrence between the polariton branches after the
vacuum contribution has been removed through a post-selection process. We find
that for entanglement to be observed it is necessary that: i) the photonic
decay rate has to be at least an order of magnitude larger than the rate at
which the source is emitting light (we find that, i.e.,
$\Gamma_{a}/\gamma_{\sigma}\approx 30$ is typically enough) and ii) the
polariton light-matter coupling should be of the order of magnitude of the
decay of our source. The first condition guarantees that there are at most two
polaritons in the system at any given time, whereas the latter prevents the
excitations to be localized in a single polariton branch, thus preserving
entanglement. The polariton quantum state that yields the maximum concurrence
is shown in Fig. 4(b) before and in Fig. 4(c) after a post-selection process
removing the vacuum contribution. Here, we find that the quantum state of the
polaritons has a 0.997 fidelity with a superposition between vacuum and the
Bell state
$\left|\\!\left|\Psi^{-}\right\rangle\\!\right\rangle=\frac{1}{\sqrt{2}}\left(\left|\\!\left|0,1\right\rangle\\!\right\rangle-\left|\\!\left|1,0\right\rangle\\!\right\rangle\right)$,
where we have used the notation
$\left|\\!\left|m,n\right\rangle\\!\right\rangle$ to label a Fock state of
polaritons with $m$ and $n$ particles in the lower and upper branch,
respectively. These remarkable results show the power of our source of
entangled photons, while providing further evidence supporting the observation
that polaritons sustain entanglement 56, which makes them an attractive
platform to perform, e.g., quantum communication and cryptography. These
applications, however, lay outside the scope of this Article and will be left
for future research.

## Discussion and conclusions

We have described a novel source of entangled photon pairs, based on the
Mollow triplet regime of resonance fluorescence. In fact, we showed that when
the laser driving the two-level system (2LS) becomes detuned from the natural
frequency of the 2LS, the lateral peaks of the triplet become the dominant
feature of the emission spectrum. Thus, using the theory of frequency resolved
correlations 13 we showed that, when one focuses on the photons emitted from
such lateral peaks, one finds that their emission is heralded. Considering the
dressed-atom picture, we showed that these photons, emitted at frequencies
$\omega_{+}$ and $\omega_{-}$ (cf. the scheme in Fig. 1), are the result of
two atomic transitions that change the quantum state of the 2LS, and therefore
cannot take place twice in a row. Instead, they alternate in such a way that
the emission of a photon $\omega_{+}$ favours the emission of a photon
$\omega_{-}$ (note that the opposite order can be achieved if the laser is
blue- instead of red-shifted with respect to the 2LS). We have used a quantum
Monte Carlo experiment to show that, taking the laser out of resonance from
the 2LS, not only leads to the heralding behaviour, but also suppresses
emissions consisting of more than one photons; namely, our system operates as
a source of heralded single photons.

Analysing the quantum correlations between the single photons emitted from the
lateral peaks of the detuned Mollow triplet, we find that they are emitted in
a superposition of vacuum and the Bell state $\left|\Phi^{-}\right\rangle$.
Furthermore, we have obtained the volume in the parameter space for which the
photons violate the Cauchy-Schwarz inequality. Thus, we observe that for the
photons to display quantum energy-time entanglement, the linewidth of the
detector—or, more generally, the optical target of the emission—has to be
narrower than the separation between the lateral peaks of the triplet.
Otherwise, the degree of indistinguishability of the photons emitted from the
two peaks decreases, which in turn leads to the loss of entanglement.
Furthermore, we showed the relation between the emission rate of our source
and the degree of entanglement, as measured through the logarithmic
negativity. Thus, our manuscript can be used as a road map for the
experimental implementation of a source of entangled photons based on
resonance fluorescence.

Lastly, to demonstrate the power of our source, we used it to excite an
ubiquitous quantum system, commonly recurring in condensed-matter physics,
namely a pair of coupled harmonic oscillators embedded in a dissipative
environment. Thus, we use the example of exciton-polaritons (although our
results are also applicable to systems composed or containing phonons,
plasmons, bosonic nanoparticles, and photonics in general) and showed that our
source is capable to inject entangled particles into the polariton system, and
that, in turn, they are able to maintain such quantum correlations in spite of
the decoherence introduced by spontaneous decay. We note that our model for
polaritons has already been used to explain experimental results 56, and while
further elements of decoherence may be relevant in the analysis of the
entanglement in polaritons, their consideration would take us far away the
scope of the present paper and therefore are left for future references.

## References

  * Kimble and Mandel 1976 H. J. Kimble and L. Mandel, Theory of Resonance Fluorescence, Phys. Rev. A 13, 2123 (1976).
  * Allen and Eberly 1987 L. Allen and J. H. Eberly,  _Optical Resonance and Two-Level Atoms_ (Dover, 1987).
  * Kimble  _et al._ 1977 H. J. Kimble, M. Dagenais, and L. Mandel, Photon Antibunching in Resonance Fluorescence, Phys. Rev. Lett. 39, 691 (1977).
  * Chang  _et al._ 2007 D. E. Chang, A. S. Sørensen, E. A. Demler, and M. Lukin, A Single-Photon Transistor Using Nanoscale Optical Transistor, Nat. Phys. 3, 76 (2007).
  * Senellart  _et al._ 2017 P. Senellart, G. Solomon, and A. White, High-Performance Semiconductor Quantum-Dot Single-Photon Sources, Nat. Nanotech. 12, 1026 (2017).
  * Sinha  _et al._ 2019 U. Sinha, S. N. Sahoo, A. Singh, K. Joarder, R. Chatterjee, and S. Chakraborti, Single-Photon Sources, Optics & Photonics News 30, 32 (2019).
  * Fattal  _et al._ 2004 D. Fattal, E. Diamanti, K. Inoue, and Y. Yamamoto, Quantum Teleportation with a Quantum Dot Single Photon Source, Phys. Rev. Lett. 92, 037904 (2004).
  * Beveratos  _et al._ 2002 A. Beveratos, R. Brouri, T. Gacoin, A. Villing, J.-P. Poizat, and P. Grangier, Single Photon Quantum Cryptography, Phys. Rev. Lett. 89, 187901 (2002).
  * Dousse  _et al._ 2010 A. Dousse, J. Suffczyński, A. Beveratos, O. Krebs, A. Lemaître, I. Sagnes, J. Bloch, P. Voisin, and P. Senellart, Ultrabright Source of Entangled Photon Pairs, Nature 466, 217 (2010).
  * Gazzano  _et al._ 2012 O. Gazzano, M. de Vasconcellos, C. Arnold, A. Nowak, E. Galopin, I. Sagnes, L. Lanco, A. Lemaitre, and P. Senellart, Bright Solid-State Sources of Indistinguishable Single Photons, Nat. Comm. 4, 1425 (2012).
  * Loredo  _et al._ 2016 J. C. Loredo, N. A. Zakaria, N. Somaschi, C. Anton, L. D. Santis, V. Giesz, T. Grange, M. A. Broome, O. Gazzano, G. Coppola, I. Sagnes, A. Lemaitre, A. Auffeves, P. Senellart, M. P. Almeida, and A. G. White, Scalable Performance in Solid-State Single-Photon Sources, Optica 3, 433 (2016).
  * Somaschi  _et al._ 2016 N. Somaschi, V. Giesz, L. D. Santis, J. C. Loredo, M. P. Almeida, G. Hornecker, S. L. Portalupi, T. Grange, C. Antón, J. Demory, C. Gómez, I. Sagnes, N. D. Lanzillotti-Kimura, A. Lemaître, A. Auffeves, A. G. White, L. Lanco, and P. Senellart, Near-Optimal Single-Photon Sources in the Solid State, Nat. Photon. 10, 340 (2016).
  * del Valle  _et al._ 2012 E. del Valle, A. González-Tudela, F. P. Laussy, C. Tejedor, and M. J. Hartmann, Theory of Frequency-Filtered and Time-Resolved N-Photon Correlations, Phys. Rev. Lett. 109, 183601 (2012).
  * González-Tudela  _et al._ 2013 A. González-Tudela, F. P. Laussy, C. Tejedor, M. J. Hartmann, and E. del Valle, Two-Photon Spectra of Quantum Emitters, New J. Phys. 15, 033036 (2013).
  * Peiris  _et al._ 2015 M. Peiris, B. Petrak, K. Konthasinghe, Y. Yu, Z. C. Niu, and A. Muller, Two-Color Photon Correlations of the Light Scattered by a Quantum Dot, Phys. Rev. B 91, 195125 (2015).
  * López Carreño  _et al._ 2017 J. C. López Carreño, E. del Valle, and F. P. Laussy, Photon Correlations from the Mollow Triplet, Laser Photon. Rev. 11, 1700090 (2017).
  * Sánchez Muñoz  _et al._ 2014a C. Sánchez Muñoz, E. del Valle, C. Tejedor, and F. P. Laussy, Violation of Classical Inequalities by Photon Frequency Filtering, Phys. Rev. A 90, 052111 (2014a).
  * Díaz-Camacho  _et al._ 2021 G. Díaz-Camacho, E. Zubizarreta Casalengua, J. C. López Carreño, S. Khalid, C. Tejedor, E. del Valle, and F. P. Laussy, Multiphoton Emission 10.48550/ARXIV.2109.12049 (2021).
  * Zubizarreta Casalengua  _et al._ 2022 E. Z. Zubizarreta Casalengua, E. del Valle, and F. P. Laussy, Two-photon emission in detuned resonance fluorescence 10.48550/ARXIV.2210.03733 (2022).
  * Masters  _et al._ 2022 L. Masters, X. Hu, M. Cordier, G. Maron, L. Pache, A. Rauschenbeutel, M. Schemmer, and J. Volz, Will a single two-level atom simultaneously scatter two photons? 10.48550/ARXIV.2209.02547 (2022).
  * Mollow 1969 B. R. Mollow, Power Spectrum of Light Scattered by Two-Level Systems, Phys. Rev. 188, 1969 (1969).
  * Sánchez Muñoz  _et al._ 2014b C. Sánchez Muñoz, E. del Valle, A. González-Tudela, K. Müller, S. Lichtmannecker, M. Kaniber, C. Tejedor, J. J. Finley, and F. P. Laussy, Emitters of N-photon bundles, Nature Photon 8, 550 (2014b).
  * López Carreño  _et al._ 2015 J. C. López Carreño, C. Sánchez Muñoz, D. Sanvitto, E. del Valle, and F. P. Laussy, Exciting Polaritons with Quantum Light, Phys. Rev. Lett. 115, 196402 (2015).
  * López Carreño and Laussy 2016 J. C. López Carreño and F. P. Laussy, Excitation with Quantum Light. I. Exciting a Harmonic Oscillator, Phys. Rev. A 94, 063825 (2016).
  * López Carreño  _et al._ 2016 J. C. López Carreño, C. Sánchez Muñoz, E. del Valle, and F. P. Laussy, Excitation with Quantum Light. II. Exciting a Two-Level System, Phys. Rev. A 94, 063826 (2016).
  * Grünwald and Vogel 2010 P. Grünwald and W. Vogel, Entanglement in Atomic Resonance Fluorescence, Phys. Rev. Lett. 104, 233602 (2010).
  * del Valle 2013 E. del Valle, Distilling One, Two and Entangled Pairs of Photons from a Quantum Dot with Cavity QED Effects and Spectral Filtering, New J. Phys. 15, 025019 (2013).
  * Akopian  _et al._ 2006 N. Akopian, N. H. Lindner, E. Poem, Y. Berlatzky, J. Avron, D. Gershoni, B. D. Gerardot, and P. M. Petroff, Entangled Photon Pairs from Semiconductor Quantum Dots, Phys. Rev. Lett. 96, 130501 (2006).
  * Stevenson  _et al._ 2006 R. M. Stevenson, R. J. Young, P. Atkinson, K. Cooper, D. A. Ritchie, and A. J. Shields, A Semiconductor Source of Triggered Entangled Photon Pairs, Nature 439, 179 (2006).
  * Schumacher  _et al._ 2012 S. Schumacher, J. Förstner, A. Zrenner, M. Florian, C. Gies, P. Gartner, and F. Jahnke, Cavity-Assisted Emission of Polarization-Entangled Photons from Biexcitons in Quantum Dots with Fine-Structure Splitting, Opt. Express 20, 5335 (2012).
  * Johne  _et al._ 2008 R. Johne, N. A. Gippius, G. Pavlovic, D. D. Solnyshkov, I. A. Shelykh, and G. Malpuech, Entangled Photon Pairs Produced by a Quantum Dot Strongly Coupled to a Microcavity, Phys. Rev. Lett. 100, 240404 (2008).
  * Pathak and Hughes 2009 P. K. Pathak and S. Hughes, Generation of Entangled Photon Pairs from a Single Quantum Dot Embedded in a Planar Photonic-Crystal Cavity, Phys. Rev. B 79, 205416 (2009).
  * Konthasinghe  _et al._ 2012 K. Konthasinghe, J. Walker, M. Peiris, C. K. Shih, Y. Yu, M. F. Li, J. F. He, L. J. Wang, H. Q. Ni, Z. C. Niu, and A. Muller, Coherent versus Incoherent Light Scattering from a Quantum Dot, Phys. Rev. B 10.1103/PhysRevB.85.235315 (2012).
  * Ulhaq  _et al._ 2012 A. Ulhaq, S. Weiler, S. M. Ulrich, R. Roßbach, M. Jetter, and P. Michler, Cascaded Single-Photon Emission from the Mollow Triplet Sidebands of a Quantum Dot, Nat. Photon. 6, 238 (2012).
  * Muller  _et al._ 2007 A. Muller, E. B. Flagg, P. Bianucci, X. Y. Wang, D. G. Deppe, W. Ma, J. Zhang, G. J. Salamo, M. Xiao, and C. K. Shih, Resonance Fluorescence from a Coherently Driven Semiconductor Quantum Dot in a Cavity, Phys. Rev. Lett. 99, 187402 (2007).
  * Vamivakas  _et al._ 2009 A. N. Vamivakas, Y. Zhao, C.-Y. Lu, and M. Atatüre, Spin-resolved quantum-dot resonance fluorescence, Nat. Phys. 5, 198 (2009).
  * Flagg  _et al._ 2009 E. B. Flagg, A. Muller, J. W. Robertson, S. Founta, D. G. Deppe, M. Xiao, W. Ma, G. J. Salamo, and C. K. Shih, Resonantly Driven Coherent Oscillations in a Solid-State Quantum Emitter, Nat. Phys. 5, 203 (2009).
  * Ates  _et al._ 2009 S. Ates, S. M. Ulrich, S. Reitzenstein, A. Löffler, A. Forchel, and P. Michler, Post-Selected Indistinguishable Photons from the Resonance Fluorescence of a Single Quantum Dot in a Microcavity, Phys. Rev. Lett. 103, 167402 (2009).
  * Wrigge  _et al._ 2008 G. Wrigge, I. Gerhardt, J. Hwang, G. Zumofen, and V. Sandoghdar, Efficient Coupling of Photons to a Single Molecule and the Observation of Its Resonance Fluorescence, Nat. Phys. 4, 60 (2008).
  * Ortiz-Gutiérrez  _et al._ 2019 L. Ortiz-Gutiérrez, R. C. Teixeira, A. Eloy, D. Ferreira da Silva, R. Kaiser, R. Bachelard, and M. Fouché, Mollow triplet in cold atoms, New J. Phys. 21, 093019 (2019).
  * Ng  _et al._ 2022 B. L. Ng, C. H. Chow, and C. Kurtsiefer, Observation of the Mollow triplet from an optically confined single atom, Phys. Rev. A 106, 063719 (2022).
  * Cui  _et al._ 2021 C. Cui, L. Zhang, and L. Fan, Photonic analog of Mollow triplet with on-chip photon-pair generation in dressed modes, Opt. Lett. 46, 4753 (2021).
  * Astafiev  _et al._ 2010 O. Astafiev, A. M. Zagoskin, A. A. A. Jr, Y. A. Pashkin, T. Yamamoto, K. Inomata, Y. Nakamura, and J. S. Tsai, Resonance Fluorescence of a Single Artificial Atom, Science 327, 840 (2010).
  * van Loo  _et al._ 2013 A. F. van Loo, A. Fedorov, K. Lalumière, B. C. Sanders, A. Blais, and A. Wallraff, Photon-Mediated Interactions Between Distant Artificial Atoms, Science 342, 1494 (2013).
  * Toyli  _et al._ 2016 D. M. Toyli, A. W. Eddins, S. Boutin, S. Puri, D. Hover, V. Bolkhovsky, W. D. Oliver, A. Blais, and I. Siddiqi, Resonance Fluorescence from an Artificial Atom in Squeezed Vacuum, Phys. Rev. X 6, 031004 (2016).
  * Bouchiat  _et al._ 1998 V. Bouchiat, D. Vion, P. Joyez, D. Esteve, and M. H. Devoret, Quantum Coherence with a Single Copper Pair, Phys. Scr. 1998, 165 (1998).
  * Nakamura  _et al._ 1999 Y. Nakamura, Y. A. Pashkin, and J. S. Tsai, Coherent Control of Macroscopic Quantum States in a Single-Cooper-Pair Box, Nature 398, 786 (1999).
  * Mooji  _et al._ 1999 J. E. Mooji, T. P. Orlando, L. Levitov, L. Tian, H. van der Wal, and S. Lloyd, Josephson Persistent-Current Qubit, Science 285, 1036 (1999).
  * Vion  _et al._ 2002 D. Vion, A. Aassime, A. Cottet, P. Joyez, H. Pothier, C. Urbina, D. Esteve, and M. H. Devoret, Manipulating the Quantum State of an Electrical Circuit, Science 296, 886 (2002).
  * Martinis  _et al._ 2002 J. M. Martinis, S. Nam, J. Aumentado, and C. Urbina, Rabi Oscillations in a Large Josephson-Juction Qubit, Phys. Rev. Lett. 89, 117901 (2002).
  * Koch  _et al._ 2007 J. Koch, T. M. Yu, J. M. Gambetta, A. A. Houck, D. I. Schuster, J. Majer, A. Blais, M. H. Devoret, S. M. Girvin, and R. J. Schoelkopf, Charge-insensitive qubit design derived from the Cooper pair box, Phys. Rev. A 76, 042319 (2007).
  * Schreier  _et al._ 2008 J. A. Schreier, A. A. Houck, J. Koch, D. I. Schuster, B. R. Johnson, J. M. Chow, J. M. Gambetta, J. Majer, L. Frunzio, M. H. Devoret, S. M. Girvin, and R. J. Schoelkopf, Suppressing charge noise decoherence in superconducting charge qubits, Phys. Rev. B 77, 180502 (2008).
  * Bányai and Koch 1993 L. Bányai and S. W. Koch,  _Semiconductor Quantum Dots_ , World Scientific Series on Atomic, Molecular, and Optical Physics No. vol. 2 (World Scientific, Singapore ; River Edge, NJ, 1993).
  * O’Brien  _et al._ 2009 J. L. O’Brien, A. Furusawa, and J. Vučković, Photonic Quantum Technologies, Nat. Phys. 3, 687 (2009).
  * Lodahl  _et al._ 2015 P. Lodahl, S. Mahmoodian, and S. Stobbe, Interfacing Single Photons and Single Quantum Dots with Photinc Nanostructures, Rev. Mod. Phys. 87, 347 (2015).
  * Cuevas  _et al._ 2018 Á. Cuevas, J. C. López Carreño, B. Silva, M. D. Giorgi, D. G. Suárez-Forero, C. Sánchez Muñoz, A. Fieramosca, F. Cardano, L. Marrucci, V. Tasco, G. Biasiol, E. del Valle, L. Dominici, D. Ballarini, G. Gigli, P. Mataloni, F. P. Laussy, F. Sciarrino, and D. Sanvitto, First Observation of the Quantized Exciton-Polariton Field and Effect of Interactions on a Single Polariton, Sci. Adv. 4, eaao6814 (2018).
  * Weisbuch  _et al._ 1992 C. Weisbuch, M. Nishioka, A. Ishikawa, and Y. Arakawa, Observation of the Coupled Exciton-Photon Mode Splitting in a Semiconductor Quantum Microcavity, Phys. Rev. Lett. 69, 3314 (1992).
  * Gardiner 1993 C. W. Gardiner, Driving a Quantum System with the Output Field from Another Driven Quantum System, Phys. Rev. Lett. 70, 2269 (1993).
  * Carmichael 1993 H. J. Carmichael, Quantum Trajectory Theory for Cascaded Open Systems, Phys. Rev. Lett. 70, 2273 (1993).
  * López Carreño  _et al._ 2018 J. C. López Carreño, E. del Valle, and F. P. Laussy, Frequency-Resolved Monte Carlo, Sci. Rep. 8, 6975 (2018).
  * López Carreño  _et al._ 2022 J. C. López Carreño, E. Zubizarreta Casalengua, B. Silva, E. del Valle, and F. P. Laussy, Loss of antibunching, Phys. Rev. A 105, 023724 (2022).
  * Gühne and Tóth 2009 O. Gühne and G. Tóth, Entanglement detection, Phys. Rep. 474, 1 (2009).
  * Duan  _et al._ 2000 L.-M. Duan, G. Giedke, J. I. Cirac, and P. Zoller, Inseparability Criterion for Continuous Variable Systems, Phys. Rev. Lett. 84, 2722 (2000).
  * Simon 2000 R. Simon, Peres-Horodecki Separability Criterion for Continuous Variable Systems, Phys. Rev. Lett. 84, 2726 (2000).
  * Vidal and Werner 2002 G. Vidal and R. F. Werner, Computable Measure of Entanglement, Phys. Rev. A 65, 032314 (2002).
  * Plenio 2005 M. B. Plenio, Logarithmic Negativity: A Full Entanglement Monotone That is not Convex, Phys. Rev. Lett. 95, 090503 (2005).
  * Horodecki  _et al._ 2009 R. Horodecki, P. Horodecki, M. Horodecki, and K. Horodecki, Quantum entanglement, Rev. Mod. Phys. 81, 865 (2009).
  * Plumhof  _et al._ 2014 J. D. Plumhof, T. Stöferle, L. Mai, U. Scherf, and R. F. Mahrt, Room-Temperature Bose–Einstein Condensation of Cavity Exciton-Polaritons in a Polymer, Nat. Mater. 13, 247 (2014).
  * Lerario  _et al._ 2017 G. Lerario, A. Fieramosca, F. Barachati, D. Ballarini, K. S. Daskalakis, L. Dominici, M. D. Giorgi, S. A. Maier, G. Gigli, S. Kená-Cohen, and D. Sanvitto, Room-Temperature Superfluidity in a Polariton Condensate, Nat. Phys. 13, 837 (2017).
  * Cohen-Tannoudji and Reynaud 1977 C. N. Cohen-Tannoudji and S. Reynaud, Dressed-Atom Description of Resonance Fluorescence and Absorption Spectra of a Multi-Level Atom in an Intense Laser Beam, J. Phys. B.: At. Mol. Phys. 10, 345 (1977).

## Methods

Dynamics of the source
We describe resonance fluorescence as a two-level system (2LS) with natural
frequency $\omega_{\sigma}$ driven by a laser with intensity $\Omega$.
Formally, this is described by the Hamiltonian (we take $\hbar=1$ along the
paper)

$H_{\sigma}=\Delta{\sigma^{\dagger}}\sigma+\Omega({\sigma^{\dagger}}+\sigma)\,,$
(2)

where $\Delta=(\omega_{\sigma}-\omega_{\mathrm{L}})$ is the detuning between
the 2LS and the laser of frequency $\omega_{\mathrm{L}}$, and
${\sigma^{\dagger}}$ ($\sigma$) is the creation (annihilation) operator of the
2LS, which follow the pseudo-spin algebra. The dissipation of the system is
taken into account through the master equation

$\partial_{t}\rho=i[\rho,H_{\sigma}]+\frac{\gamma_{\sigma}}{2}\mathcal{L}_{\sigma}(\rho)\,,$
(3)

where $H_{\sigma}$ is the Hamiltonian in Eq. (2) and
$\mathcal{L}_{\sigma}(\rho)\equiv
2\sigma\rho{\sigma^{\dagger}}-\rho{\sigma^{\dagger}}\sigma-{\sigma^{\dagger}}\sigma\rho$.
When the laser drives the 2LS with a large intensity, the system enters into
the so-called _Mollow_ regime 21, which is characterized by its emission
spectrum in the shape of a triplet, with a central line flanked by a pair of
symmetric peaks. An archetypal spectrum of a driven 2LS with $\Delta=0$ is
shown as a dashed line in Fig. 1(a). The origin of the three peaks has a
natural explanation in the context of the “dressed atom” picture 70, which
relates each of the four possible transitions between consecutive energy
manifolds [shown in the inset of panel (a)] to the emission peaks. Keeping the
intensity of the excitation constant and taking the driving laser out of
resonance from the 2LS, the triplet splits further and the satellite peaks
become the dominant feature of the spectrum as the central peak loses its
intensity. This is shown in solid lines in Fig. 1(a). The processes that yield
the emission spectrum can be obtained through the diagonalization of the
Liouvillian of the system 25. Thus, rewriting Eq. (3) as
$\partial_{t}\rho=-M\rho$, one can find the energy of the transitions as the
imaginary part of the eigenvalues of the matrix $M$. Figure 1(c) shows in
solid lines the three energy lines available for the emission of the 2LS
driven out of resonance. For comparison, we also show in dashed lines the
energies that unfold when the excitation is resonant [and which give rise to
the spectrum shown in dashed lines in panel (a)]. The main distinction between
these two cases is that in the detuned case the lines are always splitted,
even in the limit when $\Omega/\gamma_{\sigma}\rightarrow 0$. In the opposite
regime, when the intensity of the excitation dominates over the detuning,
i.e., when $\Omega\gg\Delta$, the splitted lines coincide again. For
intermediate intensities, the energy lines are approximately given by
$\omega_{\pm}=\omega_{\mathrm{L}}\pm\sqrt{4\Omega^{2}+\Delta^{2}}$ (the exact
expression that takes into account the dissipation is given in Ref. 16). The
energies $\omega_{\pm}$ are associated to the transitions
$\left|\pm\right\rangle\rightarrow\left|\mp\right\rangle$ (shown in Fig. 1(c)
in blue and red); namely, quantum jumps that change the quantum state of the
2LS. These types of transitions dominate the dynamics of the driven 2LS when
$\Delta\gtrsim\Omega$ and, because they change the quantum state of the 2LS,
the same transition cannot take place twice in a row. Instead, they take place
one after the other and yield a scheme of photon heralding, as we showed in
the results section.

Excitation with quantum light
The theory of cascaded systems 58, 59 allows us to describe the excitation of
an optical target with the emission from our source. Thus, the master equation
describing our system is upgraded to (cf. section I of the Supplemental
Material for details of the derivation)

$\partial_{t}\rho=i[\rho,H_{\sigma}+H_{t}]+\frac{\gamma_{\sigma}}{2}\mathcal{L}_{\sigma}(\rho)+\frac{\Gamma_{1}}{2}\mathcal{L}_{a_{1}}(\rho)+\frac{\Gamma_{2}}{2}\mathcal{L}_{a_{2}}(\rho)-{}\\\
{}-\sqrt{\gamma_{\sigma}\Gamma_{1}/2}\left\\{[{a_{1}^{\dagger}},\sigma\rho]+[\rho{\sigma^{\dagger}},a_{1}]\right\\}-{}\\\
{}-\sqrt{\gamma_{\sigma}\Gamma_{2}/2}\left\\{[{a_{2}^{\dagger}},\sigma\rho]+[\rho{\sigma^{\dagger}},a_{2}]\right\\}\,,$
(4)

where $H_{\sigma}$ is the Hamiltonian in Eq. (2) and $H_{t}$ is the
Hamiltonian describing the internal degrees of freedom of the optical target.
Note that the terms in the bottom two lines of Eq. (4) are responsible for the
unidirectional coupling between the source of light and the optical targets.
For the case of the detectors used to measure the photon correlations in the
first part of the Article, we set $H_{t}=H_{d}$ with the latter defined as
$H_{d}=(\omega_{1}-\omega_{\mathrm{L}}){a_{1}^{\dagger}}a_{1}+(\omega_{2}-\omega_{\mathrm{L}}){a_{2}^{\dagger}}a_{2}$,
which takes into account the free energy of the detectors, and we let
$\Gamma_{1}=\Gamma_{2}=\Gamma$ to be the rate at which the detectors decay.
The terms in the second line of Eq. (4) are responsible for the unidirectional
coupling between the source of light and the detectors. Conversely, for the
case of the excitation of polaritons, we use the convention $a_{1}\rightarrow
a$, $a_{2}\rightarrow b$, and we replace $H_{t}$ with the polariton
Hamiltonian
$H_{p}=(\omega_{a}-\omega_{\mathrm{L}}){a^{\dagger}}a+(\omega_{b}-\omega_{\mathrm{L}}){b^{\dagger}}b+g({a^{\dagger}}b+{b^{\dagger}}a)$,
where we take into account a photon with energy $\omega_{a}$ and an exciton
with energy $\omega_{b}$ coupled with a strength $g$. Finally, we introduced
$\Gamma_{1}=\Gamma_{a}$ and $\Gamma_{2}=\Gamma_{b}$, the decay rates of the
photon and exciton modes of the polaritons.

## Data availability

The data that support the plots within this paper and other findings of this
study are available from the corresponding Author upon reasonable request.

## Code availability

The various codes used for modelling the data are available from the
corresponding Author upon reasonable request.

## Acknowledgements

J.C.L.C. was supported by the Polish National Agency for Academic Exchange
(NAWA) under project TULIP with number PPN/ULM/2020/1/00235 and from the
Polish National Science Center (NCN) “Sonatina” project CARAMEL with number
2021/40/C/ST2/00155. M.S. was supported by the European Union’s Horizon 2020
research and innovation programme under the Marie Skłodowska-Curie project
“AppQInfo” No. 956071, the National Science Centre “Sonata Bis” project No.
2019/34/E/ST2/00273, and the QuantERA II Programme that has received funding
from the European Union’s Horizon 2020 research and innovation programme under
Grant Agreement No 101017733, project “PhoMemtor” No. 2021/03/Y/ST2/00177.

## Author contributions

J.C.L.C. proposed the idea. J.C.L.C. and S.B.F. developed the theoretical
formalism and the conceptual tools, and then implemented the theoretical
methods and analysed the data. J.C.L.C. and M.S. contributed material,
analysis tools and expertise. J.C.L.C. wrote the main paper, the Supplemental
Material, and supervised the research. All authors discussed the results and
its implications and commented on the manuscript.

## Competing interests

The authors declare no competing interests.