# The Interplay between Disorder, Local Relaxation and Collective Behaviors
for an ensemble of emitters outside vs inside cavity

Zeyu Zhou  Department of Chemistry, University of Pennsylvania, 231 South 34th
Street, Philadelphia, Pennsylvania 19104, United States  Hsing-Ta Chen
Department of Chemistry, University of Pennsylvania, 231 South 34th Street,
Philadelphia, Pennsylvania 19104, United States  Department of Chemistry and
Biochemistry, University of Notre Dame, 251 Nieuwland Science Hall, Notre
Dame, Indiana 46556, USA  Joseph E. Subotnik  Department of Chemistry,
University of Pennsylvania, 231 South 34th Street, Philadelphia, Pennsylvania
19104, United States  Abraham Nitzan  Department of Chemistry, University of
Pennsylvania, 231 South 34th Street, Philadelphia, Pennsylvania 19104, United
States  Department of Chemistry, Tel Aviv University, Tel Aviv 69978, Israel

###### Abstract

The interplay between collective optical response and molecular static and
dynamic disorder is studied using simple effective Hamiltonians for an
ensemble of two-level emitters inside and outside a single-mode cavity. We
model environmental disorder by randomly modulating the molecular transition
frequencies and the coupling between the emitters and the electromagnetic
field. We also consider effects of intermolecular interactions and
orientational disorder. We investigate how these effects lead to new features
in the steady state absorption (outside the cavity), transmission spectra
(inside the cavity) and the yield of local molecular processes such as a
unimolecular reaction. Outside the cavity, the collective behavior is
manifested in the linewidth of the steady state absorption, the emission
spectrum, and the local chemical yield. Inside the cavity, however, the
collective behavior primarily determines the Rabi splitting. Effects of
intermolecular interactions under orientational disorder are also studied. For
the most part, for all types of disorder, if we increase disorder, we find a
reduction in the collective nature of the molecular response (smaller
effective $N$), and therefore, the Rabi splitting contraction occurs with
orientational disorder. Moreover, we find that static disorder is more
destructive to collective behavior than dynamic disorder.

##  I Introduction

Understanding the interplay between collective molecular optical response and
molecular disorder has drawn a great deal of interest in recent years.
Collective optical response could be seen in the time domain, for instance,
Dicke superradiance emission[1, 2, 3], superfluorescence[4, 5],
superabsorption[6] have been discussed extensively in the literature, or in
the frequency domain, e.g. in optical cavities where collective response is
manifested in the observed Rabi splitting that characterizes the strong light-
matter coupling regime in such a system. Near metal interface and in optical
cavities, the response of an individual molecule is affected also by the
confined character of the cavity as demonstrated by the Purcell effect as well
as different realizations of surface enhanced spectroscopies. [7, 8, 9, 10,
11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26] Collective
response is characterized by a non-trivial dependence of experimental
observables on the number $N$ of involved molecules; for example, one finds a
rate linear in $N$ for superradiance emission and a Rabi splitting $\sqrt{N}$
in the optical response of molecular optical cavities. The dependence of these
behaviors on local disorder has been the subject of several recent
studies.[27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43,
44, 45, 46, 47, 48, 49, 50, 51, 52, 53]

In this paper, we discuss the effect of different types of disorder on the
collective optical response of molecular systems in and out of the cavity
using two simple models: site energy disorder and coupling disorder. We also
discuss how the interplay between collective optical response (outside vs
inside a single-mode cavity) and local molecular processes might affect the
yield of molecular photoprocesses such as photochemical reactions.

##  II Model and Calculation Method

###  II.1 A molecular ensemble in the single exciton subspace

Consider an ensemble of 2-level molecules whose spatial extent is assumed for
shorter than the wavelength of its resonant absorption, subjected to pumping
by a near resonant classical field and further exhibiting collective emission.
In addition, each molecular excited state can decay due to interaction with
its local environment or a local chemical process. The time evolution of this
system can be described by the following Hamiltonian:

$\displaystyle\hat{H}=$
$\displaystyle\hat{H}_{0}+\hat{V}_{\text{pump}}-i\Gamma_{rad}\hat{Q}-i\Gamma_{loc}\hat{C}$
(1) $\displaystyle\hat{H}_{0}=$ $\displaystyle
E_{0}\left|0\right\rangle\left\langle
0\right|+\sum_{m}E_{m}\left|m\right\rangle\left\langle m\right|$ (2)
$\displaystyle\hat{V}_{\text{pump}}=$ $\displaystyle\sum_{m}V_{0m}\cos(\omega
t)\left|0\right\rangle\left\langle m\right|+V_{m0}\cos(\omega
t)\left|m\right\rangle\left\langle 0\right|$ (3) $\displaystyle\hat{Q}=$
$\displaystyle\sum_{m,m^{\prime}}\left|m\right\rangle\left\langle
m^{\prime}\right|$ (4) $\displaystyle\hat{C}=$
$\displaystyle\sum_{m}\left|m\right\rangle\left\langle m\right|$ (5)

In eq (1)-(5), the molecular states are
$\left\\{\left|0\right\rangle,\left|m\right\rangle\right\\}$:
$\left|0\right\rangle=\prod_{j}\left|g_{j}\right\rangle$ is the total ground
state and $\left|m\right\rangle=\left|x_{m}\right\rangle\prod_{j\neq
m}\left|g_{j}\right\rangle$ is a singly-excited molecular state. Note that we
are working in a single-exciton basis, so that our model can describe only
weak excitations where many exciton correlations can be disregarded. The
radiative emission of the molecules is treated implicitly using the effective
non-Hermitian operator $-i\Gamma_{rad}\hat{Q}$ that arises when one considers
one idealized wide band radiative continuum (with density of state
$\rho_{dos}$) coupled identically to all molecular emitters
($V_{j}=V_{rad},\forall j$), leading to the purely imaginary self-energy

$\displaystyle\Gamma_{rad}=2\pi\rho_{dos}|V_{rad}|^{2}$ (6)

If we denote the wavefunction of the molecular system by $c_{m}(t)$, then
under radiative pumping of the molecules, the outgoing radiative flux at any
given time $t$ is

$\displaystyle
J_{rad}(t)=\sum_{m,m^{\prime}}\Gamma_{rad}c_{m}(t)c^{*}_{m^{\prime}}(t)$ (7)

Finally, $-i\Gamma_{loc}\hat{C}$ represents decay through a local channel
involving a single molecule process. In contrast to the radiative decay, this
channel describes a process defined by the state of each individual molecule
(and does not depend on intermolecular coherence). The outgoing local flux at
any given time $t$ is

$\displaystyle J_{loc}(t)=\sum_{m}\Gamma_{loc}|c_{m}(t)|^{2}$ (8)

Note that in the Hamiltonian (eq 1), pumping is done by a continuous wave (CW)
field at the driving frequency $\omega$. Therefore, if we plot total steady
state outgoing flux ($J_{rad}+J_{loc}$) as a function of this frequency
$\omega$, we obtain the steady state absorption spectrum.

###  II.2 A molecular ensemble in an optical cavity

To describe the same molecular ensemble inside a single-mode optical cavity,
we modify the Hamiltonian (eq (1)-(5)) in several ways.[54] First, the
external pumping does not directly couple to the emitters, but instead to the
cavity mode $\left|v\right\rangle$. Second, the molecules do not couple
directly to the far EM field, and are assumed to couple only to the cavity
mode. Hence, the term $-i\Gamma_{rad}\hat{Q}$ is replaced by the coupling term
(eq 11) below). Finally, the radiative ($\Gamma_{cav}^{(R)}$) and nonradiative
($\Gamma_{cav}^{(NR)}$) decay rates of the cavity mode to the far field and to
damping in the cavity mirrors, respectively, are taken into account. The
corresponding Hamiltonian is

$\displaystyle\hat{H}=$ $\displaystyle\hat{H}_{0}+\hat{V}_{\text{cav-
mol}}+\hat{V}_{\text{pump}}-i\Gamma_{\text{loc}}\hat{C}-i\Gamma_{\text{cav}}\left|v\right\rangle\left\langle
v\right|$ (9) $\displaystyle\hat{H}_{0}=$ $\displaystyle
E_{v}\left|v\right\rangle\left\langle
v\right|+E_{0}\left|0\right\rangle\left\langle
0\right|+\sum_{m}E_{m}\left|m\right\rangle\left\langle m\right|$ (10)
$\displaystyle\hat{V}_{\text{cav-mol}}=$
$\displaystyle\sum_{m}V_{vm}\left|v\right\rangle\left\langle
m\right|+V_{mv}\left|m\right\rangle\left\langle v\right|$ (11)
$\displaystyle\hat{V}_{\text{pump}}=$ $\displaystyle V_{0v}\cos(\omega
t)\left|0\right\rangle\left\langle v\right|+V_{v0}\cos(\omega
t)\left|v\right\rangle\left\langle 0\right|$ (12) $\displaystyle\hat{Q}=$
$\displaystyle\sum_{m,m^{\prime}}\left|m\right\rangle\left\langle
m^{\prime}\right|$ (13)

