# Scaling law for Kasha’s rule in photoexcited molecular aggregates

R. Holzinger  Institut für Theoretische Physik, Universität Innsbruck,
Technikerstrasse 21a, A-6020 Innsbruck, Austria  N. S. Bassler  Max Planck
Institute for the Science of Light, Staudtstraße 2, D-91058 Erlangen, Germany
Department of Physics, Friedrich-Alexander-Universität Erlangen-Nürnberg,
Staudtstraße 7, D-91058 Erlangen, Germany  H. Ritsch  Institut für
Theoretische Physik, Universität Innsbruck, Technikerstrasse 21a, A-6020
Innsbruck, Austria  C. Genes  Max Planck Institute for the Science of Light,
Staudtstraße 2, D-91058 Erlangen, Germany  Department of Physics, Friedrich-
Alexander-Universität Erlangen-Nürnberg, Staudtstraße 7, D-91058 Erlangen,
Germany

###### Abstract

We study the photophysics of molecular aggregates from a quantum optics
perspective, with emphasis on deriving scaling laws for the fast non-radiative
relaxation of collective electronic excitations, referred to as Kasha’s rule.
Aggregates exhibit an energetically broad manifold of collective states with
delocalized electronic excitations originating from near field dipole-dipole
exchanges between neighboring monomers. Photo-excitation at optical
wavelengths, much larger than the monomer-monomer average separation,
addresses almost exclusively symmetric collective states, which for an
arrangement known as H-aggregate, show an upward hypsochromic shift. The
extremely fast subsequent non-radiative relaxation via intramolecular
vibrational modes populates lower energy, subradiant states, resulting in an
effective inhibition of fluorescence. Our analytical treatment allows for the
derivation of an approximate scaling law of this relaxation process, linear in
the number of available low energy vibrational modes and directly proportional
to the dipole-dipole interaction strength between neighbouring monomers.

Molecular aggregates Hestand and Spano (2018); Saikin _et al._ (2013); Ma _et
al._ (2021) are self-ordered arrangements of monomers showing strong
collective optical transition dipole strengths. Owing to the dense packing of
monomers, with monomer-monomer separations at the order of a nanometer and
total aggregate lengths of tens or hundred of nanometers, thus much below an
optical wavelength and despite inhomogeneous broadening and separation
disorder, they exhibit delocalized excitons Abramavicius _et al._ (2009). This
allows strong coupling to external light modes, resulting in collectively
modified fluorescence rates. Following the discovery of J- and H-aggregates in
the $1930$s by Scheibe Scheibe _et al._ (1937) and Jelley Jelley (1936) their
standard understanding is based on the original approach introduced by Kasha
in the 1960s Kasha (1963). Currently, J-aggregates are widely employed in
light-matter coupling experiments aiming at the modification of material
properties via the manipulation of the electromagnetic vacuum mode density
around electronic resonances Garcia-Vidal _et al._ (2012). Numerous approaches
to understand their aggregation behavior and the subsequent response to
external illumination have been taken, as thoroughly reviewed in Ref. Hestand
and Spano (2018).
There has been recently a new surge in interest to theoretically understand
the behavior of organic molecular systems under the action of confined vacuum
fields, i.e. within the context of cavity quantum electrodynamics Berman
(1994); Haroche and Kleppner (1989); Walther _et al._ (2006); Haroche and
Raimond (2013). Among many approaches taken, a purely quantum optics formalism
has been introduced, which is based on a simple model for electron-vibron
interactions contained within the formalism of the Holstein Hamiltonian and
which can incorporate loss of photons or vibrations via open system approach
following either a master equation of a quantum Langevin equations approaches
Reitz _et al._ (2019); Holzinger _et al._ (2022).

Figure 1: Illustration of the Kasha process in an H-aggregate for a dimer
(left) and larger aggregate (right). The population migrates from the
optically addressed state $\ket{E}$ (superradiant, therefore color coded in
red) towards the lower and subradiant collective states $\ket{A}$ (color coded
towards blue) via a multitude of vibrational relaxation channels (indicated as
green arrows). The energy splitting in the large aggregate is four times the
nearest neighbor dipole-dipole coupling strength $\Omega$ and strongly depends
on the aggregate separation at the level of nanometers. A resulting effective
relaxation rate from the high energy optically addressed state to the bottom
of the energy manifold results, on the femtosecond timescale.