where, $v$ is the excited state of the cavity mode (and the ground state
$\left|0\right\rangle$ refers to both the cavity and the molecular systems).
Similar to many previous works, we consider only a single cavity mode and in
correspondence with the single-exciton model, we include only the ground and
first excited states of this mode.

Here, the cavity photon leaking rate
($\Gamma_{\text{cav}}=\Gamma_{\text{cav}}^{(R)}+\Gamma_{\text{cav}}^{(NR)}$)
determines the quality factor $\hbar\omega_{\text{cav}}/\Gamma_{\text{cav}}$.
In addition to the local relaxation channel (which is the same as in eq (8)),
the outgoing radiative flux through the cavity is

$\displaystyle J_{cav}(t)=\Gamma_{cav}^{(R)}|c_{v}(t)|^{2}$ (14)

This flux represents the combined transmission spectrum.

###  II.3 Disorder and intermolecular couplings

Disorder can be implemented in the Hamiltonian (1) and (9) by assuming that
either the molecular transition energies $E_{m}$ or the molecular couplings
with its radiative environment ($V_{m0}$ in eq 3 or $V_{mv}$ in eq 11) have a
random component. This random component can be constant in time (static
disorder) or time-dependent (dynamic disorder)

####  II.3.1 Energy disorder

If we focus first on molecular transition energy disorder, the static disorder
limit is described by including a random component in the individual molecular
frequencies ($\omega_{m}=E_{m}/\hbar$)

$\omega_{m}=\omega_{0}+\Omega_{m}.$ (15)

In the calculations reported below, the random component $\Omega_{m}$ is
sampled from a Gaussian distribution

$P[\Omega_{m}]=\frac{1}{\sqrt{2\pi\sigma^{2}}}\exp{\left(-\Omega_{m}^{2}/2\sigma^{2}\right)}$
(16)

where different molecules are assumed to be uncorrelated,
$\langle\Omega_{m}\Omega_{m^{\prime}}\rangle=\delta_{mm^{\prime}}\sigma^{2}$
and here the variance $\sigma^{2}$ indicates the disorder strength. This
static disorder model may be viewed as the static limit of Kubo’s stochastic
modulation model where

$\omega_{m}=\omega_{0}+\Omega_{m}(t)$ (17)

and $\Omega_{m}(t)$ are stochastic random variables that change over time. We
choose $\Omega_{m}(t)$ to be a Gaussian stochastic variable that satisfies
$\langle\Omega_{m}(t)\rangle=0$ and

$\langle\Omega_{m}(t_{1})\Omega_{m^{\prime}}(t_{2})\rangle=\delta_{mm^{\prime}}\sigma^{2}e^{-|t_{1}-t_{2}|/\tau_{c}}$
(18)

where $\delta_{mm^{\prime}}$ is a Kronecker delta function. The Gaussian
stochastic variables are characterized by the amplitude
$\sigma=\langle\delta\Omega^{2}\rangle^{1/2}$ and correlation lifetime
$\tau_{c}$ of local energy fluctuations. For a detailed algorithm to generate
such a stochastic process, see Appendix A.

####  II.3.2 Orientational Disorder

The most obvious source of coupling disorder is the inherent orientational
disorder that characterizes most molecular liquids, and here we focus on this
aspect of our systems. The simplest way to model such disorder is to assign an
angle $\theta_{m}$ for each emitter $m$, where $\theta_{m}$ is oriented with
respect to the cavity mode (see Fig. 1). We assume that the field of the
cavity mode that may be excited by the external pumping is polarized in the
$z$ direction and denote $\theta_{m}$ as the angle between the $z$ axis and
the molecular transition dipole direction.

Figure 1: Schematic diagram of the molecular emitters (blue arrows) confined
in an optical cavity (grey walls) and pumped by the incoming CW field (black
line). The field of the cavity mode is assumed to be in the $z$ direction. The
transition dipole moments of the emitters have different angles $\theta_{m}$
with respect to the $z$ axis.

This disorder leads to three modifications of the model Hamiltonian: First,
the cavity-molecule coupling becomes

$\hat{V}_{\text{cav-
mol}}=\sum_{m}V_{vm}\cos(\theta_{m})\left|v\right\rangle\left\langle
m\right|+V_{mv}\cos(\theta_{m})\left|m\right\rangle\left\langle v\right|$ (19)

Similarly, for the case outside the cavity, the couplings to the external
driving EM field also depend on the emitters orientation

$\hat{V}_{\text{pump}}=\sum_{m}V_{0m}\cos(\theta_{m})\cos(\omega
t)\left|0\right\rangle\left\langle m\right|+V_{m0}\cos(\theta_{m})\cos(\omega
t)\left|m\right\rangle\left\langle 0\right|$ (20)

Also, outside the cavity, the elements of the non-Hermitian $\hat{Q}$ matrix
depend on the relative angle $\cos(\theta_{m}-\theta_{m^{\prime}})$[55]

$\hat{Q}=\sum_{m,m^{\prime}}\cos(\theta_{m}-\theta_{m^{\prime}})\left|m\right\rangle\left\langle
m^{\prime}\right|$ (21)

The disorder of $\\{\theta_{m}\\}$ can be either static or dynamic as
described in section II.3. At room temperature, the characteristic timescale
of the rotational motion in solution is of order $10ps$, which implies that in
most situations including light-matter strong coupling regime in molecular
system, orientational disorder may be assumed static. Nevertheless, for
completeness, we will consider below a full range of dynamic disorder. In the
static case, we choose a time-independent $\\{\theta_{m}\\}$ according to
random orientation angles sampled from a Gaussian distribution

$P[\theta_{m}]=\frac{1}{\sqrt{2\pi\sigma^{2}}}\exp{\left(-\theta_{m}^{2}/2\sigma^{2}\right)}$
(22)

This again is the limit $\tau_{c}\rightarrow\infty$ of a process where
orientations are dynamically modulated according to
$\langle\theta_{m}(t)\rangle=0$ and

$\langle\theta_{m}(t_{1})\theta_{m^{\prime}}(t_{2})\rangle=\delta_{mm^{\prime}}\sigma^{2}e^{-|t_{1}-t_{2}|/\tau_{c}}$
(23)

####  II.3.3 Orientational disorder with intermolecular coupling

To complete our investigation of linear response of molecules outside and
inside optical cavities as expressed by the absorption and transmission
signals, we will also consider models with intermolecular couplings (dipole-
dipole couplings), again focusing on static and dynamic orientational
disorder. The intermolecular couplings is modeled by a dipole-dipole coupling:

$\displaystyle\hat{V}_{\text{dip}}=$
$\displaystyle\sum_{m=k+1}^{N}\sum_{k=1}^{n}W_{k}f(\theta_{m},\theta_{m-k})\left|m\right\rangle\left\langle
m-k\right|+\sum_{m=k+1}^{N}\sum_{k=1}^{n}f(\theta_{m-k},\theta_{m})W_{k}\left|m-k\right\rangle\left\langle
m\right|$ (24)

Here, $W_{k}=W/k^{3}$ and in the calculations reported below we have used
$W=0.1$,111The choice of intermolecular coupling is $0.1$ in unit of Rabi
splitting. Rabi splitting is approximately on the order of $10meV$, and if the
distance between dipoles are $2nm$, the dipole moments are on the order of
$10Debye$. unless otherwise specified. In principle, this Hamiltonian
describes a chain of emitters and the interaction between them decays as one
over their distance cubed, $1/r^{3}$, as sketched in Fig 1.

The function $f(\theta_{m},\theta_{m-k})$ in eq 24 reflects the angular
dependence of the dipole-dipole coupling between emitters and is given by

$\displaystyle
f(\theta_{m},\theta_{m-k})=\cos(\theta_{m}-\theta_{m-k})-3\cos\theta_{m}\cos\theta_{m-k}$
(25)

Again, these angles are defined to be relative to the polarization of the
incoming driving electric field.

###  II.4 The Steady State flux

In order to simulate cavity-molecular dynamics, we propagate the Schrödinger
equation for the Hamiltonians above under steady state boundary conditions.
The wavefunction for the driven, open quantum system can be written as
$\left|\Psi\right\rangle=c_{0}\left|0\right\rangle+c_{v}\left|v\right\rangle+\sum_{m}c_{m}\left|m\right\rangle$.
This implies propagating the Schrödinger equation
$i\hbar\frac{d}{dt}\left|\Psi\right\rangle=\hat{H}\left|\Psi\right\rangle$
under the Hamiltonians of eqs 1 and 9 keeping $|c_{0}|^{2}=1$. The long-time
evolution yields the steady state forms of the coefficients $c_{m}(t)$ and
$c_{v}(t)$ that are used to evaluate the steady state fluxes $J_{rad}$ (eq 7),
$J_{loc}$ (eq 8), and $J_{cav}$ (eq 14). Again, the total steady state flux
$J_{rad}+J_{loc}$ as a function of the incoming driving frequency yields the
steady state absorption spectrum for emitters outside optical cavities; while
inside the cavity, the steady state flux $J_{cav}$ gives
reflection/transimission spectrum and $J_{loc}$ represents the local
relaxation flux.

With this single-exciton model, there are many parameters we can adjust to
account for different physical conditions. We will discuss these different
possibilities in the following sections.

##  III Results

###  III.1 Energy disorder

In Fig 2, we plot the absorption spectrum for an ensemble of emitters outside
the cavity. In the limit of zero disorder, the linewidth of the absorption
spectrum is associated with the decay rate for the radiative and local non-
radiative channels. For the static disorder results (a), as the disorder
strength $\sigma$ increases, the absorption lineshape become lower and
broader. This feature reflects inhomogeneous broadening and the static
uncertainty in the molecular excitation energy. Under dynamic disorder (b), as
the correlation time $\tau_{c}$ decreases, the absorption linewidth is
smaller, which reflects motional narrowing. We can also predict the absorption
lineshape in Fig 2 by performing a Fourier transform of $\phi(t)$222This
formula is by modifying equation (7.104) in ref. 62 by taking both decay rates
$N\Gamma_{rad}+\Gamma_{loc}$ into account according to Voigt theorem.

$\displaystyle\phi(t)=\exp[-\sigma^{2}\tau_{c}^{2}(\frac{t}{\tau_{c}}-1+e^{-\frac{t}{\tau_{c}}})-(N\Gamma_{rad}+\Gamma_{loc})t].$
(26)

For the case of static disorder ($\tau_{c}\rightarrow\infty$), we can simplify
the formula by expanding the exponential of $t/\tau_{c}$

$\displaystyle\phi(t)\rightarrow\exp(-\sigma^{2}t^{2}/2-(N\Gamma_{rad}+\Gamma_{loc})t)$
(27)

which, after a Fourier transform, yields a convolution between the gaussian
associated with static disorder and the lorentzian associated with the
radiative and non-radiative decay rates. For more details and a comparison of
analytic versus numerical results, see Appendix B.

Figure 2: The steady state absorption spectrum for a molecular ensemble
outside the cavity with (a) different static energetic disorder strengths and
(b) different correlation time $\tau_{c}$ but constant disorder strength
$\sigma=0.2$. The molecular ensemble includes $N=40$ molecules. If there is no
disorder (panel (a), $\sigma=0$, blue line), the linewidth of the peak is
$2N\Gamma_{rad}+2\Gamma_{loc}$. As the static disorder strength increases, we
obtain a broader inhomogeneous peak. As the correlation time decreases under
fixed disorder strength (panel (b), $\sigma=0.2$), the peak gets narrower,
reflecting motional narrowing. As shown in eqs 26 and 27, as well as Appendix
B, these features can be explained analytically.

Next, consider the case where the same ensemble of emitters interacts with a
single-mode cavity under steady state pumping. As discussed in section II.2,
we calculate two observables: (1) the steady state flux $J_{cav}$ through the
cavity leakage channel (which gives the combined transmission spectrum) and
the steady state flux $J_{loc}$ through the local relaxation channel of
individual emitters. The results in Fig. 3 show that static disorder leads to
an inhomogeneous broadening of the two polariton peaks that are present in
both of the relevant channels. Interestingly, Fig 4, in which we suppress the
homogeneous contributions (i.e. when non-Hermitian parts $\Gamma_{rad}$,
$\Gamma_{loc}$, $\Gamma_{cav}$ are small) to the lineshape so that the width
is dominated by the inhomogeneous contribution, shows that this inhomogeneous
broadening is much weaker inside a cavity relative to outside a cavity . This
effect has already been seen in reference 58.

Returning to Fig. 3, we note that, for increasing $\sigma$, the baseline
between the two polariton peaks (around driving frequency $\omega_{d}=0$) for
the local relaxation channel becomes larger, whereas in the transmission
spectrum, the baseline around driving frequency $\omega_{d}=0$ remains
approximately zero. This different behavior arises because the dark modes
cannot contribute to the transmission spectrum but they do contribute to the
local relaxation channel. Another observation in Fig. 3 is the increase of the
effective Rabi splitting in both figures (a) and (b), as the energetic static
disorder strength increases. This increase can be captured approximately by
the following formula

$\Omega_{R}(\sigma)/\Omega_{R}(\sigma=0)=1+2\sigma^{2}$ (28)

One can obtain such a formula by using second-order perturbation theory for
the polariton states. For a simple derivation, see Appendix C.

Lastly, for the case with dynamic disorder, as shown in Fig. 5, as one might
expect, we find motional narrowing for both channels as the correlation time
$\tau_{c}$ decreases (similar to Fig. 2). Moreover, we observe the same
behavior as in Fig. 3: in between the two polariton peaks (around driving
frequency $\omega_{d}=0$), the presence of dark modes leads to a finite steady
state flux for the local relaxation channel.

Figure 3: Steady state (a) transmission spectrum and (b) local relaxation flux
spectra, as calculated for molecular ensembles with different static disorder
strengths $\sigma$ interacting with a single-mode cavity. The molecular
ensemble and the cavity mode form two polariton peaks in both figures. As the
amplitude of static disorder increases, the split peaks become broadened;
however, the baseline in between the two peaks (around driving frequency
$\omega_{d}=0$) grows larger only for the local signal. This growth arises
because the dark modes contribute only to the local relaxation flux but not to
the transmission spectrum. Lastly, note that as the energetic static disorder
increases, the effective Rabi splitting increases. This increase can be
captured approximately by perturbation theory (eq 28). Figure 4: Steady state
(a) absorption spectra for outside the cavity, (b) transmission and (c) local
relaxation flux spectra interacting with a single-mode cavity, as calculated
for molecular ensembles with different static disorder strengths $\sigma$. The
black dotted lines in panel (a) are analytical results calculated in Appendix
B. Here, the contributions to homogeneous broadening, namely the non-Hermitian
parts of the Hamiltonian matrices are taken
$\Gamma_{rad}=0.005/N,\Gamma_{cav}=0.005,\Gamma_{loc}=0.005$, which are $10$
times smaller than in Fig 3, so (except from the case without disorder) the
linewidth outside the cavity is dominated by the inhomogeneous broadening. In
agreement with earlier observations[58], the inhomogeneous broadening observed
outside the cavity is not manifested in the linewidth of the polariton peaks
inside the cavity which are therefore much narrower as seen in panel (b) and
(c).  Figure 5: Steady state (a) transmission spectra and (b) local relaxation
for the molecular ensemble inside a single-mode cavity with dynamic disorder.
We fix the disorder strength ($\sigma=0.2$) and scan the correlation time
$\tau_{c}=100,10,1,0.1$. The molecular ensemble and the cavity mode forms two
polariton peaks in both figures. Note that the long correlation limit (blue
lines) recovers the results with static disorder in Fig. 3 (Yellow lines).
Note also that similar to Fig. 3, the dark modes contribute only to the local
relaxation channel. As we decrease $\tau_{c}$, the two peaks become narrower
in both figures, and the dark modes contribute less to the local relaxation
channel.

###  III.2 Orientational disorder

In this subsection, we present results for angular disorder as defined in
section II.3.2. Note that the way disorder is defined in this case implies the
average contribution is nonzero ($\langle\cos\theta_{m}\rangle\neq 0$), and
therefore, some of the results shown below reflect this situation. As shown in
Fig 6, for the static angular disorder case, as the disorder strength $\sigma$
increases, both the height and the linewidth of the absorption peak decrease
(see eq 22). This decrease arises because each dipole in the ensemble of
emitters couples to the driving field differently, resulting in less
collective behavior. For the dynamic angular disorder case, in agreement with
the results in section III.1, as we decrease the correlation time of the
angular disorder, we find motional narrowing.

Figure 6: Steady state absorption spectra for the ensemble of emitters outside
the cavity. (a) static angular disorder; (b) dynamic angular disorder. As
shown in the static disorder figure, when the disorder strength $\sigma$
increases, the absorption peaks decrease. More importantly, the linewidth
decreases because the emitters act less collectively. For the case of dynamic
disorder, similar to Fig. 2, we observe motional narrowing.

Next, we consider this molecular ensemble with angular disorder inside the
cavity. The results are shown in Fig. 7. Unlike the results in Fig. 3, (where
we investigated disorder in the excitation energies), here the Rabi splitting
consistently decreases as we increase the static angular disorder strength and
lose collectivity of the response. Clearly, when interpreting the spectra from
cavity, one must consider geometry and not just the energy levels of the
molecules in the cavity. In principle, in a cavity, one might expect to have
both energy and angular disorder and so, as the temperature rises, the Rabi
splitting should decrease. Unfortunately, it remains unclear which effect
dominates. Note that the increase of Rabi splitting with increasing static
energy disorder does not suggest that more emitters are behaving collectively.
As shown both in Fig. 3 and Fig. 7, the maximal radiative steady state flux
decreases as static disorder strength $\sigma$ increases. Lastly, as in Fig.
3, the dark modes contribute to the rise of baseline between the two polariton
peaks only for the local relaxation channel.