Within such a quantum optics approach, rate equations for the dynamics of a
molecular aggregate upon photo-excitation with light at visible wavelengths
$\lambda$ can be derived. As pointed out in Ref. Hestand and Spano (2018), the
dynamics of the exciton in the aggregate can be understood in terms of bright
(or symmetric) and dark (or asymmetric) states. In the simplest example, a
molecular dimer in the side-by-side H-aggregate configuration (see Fig.2,
exhibits a symmetric superposition $\ket{S}=(\ket{eg}+\ket{ge})/\sqrt{2}$
sitting at energy $\omega_{0}+\Omega$. Here $\ket{eg}$ denotes that the first
monomer is in the excited electronic state while the second one is in the
ground state and $\Omega$ is the nearest neighbor dipole-dipole interaction
scaling with the monomer separation as $d^{-3}$. As $d$ is typically much
smaller than $\lambda$, the asymmetric state, at energy $\omega_{0}-\Omega$,
is decoupled from external light and instead is reached indirectly from the
symmetric state $\ket{S}$ via redistribution of energy into molecular
vibrations, consistent with the standard understanding of a Kasha process.
This can be easily generalized to a mesoscopic sized aggregate, where an
additional large number of dark states fill the $4\Omega$ frequency split
between the photoexcited symmetric state and the lowest energy state. The
Kasha incoherent loss of energy is mediated by the multitude of vibrational
relaxation pathways. For a molecule with $\mathcal{N}$ atoms a large number of
the order $3\mathcal{N}-5$ vibrational modes are available; some of them, up
to some index $n_{\text{max}}$, find themselves within the interval $4\Omega$
thus providing resonant transfer between the set of collective electronic
states.

We show that an approach based on rate equations allows for an approximate
expression for the rate of a Kasha process as

$\kappa_{\mathcal{S}}\approx\frac{4s\Omega}{3}\frac{(n_{\text{max}}+1)(2n_{\text{max}}+1)}{n_{\text{max}}}.$
(1)

The simplicity of the result lies in the small number of parameters involved:
the average estimated Huang-Rhys factor
$s=\sum_{j=1}^{n_{\text{max}}}s_{j}/{n_{\text{max}}}$, the monomer separation
$d$ (appearing via $\Omega$) and the estimated number of intramolecular
vibrations $n_{\text{max}}$ of relevant frequencies. As an example, we
consider the aqueous dye molecule Cresyl Violet which aggregates in an
H-configuration, featuring more then $30$ vibrational modes in the range $\sim
5-50$ THz with a zero-phonon transition wavelength $\lambda_{0}=590$ nm and a
monomer fluorescence linewidth $\gamma_{0}=454$ MHz Jafari _et al._ (2007).
The nearest neighbour dipole-dipole coupling is estimated at around
$\Omega=23$ THz and the vibrational spectrum is extracted from Ref. Barclay
_et al._ (2022). We estimate around $n_{\mathrm{max}}=12$ vibrational modes to
have a non-negligible Huang-Rhys factor larger than $0.01$ and to fall into
the range of the electronic energy band $[0,4\Omega]$; this leads to an
average estimated Huang-Rhys factor $s=0.09$. The estimated timescale for the
Kasha process is then around $\kappa_{\mathcal{S}}^{-1}\approx 13$ fs. Similar
estimates for Rhodamine $800$ Christensson _et al._ (2010) and
Bacteriochlorophyll a (BChla) Rätsep _et al._ (2011) give values for the Kasha
process timescale of order $30$ fs and $34$ fs, respectively.

Figure 2: (a) Diagramatic description of electron-vibron interactions within
each monomer. A single optically addressable electronic transition operator
$\sigma$ is coupled to $n$ vibrational mode operators $b_{m}$ with frequencies
$\nu_{m}$, linewidths $\Gamma_{m}$ and Huang-Rhys factors $s_{m}$, where $m$
runs from $1$ to $n$. The vibrational spectral density $J(\nu)$ has a number
vibrational frequencies $\nu_{m}$ and associated relaxation rates $\Gamma_{m}$
which dictate their linewidths. (b) Ground and excited electronic pootential
energy landscapes along a single nuclear coordinate under the harmonic
assumption. The vertical lines illustrate the Franck-Condon principle where
substantial transition rates only are possible when a match between the
symmetry of the nuclear wavepacket in the ground and excited electronic
states.

Photophysics of molecular aggregates \- The minimal model we consider involves
a one dimensional molecular chain of identical monomers with a typical inter-
monomer distance $d$ in the nm range. As photo-excitation with wavelength
$2\pi/k_{0}$ (with $k_{0}$ the wavevector) is in the $\mu$m range, the
condition $k_{0}d\ll 1$ is fulfilled meaning that the chain is uniformly
illuminated. Each monomer is assumed to undergo a single electronic transition
which in turn is coupled to a $n$ vibrational modes (see Fig. 2). Each
vibrational mode has a frequency $\nu_{m}$ and relaxation rate $\Gamma_{m}$,
where $m$ runs from $1$ to $n$. The monomer can undergo spontaneous emission
at rate $\gamma_{0}$, owing to the coupling to the electromagnetic
environment. For each monomer $j$, the electronic transition is at frequency
splitting $\omega_{0}$ ($\hbar=1$) and is described by the collapse operator
$\sigma_{j}=\ket{g}_{j}\bra{e}_{j}$. The vibrational degrees of freedom are
described by bosonic operators $b_{jm}$ satisfying the commutation relations
$\left[b_{jm},b_{j^{\prime}m^{\prime}}^{\dagger}\right]=\delta_{jj^{\prime}}\delta_{mm^{\prime}}$.
The vibronic couplings are illustrated in Fig. 2(a) as links between the
electronic and vibration operators with magnitude characterized by the Huang-
Rhys factors $s_{m}$. The electronic and vibrational degrees of freedom are
subject to loss quantified by the spontaneous emission rate $\gamma_{0}$ and
by the vibrational relaxation rates $\Gamma_{m}$, respectively. The
Hamiltonian for all $\mathcal{N}$ monomers is obtained as a sum over each
particle’s Hamiltonian

$h^{(j)}=\bar{\omega}_{0}\sigma^{\dagger}_{j}\sigma_{j}+\sum_{m=1}^{n}\nu_{m}\Big{(}b_{jm}^{\dagger}b_{jm}-\sqrt{s_{m}}\sigma_{j}^{\dagger}\sigma_{j}(b_{jm}^{\dagger}+b_{jm})\Big{)}.$
(2)

Notice that the bare frequency is Stokes shifted
$\bar{\omega}_{0}=\omega_{0}+\sum_{m}s_{m}\nu_{m}$. This shift will later be
eliminated after a polaron transformation (see SI for more details). The model
and its validity and relevance is extensively discussed in Ref. Reitz _et al._
(2019).

Dipole-dipole exchanges at rates $\Omega_{jj^{\prime}}$ have a strong imprint
at nm distances, owing to their scaling with the inverse cube of the particle
separation $|\vec{r}_{j}-\vec{r}_{j^{\prime}}|^{-3}$ in the near field region
Reitz _et al._ (2022). This can be listed in the Hamiltonian as
$\mathcal{H}_{\text{d-d}}=\textstyle\sum_{j\neq
j^{\prime}}\Omega_{jj^{\prime}}\sigma^{\dagger}_{j}\sigma_{j^{\prime}}$ and
describes an excitation transfer between pairs of monomers via a virtual
photon exchange. The coherent exchange is mediated by the dipole-dipole
frequency shifts $\Omega_{jj^{\prime}}$, which in units of the optical
emission rate $\gamma_{0}$ are given by
$\Omega_{jj^{\prime}}/\gamma_{0}=-3\pi/k_{0}\ \vec{\mu}^{*}\cdot\mathrm{Re}\
\boldsymbol{G}(\vec{r}_{j},\vec{r}_{j^{\prime}},\omega_{0})\cdot\vec{{\mu}}$,
namely, proportional to the real part of the Green’s tensor in free space
(shown in the Supporting Information). In the following we will consider the
particular case of side-by-side arrangement, where all transition dipoles
$\vec{\mu}$ are parallel to each other and perpendicular to the chain
direction.

The model can be better tackled in a collective basis Holzinger _et al._
(2022) resulting from the diagonalization of the first excitation subspace
including dipole-dipole interactions. This can be done either in the discrete
space where states are indexed from $1$ to $\mathcal{N}$ (see Fig. 3(a)) or
alternatively, with the quasi-momentum index $q=2\pi k/(\mathcal{N}d)$
obtained by a re-scaling of the index of the mode $k=\pm
1,...,\pm(\mathcal{N}-1)/2$ (see Fig. 3(b)). This is based on the observation
that, for aggregates with $\mathcal{N}$ of the order of hundred monomers, the
mesoscopic limit can be assumed such that the system can be considered
translationally invariant and periodic boundary conditions can be invoked. A
single symmetric mode can be distinguished with state vector obtained by the
application of the symmetric operator
$\mathcal{S}^{\dagger}=\sum_{j}\sigma_{j}^{\dagger}/\sqrt{\mathcal{N}}$ to the
collective ground state $\ket{G}$. To a good approximation the system can be
considered to be in the Dicke limit where a single superradiant emission rate
roughly estimated by $\gamma_{\mathcal{S}}=\mathcal{N}\gamma_{0}$
characterizes this ’bright’ state. In addition, the other orthogonal
$\mathcal{N}-1$ asymmetric states are obtained via the application of
asymmetric operators
$\mathcal{A}_{q}=\sum_{j=1}^{\mathcal{N}}{e^{iqjd}}\sigma_{j}/\sqrt{\mathcal{N}}$
to $\ket{G}$. These states are non-radiative at such deep subwavelength
molecular separations and we dub them therefore as ’dark’.

The collective excitations are eigenstates of the dipole-dipole interaction
Hamiltonian
$\mathcal{H}_{\mathrm{d-d}}\mathcal{A}_{q}|g\rangle^{\otimes\mathcal{N}}=\Omega_{q}\mathcal{A}_{q}|g\rangle^{\otimes\mathcal{N}}$
where by definition we fix the symmetric shift
$\Omega_{\mathcal{S}}\equiv\Omega_{q=0}$. With periodic boundary conditions
imposed, in the mesoscopic limit, one can derive the collective shifts as
$\Omega_{q}=2\sum_{j}\Omega_{1j}\cos(q(j-1)d)$ Asenjo-Garcia _et al._ (2017);
Holzinger _et al._ (2022). Further simplifications occur by considering the
nearest neighbour approximation and the collective eigenenergies become
$\Omega_{q}=2\Omega\cos(qd)$, where the nearest-neighbour coupling is simply
denoted by $\Omega\equiv\Omega_{12}$. Closely following the procedure
introduced in Ref. Holzinger _et al._ (2022), one can analyze the coupling
between states of different symmetries via electron-vibron couplings by an
additional transformation to a collective basis for the vibrational degrees of
freedom as well. This is done by introducing collective vibrational modes
$Q_{q}^{(m)}=\sum_{j=1}^{\mathcal{N}}(b_{m,j}+b_{m,j}^{\dagger})e^{iqjd}/\sqrt{\mathcal{N}},$
with the momentum quadratures satisfying $[Q^{(m)}_{q},P^{(m)}_{q}]=i$ and $m$
labels the vibrational mode running from $1$ to $n$.

Dynamics of the Kasha process in H-aggregates \- Owing to small dimension of a
molecular aggregate compared to an optical wavelength, under photo-excitation
only the symmetric collective mode is directly accessible. This immediately
leads to a rescaling of the Rabi driving by a factor of $\sqrt{\mathcal{N}}$:
this can be seen equivalently as an increase of the oscillator strength by
$\sqrt{\mathcal{N}}$, thus rendering aggregates of any kind as good candidates
for strong light-matter coupling. In addition, the particularity of the
subsequent aggregate electronic dynamics lies within the shape of the energy
band. For example, J-aggregates present an energy band where the symmetric
mode lies at the bottom of the band thus leading to a bathochromic frequency
shift (to the left of the bare monomer frequency) and subsequently shows not
only an enhanced absorption cross section but also enhanced fluorescence at a
superradiant rate. In contrast, H-aggregates have the symmetric state located
at the top of the energy band corresponding to a hypsochromic shift. Most
importantly, quick dynamics follows the optical excitation involving the
relaxation of the collective state towards low energy dark states. This
dynamics is illustrated in Fig. 3 in both the discrete and continuous cases.
In Fig. 3 a numerical simulation for $\mathcal{N}=100$ at a inter-molecular
spacing $k_{0}d=0.0126$. We have considered a large n spanning the range of
frequencies from $0$ to $4\Omega$, with $\Gamma_{m}=\nu_{m}/10$ and identical
$s_{m}=0.01$.

Figure 3: Dynamics of the collective state populations upon photoexcitation of
the symmetric state in H-aggregates. (a) Illustration of the population
migration in time from the optically excited state (superradiant and therefore
color coded in red) towards the lower and subradiant collective states (color
coded towards blue) via a multitude of vibrational relaxation channels. (b)
Exact numerical time dynamics for a chain of $\mathcal{N}=100$ monomers. In
such a mesoscopic limit, the analysis can be performed on a band diagram where
the symmetric state, photoexcitable, sits at $q=0$. The numerical simulation
shows the time dynamics of populations on a time scale normalized to
$\Omega^{-1}$.

Such time dynamics can be analytically tackled at the level of rate equations.
We largely follow the derivation in Ref. Holzinger _et al._ (2022) which we
generalize here to incorporate the crucial aspect that many vibrational modes
have to be taken into account. Under the assumption that the vibrational
relaxation rates are fast compared to the coherent couplings and radiative
loss rates, a set of rate equations for the populations of the symmetric state
$p_{\mathcal{S}}=\langle\mathcal{S}^{\dagger}\mathcal{S}\rangle$ and all dark
states $p_{q}=\langle\mathcal{A}_{q}^{\dagger}\mathcal{A}_{q}\rangle$ can be
derived (see SI). With the definitions $\kappa_{\mathcal{S}}$ \- total loss
rate of the symmetric state, $\kappa_{q}$ \- loss rate for dark state $q$,
$\kappa_{q\rightarrow\mathcal{S}}$ \- incoherent repopulation rate from the
dark to the bright state and $\kappa_{q^{\prime}\rightarrow q}$ \- incoherent
rate for redistribution of energy within the dark state manifold, one can
write

$\displaystyle\dot{p}_{\mathcal{S}}$
$\displaystyle=-(\gamma_{\mathcal{S}}+\kappa_{\mathcal{S}})p_{\mathcal{S}}+\sum_{q\neq
0}\kappa_{q\rightarrow\mathcal{S}}p_{q},$ (3a) $\displaystyle\dot{p}_{q}$
$\displaystyle=-\kappa_{q}p_{q}+\sum_{q^{\prime}\neq
q}\kappa_{q^{\prime}\rightarrow q}p_{q^{\prime}}.$ (3b)

The rate equations show that the symmetric state energy spills into the whole
dark state manifold via rate $\kappa_{\mathcal{S}}$ and in addition, higher
energy dark states spill into the lower energy ones via $\kappa_{q}$. The
quasi-unidirectionality of the process is ensured by the fact that, in this
perturbative treatment, the coherent coupling between states is followed by
quick vibrational relaxation, making the reverse process, governed by rates
$\kappa_{q\rightarrow\mathcal{S}}$ and $\kappa_{q^{\prime}\rightarrow q}$ from
lower energy state to higher ones, very unlikely (as shown in Ref. Holzinger
_et al._ (2022) for the two monomer case). Mathematically, the condition is
$\sqrt{s_{m}}\nu_{m}/\sqrt{\mathcal{N}}\ll\Gamma$. Analytically, one can get
an expression for the transfer rate between the symmetric mode and any dark
mode $q$ mediated by vibrational mode $m$ as

$\kappa^{(m)}_{\mathcal{S}\rightarrow
q}=\frac{2s_{m}\nu_{m}^{2}(\Gamma_{m}+\gamma_{\mathcal{S}})/{\mathcal{N}}}{(\Gamma_{m}+\gamma_{\mathcal{S}})^{2}+4(\Omega_{\mathcal{S}}-\Omega_{q}-\nu_{m})^{2}}$
(4)

and construct the total rate to all states spanned by the index $q$ by summing
$\kappa_{\mathcal{S}}=\sum_{m=1}^{n_{\text{max}}}\sum_{q}\kappa^{(m)}_{\mathcal{S}\rightarrow
q}$ over all vibrations up to an index $n_{\text{max}}$ within the frequency
interval covered by $4\Omega$ where the electronic collective states are
positioned in energy. Moreover, we will consider the standard underdamped
harmonic oscillator model for the molecular vibrations, i.e. the dissipation
rate is much smaller than the resonance frequency for any mode
$\Gamma_{m}\ll\nu_{m}$. Equivalently, one can state that the quality factor of
any vibrational mode is much larger than unity $\nu_{m}/\Gamma_{m}\gg 1$.

Figure 4: (a) Only a limited number of vibrational modes, from $1$ to
$n_{\mathrm{max}}$ can efficiently mediate resonant transfer between the
symmetric state and the dark manifold. These modes are within the $4\Omega$
bandwidth. (b) Linear scaling of the transfer rate $\kappa_{\mathcal{S}}$ as a
function of the number of vibrational modes, with the same parameters as in
(c) and $\mathcal{N}=20$ molecules at $k_{0}d=0.0126$ separation. The scaling
law from Eq. (5) provides a very good fit to the rate equations in Eq. (4).
(c) Time dynamics of $20$ molecules initialized in the symmetric state.
Numerical results in the single excitation manifold show excellent agreement
with an exponential decay given by $e^{-\kappa_{\mathcal{S}}t}$, and governed
by the analytical formula in Eq. (5) (Same parameters as in (d)-(e)). (d)
Numerical results for $\kappa_{\mathcal{S}}$ plotted against increasing
$n_{\text{max}}$ for a randomly drawn set of vibrational frequencies
$\nu_{1},...\nu_{n_{\text{max}}}$ in the window $0$ to $4\Omega$ and Huang-
Rhys factors $s_{1},...s_{n_{\text{max}}}$ in the interval $[0,0.2]$. The fit
is performed with the analytical scaling in Eq. (5) for evenly spaced
vibrational frequencies and Huang-Rhys factor at the level of the distribution
average $s=0.1$. (e) Further comparison of the analytical scaling from Eq. (5)
with an average over $200$ random realizations shows almost perfect agreement.
(f) Decay of the initial symmetric state population under the influence of
random static frequency disorder with fluctuation $\delta\omega$ around
$\omega_{0}$. Further parameters in all plots: $\mathcal{N}=20$ molecules,
$k_{0}d=0.0126$, $\Gamma_{m}=\nu_{m}/10$, $\Omega\approx 3.759\times
10^{5}\gamma_{0}$.

In order to further proceed with analytical estimates, let us first make some
comments regarding typical timescales. Given that vibrational relaxation is in
the order of tens to hundreds of GHz while spontaneous emission is in the
range of tens of MHz, a very quick non-radiative path from the symmetric to
low energy asymmetric states can be achieved on ps timescales. For monomer
separations in the nm range, expected near field shifts $\Omega$ in the range
of $1$ THz to tens of THz are expected. This means that only a few, low
energy, molecular vibrations can fit in the window of $4\Omega$ (see Fig.
4(a)) and aid the relaxation process. This allows us to derive an approximate
scaling law for $\kappa_{\mathcal{S}}$ as a function of a given number of
vibrational modes $n_{\text{max}}$ that can efficiently mediate the relaxation
of the symmetric state into the dark state manifold. We proceed by first
consider a given vibrational mode $m$ and asking for the condition that this
mode can transfer excitation from the symmetric mode to any of the dark
states. We notice that with the condition that $\Gamma_{m}\ll\nu_{m}$ (it is
also implied that $\gamma_{\mathcal{S}}\ll\nu_{m}$ even for large
$\mathcal{N}$) the resonance condition requires that the states to which
resonant transfer can take place are only in the vicinity of the mode $q$
fulfilling $\Omega_{q}=\Omega_{\mathcal{S}}-\nu_{m}$. Of course, in this case
one can immediately observe that any modes with $\nu_{m}>4\Omega$ cannot take
part in this transfer. Assuming a constant density of all $\mathcal{N}$
collective states spread within the interval $4\Omega$, we can then estimate
that a number of approximately $\epsilon_{m}=\Gamma_{m}\mathcal{N}/(4\Omega)$
states fall close to the resonance $\Omega_{q}=\Omega_{\mathcal{S}}-\nu_{m}$,
i.e. within the linewidth $\Gamma_{m}$. Summing over all these contributions
gives the total rate for all transitions mediated by mode $m$ to states close
to $q$. Summing over all possible relaxation paths that participate in the
transfer giving thus an estimate for the total rate

$\kappa_{\mathcal{S}}=\sum_{q}\kappa^{(m)}_{\mathcal{S}\rightarrow
q}\approx\sum_{m=1}^{n_{\text{max}}}\frac{\Gamma_{m}}{2\Omega}\frac{s_{m}\nu_{m}^{2}}{(\Gamma_{m}+\gamma_{\mathcal{S}})}.$
(5)

In order to estimate this quantity, knowledge of the particular nature of the
monomer’s frequencies, Huang-Rhys factors and vibrational relaxation rates is
necessary. However, while later we will numerically investigate random
distributions of frequencies and Huang-Rhys factors, we aim first at deriving
a simple scaling law. To this end we will proceed by making some simplifying
assumptions, among which the first is that the vibrational spectrum is equally
spaced in the interval from $0$ up to $4\Omega$. We denote the frequency of
mode $m$ by $\nu_{m}=m4\Omega/n_{\text{max}}$. Let us also consider that all
Huang-Rhys factors are equal to $s$ (later we compare with a randomized
distribution with an average $s$). Moreover, we neglect the contribution of
$\gamma_{\mathcal{S}}$ as it is much smaller than $\Gamma_{m}$ (this should be
typically very well fulfilled as $\gamma_{0}/\Gamma_{m}$ is expected to be
around $10^{-5}$). Summing over all vibrational modes within the interval of
$4\Omega$ gives us an approximated scaling law as the expression already
introduced in Eq. (1). The result shows independence of the total number of
monomers $\mathcal{N}$ and a quasi linear dependence on the total number of
available low frequency vibrational modes which can resonantly participate in
Kasha’s relaxation process from the high energy symmetric state to the bottom
of the dark state manifold. Notice that the predicted timescale is dictated by
the nearest neighbor dipole-dipole coupling strength $\Omega$ which in turn
depends on the inverse cube of the monomer-monomer separation. In a first
step, we can estimate that the analytical scaling is in very good agreement
with the results of the rate equations, as seen in Fig. 4(b).

Dye molecule | d(nm) | $|\vec{\mu}|(\mathrm{D})$ |  $\gamma_{0}$(MHz) | $n_{\mathrm{max}}$ | $s$ |  $\kappa_{\mathcal{S}}^{-1}$(fs)
---|---|---|---|---|---|---
Cresyl Violet |  $\sim$ 2 | 1.8 | 455 | 12 | 0.09 | 13
Rhodamine 800 |  $\sim$ 2 | 1.8 | 525 | 16 | 0.067 | 30
BChl a |  $\sim$ 1 | 6.4 | 335 | 25 | 0.07 | 34
Table 1: Estimates of the timescale $\kappa_{\mathcal{S}}^{-1}$ (femtoseconds)
of the Kasha process according to Eq. (1). The specific dyes considered are
Cresyl Violet (CV) with $\Omega_{0}\approx 23$ THz Jafari _et al._ (2007);
Barclay _et al._ (2022), Rhodamine 800 with $\Omega_{0}\approx 8$ THz
Christensson _et al._ (2010) and Bacteriochlorophyll a (BChla) with
$\Omega_{0}\approx 6$ THz, measured at room temperature Rätsep _et al._
(2011).

The validity of the analytically obtained result can be easily checked against
numerical simulations. In a first step, we compare the numerical results with
the analytical scaling of Eq. (1), for an equidistant vibrational spectrum, in
Fig. 3(c). A very good agreement is obtained showing that the time evolution
of the symmetric state is well reproduced by an exponential decay following
$e^{-\kappa_{\mathcal{S}}t}$. In the next step we pick a set of randomly drawn
vibrational frequencies $\nu_{1},...\nu_{n_{\text{max}}}$ in the window $0$ to
$4\Omega$ and Huang-Rhys factors $s_{1},...s_{n_{\text{max}}}$ in the interval
$[0,0.2]$ (with an average $s=0.1$). The results plotted in Fig. 4(a) show
that a good fit is obtained with the simplified result of Eq. (1) which
assumes evenly spaced vibrational frequencies and a Huang-Rhys factor at the
level of the distribution average $s=0.1$. Furthermore, an average of the
numerical results over $200$ random realizations predicts an excellent
agreement to the linear fit predicted by Eq. (1).
As aggregates are immersed in solvents and usually under room temperature
conditions, they are expected to present a large inhomogeneous broadening, at
the level of THz. Moreover, monomers are not placed in an exactly equidistant
chain resulting in further shifts in molecular energies. We consider a
distribution of the $\mathcal{N}$ monomer frequencies around $\omega_{0}$ such
that the frequency of each monomer becomes $\omega_{0}+\delta_{j}$ where
$\delta_{j}$ is randomly drawn from a distribution of width $\delta\omega$. In
the collective basis, the symmetric state couples to any asymmetric state and
acquires a shift as well Sommer _et al._ (2021). The coupling of $\mathcal{S}$
to a state $q$ mediated by disorder is simply given by the Fourier transform
of the distribution
$\delta_{q}=\frac{1}{\sqrt{\mathcal{N}}}\sum_{j=1}^{\mathcal{N}}{e^{iqjd}}\delta_{j},$
while the shift of the collective state
$\delta_{\mathcal{S}}=1/\sqrt{\mathcal{N}}\sum_{j=1}^{\mathcal{N}}\delta_{j}$
is the average of the distribution, thus close to zero. According to Ref.
Sommer _et al._ (2021), disorder induced couplings introduce an additional
loss channel, thus slightly increasing the Kasha rate. This is indeed
consistent with numerical simulations shown in Fig. 4(f) where the dynamics of
the symmetric state without disorder and with considerable disorder at
$\delta\omega=\Omega$ and $\delta\omega=2\Omega$ are compared. The upshot is
that the analytically derived loss rate holds well even for considerable
disorder levels.

Finally, let us provide some estimates of the expected timescale for such a
Kasha process in the aqueous dye molecule Cresyl Violet, Rhodamine $800$ and
Bacteriochlorophyll a (BChla). The results are listed in Table. 1 with
parameters estimated from Refs. Jafari _et al._ (2007); Barclay _et al._
(2022); Christensson _et al._ (2010); Rätsep _et al._ (2011). The upshot is
that analytically estimated values for the Kasha timescale run in the tens of
femtoseconds.

Conclusions and outlook \- In conclusion, we have investigated the time
dynamics following the photoexcitation of molecular aggregates in the side-by-
side configuration. These are a perfect showcase of cooperative phenomena as
their photophysics is naturally characterized by coherent and incoherent
effects brought on by the positioning of individual monomers in the near field
of each other. Effects widely explored in quantum optics, such as Dicke
superradiance and subradiance, naturally occur in the theoretical description
of such compounds, albeit in the presence of more complex, additional
interactions between electrons and a vast number of molecular vibrations. We
have provided a theoretical approach to the dynamics of collective electronic
states making the connection, at the analytical and numerical level, with the
physical mechanism introduced long ago by Kasha, stating that photon emission
(fluorescence or phosphorescence) occurs in appreciable yield only from the
lowest excited state of a given multiplicity. Our analytical conclusions
predict that the Kasha loss rate from the symmetric, high energy, optically
addressable state to the lower energy states of an H-aggregate is roughly
independent of the number of monomers but strongly dependent on the number of
low energy vibrational modes which can be excited and then dissipate the
accumulated energy afterwards. Further investigations will focus on aspects
such as the role of quantum coherence in such systems. An important direction
is the application of the methods presented in this manuscript to
photosynthetic systems under various conditions of illumination, ranging from
spatially and time coherent laser light to spatially and time incoherent light
sources.

##  I Methods

Numerical simulations are performed in the single excitation subspace. We
follow the time evolution of the system assuming unit population of the
symmetric state at the initial time $t=0$. As $\Gamma_{m}$ is larger than the
coherent rates, any excitation of a vibration is followed by quick relaxation,
justifying the assumption that any double excitation can be neglected. The
basis set is picked as a tensor product $\ket{j}\otimes\ket{j^{\prime}}^{(m)}$
where by definition $\ket{j}=\ket{g,g,...e_{j},...}$ \- only emitter $j$
excited electronically and
$\ket{j^{\prime}}=\ket{0,0,...1_{j^{\prime}},...}^{(m)}$ \- only mode $m$ in
emitter $j^{\prime}$ has one vibrational excitation. Instead of solving the
master equation directly, one can use the quantum jump formalism to evaluate
single stochastic quantum trajectories using the Monte Carlo wave function
method (MCWF). The advantage is that instead of describing the state of the
quantum system by a density matrix of size $\mathcal{N}^{4}\times n^{2}$ the
stochastic method only requires state vectors of size $\mathcal{N}^{2}\times
n$. This is somewhat counteracted by the stochastic nature of the formalism
which makes it necessary to repeat the simulation until the wanted accuracy is
reached. However, in most scenarios, especially for higher dimensional quantum
systems, the necessary number of repetitions is much smaller than the system
size $\mathcal{N}^{2}\times n$ and therefore using the MCWF method is
advantageous (see SI). For the simulation we write

$\displaystyle|\Psi\rangle=\sum_{j,j^{\prime}=1}^{\mathcal{N}}\sum_{m=1}^{n}\alpha_{jj^{\prime}}^{(m)}|g,g,...e_{j},...\rangle\otimes|0,0,...1_{j^{\prime}},...\rangle^{(m)},$
(6)

with coefficients $\alpha_{jj^{\prime}}^{(m)}$ and where the first part refers
to the electronic excitation of molecule $j$ and the second part to the
excitation of the $m$-th vibrational mode of molecule $j^{\prime}$.

##  II Acknowledgments

We acknowledge financial support from the Max Planck Society and the Deutsche
Forschungsgemeinschaft (DFG, German Research Foundation) – Project-ID
429529648 – TRR 306 QuCoLiMa (“Quantum Cooperativity of Light and Matter”). R.
H. acknowledges funding from the Austrian Science Fund (FWF) doctoral college
DK-ALM W1259-N27. We acknowledge fruitful discussions with Michael Reitz and
Johannes Feist. The numerical simulations were performed with the open-source
framework QuantumOptics.jl Krämer _et al._ (2018).

## References

  * Hestand and Spano (2018) N. J. Hestand and F. C. Spano, “Expanded theory of H- and J-molecular aggregates: The effects of vibronic coupling and intermolecular charge transfer,” Chemical Reviews 118, 7069–7163 (2018).
  * Saikin _et al._ (2013) S. K. Saikin, A. Eisfeld, S. Valleau, and A. Aspuru-Guzik, “Photonics meets excitonics: natural and artificial molecular aggregates,” Nanophotonics 2, 21–38 (2013).
  * Ma _et al._ (2021) S. Ma, S. Du, G. Pan, S. Dai, B. Xu, and W. Tian, “Organic molecular aggregates: From aggregation structure to emission property,” Aggregate 2, e96 (2021).
  * Abramavicius _et al._ (2009) D. Abramavicius, B. Palmieri, D. V. Voronine, F. Sanda, and S. Mukamel, “Coherent multidimensional optical spectroscopy of excitons in molecular aggregates; quasiparticle versus supermolecule perspectives,” Chem. Rev. 109, 2350 (2009).
  * Scheibe _et al._ (1937) G. Scheibe, L. Kandler, and H. Ecker, “Polymerisation und polymere Adsorption als Ursache neuartiger Absorptionsbanden von organischen Farbstoffen,” Naturwissenschaften 25, 75 (1937).
  * Jelley (1936) E. E. Jelley, “Spectral absorption and fluorescence of dyes in the molecular state,” Nature 138, 1009 (1936).
  * Kasha (1963) M. Kasha, “Energy transfer mechanisms and the molecular exciton model for molecular aggregates(1, 2). 1963,” Radiat. Res. 178, AV27–34 (1963).
  * Garcia-Vidal _et al._ (2012) F. J. Garcia-Vidal, C. Ciuti, and T. W. Ebbesen, “Manipulating matter by strong coupling to vacuum fields,” Science 373, 6551 (2012).
  * Berman (1994) P. Berman, _Cavity Quantum Electrodynamics_ , Advances in atomic, molecular, and optical physics (Academic Press, 1994).
  * Haroche and Kleppner (1989) S. Haroche and D. Kleppner, “Cavity Quantum Electrodynamics,” Phys. Today 42, 24 (1989).
  * Walther _et al._ (2006) H. Walther, B. T. Varcoe, B. Englert, and T. Becker, “Cavity Quantum Electrodynamics,” Rep. Prog. Phys. 69, 1325 (2006).
  * Haroche and Raimond (2013) S. Haroche and J.-M. Raimond, _Exploring the Quantum: Atoms, Cavities, and Photons_ (Oxford University Press, 2013).
  * Reitz _et al._ (2019) M. Reitz, C. Sommer, and C. Genes, “Langevin approach to quantum optics with molecules,” Phys. Rev. Lett. 122, 203602 (2019).
  * Holzinger _et al._ (2022) R. Holzinger, S. A. Oh, M. Reitz, H. Ritsch, and C. Genes, “Cooperative subwavelength molecular quantum emitter arrays,” Phys. Rev. Res. 4, 033116 (2022).
  * Jafari _et al._ (2007) A. Jafari, A. Ghanadzadeh, H. Tajalli, M. Yeganeh, and M. Moghadam, “Electronic absorption spectra of cresyl violet acetate in anisotropic and isotropic solvents,” Spectrochimica Acta Part A: Molecular and Biomolecular Spectroscopy 66, 717–725 (2007).
  * Barclay _et al._ (2022) M. S. Barclay, J. S. Huff, R. D. Pensack, P. H. Davis, W. B. Knowlton, B. Yurke, J. C. Dean, P. C. Arpin, and D. B. Turner, “Characterizing mode anharmonicity and huang–rhys factors using models of femtosecond coherence spectra,” The Journal of Physical Chemistry Letters 13, 5413–5423 (2022).
  * Christensson _et al._ (2010) N. Christensson, B. Dietzek, A. Yartsev, and T. Pullerits, “Electronic photon echo spectroscopy and vibrations,” Vibrational Spectroscopy 53, 2–5 (2010).
  * Rätsep _et al._ (2011) M. Rätsep, Z.-L. Cai, J. R. Reimers, and A. Freiberg, “Demonstration and interpretation of significant asymmetry in the low-resolution and high-resolution qy fluorescence and absorption spectra of bacteriochlorophyll a,” The Journal of Chemical Physics 134 (2011).
  * Reitz _et al._ (2022) M. Reitz, C. Sommer, and C. Genes, “Cooperative quantum phenomena in light-matter platforms,” PRX Quantum 3, 010201 (2022).
  * Asenjo-Garcia _et al._ (2017) A. Asenjo-Garcia, M. Moreno-Cardoner, A. Albrecht, H. J. Kimble, and D. E. Chang, “Exponential improvement in photon storage fidelities using subradiance and “selective radiance” in atomic arrays,” Phys. Rev. X 7, 031024 (2017).
  * Sommer _et al._ (2021) C. Sommer, M. Reitz, F. Mineo, and C. Genes, “Molecular polaritonics in dense mesoscopic disordered ensembles,” Phys. Rev. Res. 3, 033141 (2021).
  * Krämer _et al._ (2018) S. Krämer, D. Plankensteiner, L. Ostermann, and H. Ritsch, “Quantumoptics.jl: A Julia framework for simulating open quantum systems,” Computer Physics Communications 227, 109 (2018).
  * Hu _et al._ (2015) D. Hu, S.-Y. Huang, J.-Q. Liao, L. Tian, and H.-S. Goan, “Quantum coherence in ultrastrong optomechanics,” Phys. Rev. A 91, 013812 (2015).

##  Appendix A Vibronic coupling

Let us justify the form of the Holstein Hamiltonian in Eq. (2) by following a
first-principle derivation for a single nuclear coordinate $R$ of effective
mass $\mu$. We assume that, along the nuclear coordinate, the equilibria for
ground (coordinate $R_{g}$ , state vector $|g\rangle$) and excited (coordinate
$R_{e}$ and state vector $|e\rangle$) electronic orbitals are different.
Assuming equilibrium positions $R_{g}$ and $R_{e}$ for the potential surfaces
of electronic ground and excited states, one can write the total molecular
Hamiltonian describing both electronic and vibrational dynamics as

$\displaystyle\mathcal{H}_{\text{mol}}=\left[\omega_{0}+\frac{\hat{P}^{2}}{2\mu}+\frac{1}{2}\mu\nu^{2}\left(\hat{R}-R_{e}\right)^{2}\right]\sigma^{\dagger}\sigma+\left[\frac{\hat{P}^{2}}{2\mu}+\frac{1}{2}\mu\nu^{2}\left(\hat{R}-R_{g}\right)^{2}\right]\sigma\sigma^{\dagger},$
(7)

where $\mu$ is the reduced mass of the vibrational mode and
$\sigma=|g\rangle\langle e|$. The kinetic and potential energies are written
in terms of the position $\hat{Q}$ and momentum operator $\hat{P}$ describing
the nuclear coordinate under consideration, with commutation
$[\hat{Q},\hat{P}]=i$. Introducing oscillations around the equilibria
$\hat{Q}=\hat{R}-R_{\text{g}}$ and subsequently
$\hat{R}-R_{\text{e}}=\hat{Q}+R_{\text{g}}-R_{\text{e}}=:\hat{Q}-R_{\text{ge}}$
we obtain

$\mathcal{H}_{\text{mol}}=\frac{\hat{P}^{2}}{2\mu}+\frac{1}{2}\mu\nu^{2}\hat{Q}^{2}+\omega_{0}\sigma^{\dagger}\sigma-\mu\nu^{2}\hat{Q}R_{\text{ge}}\sigma^{\dagger}\sigma+\frac{1}{2}\mu\nu^{2}R_{\text{ge}}^{2}\sigma^{\dagger}\sigma.$
(8)

We can now rewrite the momentum and position operators in terms of bosonic
operators $\hat{Q}=q_{\text{zpm}}(b^{\dagger}+b)$,
$\hat{P}=ip_{\text{zpm}}(b^{\dagger}-b)$. The bosonic operators satisfy the
usual commutation relation $[{b},{b^{\dagger}}]=1$ and the zero-point motion
displacement and momentum are defined as $q_{\text{zpm}}=1/\sqrt{2\mu\nu}$ and
$p_{\text{zpm}}=\sqrt{\mu\nu/2}$. Reexpressing the terms above yields the
Holstein Hamiltonian

$\displaystyle\mathcal{H}_{\text{mol}}=(\omega_{0}+s\nu){\sigma}^{\dagger}{\sigma}+\nu{b}^{\dagger}{b}-\sqrt{s}\nu({b}^{\dagger}+{b}){\sigma}^{\dagger}{\sigma}.$
(9)

The dimensionless vibronic coupling strength $s$ is given by $\sqrt{s}=\mu\nu
R_{\text{ge}}q_{\text{zpm}}$ ($s$ is called the Huang-Rhys factor and is
typically on the order of $\sim 0.01-1$).

##  Appendix B Radiative and vibrational loss

In a master equation formulation for the system density operator $\rho$
written as $\partial_{t}\rho=i[\rho,\mathcal{H}]+\mathcal{L}[\rho]$ loss can
be included via the Lindblad superoperator
$\mathcal{L}_{\gamma}[\rho]=\gamma_{\mathcal{O}}/2\left[2\mathcal{O}\rho(t)\mathcal{O}^{\dagger}-\mathcal{O}^{\dagger}\mathcal{O}\rho(t)-\rho(t)\mathcal{O}^{\dagger}\mathcal{O}\right],$
describing decay at generic rate $\gamma_{\mathcal{O}}$ through a single
channel with a generic collapse operator $\mathcal{O}$. For vibrational loss,
the collapse rate for each mode $m$ is $\Gamma_{m}$ and the corresponding
collapse operator is $b_{jm}-\sqrt{s_{m}}\sigma_{j}^{\dagger}\sigma_{j}$. This
form for the collapse operator is derived in analogy to the dissipative
physics of optomechanical systems in the ultrastrong coupling regime Hu _et
al._ (2015). The radiative loss is not in diagonal Lindblad form but achieves
the following expression
$\mathcal{L}_{e}[\rho]=\sum_{j,j^{\prime}}\gamma_{jj^{\prime}}/2\left[2\sigma_{j}\rho\sigma^{\dagger}_{j^{\prime}}-\sigma^{\dagger}_{j}\sigma_{j^{\prime}}\rho-\rho\sigma^{\dagger}_{j}\sigma_{j^{\prime}}\right]$.
This form can be diagonalized and it shows the emergence of $\mathcal{N}$
independent decay channels, each corresponding to some collective electronic
superposition state Reitz _et al._ (2022). At very small separation, deep into
the subwavelength regime, the fully symmetric superposition decays at a
superradiant rate roughly equal to $\mathcal{N}\gamma_{0}$ while all other
states have vanishingly small decay rates (which we will assume in the
following to be exactly zero). This is by no means a limitation of our
treatment as one can easily generalize this to the case of non-zero decay
rates of the dark manifold Holzinger _et al._ (2022); Reitz _et al._ (2022).

##  Appendix C Hamiltonian in the collective basis

The Hamiltonian coupling the symmetric state to the dark state manifold is
then given by

$\displaystyle\mathcal{H}_{\text{int}}^{\mathcal{S}\mathcal{A}}=-\sum_{m=1}^{n}\sum_{q\neq
0}\frac{\sqrt{s_{m}}\nu_{m}}{\sqrt{\mathcal{N}}}\Big{(}Q_{q}^{(m)}\mathcal{S}^{\dagger}\mathcal{A}_{q}+\mathrm{h.c.}\Big{)},$
(10)

via collective vibrations. This coupling is responsible for funneling
population into the long lived dark state manifold after the initial driving
of the fully symmetric state under uniform illumination. This mechanism is
fundamental to understand the dynamics associated with Kasha’s rule. In
addition, within the dark state manifold an all-to-all coupling Hamiltonian
acts with the following form

$\displaystyle\mathcal{H}_{\text{int}}^{\mathcal{A}\mathcal{A}}=-\sum_{m=1}^{n}\sum_{q\neq
q^{\prime}}\frac{\sqrt{s_{m}}\nu_{m}}{\sqrt{\mathcal{N}}}\Big{(}Q_{q-q^{\prime}}^{(m)}\mathcal{A}_{q}^{\dagger}\mathcal{A}_{q^{\prime}}+\mathrm{h.c.}\Big{)},$
(11)

and the sum implies that $q,q^{\prime}\neq 0$. This indicates that a
redistribution of energy takes place within the whole manifold of dark states.
After transforming the system Hamiltonian, the energies of the collective
states are shifted by the contribution of the symmetric vibrational mode
$-\sum_{m}\sqrt{s_{m}}\nu_{m}Q^{(m)}_{0}/\sqrt{\mathcal{N}}$. The energy
shifts can be removed by the collective polaron transformation
$U=\prod_{q}\prod_{m=1}^{n}e^{i\sqrt{s_{m}}/\sqrt{\mathcal{N}}P_{0}^{(m)}\mathcal{A}_{q}^{\dagger}\mathcal{A}_{q}}$
which leads to a renormalization of the collective state energies as
$\bar{\omega}_{q}=\omega_{0}+\sum_{m}s_{m}\nu_{m}/2+\Omega_{q}$.

##  Appendix D Vacuum mediated dipole-dipole coupling rates

The vacuum mediated dipole-dipole interactions for an electronic transition at
wavelength $\lambda_{0}$ (corresponding wave vector $k_{0}=2\pi/\lambda_{0}$)
between identical pairs of emitters separated by $r_{ij}$ is given in terms of
the free-space electromagnetic Green’s tensor
$\boldsymbol{G}(\vec{r}_{i}-\vec{r}_{j},\omega_{0})\equiv\boldsymbol{G}(\vec{r}_{ij},\omega_{0})$,
with $\vec{r}_{ij}=\vec{r}_{i}-\vec{r}_{j}$, which reads

$\boldsymbol{G}(\vec{r},\omega_{0})=\frac{e^{ik_{0}r}}{4\pi
k_{0}^{2}r^{3}}\Big{[}(k_{0}^{2}r^{2}+ik_{0}r-1)\mathbb{I}+(-k_{0}^{2}r^{2}-3ik_{0}r+3)\frac{\vec{r}\otimes\vec{r}}{r^{2}}\Big{]},$
(12)

where $r=|\vec{r}|$. The Green’s function $\boldsymbol{G}_{\alpha\beta}$ is a
tensor quantity, with $\\{\alpha,\beta\\}=\\{x,y,z\\}$ which is determined by
the polarization direction of the dipoles. In order to obtain the dipole-
dipole couplings for H-aggregates we chose linear polarization in z-direction
for all molecules, namely we take
$\vec{\mu}_{z}^{*}\cdot\boldsymbol{G}\cdot\vec{\mu}_{z}$.

##  Appendix E Deriving rate equations

Starting from the Holstein Hamiltonian for $\mathcal{N}$ identical molecules
with $n$ vibrational modes each, the Heisenberg equations for the collective
electronic modes are given by

$\displaystyle\dot{\mathcal{S}}$
$\displaystyle=-i\Big{(}\Omega_{\mathcal{S}}-\frac{\gamma_{\mathcal{S}}}{2}\Big{)}\mathcal{S}+\frac{i\sqrt{s}_{m}\nu_{m}}{\sqrt{\mathcal{N}}}\sum_{m=1}^{n}\sum_{q}Q_{q}^{(m)}\mathcal{A}_{q}+\mathrm{noise},$
(13a) $\displaystyle\dot{\mathcal{A}_{q}}$
$\displaystyle=-i\Omega_{q}\mathcal{A}_{q}+\sum_{m=1}^{n}\frac{i\sqrt{s_{m}}\nu_{m}}{\sqrt{\mathcal{N}}}\Big{(}{Q^{(m)}_{q}}^{\dagger}\mathcal{S}+\sum_{q^{\prime}\neq
q}Q^{(m)}_{q-q^{\prime}}\mathcal{A}_{q^{\prime}}\Big{)}+\mathrm{noise}.$ (13b)

The collective noise terms will be neglected from now on as they do not
contribute to the transfer process.

To calculate the transfer rate from the symmetric state to the antisymmetric
states we assume some initial population in the symmetric state and no
population in the antisymmetric states, additionally we assume that the
symmetric state decays independently and formally integrate

$\displaystyle\mathcal{S}(t)$
$\displaystyle=S(0)e^{-(i\Omega_{\mathcal{S}}+\gamma_{\mathcal{S}}/2)t},$
(14a) $\displaystyle\mathcal{A}_{q}(t)$
$\displaystyle=A_{q}(0)e^{-i\Omega_{q}t}+\sum_{m=1}^{n}\frac{i\sqrt{s_{m}}\nu_{m}}{\sqrt{\mathcal{N}}}\int^{t}_{0}dt^{\prime}e^{-i\Omega_{q}(t-t^{\prime})}\Big{(}Q_{q}^{(m)}(t^{\prime})\mathcal{S}(t^{\prime})+\sum_{q^{\prime}\neq
q}Q^{(m)}_{q-q^{\prime}}(t^{\prime})\mathcal{A}_{q^{\prime}}(t^{\prime})\Big{)},$
(14b)

and for the expectation value of the populations we get

$\displaystyle\dot{\langle\mathcal{S}^{\dagger}\mathcal{S}\rangle}$
$\displaystyle=-\gamma_{\mathcal{S}}{\langle\mathcal{S}^{\dagger}\mathcal{S}\rangle}-\sum_{m=1}^{n}\sum_{q}\frac{2\sqrt{s_{m}}\nu_{m}}{\sqrt{\mathcal{N}}}\,\mathrm{Im}\langle\mathcal{S}^{\dagger}\mathcal{A}_{q}Q_{q}^{(m)}\rangle,$
(15a)
$\displaystyle\dot{\langle\mathcal{A}_{q}^{\dagger}\mathcal{A}_{q}\rangle}$
$\displaystyle=-\sum_{m=1}^{n}\frac{2\sqrt{s_{m}}\nu_{m}}{\sqrt{\mathcal{N}}}\,\mathrm{Im}\Big{(}\langle\mathcal{A}_{q}^{\dagger}\mathcal{S}Q_{q}^{(m)}\rangle+\sum_{q^{\prime}\neq
q}\langle\mathcal{A}_{q}^{\dagger}\mathcal{A}_{q^{\prime}}Q_{q^{\prime}-q}^{(m)}\rangle\Big{)}.$
(15b)

Therefore the terms
$-2\sqrt{s}_{m}\nu_{m}/\sqrt{\mathcal{N}}\,\mathrm{Im}\langle\mathcal{A}_{q}^{\dagger}\mathcal{S}Q_{q}^{(m)}\rangle$
will be responsible for population transfer from the symmetric to the
antisymmetric state with quasi-momentum $q$ at a rate
$\kappa^{(m)}_{\mathcal{S}\rightarrow{q}}$. We can calculate the rates
explicitly up to order $\mathcal{O}(s_{m}\nu_{m}^{2})$ and assuming that
correlations between vibronic and electronic operators factorize.

$\displaystyle-2\sqrt{s}_{m}\nu_{m}/\sqrt{\mathcal{N}}\langle\mathcal{A}_{q}^{\dagger}\mathcal{S}Q_{q}^{(m)}\rangle$
$\displaystyle=-is_{m}\nu_{m}^{2}\int_{0}^{t}dt^{\prime}e^{-\Omega_{q}(t-t^{\prime})}\langle
Q_{q}^{(m)}(t^{\prime})Q_{q}^{(m)}(t)\rangle\langle\mathcal{S}^{\dagger}(0)\mathcal{S}(0)\rangle
e^{-\epsilon_{\mathcal{S}}t^{\prime}}e^{-\epsilon^{*}_{\mathcal{S}}t}$
$\displaystyle=-is_{m}\nu_{m}^{2}\langle\mathcal{S}^{\dagger}(0)\mathcal{S}(0)\rangle\frac{e^{-\gamma_{\mathcal{S}}t}-e^{-((\Gamma_{m}+\gamma_{\mathcal{S}})/2+i(\Omega_{\mathcal{S}}-\Omega_{q}-\nu_{m}))t}}{(\Gamma_{m}+\gamma_{\mathcal{S}})/2+i(\Omega_{\mathcal{S}}-\Omega_{q}-\nu_{m})},$
(16)

where we defined
$\epsilon_{\mathcal{S}}=-(\gamma_{\mathcal{S}}/2-i\Omega_{\mathcal{S}})$ and
used the fact that different vibrational modes are uncorrelated at all times,
i.e. $\langle Q^{(m^{\prime})}_{q}(t^{\prime})Q^{(m)}_{q}(t)\rangle=0$ for
$m^{\prime}\neq m$. The correlations for $Q_{q}^{(m)}$ are evaluated assuming
free evolution of the vibrations (to lowest order) and zero temperature for
the vibrational modes:

$\langle
Q_{q}^{(m)}(t^{\prime})Q_{q}^{(m)}(t)\rangle=\frac{1}{\mathcal{N}}\sum_{j=1}^{\mathcal{N}}\langle
b_{jm}(t^{\prime})b_{jm}^{\dagger}(t)\rangle=e^{-(\Gamma_{m}/2-i\nu_{m})(t-t^{\prime})}.$
(17)

The transfer rate can be written as

$\kappa^{(m)}_{\mathcal{S}\rightarrow
q}=\frac{2s_{m}\nu_{m}^{2}(\Gamma_{m}+\gamma_{\mathcal{S}})/{\mathcal{N}}}{(\Gamma_{m}+\gamma_{\mathcal{S}})^{2}+4(\Omega_{\mathcal{S}}-\Omega_{q}-\nu_{m})^{2}},$
(18)

given fast vibrational relaxation rates $\Gamma_{m}\gg\gamma_{\mathcal{S}}$
compared to the electronic decay rates.

##  Appendix F Single excitation subspace

The numerical diagonalization and subsequent time dynamics are evaluated in
the single-excitation sector for both the electronic and vibrational degrees
of freedom. This allows to rewrite the effective Hamiltonian in non-hermitian
form as ($\hbar=1$)

$\displaystyle\mathcal{H}_{\mathrm{eff}}=\sum_{j=1}^{\mathcal{N}}\Big{(}h^{(j)}+\sum_{j^{\prime}=1}^{\mathcal{N}}\Big{(}\Omega_{jj^{\prime}}-i\frac{\gamma_{jj^{\prime}}}{2}\Big{)}\sigma^{\dagger}_{j}\sigma_{j^{\prime}}-\frac{i}{2}\sum_{m=1}^{n}\Gamma_{m}\mathcal{O}^{\dagger}_{jm}\mathcal{O}_{jm}\Big{)},$
(19)

where $h^{(j)}$ is defined in Eq. (2) and
$\mathcal{O}_{jm}=b_{jm}-\sqrt{s_{m}}{\sigma}^{\dagger}_{j}{\sigma}_{j}$. The
dynamics of the electron-vibron density matrix $\rho$ can be described by a
von Neumann equation of the form

$\displaystyle
i\frac{\mathrm{d}}{\mathrm{dt}}\rho(t)=[\mathcal{H}_{\mathrm{eff}}\rho-\rho\mathcal{H}_{\mathrm{eff}}],$
(20)

and the expectation value of observable $\mathcal{O}$ becomes
$\bar{\mathcal{O}}=\mathrm{tr}(\rho{\mathcal{O}})$. However, instead of
solving the von Neumann equation directly, one can use the quantum jump
formalism to evaluate single stochastic quantum trajectories using the Monte
Carlo wave function method (MCWF). For large numbers of trajectories, the
statistical average then approximates the result of the Master equation. The
huge advantage is that instead of describing the state of the quantum system
by a density matrix of size $\mathcal{N}^{4}\times n^{2}$ these trajectories
work in terms of state vectors of size $\mathcal{N}^{2}\times n$. This is
somewhat counteracted by the stochastic nature of the formalism which makes it
necessary to repeat the simulation until the wanted accuracy is reached. It
turns out, however, that for many cases, especially for high dimensional
quantum systems, the necessary number of repetitions is much smaller than the
system size $\mathcal{N}^{2}\times n$ and therefore using the MCWF method is
advantageous.

The system size stems from the fact, that in the single excitation subspace
for both electronic and vibronic modes a general state vector can be written
as

$\displaystyle|\Psi\rangle=\sum_{j=1}^{\mathcal{N}}\alpha_{jj^{\prime}}^{(m)}|g,g,...e_{j},...\rangle\otimes\sum_{m=1}^{n}\sum_{j^{\prime}=1}^{\mathcal{N}}|0,0,...1_{j^{\prime}},...\rangle^{(m)},$
(21)

with coefficients $\alpha_{jj^{\prime}}^{(m)}$ and where the first part refers
to the electronic excitation of molecule $j$ and the second part to the
excitation of the $m$-th vibrational mode of molecule $j^{\prime}$. Thus, the
single excitation assumption substantially reduces the Hilbert space dimension
from $2^{\mathcal{N}}\times n_{\mathrm{cut}}^{\mathcal{N}}\times n$ to
$\mathcal{N}^{2}\times n$ (where $n_{\mathrm{cut}}$ is the cut-off of the Fock
space dimension for the vibrational modes), allowing the simulation of
mesoscopic numbers of molecules.


# Recognising Known Configurations of Garments For Dual-Arm Robotic Flattening

Li Duan  School of Computing Science
University of Glasgow
Glasgow, United Kingdom
<EMAIL_ADDRESS>Gerardo Aragon-Camarasa  School of Computing
Science
University of Glasgow
Glasgow, United Kingdom
<EMAIL_ADDRESS>School of Computing Science
University of Glasgow
Glasgow, United Kingdom
<EMAIL_ADDRESS>School of Computing Science
University of Glasgow
Glasgow, United Kingdom
<EMAIL_ADDRESS>
###### Abstract

Robotic deformable-object manipulation is a challenge in the robotic industry
because deformable objects have complicated and various object states.
Predicting those object states and updating manipulation planning is time-
consuming and computationally expensive. In this paper, we propose learning
known configurations of garments to allow a robot to recognise garment states
and choose a pre-designed manipulation plan for garment flattening.

###### Index Terms:

pattern recognition; robot vision; depth images; deformable objects

Recognising Known Configurations of Garments For Dual-Arm Robotic Flattening

Li Duan Gerardo Aragon-Camarasa

##  I Introduction

Robots have been widely applied in many aspects of life, for example in
robotic assistance for humans, life-threatening rescue, and manufacturing.
Among these applications, robotic deformable object manipulation plays a
critical part. Autonomous garment laundry and sorting, soft-object
manipulation in life-threatening rescue and manufacturing, and robotic
threading in the textile industry require an efficient and robust robotic
deformable-object manipulation strategy.

However, robotic manipulation of deformable objects remains an open problem in
robotics because deformable objects can take unpredictable object states
(crumpled, stretched, or bent configurations). That is, their states
frequently vary when they are being manipulated and can not be predicted
precisely and timely. Deformable objects are easily slipped from robots and
take a long time to be manipulated because of their properties. Therefore,
current research efforts focus on reducing object-state or robotic action-
state complexity.

Current approaches that reduce deformable object complexity for robotic
manipulation can be divided into two categories: those that propose to find
simplified models to represent object states or action states first [1, 2, 3,
4, 5, 6], and those that rely on model-free reinforcement learning solutions
where the agent learns to model the object while interacting with it [7, 8,
9]. When deformable object models are simplified, it is possible to reduce the
computational costs, making it possible to predict the object states of
manipulated deformable objects close to real-time. An efficient prediction of
future object states enables a robot to update its manipulation plan on-the-
fly, which is critical to manipulating deformable objects in a real-time
setting. However, current research is geared towards simple geometries such as
ropes, towels or cube-shaped objects. When a reinforcement learning solution
is used, data-driven reinforcement learning agent learns manipulation policies
using a reward-based mechanism. Robots are rewarded if their manipulation
policies produce actions that contribute to a step closer to the object’s goal
state. However, most of the reinforcement learning approaches rely on
simulated environments of which the physics of deformable objects are
simplified.

In this paper, we propose taking advantage of gravity to model and learn the
known configurations of deformable objects (garments in this paper), which can
be later used to choose a pre-designed manipulation plan based on the
recognised known configurations to flatten these objects. Our approach does
not need to find simplified models that provide good representations of object
states for deformable objects. Similarly, our approach does not require manual
action labelling compared to reinforcement-learning approaches. In this paper,
we focus on recognising the known configurations of garments and develop a
pipeline to flatten garments, see Figure 1.

In summary, our main contributions in this paper are two-fold: we propose
learning the known configurations of garments, which is the foundation for a
full pipeline of an effective robotic garment flattening pipeline. Compared
with state-of-art, our approach has a higher accuracy (89% versus 73% in [6,
5]) and we captured a known configurations database which comprises RGB and
depth images of garments

We introduce the related work in section II, and our methodology in section
III. Our experiment results are in section IV and a discussion in section V.
Finally, we conclude this paper and suggest future work in section VI.

##  II Related Work

Research on robotic deformable-object manipulation is divided into two
categories: model-based approaches and model-free approaches. Model-based
approaches focus on engineering methods, simplifying models representing
object and action states or applying dynamics modelling rather than state-
space modelling. On the other hand, model-free approaches concentrate on
devising data-driven methods mainly relating to reinforcement or imitation
learning. Robots learn to manipulate deformable objects by applying
manipulation policies, usually trained on simulated environments.

The main two approaches of models-based approaches are object-state modelling
and action-state modelling. Lin et al. [10] have proposed to use finite
element methods (FEM) to model the dynamics of bottles, where they defined a
‘liftability ’ test. A bottle can only be lifted if its deformation passes the
‘liftability’ test, so the bottle can be stably manipulated without slipping
from the robot’s gripper. A mass-spring system has been applied in [4] to
grasp deformable objects using a multi-fingered robotic hand. The Nvidia PhysX
simulator has also been utilised in [11] to teach a robot for assisted human
dressing. Miller et al. [12] have proposed a robotic laundry folding approach
by introducing a quasi-static cloth model that represents states of garments
and an algorithm that outputs a motion planning based on a 2D cloth polygon
generated from the quasi-static cloth model and a desired sequence of folds.
Yamakawa et al. [13] have also proposed learning a simplified model to
represent the state space of ropes and derive motion planning from the inverse
dynamics of rope states. McConachie et al.[14] have suggested learning reduced
state spaces by using a classifier to bias a planner away from state-action
pairs that are not feasible under general robotic working environments.

The problem with model-based approaches is that due to the high dimensionality
of the deformable object and action spaces, only objects with simple
geometries (such as towels, ropes, and sponges) can be modelled and
manipulated. Objects with complex geometries such as garments can not be used
in these approaches because modelling and manipulation are time-consuming and
computationally costly. An alternative to simplify model representations is to
convert unknown object states to known object states. Li et al. [6] [5] have
proposed to estimate the pose of hanging garments and apply manipulation plans
based on the recognised garments’ poses. Poses only relate to grasping points
rather than initial object states of garments, effectively converting garments
from unknown object states to known object states. But their approach requires
re-grasping and rotating garments to recognise the poses, adding a time
overhead for executing a garment manipulation task.

Data-driven approaches have been recently proposed, namely imitation and
reinforcement learning. That is, Pignat et al. [7] have proposed to encode
sensory information and motor commands as a joint distribution in a hidden
semi-Markov model, of which parameters are learned from a set of human
demonstrations. Each model set represents a sensorimotor pattern whose
sequencing can produce complex behaviours. Balaguer et al. [9] have presented
a combined approach of imitation and reinforcement learning, where human
demonstrations are used to reduce the search space of the reinforcement
learning agent and allow the solution to converge quickly. Tsumine et al. [15]
have proposed a deep reinforcement learning for robotic cloth manipulation
with a smooth policy update. Their networks combine the nature of smooth
policy updates in value-function-based reinforcement learning with automatic
feature extraction from high-dimensional observations in deep neural networks
to enhance the sample efficiency and learning stability with fewer samples.
Twardon et al. propose restricting robotic action-space modelling to the
actions that are safe for robotic hands and arms, garments and the head when a
robot is trained to put a knit cap on a Styrofoam head. A direct policy-search
algorithm finds appropriate trajectories in the restricted action space in
order to enable the robotic knit-cap operation. Their approach has
demonstrated an effective manipulation planning under a restricted action
space. McConachie et al. [16] have proposed formulating a deformable object
manipulation task as a multi-armed bandit problem, with each arm representing
a model of the deformable object. Matas et al. [17] use a tailored version of
the deep deterministic policy gradients (DDPG [18]) for deformable object
manipulation, where they have learned a policy from simulation data while
testing the policy on real data. The limitation of imitation and reinforcement
learning approaches is that they need several training iterations to learn
tasks. Matas et al. [17] trained a robot to fold a towel and drape a piece of
cloth through a hanger by reinforcement learning. The robot is trained in a
simulated environment for around 250 iterations. Training such a robot in real
environments is difficult because some actions need to be validated in a
simulated environment. Hazardous actions (actions that may harm robots or
human beings) are easy to handle in simulated environments but need to be
avoided in real environments.

##  III Methodology

As described in section I, state-of-art approaches conduct manipulations on
tables or platforms, where deformable objects have complicated and various
configurations (or object states) and result in computationally costly updates
to find a manipulation plan. In contrast, we consider taking advantage of
gravity to control the variety of object configurations and reduce the
complexity of deformable object configurations.

In this paper, we start by assuming that garments of the same categories (e.g.
jeans, towels, tshirts, etc.) lying on a table have different configurations.
If the garments are grasped from similar grasping points, they will have
similar configurations when the robot picks them up, which we call known
configurations in this paper. For example two towels have different crumpled
configurations on a table (i.e. starting configurations), but they have
similar known configurations after a robot grasps them. These known
configurations only depend on the grasping points because of gravity. That is,
we convert complex starting configurations to simple known configurations from
which the robot can follow pre-defined manipulation plans to flatten them.
Recognising the known configurations of garments is crucial for the robotic
garment flattening pipeline shown in Figure 1.

After the garment’s known configuration is obtained after grasping them from
the table, we find the second grasping points based on how gravity effect on
the garment. For example, as shown in Figure 1 step 1 in the manipulation
routine box, the second grasping point is the lowest point of the towel.
Different known configurations have different locations of the second grasping
points, so recognising the garment’s known configurations are critical to
localise the second grasping point. After the robot finds the second grasping
point, the next step (step two in Figure 1) is to find the third grasping
point, which is the opposing ending corner of the towel. Then, the robot
stretches the towel from the two grasping points in step three in Figure 1.
Because of the stretching and gravity, the towel is in flattened state;
therefore, the final step consists of placing the towel on the table by
sliding it on the edge of the table (shown in steps four and five in Figure
1).

Figure 1: Known Configurations Pipeline: a towel with a crumpled configuration
(starting configuration) is grasped by a robot from a grasping point (the
first grasping point as described in Section III), giving rise to a known
configuration. The known configuration is recognised, and a pre-designed five-
step manipulation plan (described in section III) can be chosen based on the
recognised known configuration. The robot manipulates the towel with the
chosen manipulation plan.

To recognise garment’s known configurations, we propose to use an off-the-
shelf deep neural network (named KCNet in this paper) based on ResNet-18[19].
KCNet is mathematically expressed as: $O=F(C(I))$, where $O$ is the output of
the KCNet, which is the recognised known configurations (represented as
classification class), $I$ is the input image that captures known
configurations, $C$ are convolutional layers (a ResNet-18[19]), and $F$ are
fully connected layers. We use a negative log-likelihood loss (NLLLoss) to
train this network:
$L(\theta)=\sum_{i=1}^{n}(y_{i}\log{\hat{y}_{\theta,i}}+(1-y_{i})\log(1-\hat{y}_{\theta,i}))$,
where $\theta$ is the weight parameters of the network, $y_{i}$ is the ground
truth probability that the $i$th data point is positive, and
$\hat{y}_{\theta,i}$ is the predicted probability that the $i$th data point is
positive. Figure 1 shows how KCNet works in the known configurations pipeline.

As stated in section III, a known configuration highly depends on the grasping
points, therefore we label known configurations in terms to its corresponding
grasp points. Known configurations should be recognisable by KCNet, so they
should be distinguishable from each other. Therefore, a known configuration
represents all grasping points in a segmented area of a garment, rather than a
single grasping point. We proposed to discritise garments into different
segments and define a grasping point to represent one segment. Our experiments
found that dividing a garment into ten segments can provide the most known
configurations that are recognisable. In each segment, we locate its
corresponding grasping point as the centre of the area.

State-of-art focuses on real-time policy updating in robotic deformable object
manipulation. A robot will update its manipulation plan after observing a new
object state. The updates are computationally intensive and time-consuming. In
contrast, we propose designing manipulation plans in advance to avoid updating
the manipulation plan to speed up the task execution.

A manipulation plan comprises sequences of 6D poses and the the robot’s
gripper state (i.e. open or close). As described in section III, our
manipulation strategy include finding the second and third grasping points,
stretching to flatten garments and lifting garments down to the table. Our
manipulation plan design is thus based on these strategy and ensures that each
step requires the least actions. We defined these sequences in CSV files to
transfer manipulation commands to the robot after a known configuration is
recognised.

##  IV Experiments And Results

###  IV-A Experimental Methodology

Our known configurations pipeline consists of three stages: the robot
recognises the known configuration of a grasped garment, based on the
recognised known configuration, the robot chooses a pre-defined manipulation
plan, and then the robot flattens the garment.

To recognise garments’ known configurations, we captured a novel dataset
comprising depth and RGB images of garments from five categories: jeans,
shirts, sweaters, towels and tshirts. There are four garment instances in each
category, and as described in section III, each garment instance has ten
grasping segments represented by ten grasping positions. We collected 100
images for each position, where each image captures the known configuration of
a garment grasped from a specific position. We captured a total of 19,269
depth and RGB images, respectively, with a resolution of $256\times 256$
pixels. Figure 2 shows examples of images in this dataset. For each garment
category, we pre-defined 10 manipulation plans for each the 10 segments as
described in section III. Therefore, we have 50 manipulation plans.

With this dataset, we trained a KCNet which is a classification network that
consists of a pre-trained ResNet 18 structure and a fully connected network.
The detailed implementation of KCNet can be found at
https://liduanatglasgow.github.io/known_configurations/. The learning rate is
set to $10^{-3}$ and we allow it to decay during the training with a step
scheduler (a decay factor of 0.1 and a step size of 8).

We have implemented a k-fold cross-validation approach (k-fold CVA) to train
and test KCNet, rather than the traditional approach of train-validate-test
splits. The $k$ value in our experiment is set to four, which means that our
database is split into four groups for the k-fold CVA training and testing
sessions. There are four garment instances in a group. Three groups are
assigned as training groups for each session and one group as a testing group.
We ensured that the garment instances in the testing group were ‘unseen’ by
KCNet. We iterate the testing group to include all garment instances in a
category, and we averaged the classification accuracies as the classification
accuracy for the network.

The robot in this experiment is a Baxter dual-arm robot with a table at the
front to place the garments. We use a computer with Ubuntu 16.04 and an NVIDIA
1080 Ti GPU to train the KCNet. The robot is controlled by the Robot Operating
System, and an Xtion camera captures the images. In section IV-B, we
demonstrate an example of a five-step manipulation plan based on the
recognised known configurations of a towel for flattening the towel. We will
improve pre-designed manipulation plans in our future research.

Figure 2: Database Examples: There are five garment categories in our
database: jeans, shirts, sweaters, towels and tshirts. Our database contains
masked depth and RGB images of four garment instances in each category, which
have five categories (Top: original images, Bottom: masked depth images).

###  IV-B KCNet Validation Results and Manipulation Demonstration

We compared the performance of KCNets trained on RGB, depth and RGBD images.
Table I shows our results for each category in each k-fold cross-validation
described in section IV-A. Table I shows classification accuracies that
represent the percentages of known configurations recognised by the KCNet.

We found that training a KCNet on depth images (89%) outperforms a KCNet on
RGB images (73%). Depth images capture structural and physical characteristics
of garments (as found in [20]), and therefore, enables a KCNet to recognise
the known configurations of unseen garment instances. The RGB images capture
texture characteristics that depend on the lighting condition and shadows
which hinders the ability of KCNet to predict known configurations. Meanwhile,
We can observe an increase from 73% to 78.5% for the KCNet trained on RGBD
images but still a significant gap behind the KCNet trained only on depth
images (89.0%). Depth information facilitates a KCNet to recognise the known
configurations of unseen garments, while RGB information may affect the
effectiveness of depth information.

Based on this ablation study, we chose the KCNet trained on the depth images
to conduct the manipulation experiments. Figure 3 shows an example of
manipulations in our experiments. As described in section III, the robot
firstly recognises the known configurations of the garments, then chooses a
manipulation plan based on the recognised known configurations, and finally
flattens the garments with the chosen manipulation plan. We can observe that
the towel shown in Figure 3 are correctly flattened from its crumpled
configuration. A robotic manipulation example can be found on our project
website https://liduanatglasgow.github.io/known_configurations/.

TABLE I: k-fold cross validation experiment results: comparison between depth, RGB and RGBD images (unit: %) Depth | 1 | 2 | 3 | 4 | RGB | 1 | 2 | 3 | 4 | RGBD | 1 | 2 | 3 | 4
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
towel | 94.4 | 96.4 | 93.1 | 86.2 | towel | 67.0 | 55.9 | 74.1 | 64.5 | towel | 71.8 | 79.7 | 85.2 | 76.7
tshirt | 86.8 | 87.2 | 96.3 | 94.7 | tshirt | 70.8 | 71.2 | 72.1 | 76.9 | tshirt | 74.6 | 68.3 | 84.4 | 83.9
shirt | 78.2 | 80.3 | 75.9 | 92.5 | shirt | 87.9 | 58.4 | 75.4 | 91.4 | shirt | 87.9 | 58.4 | 68.8 | 91.9
sweater | 78.4 | 85.4 | 87.0 | 86.2 | sweater | 54.6 | 42.0 | 68.1 | 85.2 | sweater | 52.7 | 51.6 | 74.7 | 77.0
jean | 99.3 | 95.8 | 95.3 | 99.1 | jean | 76.4 | 87.8 | 87.0 | 98.8 | jean | 87.5 | 92.0 | 97.2 | 99.5
average | 87.0 | 89.0 | 89.0 | 92.0 | average | 73.0 | 62.0 | 75.0 | 83.0 | average | 76.0 | 70.0 | 82.0 | 86.0
AVERAGE: 89.0 | AVERAGE: 73.0 | AVERAGE: 78.5
Figure 3: An Example of Manipulations: an example are demonstrated: firstly,
the robot grasps the towel from the table to recognise its known
configurations. Then a manipulation plan is chosen. Finally, the robot
flattens the towel with the chosen manipulation plan.

##  V Discussion

redIn our experiments, we proposed a robotic garment manipulation by
recognising garments’ known configurations and choosing a pre-designed
manipulation plans to flatten garments. This paper introduces the first part
of this pipeline: how to recognise the ‘known configurations’ of garments. The
known configurations can be recognised because we took advantage of gravity
and reduced deformable objects’ complexity. Various and complicated initial
configurations are converted into simple known configurations by taking
advantage of gravity after the robot has grasped the garment.

##  VI Conclusion And Future Work

We proposed an effective robotic garment flattening approach by recognising
the known configurations of garments and choosing a pre-designed manipulation
plan to flatten them. Our approach also features training KCNets on depth
images of real garments, resulting in higher recognition accuracy compared
with previous work [6].

In our experiments, however, a robot can only recognise the known
configurations of the garments with shapes from the five categories, while it
is unable to recognise those of the garments with an unknown category. Our
future work aims to devise a continual learning framework so that a robot
learns unknown shapes when it deals with those shapes. Additionally, we
believe that a robot can benefit from prior knowledge of shapes [21] or
physics properties [20] of garments to recognise known configurations, so we
propose learning that prior knowledge during the robot’s grasping garments and
improving recognition accuracy on known configurations.

##  VII Acknowledgement

We thank the members of the Computer Vision and Automation System group from
the University of Glasgow for discussing the experimental results and
providing experiment suggestions.

## References

  * [1] W. Yan, A. Vangipuram, P. Abbeel, and L. Pinto, “Learning predictive representations for deformable objects using contrastive estimation,” _arXiv preprint arXiv:2003.05436_ , 2020.
  * [2] Y. Li, Y. Yue, D. Xu, E. Grinspun, and P. K. Allen, “Folding deformable objects using predictive simulation and trajectory optimization,” in _2015 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_. IEEE, 2015, pp. 6000–6006.
  * [3] D. Agrawal, M. Choy, H. Va Leong, and A. K. Singh, “Maya: A simulation platform for distributed shared memories,” in _Proceedings of the Eighth Workshop on Parallel and Distributed Simulation_ , ser. PADS ’94. New York, NY, USA: Association for Computing Machinery, 1994, p. 151–155. [Online]. Available: https://doi.org/10.1145/182478.182583
  * [4] L. Zaidi, J. A. Corrales, B. C. Bouzgarrou, Y. Mezouar, and L. Sabourin, “Model-based strategy for grasping 3d deformable objects using a multi-fingered robotic hand,” _Robotics and Autonomous Systems_ , vol. 95, pp. 196–206, 2017.
  * [5] Y. Li, D. Xu, Y. Yue, Y. Wang, S.-F. Chang, E. Grinspun, and P. K. Allen, “Regrasping and unfolding of garments using predictive thin shell modeling,” in _2015 IEEE International Conference on Robotics and Automation (ICRA)_ , 2015, pp. 1382–1388.
  * [6] Y. Li, Y. Wang, Y. Yue, D. Xu, M. Case, S.-F. Chang, E. Grinspun, and P. K. Allen, “Model-driven feedforward prediction for manipulation of deformable objects,” _IEEE Transactions on Automation Science and Engineering_ , vol. 15, no. 4, pp. 1621–1638, 2018.
  * [7] E. Pignat and S. Calinon, “Learning adaptive dressing assistance from human demonstration,” _Robotics and Autonomous Systems_ , vol. 93, pp. 61–75, 2017\.
  * [8] J. Schulman, J. Ho, C. Lee, and P. Abbeel, “Learning from demonstrations through the use of non-rigid registration,” in _Robotics Research_. Springer, 2016, pp. 339–354.
  * [9] B. Balaguer and S. Carpin, “Combining imitation and reinforcement learning to fold deformable planar objects,” in _2011 IEEE/RSJ International Conference on Intelligent Robots and Systems_. IEEE, 2011, pp. 1405–1412.
  * [10] H. Lin, F. Guo, F. Wang, and Y.-B. Jia, “Picking up a soft 3d object by “feeling” the grip,” _The International Journal of Robotics Research_ , vol. 34, no. 11, pp. 1361–1384, 2015.
  * [11] W. Yu, A. Kapusta, J. Tan, C. C. Kemp, G. Turk, and C. K. Liu, “Haptic simulation for robot-assisted dressing,” in _2017 IEEE international conference on robotics and automation (ICRA)_. IEEE, 2017, pp. 6044–6051.
  * [12] S. Miller, J. Van Den Berg, M. Fritz, T. Darrell, K. Goldberg, and P. Abbeel, “A geometric approach to robotic laundry folding,” _The International Journal of Robotics Research_ , vol. 31, no. 2, pp. 249–267, 2012.
  * [13] Y. Yamakawa, A. Namiki, and M. Ishikawa, “Simple model and deformation control of a flexible rope using constant, high-speed motion of a robot arm,” in _2012 IEEE International Conference on Robotics and Automation_. IEEE, 2012, pp. 2249–2254.
  * [14] D. McConachie, T. Power, P. Mitrano, and D. Berenson, “Learning when to trust a dynamics model for planning in reduced state spaces,” _IEEE Robotics and Automation Letters_ , vol. 5, no. 2, pp. 3540–3547, 2020.
  * [15] Y. Tsurumine, Y. Cui, E. Uchibe, and T. Matsubara, “Deep reinforcement learning with smooth policy update: Application to robotic cloth manipulation,” _Robotics and Autonomous Systems_ , vol. 112, pp. 72–83, 2019\.
  * [16] D. McConachie and D. Berenson, “Estimating model utility for deformable object manipulation using multiarmed bandit methods,” _IEEE Transactions on Automation Science and Engineering_ , vol. 15, no. 3, pp. 967–979, 2018.
  * [17] J. Matas, S. James, and A. J. Davison, “Sim-to-real reinforcement learning for deformable object manipulation,” in _Conference on Robot Learning_. PMLR, 2018, pp. 734–743.
  * [18] T. P. Lillicrap, J. J. Hunt, A. Pritzel, N. Heess, T. Erez, Y. Tassa, D. Silver, and D. Wierstra, “Continuous control with deep reinforcement learning,” _arXiv preprint arXiv:1509.02971_ , 2015.
  * [19] K. He, X. Zhang, S. Ren, and J. Sun, “Deep residual learning for image recognition,” in _Proceedings of the IEEE conference on computer vision and pattern recognition_ , 2016, pp. 770–778.
  * [20] L. Duan, L. Boyd, and G. Aragon-Camarasa, “Learning physics properties of fabrics and garments with a physics similarity neural network,” 2021. [Online]. Available: https://arxiv.org/abs/2112.10727
  * [21] L. Duan and G. Aragon-Camarasa, “Garnet: A continuous robot vision approach for predicting shapes and visually perceived weights of garments,” 2022. [Online]. Available: https://arxiv.org/abs/2109.07831