Figure 7: Steady state transmission spectrum and local relaxation flux for
different static angular disorder strength ($\sigma=0,0.2\pi,0.4\pi,0.6\pi$)
inside a cavity. As shown in the figures, when the disorder strength $\sigma$
increases, the absorption peak height decreases and the dark modes contribute
only to the local relaxation flux but not to the transmission spectrum
(similar to Fig. 3). However, unlike Fig. 3, disorder in the angular momentum
leads to a decrease in the Rabi splitting.

The results for the dynamic angular disorder case are slightly more non-
intuitive. In Fig. 8, we fix the disorder strength ($\sigma=0.4\pi$) and
calculate transmission and local relaxation fluxes for different correlation
time $\tau_{c}$. In the long correlation time limit (blue line), we recover
the static disorder result. Two features are worth noting in Fig. 8. First, as
$\tau_{c}$ becomes smaller, we find that, as expected, the polaritonic spectra
undergo motional narrowing. However, reducing the correlation time is not
equivalent to reducing disorder insofar as the fact that the Rabi splitting in
Fig. 8 decreases as $\tau_{c}\rightarrow 0$ (whereas the Rabi splitting
increases as $\sigma\rightarrow 0$ in Fig. 7). For a simple explanation of
this Rabi splitting contraction, consider the two limiting scenarios, (i) the
static limit ($\tau_{c}=\infty$) and (ii) the no correlation limit
($\tau_{c}=0$). For the static limit ($\tau_{c}=\infty$), the averaged Rabi
splitting is calculated by estimating the eigenvalues $E_{LP/UP}$.

$\displaystyle
E_{LP/UP}=\sqrt{N|V_{mv}|^{2}\langle\cos^{2}(\theta_{m})\rangle}$ (29)

Here, $\langle\rangle$ stands for ensemble average. For the no correlation
limit ($\tau_{c}=0$), a good estimate is

$\displaystyle E_{LP/UP}=\sqrt{N|V_{mv}|^{2}}\langle\cos(\theta_{m})\rangle$
(30)

However, we are unable to predict what Rabi splittings are for finite
$\tau_{c}$. For detailed derivations of these expressions and an analysis of
how they match up with numerical results, see Appendix D. Note that in both
limits, the Rabi splittings are smaller than the case with no disorder,
reflecting the fact that the average coupling in this model is smaller than
$V_{mv}$.

Another non-intuitive result in Fig. 8 pertains to the height of the peaks
themselves. Here, we see that the height behaves in a non-monotonic fashion as
a function of $\tau_{c}$; starting from static disorder and decreasing
$\tau_{c}$, we find that the peak height and integral decreases and then
increases. Qualitatively, we believe this non-monotonic behavior can be
understood by recognizing that with moderate correlation times, introducing
disorder and slowly distorting the coherences removes the possibility of $(i)$
seeing individual configurations, each absorbing for a long time or $(ii)$
seeing many configurations averaged to one static configuration. As a result,
there is less absorption at intermediate correlation times. Pursuing a more
rigorous analysis of this effect will certainly be an important research
direction in the future. Finally, just as in Fig. 7, the dark modes contribute
to the local relaxation, and as the correlation time decreases, the Rabi
splitting decreases.

Figure 8: Steady state absorption spectra for dynamic angular disorder. We fix
the angular disorder strength ($\sigma=0.4\pi$) and scan over correlation time
$\tau_{c}=100,10,1,0.1$. We also plot the result with no disorder for
reference. As shown in both figures, when the correlation time decreases, the
collective effect diminishes, leading to a smaller Rabi splitting. Similar to
Fig. 7, the dark modes contribute only to the local relaxation signals.

###  III.3 Intermolecular coupling and orientational disorder

Finally, let us address the presence of how intermolecular coupling affects
orientational disorder. Let $\left|B\right\rangle$ be the bright state of the
system outside the cavity

$\displaystyle\left|B\right\rangle$
$\displaystyle=\frac{1}{\sum_{j}\cos^{2}\theta_{j}}(\cos\theta_{1},\cos\theta_{2},...)^{T}.$
(31)

Before presenting the absorption spectra, we start by showing the density of
states of the system outside/inside a cavity, weighted by the brightness of
the states in Fig 9. For the case outside the cavity, we define the brightness
$W_{j}$ of each eigenstate $\phi_{j}$ of $\hat{H}_{sub}$ (keeping only the
real part of $m$ and $v$ elements, see eqs 36-38 in Appendix C) to be:

$\displaystyle w_{j}$ $\displaystyle=|\left\langle
B\middle|\phi_{j}\right\rangle|^{2}$ (32)

For the case inside the cavity, the bright mode must contain the cavity mode
and thus the following definition makes more sense:

$\displaystyle w_{j}^{\text{cav}}$ $\displaystyle=|\left\langle
v\middle|\phi_{j}\right\rangle|^{2}$ (33)

where $\left|v\right\rangle$ is the bare cavity state (see eqs 9-13).

Figure 9: Histogram of the brightness ($w_{j}$ in eq 32 and
$w_{j}^{\text{cav}}$ in eq 33) for the Hamiltonian eigenstates that arise with
intermolecular coupling (see eqs 24 and 25) as averaged over realizations of
$\theta_{m}$ following the distribution in eqs 22 and 23. These brightness
factors are relevant to pumping the emitters outside/inside the cavity. We
consider that the dipoles are fully aligned in the limit of zero disorder and
thus, intermolecular couplings are negative (eq 25 equals $-2$). The
eigenvalue of the bright state (symmetric superposition state) becomes
negative. Hence, for the case outside the cavity, the brightest states have
energy $E\approx-0.3\ a.u.$, which is lower than the energy of the independent
emitters (which is set to be $0$); for the case inside the cavity, the two
brightest states yield lower and upper polariton peaks at around $E=-0.6,0.3\
a.u.$ respectively. The upper polariton peak is brighter than the lower
polariton because the former has an energy eigenvalue ($0.3\ a.u.$) closer to
the energy of the bare cavity mode ($0\ a.u.$). Again, the center of the two
polariton peaks ($(E_{LP}+E_{UP})/2\approx-0.15\ a.u.$) is lower than the
energy of the independent emitters/bare cavity mode.

As shown in Fig 9, for the case outside the cavity, the brightest eigenstates
have negative energy eigenvalues (relative to the emitters energy). For the
case inside the cavity, there are two regions of bright states corresponding
to the two polariton peaks; the center of the two peaks is shifted towards the
negative direction, just as for the case outside the cavity. The upper
polariton eigenstates are brighter than the lower polariton because they are
closer to the bare cavity mode energy.

Figure 10: Steady state absorption for emitters outside cavity with
intermolecular couplings and static/dynamic orientational disorder. (a) static
orientational disorder; (b) dynamic orientational disorder. For the static
disorder case, the absorption peaks height decreases as the disorder strength
increases. As predicted in Fig. 9, intermolecular coupling reduces the energy
of the system and the peaks appear to have energy lower than that of the
individual emitters. For the dynamic disorder case, the peak gets narrower as
we decrease the correlation time. Figure 11: Steady state (a) transmission and
(b)local relaxation for an ensemble of emitters inside a single-mode cavity
with intermolecular couplings for different static orientational disorder
strength. As predicted in Fig. 9, the upper polariton contributes more
strongly than does the lower polariton peak to the transmission spectrum.
Intermolecular couplings stabilize the system and thus the center of the two
polariton peaks corresponds to an energy lower than that of an individual
emitter/bare cavity mode. As the disorder strength increases, the Rabi
splitting decreases as the collectiveness is destroyed by the disorder. The
difference in height between the upper and lower polariton peaks is small for
the local relaxation channel because the signals are not dependent on how
strong is the cavity mode contribution to the polariton (but rather are
dictated by the dynamics of the individual emitters). Figure 12: Steady state
(a) transmission and (b) local relaxation flux for an ensemble of emitters
inside the cavity with intermolecular couplings and dynamic disorder. Both
figures are for the indicated correlation times but a fixed disorder strength
($\Delta=0.4\pi$). Again, the long correlation time limit (blue lines) is
equivalent to static disorder results (yellow lines in Fig. 11). Both the Rabi
splitting and the linewidth decreases as the correlation time decreases
similar to Fig. 8. The center of the two polariton peaks shift to lower energy
by the intermolecular coupling. Same as in Fig. 11, the difference in height
between the upper and lower polariton peaks is more significant in the
transmission spectrum than in the local relaxation signal.

With these structures in mind, in Fig 10, we plot the absorption spectra for
different orientational disorder strengths. When there is no orientational
disorder, because there is finite negative intermolecular coupling, the peaks
are red shifted, as the brightest eigenstate corresponds to the smallest
eigenvalue. As we increase the disorder strength, the total absorption
decreases. For the dynamic orientational disorder case, similar to Fig. 6, as
the correlation time decreases, the absorption linewidth decreases
(corresponding to motional narrowing). The scale of the red shift is also
reduced as we increase modulation speed.

In Fig. 11, we plot the transmission and local relaxation flux for an ensemble
of emitters placed inside the cavity. First, one might guess from what was
shown in Fig. 9, the upper polariton peak height is greater than that of the
lower polariton; the center of the two polariton peaks has a lower energy than
that of the bare cavity mode. Second, as the static disorder strength
increases, the effective Rabi splitting decreases because the disorder reduces
the collectiveness of the emitters. Third, as we compare the two channels,
because the transmission spectra relies on the cavity mode and the local
relaxation relies on the emitter states, the difference in intensity between
the upper and lower polaritons for transmission spectrum is greater than that
for the local relaxation signals. Lastly, the dark modes contribute only to
the local relaxation signals, as was also seen in Fig. 7.

Finally, we show the steady state transmission and local relaxation flux
inside a single-mode cavity with intermolecular coupling and dynamic
orientational disorder in Fig. 12. As the correlation time decreases, the
effective Rabi splitting and the linewidth of each polariton peak decreases,
just as in Fig. 8. The difference between the upper and lower polariton peak
height is less significant in local relaxation flux than the transmission
spectrum, similar to what was found in Fig. 11. All of the observations are
consistent with Section III.2.

##  IV cavity effect on local relaxation yield

The models defined by eqs 1-13 are minimal models that allow considerations of
three relaxation channels for the radiation absorbed from the driving field: a
radiative channel that combines transmission and reflection, nonradiative
damping of the cavity mode (heat production in the cavity walls) and reactive
relaxation. The latter represents a reaction that occurs following the
molecular excitation. It is of interest to examine the effect of the cavity
environment on the yield of the latter channel.

Outside the cavity, the system undergoes collective superradiance at a rate
$N\Gamma_{rad}$ and local relaxation at a rate $\Gamma_{loc}$. The yield of
the individual molecular relaxation (that we assume to be a reactive process)
is

$\displaystyle\text{Yield}=\frac{\Gamma_{loc}}{N\Gamma_{rad}+\Gamma_{loc}}$
(34)

By comparison, inside the cavity, the energy will leak through the cavity mode
at rate $\Gamma_{cav}$ with the same local relaxation rate as above,
$\Gamma_{loc}$. It is clear that the local yield is then

$\displaystyle\text{Yield}=\frac{\Gamma_{loc}}{\Gamma_{cav}+\Gamma_{loc}}$
(35)

Hence, without disorder, we can expect that the reactive relaxation yields
will be different outside and inside the cavity, both in terms of absolute
value and $N$ dependence. This difference is shown in Fig. 13, Note, for
simplicity, we assume $\Gamma_{cav}=10\Gamma_{rad}$. In agreement with eqs 34
and 35, the local relaxation yield depends on $N$ only when the emitters are
outside the cavity.

Consider now the effect of disorder. As discussed above and shown in Figs. 2
and 5, at reasonably high temperatures, nuclear motion and dephasing processes
lead to static and dynamic disorder, which results in a reduction of steady
state radiative fluxes. As shown in Fig. 13, we can clearly see that moderate
dynamic disorder leads to a reduction of radiative fluxes, and thus a greater
local relaxation yield both inside and outside the cavity. More importantly,
the slope of the ratio between radiative and local relaxation flux (blue-
dashed line in the subfigure (a)) is reduced in the presence of disorder
outside the cavity. In principle, this slope corresponds to the effective
number of emitters that contribute to the collective emission. For the
parameters in Fig. 13, if we focus on comparing the blue solid and dashed
lines, we find a slope of $0.06$ for moderate disorder (blue, dashed line), vs
$0.1$ for no disorder case (blue solid line), which implies that 60% molecules
are emitting collectively. That being said, for very strong dephasing and
disorder, we expect that this slope should approach zero and $N$ in eq 34
should approach unity; after all, with strong enough dephasing, we expect to
find independent spontaneous emission for emitters outside the cavity – in
contrast to the superradiant effect in eq 34. While this analysis demonstrate
a possible cavity effect on the yield of chemical reaction, a note of caution
should be added. This analysis was based on the assumption that $\Gamma_{loc}$
is identical following excitation of the bright mode outside the cavity as in
the case where a polariton is excited inside the cavity. This is not
necessarily the case because the initially prepared state differs in energy by
an amount determined by the Rabi splitting. A model where $\Gamma_{loc}$ is
sensitive to this difference was recently analyzed in ref 59

Figure 13: (a) The steady state ratio between radiative flux and local
relaxation flux outside/inside the cavity with/without energy dynamic
disorder. (b) The steady state local relaxation yield outside/inside cavity
with/without energy dynamic disorder. We keep all parameters constant and vary
only the number of emitters of the system. We assign the cavity leakage rate
to be the same as the local relaxation rate for the case of the system inside
the cavity. For the case without disorder, we choose the superradiant rate for
$N=10$ to be the same as the local relaxation rate for the case outside the
cavity. Hence, the local yield at $N=10$ is 0.5 both inside and outside the
cavity. For outside the cavity, as more emitters are collectively emitting,
the ratio between the radiative and local relaxation rate increases linearly
as $N$, and thus the local yield decreases. For inside the cavity, we assume
the cavity leakage rate remains approximately constant. Hence, the local
relaxation yield does not depend on $N$.[22] From this data, we conclude that
adding a single mode cavity can certainly change the importance of different
energy dissipation channels. Now, for the case with moderately fast energetic
dynamic disorder ($\sigma=0.2$, $\tau_{c}=1$), one clear pattern is that for
both outside/inside the cavity, the radiative channels are suppressed by the
disorder, and thus, local relaxation yields increase. Moreover, for the case
outside the cavity (blue dashed line), the slope decreases from $0.1$ (no
disorder) to $0.06$ (moderate dynamic disorder), which shows that there are
60% emitters behaving collectively. In applying this conclusion to a realistic
experiment, note that $N$ should be the effective number of emitters that
participate in the superradiant rate/Rabi splitting, which characterize the
collectiveness under the influence of disorder.

##  V Conclusion

In conclusion, we investigate the influence of energetic disorder,
orientational disorder and intermolecular couplings on the absorption spectra
for an ensemble of two-level emitters outside and inside a single-mode cavity.
As far as the influence of energetic and orientational disorder are concerned,
we recover inhomogeneous broadening for static disorder and motional narrowing
for dynamic disorder. For the case inside the cavity, there are two inherently
different observables: transmission spectrum and local relaxation flux. The
transmission spectrum is dictated primarily by properties of the cavity mode
and the local relaxation flux depends mainly on the property of the individual
emitters. These two observables are complementary to each other and together
provide a comprehensive picture of the system when we introduce different
types of disorder. Introducing intermolecular coupling is essentially
equivalent to adding finite detuning between the bare cavity mode and the
emitters eigenstates. As we compare the two observables, we confirm that dark
modes do not contribute to transmission signals which requires the state to be
bright. However, the dark modes do contribute to the local relaxation signals.

Lastly, we investigate the differences between excited state local relaxation
yields outside versus inside the cavity. It is clear that the presence of a
cavity mode will modify the yield in two ways: (1) changing the radiative
energy dissipation rate (from collective spontaneous emission to cavity
leakage), and (2) changing inherently how collectiveness translates into
physical observables (from superradiance to Rabi splitting)[22].

Looking forward, we note that, in this work, we have assumed that the rate of
independent spontaneous emission towards other polarization direction is very
slow and can be completely ignored. In truth, however, the influence of other
polarization directions cannot really be modeled with a single-mode cavity. To
better isolate possible cavity effects in chemistry, one of our next steps
will be to incorporate higher dimensions and more cavity modes inside a
realistic microcavity (likely through a Maxwell-Bloch calculation).

##  Appendix A Generation of Orenstein-Uhlenbeck process

In our computational procedure, for generating energy modulations, we employ
an Orenstein-Uhlenbeck process. To generate a sequence of values at time
$t_{n}=t_{0}+ndt$, we sample $\\{\Omega_{m}(t_{n})\\}$ in the following
process:[60]

  1. 1.

At $t_{0}$, pick $\\{\Omega_{m}(t_{0})\\}$ according to Gaussian distribution

$P[\Omega_{m}(t_{0})]=\frac{1}{\sqrt{2\pi\sigma^{2}}}\exp\Bigl{(}-\frac{\Omega_{m}(t_{0})^{2}}{2\sigma^{2}}\Bigr{)}$

  2. 2.

Calculate $r_{n}$:

$r_{n}=\exp\bigl{(}-(t_{n+1}-t_{n})/\tau_{c}\bigr{)}=\exp\bigl{(}-dt/\tau_{c}\bigr{)}$

  3. 3.

Pick $\\{\Omega_{m}(t_{n+1})\\}$ according to conditional probability

$P[\Omega_{m}(t_{n+1})|\Omega_{m}(t_{n})]=\frac{1}{\sqrt{2\pi\sigma^{2}(1-r_{n}^{2})}}\exp\Bigl{(}-\frac{(\Omega_{m}(t_{n+1})-r_{n}\Omega_{m}(t_{n}))^{2}}{2\sigma^{2}(1-r_{n}^{2})}\Bigr{)}$

  4. 4.

Go back to 2

##  Appendix B Voigt theory

In this appendix, we demonstrate that our numerical results from section III.1
for the case of a collection of molecules outside the cavity with
static/dynamic energetic disorder (Fig. 2) satisfy Voigt theory. According to
eq 26 and eq 27, we can analytically predict the Voigt linewidth for different
disorder strengths $\sigma$ and correlation time $\tau_{c}$. However, we
cannot determine the absolute peak height from eq 26 and eq 27. Nevertheless,
if we simply fit the lineshape with the numerical height, then Fig. 14
demonstrates that the anatlyic linewidths predicted by eq 27 and eq 26 (lack
dots) match perfectly with the numerical linewidths in Fig. 2.

Figure 14: Steady state absorption spectrum for a molecular ensemble with (a)
different static disorder strengths and with (b) different correlation time
$\tau_{c}$ but the same disorder strength $\sigma=0.2$, identical to Fig. 2.
The black dots are predicted by eq 27 and 26. The linewidths of numerical
results match perfectly with the prediction.

##  Appendix C Derivation of the Effective Rabi Splitting under Static
Energetic Disorder

In this appendix, we use second-order perturbation theory to derive eq 28. The
total Hamiltonian of the disordered system is

$\displaystyle\hat{H}_{sub}=$ $\displaystyle H_{0}+\Lambda$ (36)
$\displaystyle H_{0}=$ $\displaystyle E_{v}\left|v\right\rangle\left\langle
v\right|+\sum_{m=1}^{N}E_{m}\left|m\right\rangle\left\langle
m\right|+V_{mv}(\left|m\right\rangle\left\langle
v\right|+\left|v\right\rangle\left\langle m\right|)$ (37)
$\displaystyle\Lambda=$
$\displaystyle\sum_{m=1}^{N}\delta{E}_{m}\left|m\right\rangle\left\langle
m\right|$ (38)

Here $N$ is the effective number of emitters in the cavity. The modulation
$\delta E_{m}$ is treated as perturbation and has zero mean and $\sigma^{2}$
variance. We assume no detuning and thus, the Rabi splitting with no disorder
is

$\Omega_{R}(\sigma=0)=2\sqrt{N}V_{mv}$ (39)

The polariton states have energy $\pm\Omega_{R}(\sigma=0)/2$ and the
corresponding eigenvectors are

$\displaystyle\left|UP/LP\right\rangle=(1/\sqrt{2},\pm 1/\sqrt{2N},...,\pm
1/\sqrt{2N})$ (40)

All remaining eigenstates are dark states and have energy $E_{dark}=0$.
Second-order perturbation theory yields the energies of the polariton states

$\displaystyle E_{LP}=E_{LP}^{0}+\left\langle
LP\right|\Lambda\left|LP\right\rangle+\sum_{k\neq LP}\frac{|\left\langle
k\right|\Lambda\left|LP\right\rangle|^{2}}{E_{LP}^{0}-E_{k}^{0}}$ (41)

The first-order term is the mean of the static disorder.

$\displaystyle\left\langle
LP\right|\Lambda\left|LP\right\rangle=\sum_{m=1}^{N}\delta{E}_{m}/2N=0$ (42)

The second term is

$\displaystyle\sum_{k\neq LP}\frac{|\left\langle
k\right|\Lambda\left|LP\right\rangle|^{2}}{E_{LP}^{0}-E_{k}^{0}}=$
$\displaystyle-\frac{\sum_{k}|\left\langle
k\right|\Lambda\left|LP\right\rangle|^{2}-|\left\langle
LP\right|\Lambda\left|LP\right\rangle|^{2}-|\left\langle
UP\right|\Lambda\left|LP\right\rangle|^{2}/2}{\Omega_{R}(\sigma_{0})/2}$ (43)
$\displaystyle\sum_{k}|\left\langle
k\right|\Lambda\left|LP\right\rangle|^{2}=$ $\displaystyle\left\langle
LP\right|\Lambda^{2}\left|LP\right\rangle=\sum_{m=1}^{N}\delta
E_{m}^{2}/2N=\sigma^{2}/2$ (44) $\displaystyle|\left\langle
UP\right|\Lambda\left|LP\right\rangle|^{2}/2=$
$\displaystyle-\sum_{m=1}^{N}\delta E_{m}/4N=0$ (45)

Hence, the perturbed energy is

$\displaystyle E_{LP}=E_{LP}^{0}-\frac{\sigma^{2}}{\Omega_{R}(\sigma=0)}$ (46)

Last, because the detuning is zero, the spectrum is completely symmetrical.
The effective Rabi splitting is simply twice the absolute value of the lower
polariton energy, which yields eq 28.

##  Appendix D Derivation of the Effective Rabi Splitting under Static and
Dynamic Angular Disorder

In this appendix, we use a cumulant expansion to calculate the effective Rabi
splitting that is relevant with disorder as a function of the disorder
strength $\sigma$ and/or correlation time $\tau_{c}$. At this point, if we
ignore the pumping and non-hermitian decay from eqs 9 \- 13, the Hamiltonian
of the disordered system is

$\displaystyle\hat{H}=$ $\displaystyle H_{0}+\Lambda_{\theta}$ (47)
$\displaystyle H_{0}=$ $\displaystyle E_{v}\left|v\right\rangle\left\langle
v\right|+\sum_{m=1}^{N}E_{m}\left|m\right\rangle\left\langle m\right|$ (48)
$\displaystyle\Lambda_{\theta}=$
$\displaystyle\sum_{m=1}^{N}(V_{mv}(\theta_{m})\left|m\right\rangle\left\langle
v\right|+V_{vm}(\theta_{m})\left|v\right\rangle\left\langle m\right|$ (49)

As shown in eq 19, dynamic angular disorder enters as a modulation of $\theta$
in the coupling $V_{mv}(\theta_{m})=V\cos(\theta_{m})$ between each individual
emitter and the bare cavity mode. It is straightforward to obtain the
instantaneous eigenvalues of $H$ because the characteristic polynomial is

$\displaystyle\lambda^{N-1}(\lambda^{2}-\sum_{m=1}^{N}V^{2}\cos^{2}(\theta_{m}))=0$
(50)

Hence, the effective Rabi splitting is
$\Omega_{R}/\Omega_{R}(\sigma=0)=\sqrt{\langle\cos^{2}\theta_{m}\rangle}$,
where ($\theta_{m}(t)$) satisfies $\langle\theta_{m}(t)\rangle=0$ and
$\langle\theta_{m}(t_{j})\theta_{n}(t_{k})\rangle=\delta_{mn}\sigma^{2}\exp(-|t_{j}-t_{k}|/\tau_{c})$.
By cumulant expansion

$\displaystyle\langle\cos{\theta_{m}}\rangle=$
$\displaystyle\text{Re}(\exp(i\langle\theta_{m}\rangle-\frac{1}{2}\langle\delta\theta_{m}^{2}\rangle))=\exp(-\sigma^{2}/2)$
(51)

Hence, in the static limit ($\tau_{c}=\infty$),

$\displaystyle\langle\cos^{2}\theta_{m}\rangle=\Bigl{\langle}\frac{\cos
2\theta_{m}+1}{2}\Bigr{\rangle}=\frac{(e^{-2\sigma^{2}}+1)}{2}$ (52)

The effective Rabi splitting is

$\displaystyle\Omega_{R}(\sigma,\tau_{c}=\infty)/\Omega_{R}(\sigma=0)=\sqrt{(e^{-2\sigma^{2}}+1)/2}$
(53)

This analytical formula is verified in Fig. 15. The numerical results are
extracted from Fig. 7.

In the fast modulation limit ($\tau_{c}=0$), the external field sees an
averaged upper/lower polariton state (as defined in eq 40 in Appendix C),
which leads to average Rabi splitting

$\displaystyle\Omega_{R}(\sigma,\tau_{c}=0)/\Omega_{R}(\sigma=0)=\langle\cos\theta_{m}\rangle=e^{-\sigma^{2}/2}$
(54)

Unfortunately, we have not yet been able to derive an analytic form for the
effective Rabi splitting for finite $\tau_{c}$. For a general $\tau_{c}$, the
effective Rabi splitting is not determined exclusively by the instantaneous
eigenvalues (which allowed us to calculate the slow and fast limits (eqs. 53
and 54) above). That being said, in Fig. 15, we plot the effective Rabi
splitting for two different values of N (the number of emitters). Note that
the two lines are on top of one another. Thus, even though we do not have a
predictive theory for the Rabi splitting in the intermediate $\tau_{c}$
regime, we can at least be confident that the Rabi splitting does always scale
as $\sqrt{N}$ (which would seem to agree with the results in ref. 61).

Figure 15: Effective Rabi splitting vs static angular disorder strength and
dynamic angular disorder correlation time. As shown in figure (a), static
angular disorder leads to contraction of Rabi splitting. This contraction is
unlike the case with static energetic disorder, in which the effective Rabi
splitting increases as the disorder strength increases. The formula eq 53
obtained by cumulant expansion captures the contraction of effective Rabi
splitting induced by static angular disorder ($\tau_{c}\rightarrow\infty$).
For the case with dynamic angular disorder ($\sigma=0.4\pi$), as shown by the
Red dashed line in figure (b), we can only predict the fast modulation limit
$\tau_{c}=0$ as in eq 54. For the intermediate correlation time, the effective
Rabi splitting is not fully determined by the first and second moment. As
shown by the green dashed figure, the effective Rabi splitting is also
independent of $N$ in the system.

##  Appendix E Symbols, energies/rates/frequencies and values that appear in
this paper

In this Appendix, we will list all parameters in the paper.

[b]  Symbol Energy/Rate/Frequency Value $E_{v}$ Cavity mode energy 0* $E_{0}$
Total ground state energy -1* $E_{m}$ molecular singly excited state energy 0
(for no disorder system) $\hat{V}_{mv}$ coupling between cavity mode
$\left|v\right\rangle$ and molecular excited state $\left|m\right\rangle$
$0.5/\sqrt{N}$$\dagger$ $\hat{V}_{m0}$ coupling between ground state
$\left|0\right\rangle$ and molecular excited state $\left|m\right\rangle$
0.0005 $\omega$ Driving frequency [-1.2, 1.2] $\Gamma_{rad}$ Single molecule
radiative decay rate 0.05/N $\Gamma_{loc}$ Local relaxation rate 0.05
$\Gamma_{cav}(\Gamma_{cav}^{R}+\Gamma_{cav}^{NR})$ Cavity leakage rate
(including radiative and nonradiative) 0.05

Table 1: Symbols, energies/rates/frequencies, and values that appear in this
paper.

  * $\dagger$

We choose to normalize the Rabi splitting for no disorder system and thus all
parameters are in unit of $\Omega_{R}(\sigma=0)$.

  * *

The choice of energy difference between ground and cavity mode
energy/molecular singly excited state energy does not affect the results
because the calculation is under rotating wave approximation. We also set the
energy of the singly excited state energy/cavity mode energy for no disorder
system to be $0$.

## References

  * Dicke [1954] R. H. Dicke, “Coherence in spontaneous radiation processes,” Physical review 93, 99 (1954).
  * Spano and Mukamel [1989] F. C. Spano and S. Mukamel, “Superradiance in molecular aggregates,” The Journal of Chemical Physics 91, 683–700 (1989).
  * Lim  _et al._ [2004] S.-H. Lim, T. G. Bjorklund, F. C. Spano, and C. J. Bardeen, “Exciton Delocalization and Superradiance in Tetracene Thin Films and Nanoaggregates,” Physical Review Letters 92, 107402 (2004).
  * Bonifacio and Lugiato [1975] R. Bonifacio and L. Lugiato, “Cooperative radiation processes in two-level systems: Superfluorescence,” Physical Review A 11, 1507 (1975).
  * Rainò  _et al._ [2018] G. Rainò, M. A. Becker, M. I. Bodnarchuk, R. F. Mahrt, M. V. Kovalenko, and T. Stöferle, “Superfluorescence from lead halide perovskite quantum dot superlattices,” Nature 563, 671–675 (2018).
  * Raimond  _et al._ [1982] J. Raimond, P. Goy, M. Gross, C. Fabre, and S. Haroche, “Collective absorption of blackbody radiation by rydberg atoms in a cavity: an experiment on bose statistics and brownian motion,” Physical Review Letters 49, 117 (1982).
  * 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,” Physical Review Letters 69, 3314 (1992).
  * Heinzen  _et al._ [1987] D. Heinzen, J. Childs, J. Thomas, and M. Feld, “Enhanced and inhibited visible spontaneous emission by atoms in a confocal resonator,” Physical review letters 58, 1320 (1987).
  * Jhe  _et al._ [1987] W. Jhe, A. Anderson, E. A. Hinds, D. Meschede, L. Moi, and S. Haroche, “Suppression of spontaneous decay at optical frequencies: Test of vacuum-field anisotropy in confined space,” Physical review letters 58, 666 (1987).
  * De Martini  _et al._ [1987] F. De Martini, G. Innocenti, G. Jacobovitz, and P. Mataloni, “Anomalous spontaneous emission time in a microscopic optical cavity,” Physical review letters 59, 2955 (1987).
  * Bayer  _et al._ [2001] M. Bayer, T. L. Reinecke, F. Weidner, A. Larionov, A. McDonald, and A. Forchel, “Inhibition and enhancement of the spontaneous emission of quantum dots in structured microresonators,” Physical review letters 86, 3168 (2001).
  * Houdré  _et al._ [1994] R. Houdré, R. Stanley, U. Oesterle, M. Ilegems, and C. Weisbuch, “Room-temperature cavity polaritons in a semiconductor microcavity,” Physical Review B 49, 16761 (1994).
  * Gérard  _et al._ [1998] J.-M. Gérard, B. Sermage, B. Gayral, B. Legrand, E. Costard, and V. Thierry-Mieg, “Enhanced spontaneous emission by quantum boxes in a monolithic optical microcavity,” Physical review letters 81, 1110 (1998).
  * Yoshie  _et al._ [2004] T. Yoshie, A. Scherer, J. Hendrickson, G. Khitrova, H. Gibbs, G. Rupper, C. Ell, O. Shchekin, and D. Deppe, “Vacuum rabi splitting with a single quantum dot in a photonic crystal nanocavity,” Nature 432, 200–203 (2004).
  * Masiello and Schatz [2008] D. J. Masiello and G. C. Schatz, “Many-body theory of surface-enhanced raman scattering,” Physical Review A 78, 042505 (2008).
  * Spano [2010] F. C. Spano, “The Spectral Signatures of Frenkel Polarons in H- and J-Aggregates,” Accounts of Chemical Research 43, 429–439 (2010).
  * Mirsaleh-Kohan  _et al._ [2012] N. Mirsaleh-Kohan, V. Iberi, P. D. Simmons Jr, N. W. Bigelow, A. Vaschillo, M. M. Rowland, M. D. Best, S. J. Pennycook, D. J. Masiello, B. S. Guiton, _et al._ , “Single-molecule surface-enhanced raman scattering: can stem/eels image electromagnetic hot spots?” The Journal of Physical Chemistry Letters 3, 2303–2309 (2012).
  * Spano [2015] F. C. Spano, “Optical microcavities enhance the exciton coherence length and eliminate vibronic coupling in J-aggregates,” The Journal of Chemical Physics 142, 184707 (2015).
  * Herrera and Spano [2016] F. Herrera and F. C. Spano, “Cavity-Controlled Chemistry in Molecular Ensembles,” Physical Review Letters 116 (2016).
  * Herrera and Spano [2017] F. Herrera and F. C. Spano, “Absorption and photoluminescence in organic cavity QED,” Physical Review A 95 (2017).
  * Abid  _et al._ [2017] I. Abid, W. Chen, J. Yuan, A. Bohloul, S. Najmaei, C. Avendano, R. Péchou, A. Mlayah, and J. Lou, “Temperature-dependent plasmon–exciton interactions in hybrid au/mose2 nanostructures,” ACS photonics 4, 1653–1660 (2017).
  * Bisht  _et al._ [2018] A. Bisht, J. Cuadra, M. Wersall, A. Canales, T. J. Antosiewicz, and T. Shegai, “Collective strong light-matter coupling in hierarchical microcavity-plasmon-exciton systems,” Nano letters 19, 189–196 (2018).
  * Hertzog  _et al._ [2019] M. Hertzog, M. Wang, J. Mony, and K. Börjesson, “Strong light–matter interactions: a new direction within chemistry,” Chemical Society Reviews 48, 937–961 (2019).
  * Sidler  _et al._ [2020] D. Sidler, C. Schäfer, M. Ruggenthaler, and A. Rubio, “Polaritonic chemistry: Collective strong coupling implies strong local modification of chemical properties,” The journal of physical chemistry letters 12, 508–516 (2020).
  * Smith, Bhattacharya, and Masiello [2021] K. C. Smith, A. Bhattacharya, and D. J. Masiello, “Exact k-body representation of the jaynes-cummings interaction in the dressed basis: Insight into many-body phenomena with light,” Physical Review A 104, 013707 (2021).
  * Li  _et al._ [2022] T. E. Li, B. Cui, J. E. Subotnik, and A. Nitzan, “Molecular polaritonics: chemical dynamics under strong light–matter coupling,” Annual review of physical chemistry 73, 43–71 (2022).
  * Renaud, Whitley, and Stroud Jr [1977] B. Renaud, R. Whitley, and C. Stroud Jr, “Nonstationary two-level resonance fluorescence,” Journal of Physics B: Atomic and Molecular Physics (1968-1987) 10, 19 (1977).
  * Gross and Haroche [1982] M. Gross and S. Haroche, “Superradiance: An essay on the theory of collective spontaneous emission,” Physics reports 93, 301–396 (1982).
  * Celardo  _et al._ [2013] G. Celardo, A. Biella, L. Kaplan, and F. Borgonovi, “Interplay of superradiance and disorder in the anderson model,” Fortschritte der Physik 61, 250–260 (2013).
  * Biella  _et al._ [2013a] A. Biella, F. Borgonovi, R. Kaiser, and G. L. Celardo, “Subradiant hybrid states in the open 3D Anderson-Dicke model,” EPL (Europhysics Letters) 103, 57009 (2013a).
  * Biella  _et al._ [2013b] A. Biella, F. Borgonovi, R. Kaiser, and G. Celardo, “Subradiant hybrid states in the open 3d anderson-dicke model,” EPL (Europhysics Letters) 103, 57009 (2013b).
  * Delga  _et al._ [2014] A. Delga, J. Feist, J. Bravo-Abad, and F. Garcia-Vidal, “Quantum emitters near a metal nanoparticle: strong coupling and quenching,” Physical review letters 112, 253601 (2014).
  * Celardo, Giusteri, and Borgonovi [2014] G. L. Celardo, G. G. Giusteri, and F. Borgonovi, “Cooperative robustness to static disorder: Superradiance and localization in a nanoscale ring to model light-harvesting systems found in nature,” Physical Review B 90, 075113 (2014).
  * Pustovit, Urbas, and Shahbazyan [2014] V. N. Pustovit, A. M. Urbas, and T. V. Shahbazyan, “Energy transfer in plasmonic systems,” Journal of Optics 16, 114015 (2014).
  * Goban  _et al._ [2015] A. Goban, C.-L. Hung, J. Hood, S.-P. Yu, J. Muniz, O. Painter, and H. Kimble, “Superradiance for atoms trapped along a photonic crystal waveguide,” Physical review letters 115, 063601 (2015).
  * Norcia  _et al._ [2016] M. A. Norcia, M. N. Winchester, J. R. Cline, and J. K. Thompson, “Superradiance on the millihertz linewidth strontium clock transition,” Science advances 2, e1601231 (2016).
  * Asenjo-Garcia  _et al._ [2017] A. Asenjo-Garcia, M. Moreno-Cardoner, A. Albrecht, H. Kimble, and D. E. Chang, “Exponential improvement in photon storage fidelities using subradiance and “selective radiance” in atomic arrays,” Physical Review X 7, 031024 (2017).
  * Kirton and Keeling [2017] P. Kirton and J. Keeling, “Suppressing and restoring the dicke superradiance transition by dephasing and decay,” Physical review letters 118, 123602 (2017).
  * Shahbazyan, Raikh, and Vardeny [2000] T. Shahbazyan, M. Raikh, and Z. Vardeny, “Mesoscopic cooperative emission from a disordered system,” Physical Review B 61, 13266 (2000).
  * Sukharev and Nitzan [2017] M. Sukharev and A. Nitzan, “Optics of exciton-plasmon nanomaterials,” Journal of Physics: Condensed Matter 29, 443003 (2017).
  * Fauché  _et al._ [2017] P. Fauché, C. Gebhardt, M. Sukharev, and R. A. Vallée, “Plasmonic opals: observation of a collective molecular exciton mode beyond the strong coupling,” Scientific reports 7, 1–9 (2017).
  * Ribeiro  _et al._ [2018] R. F. Ribeiro, L. A. Martínez-Martínez, M. Du, J. Campos-Gonzalez-Angulo, and J. Yuen-Zhou, “Polariton chemistry: controlling molecular dynamics with optical cavities,” Chemical science 9, 6325–6339 (2018).
  * Herrera and Spano [2018] F. Herrera and F. C. Spano, “Theory of Nanoscale Organic Cavities: The Essential Role of Vibration-Photon Dressed States,” ACS Photonics 5, 65–79 (2018).
  * Gómez-Castaño  _et al._ [2019] M. Gómez-Castaño, A. Redondo-Cubero, L. Buisson, J. L. Pau, A. Mihi, S. Ravaine, R. A. Vallée, A. Nitzan, and M. Sukharev, “Energy transfer and interference by collective electromagnetic coupling,” Nano letters 19, 5790–5795 (2019).
  * Wang, Scholes, and Hsu [2019] S. Wang, G. D. Scholes, and L.-Y. Hsu, “Quantum dynamics of a molecular emitter strongly coupled with surface plasmon polaritons: A macroscopic quantum electrodynamics approach,” The Journal of Chemical Physics 151, 014105 (2019).
  * Wang, Scholes, and Hsu [2020a] S. Wang, G. D. Scholes, and L.-Y. Hsu, “Coherent-to-incoherent transition of molecular fluorescence controlled by surface plasmon polaritons,” The journal of physical chemistry letters 11, 5948–5955 (2020a).
  * Mattiotti  _et al._ [2020] F. Mattiotti, M. Kuno, F. Borgonovi, B. Jankó, and G. L. Celardo, “Thermal decoherence of superradiance in lead halide perovskite nanocrystal superlattices,” Nano Letters 20, 7382–7388 (2020).
  * Wang  _et al._ [2020] S. Wang, M.-W. Lee, Y.-T. Chuang, G. D. Scholes, and L.-Y. Hsu, “Theory of molecular emission power spectra. I. Macroscopic quantum electrodynamics formalism,” The Journal of Chemical Physics 153, 184102 (2020).
  * Spano [2020] F. C. Spano, “Exciton–phonon polaritons in organic microcavities: Testing a simple ansatz for treating a large number of chromophores,” The Journal of Chemical Physics 152, 204113 (2020).
  * Wang, Scholes, and Hsu [2020b] S. Wang, G. D. Scholes, and L.-Y. Hsu, “Coherent-to-Incoherent Transition of Molecular Fluorescence Controlled by Surface Plasmon Polaritons,” The Journal of Physical Chemistry Letters 11, 5948–5955 (2020b).
  * Lee, Chuang, and Hsu [2021] M.-W. Lee, Y.-T. Chuang, and L.-Y. Hsu, “Theory of molecular emission power spectra. II. Angle, frequency, and distance dependence of electromagnetic environment factor of a molecular emitter in plasmonic environments,” The Journal of Chemical Physics 155, 074101 (2021).
  * Wang, Chuang, and Hsu [2021] S. Wang, Y.-T. Chuang, and L.-Y. Hsu, “Simple but accurate estimation of light–matter coupling strength and optical loss for a molecular emitter coupled with photonic modes,” The Journal of Chemical Physics 155, 134117 (2021), publisher: American Institute of Physics.
  * Herrera and Litinskaya [2022] F. Herrera and M. Litinskaya, “Disordered ensembles of strongly coupled single-molecule plasmonic picocavities as nonlinear optical metamaterials,” The Journal of Chemical Physics 156, 114702 (2022).
  * Kavokin  _et al._ [2017] A. V. Kavokin, J. J. Baumberg, G. Malpuech, and F. P. Laussy,  _Microcavities_ , Vol. 21 (Oxford university press, 2017).
  * Akram, Ficek, and Swain [2000] U. Akram, Z. Ficek, and S. Swain, “Decoherence and coherent population transfer between two coupled systems,” Physical Review A 62, 013413 (2000).
  * Note [1] The choice of intermolecular coupling is $0.1$ in unit of Rabi splitting. Rabi splitting is approximately on the order of $10meV$, and if the distance between dipoles are $2nm$, the dipole moments are on the order of $10Debye$.
  * Note [2] This formula is by modifying equation (7.104) in ref. 62 by taking both decay rates $N\Gamma_{rad}+\Gamma_{loc}$ into account according to Voigt theorem.
  * Houdré, Stanley, and Ilegems [1996] R. Houdré, R. Stanley, and M. Ilegems, “Vacuum-field rabi splitting in the presence of inhomogeneous broadening: Resolution of a homogeneous linewidth in an inhomogeneously broadened system,” Physical Review A 53, 2711 (1996).
  * Cui and Nizan [2022] B. Cui and A. Nizan, “Collective response in light–matter interactions: The interplay between strong coupling and local dynamics,” The Journal of Chemical Physics 157, 114108 (2022).
  * Rybicki [1994] G. B. Rybicki, “Notes on gaussian random functions with exponential correlation functions (ornstein-uhlenbeck process),” https://www.lanl.gov/DLDSTP/fast/OU_process.pdf (1994).
  * Pérez-Sánchez  _et al._ [2022] J. B. Pérez-Sánchez, A. Koner, N. P. Stern, and J. Yuen-Zhou, “An effective single molecule model that simulates dynamics of collective strong light-matter coupling,” arXiv preprint arXiv:2209.04955 (2022).
  * Nitzan [2006] A. Nitzan,  _Chemical dynamics in condensed phases: relaxation, transfer and reactions in condensed molecular systems_ (Oxford university press, 2006).