Estimation of high dimensional Gamma convolutions through random projections

Estimation of high dimensional Gamma convolutions

A]Oskar<EMAIL_ADDRESS>
[A]Institut Camille Jordan UMR 5208, Université Claude Bernard Lyon 1, France,
and SCOR SE, e1

[Primary ]62H12
[; secondary ]60E07

Multivariate generalized Gamma convolutions
Thorin measures
Random projection
Infinite divisibility

The author was supported by SCOR SE through a CIFRE grant.

[1]
Ambrosio, LuigiL.,
Gigli, NicolaN. Savaré, GiuseppeG.
Gradient flows: in metric spaces and in the space of probability
Springer Science & Business Media.

[2]
Arbenz, PhilippP.,
Hummel, ChristophC. Mainik, GeorgG.
Copula based hierarchical risk aggregation through sample reordering.
Insurance: Mathematics and Economics

[3]
Bélisle, ClaudeC.,
Massé, Jean-ClaudeJ.-C. Ransford, ThomasT.
When is a probability measure determined by infinitely many
The Annals of Probability

[4]
Bernardini, A.A.,
Natalini, P.P. Ricci, P. E.P. E.
Multidimensional Bell Polynomials of Higher Order.
Computers & Mathematics with Applications

[5]
Bondesson, LennartL.
Generalized Gamma Convolutions and Related Classes of
Distributions and Densities.
Lecture Notes in Statistics
Springer New York, New York, NY.

[6]
Bondesson, LennartL.
On Univariate and Bivariate Generalized Gamma Convolutions.
Journal of Statistical Planning and Inference

[7]
Bondesson, LennartL.
A Class of Probability Distributions That Is Closed with
Respect to Addition as Well as Multiplication of
Independent Random Variables.
Journal of Theoretical Probability

[8]
Bondesson, LennartL. Simon, ThomasT.
Stieltjes Functions of Finite Order and Hyperbolic Monotonicity.
Transactions of the American Mathematical Society

[9]
Breger, AnnaA.,
Ehler, MartinM.,
Gräf, ManuelM. Peter, ThomasT.
Cubatures on Grassmannians: moments, dimension reduction, and related
In Compressed Sensing and its Applications

[10]
Chafaï, DD.
Random projections, marginals, and moments.
URL https://djalil. chafai. net/docs/projections. pdf. Accessed on

[11]
Cherian, KCK.
A Bivariate Correlated Gamma-Type Distribution Function.
Journal of the Indian Mathematical Society

[12]
Chizat, LenaicL.
Sparse optimization on measures with over-parameterized gradient
Mathematical Programming

[13]
Chizat, LenaicL. Bach, FrancisF.
On the global convergence of gradient descent for over-parameterized
models using optimal transport.
arXiv preprint arXiv:1805.09545.

[14]
Comon, PierreP. Jutten, ChristianC.
Handbook of Blind Source Separation: Independent component analysis and
Academic press.

[15]
Comte, FabienneF. Genon-Catalot, ValentineV.
Adaptive Laguerre Density Estimation for Mixed Poisson Models.
Electronic Journal of Statistics

[16]
Constantine, G. M.G. M. Savits, T. H.T. H.
A Multivariate Faa Di Bruno Formula with Applications.
Transactions of the American Mathematical Society

[17]
Cramér, HaraldH. Wold, HermanH.
Some theorems on distribution functions.
Journal of the London Mathematical Society

[18]
Cvijović, DjurdjeD.
New Identities for the Partial Bell Polynomials.
Applied Mathematics Letters

[19]
de Castro, YohannY.,
Gadat, SébastienS.,
Marteau, ClémentC. Maugis, CathyC.
SuperMix: Sparse Regularization for Mixtures.
The annals of Statistics

[20]
Di Nardo, EE.,
Guarino, GG. Senato, DD.
A new method for fast computing unbiased estimators of cumulants.
Statistics and Computing

[21]
Di Nardo, ElviraE.,
Guarino, GiuseppeG. Senato, DomenicoD.
A new algorithm for computing the multivariate Faa di Bruno’s
Applied mathematics and computation

[22]
Dozat, TimothyT.
Incorporating nesterov momentum into adam.

[23]
Dussap, FlorianF.
Anisotropic multivariate deconvolution using projection on the Laguerre
Journal of Statistical Planning and Inference

[24]
Ferriero, AlessandroA.
A New Method for Modeling Dependence via Extended Common Shock Type

[25]
Furman, EdwardE.,
Hackmann, DanielD. Kuznetsov, AlexeyA.
On Log-Normal Convolutions: An Analytical-Numerical
Method With Applications to Economic Capital Determination.
SSRN Electronic Journal.

[26]
Heppes, AA.
On the determination of probability distributions of more dimensions by
their projections.
Acta Mathematica Academiae Scientiarum Hungarica

[27]
James, Lancelot F.L. F.,
Roynette, BernardB. Yor, MarcM.
Generalized Gamma Convolutions, Dirichlet Means, Thorin
Measures, with Explicit Examples.
Probability Surveys

[28]
Kim, TaekyunT.,
Kim, DaeD. Jang, Gwan-WooG.-W.
On Central Complete and Incomplete Bell Polynomials I.

[29]
Kingma, Diederik PD. P. Ba, JimmyJ.
Adam: A method for stochastic optimization.
arXiv preprint arXiv:1412.6980.

[30]
Klugman, Stuart AS. A. Parsa, RahulR.
Fitting bivariate loss distributions with copulas.
Insurance: mathematics and economics

[31]
Laverny, OskarO.,
Masiello, EsterinaE.,
Maume-Deschamps, VéroniqueV. Rullière, DidierD.
Estimation of multivariate generalized gamma convolutions through
Laguerre expansions.
Electronic Journal of Statistics
5158 – 5202.

[32]
Mabon, GwennaëlleG.
Adaptive Deconvolution on the Non-Negative Real Line:
Adaptive Deconvolution on R+.
Scandinavian Journal of Statistics

[33]
Mathai, A. M.A. M.
Storage Capacity of a Dam with Gamma Type Inputs.
Annals of the Institute of Statistical Mathematics

[34]
Melman, Vadim SV. S.,
Shablya, Yuriy VY. V.,
Kruchinin, Dmitry VD. V. Shelupanov, Alexander AA. A.
Realization of a Method for Calculating Bell Polynomials Based on
Compositae of Generating Functions.
Journal of Informatics and Mathematical Sciences

[35]
Miatto, Filippo M.F. M.
Recursive Multivariate Derivatives of $e^{f(X_1,\ldots, X_n)}$ of
Arbitrary Order.
arXiv:1911.11722 [cs, math].

[36]
Miles, JustinJ.,
Furman, EdwardE. Kuznetsov, AlexeyA.
Risk Aggregation: A General Approach via the Class of Generalized
Gamma Convolutions.

[37]
Mishkov, Rumen L.R. L.
Generalization of the Formula of Faa Di Bruno for a Composite
Function with a Vector Argument.
International Journal of Mathematics and Mathematical Sciences

[38]
Moschopoulos, P. G.P. G.
The Distribution of the Sum of Independent Gamma Random Variables.
Annals of the Institute of Statistical Mathematics

[39]
Pérez-Abreu, VictorV. Stelzer, RobertR.
A Class of Infinitely Divisible Multivariate and Matrix Gamma
Distributions and Cone-Valued Generalised Gamma Convolutions.
arXiv:1201.1461 [math, stat].

[40]
Pérez-Abreu, VictorV. Stelzer, RobertR.
Infinitely Divisible Multivariate and Matrix Gamma Distributions.
Journal of Multivariate Analysis

[41]
Rényi, AlfrédA.
On projections of probability distributions.
Acta Math. Acad. Sci. Hungar

[42]
Roynette, BernardB.,
Vallois, PierreP. Yor, MarcM.
A Family of Generalized Gamma Convoluted Variables.

[43]
Smith, Kevin D.K. D.
A Tutorial on Multivariate $k$-Statistics and Their
arXiv:2005.08373 [math, stat].

[44]
Smith, Peter JP. J.
A recursive formulation of the old problem of obtaining moments from
cumulants and vice versa.
The American Statistician

[45]
Thorin, OlofO.
On the Infinite Divisibility of the Pareto Distribution.
Scandinavian Actuarial Journal

[46]
Thorin, OlofO.
On the Infinite Divisibility of the Lognormal Distribution.
Scandinavian Actuarial Journal

[47]
Wendland, HolgerH.
Scattered data approximation
Cambridge university press.


# Hybridization of terahertz phonons and magnons in disparate and spatially-
separated material specimens

M. Białek<EMAIL_ADDRESS>Institute of High Pressure Physics,
Polish Academy of Sciences, Warszawa, Poland  Y. Todorov  Laboratoire de
Physique de l’Ecole Normale Supérieure, ENS, Université PSL, CNRS, Sorbonne
Université, Université Paris Cité, Paris, France  K. Stelmaszczyk  Institute
of High Pressure Physics, Polish Academy of Sciences, Warszawa, Poland  D.
Szwagierczak  B. Synkiewicz-Musialska  J. Kulawik  Łukasiewicz Research
Network–Institute of Microelectronics and Photonics, Kraków Division, Kraków,
Poland  N. Pałka  Institute of Optoelectronics, Military University of
Technology, Warszawa, Poland  M. Potemski  LNCMI-EMFL, CNRS UPR3228, Univ.
Grenoble Alpes, Univ. Toulouse, Univ. Toulouse 3, INSA-T, Grenoble and
Toulouse, France  CENTERA, CEZAMAT, Warsaw University of Technology, 02-822
Warszawa, Poland  Institute of Experimental Physics, Faculty of Physics,
University of Warsaw, Warszawa, Poland  W. Knap  Institute of High Pressure
Physics, Polish Academy of Sciences, Warszawa, Poland  CENTERA, CEZAMAT,
Warsaw University of Technology, 02-822 Warszawa, Poland

###### Abstract

The interaction between light and matter in condensed matter excitations and
electromagnetic resonators serves as a rich playground for fundamental
research and lies at the core of photonic and quantum technologies. Herein, we
present comprehensive experimental and theoretical studies of the photon-
mediated hybridization of magnons and phonons in the terahertz (THz) range. We
demonstrate the intriguing concept of composite states formed by distinct
electric and magnetic quasiparticles strongly coupled to the same optical
cavity modes. Specifically, we explore magnons excited in a slab of an
antiferromagnetic crystal and phonons excited in a distinct specimen of an
insulating material. The crystal slabs form an optical cavity with Fabry-Pérot
oscillations in the THz range. We demonstrate hybridized phonon-magnon
polariton modes and their tunability by adjusting the distance between the
slabs, showing that hybridization persists even at separations up to several
millimeters. The experimental results are interpreted using both classical and
quantum electrodynamical models. The quantum description allows us to quantify
the degree of hybridization that is linked to a topological behavior of the
electric field phasor, in agreement with the classical electrodynamics
expectations. Importantly, the presented results refer to temperature
conditions and cavities of millimeter size, paving the way for engineering
realistic, frequency-tunable THz devices through the hybridization of electric
(phononics) and magnetic (spintronics) elementary excitations of matter.

In the regime of strong coupling between matter excitations and
electromagnetic cavities the energies of light and matter excitations are
periodically exchanged, resulting in mixed light-matter coupled quasiparticles
known as polaritons [1, 2]. Light-matter coupling was first reported with
single atoms [3], and later with single excitons bound in quantum dots [1, 4].
The light-matter coupling strength of $N$ oscillators placed in a cavity is
enhanced by a factor of $\sqrt{N}$, a phenomenon known as Dicke cooperativity
[5, 6, 7, 8, 9, 10, 11, 12]. It allows the achievement of observable couplings
of magnetization oscillations despite a tiny dipolar moment of a single spin.
The strong coupling of microwave cavity modes and magnons was observed as
early as 1962 [13], and since the 2010s it has been intensively studied [14,
15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27], because of the prospect
of using magnon-polaritons in quantum devices [28, 29, 30, 31].

Strong light-matter coupling can form hybridized states of two distant
resonators exchanging energy with the same electromagnetic cavity mode. Such
light-bound composite states were shown for different electric dipole
resonators [32, 33, 12]. The coupling of two ferromagnetic magnons mediated by
a cavity mode was recently achieved using superconducting circuits [34, 35]
and nanostripline antennas [36], and was discussed theoretically [34, 37, 38,
39]. We recently showed that magnons in distant antiferromagnets can be
coupled using Fabry-Pérot cavity modes [40].

An even more intriguing concept is the possibility of achieving composite
states of two different optically active quasiparticles, such as phonons,
magnons, collective excitations such as plasmons, or Landau level excitations,
many of which have energies in the THz range. Such couplings of magnons are
commonly observed with other excitations in the same material, such as
acoustic phonons [41, 42], optical phonons in a multiferroic [43] and in
layered magnets [44, 45, 46], excitons in a layered semiconductor [47], and
electron paramagnetic resonance [11]. However, photon-mediated coupling of
distinct excitations in different materials, was only shown in the weak
regime, i.e. without frequency matching of phonons and magnons, in a
heterostructure of an antiferromagnet and lithium niobate (LiNbO3) [48].

Here, we demonstrate frequency-matched coupling of antiferromagnetic magnons
and optical phonons hosted exclusively in two different and distant samples.
This coupling is mediated by cavity modes formed by the sample slabs
themselves.

Figure 1:  The measured system comprised NiO crystal placed on a copper mirror
and a slab of CBO. We measured reflection as a function of the gap $d$ between
the slabs and the temperature $T$ of the NiO crystal. In the top-right corner
of the figure, we show a sample time trace.

Experiments were performed close to room temperature and without external
fields. These two sample slabs were placed next to each other at a well-
controlled gap, and the slabs themselves formed the cavity, that is, the
cavity fields penetrated both materials (Fig. 1). We show phonon-magnon-
polariton modes, under frequency matching of magnon, phonon, and cavity mode.
Such hybrid modes have magnetic and electric dipole moments originating from
their constituents. We show that our experimental data can be interpreted with
a model based on classical electrodynamics, and a microscopic model. The
magnon-photon coupling remains observable when the crystal separation reaches
the millimeter range. The separate magnon and phonon modes are coupled due to
the strong coupling of each resonance to the same mode of a Fabry-Pérot
cavity. This research demonstrates hybrid resonators that could be used in
communication, quantum technologies, or for the detection of hypothetical
axion particles [49, 50], using materials hosting topological magnons [51, 52]
with the advantage of large resonator volumes and tunable detection frequency
[53].

We used a single crystal of nickel oxide (NiO) 111 cut of thickness
$h_{m}=0.26$ mm and lateral diameter of 5 mm. The antiferromagnetic resonance
in this antiferromagnetic material, the $\alpha$-mode [54], has a frequency
close to 1.0 THz at room temperature. It decreases with rising temperature,
reaching 0.7 THz at about 450 K. We chose ceramics of Cu$B_{2}$O4 (CBO) [55],
as a material that hosts a relatively narrow phonon mode at about 0.92 THz,
which can be crossed with the magnon in NiO. The oscillator strength of the
phonon in CBO is weak for an electric dipole excitation, while the oscillator
strength of the magnon in NiO is quite strong for a magnetic dipole. This
allows for comparable light-matter coupling strengths in both cases, so the
electric dipole excitation does not completely dominate the magnetic dipole
excitation. The NiO crystal was fixed on a copper mirror, the temperature of
which was measured with a K-type thermocouple. Temperature was stabilized with
a software proportional–integral–derivative controller setting electric
current applied to a stack of Peltier elements beneath the copper mirror. The
NiO crystal and CBO ceramics were placed on independent kinematic mounts,
allowing for precise parallelism control. The size $d$ of the gap between the
samples was controlled with a motorized stage that supports the NiO part of
the cavity. Thermal transfer due to convection was small for $d\gtrsim 0.2$
mm, and the phonon frequency does not depend on temperature. The experimental
setup was purged with dry air to minimize the effect of atmospheric water
absorption with a few strong lines around 1.0 THz.

We used the THz time-domain spectrometer system (Toptica). Optical pulses were
guided in single-mode zero-dispersion fibers, and a biased source antenna was
emitting pulses of electromagnetic radiation with spectrum power centered
around 1.0 THz. The THz beam was propagated from the source using a quasi-
optical setup with parabolic mirrors. We used a beam splitter made of Kapton
tape, which directs pulses reflected from the sample to the detector antenna.
We used a 200 ps delay line, which gives a 5 GHz resolution. In the presented
spectra, we extended these time-domain data to 1000 ps with zeros,
interpolating the frequency-domain resolution to 1 GHz. Spectra were collected
as a function of gap $d$ between the crystals for a fixed NiO temperature $T$.
Then, the temperature was changed with a step of 1.0 K and spectra were
collected again as a function of the gap $d$. We Fourier-transform measured
time traces of electric field, which yields complex spectra. This report shows
the standard magnitude of obtained spectra, clearly showing strong coupling
and magnons, phonons, and cavity modes.

Figure 2:  Comparison of the (a) magnitude spectrum with its (b) phase and (c)
phase differential to frequency. The latter shows typically, much narrower
features around resonances. Also, the phase jump can be positive or negative,
whereas optical and matter resonances always result in a minimum amplitude of
reflection spectra.

However, often omitted, phase carries more information of better quality than
magnitude. This is mostly because spectral features in phase are narrower than
in magnitude, and phase jumps at resonances can have positive or negative
signs. In our convention, the reflection phase $\Phi$ drops almost linearly
with rising frequency, as expected in the case of electromagnetic waves
propagating in a non-dispersive infinite medium. We remove this constant trend
by calculating phase differential to frequency $d\Phi/df$, in the unit of
rad/GHz. This treatment of phase spectra reveals rapid phase changes
characterizing resonances, both cavity modes and optically active transitions
in the solid state such as phonons or magnons. Such phase data do not require
any further normalization. In our spectra, we observe that the sign of optical
resonances is typically negative, that is, in the direction of linear change
of phase with increasing frequency. On the other hand, the signs of NiO magnon
and CBO phonon phase jumps are always positive in our convention (white and
red in 3). Such pure matter modes are barely visible due to either, the low
amplitude of the magnon, or the large width of the phonon. However, they
become much more pronounced when they couple with cavity modes. Moreover, we
observed that the sign of polariton modes is also always positive in the
strong-coupling regime. We observed, that the sign of phase jump is an
indicator of the strong coupling, and only cavity modes with positive phase
produce show strong coupling, as visible in Fig. 3b, where the magnon just
below 1.0 THz shows only weak coupling with FP modes and cavity sign is
negative (blue), while the phonon forms polaritons modes with FP modes and its
sign is positive (red).

Mode frequencies of the FP cavity are lowering, with the rising gap $d$ (Fig.
3). When they cross with the CBO phonon mode at around 0.92 THz, they form two
phonon-polariton modes, about 20 GHz apart. In Fig. 3b phase jumps of these
polariton modes are characterized by positive jumps (red color). By our
convention, cavity modes have typically negative signs represented with blue
colors (Fig. 3). However, some strong cavity modes may have positive signs,
like the one at around 1.05 THz. Magnon mode in NiO is only slightly visible
near 1.0 THz because the results were obtained at 300 K. In this case, the
interaction of the magnon with FP cavity modes is weak, resulting in only
small splittings, comparable with polariton widths.

Figure 3:  Reflection spectra: (a) amplitude and (b) phase differential to
frequency. Spectra are plotted as a function of the gap $d$ between the slabs.

The magnon frequency in NiO was tuned in the range of about 0.7–1.0 THz, using
its dependence on crystal temperature in the range of about 300–450 K. When
the magnon frequency passes the frequency of the phonon in CBO (0.92 THz),
these two matter modes would cross without any interaction, as visible in Fig.
6g ($d=305$ $\mu$m). However, when the gap between the slabs is such that one
of the cavity modes has a frequency close to the CBO phonon, we observe
phonon-magnon-polariton modes as a function of NiO temperature. At a
particular distance of $d=260$ $\mu$m (Fig. 6d), it is visible that the magnon
mode reaches the phonon frequency at about $T=360$ K. The magnon interacts
with the phonon-polaritons, making them narrower and more separated. The
narrowing is because magnons have a line width much smaller than that of
phonons, thus, owing to hybridization, the resulting phonon-magnon-polaritons
are narrower than pure phonon-polaritons. Phonon-magnon polaritons are also
slightly more split than pure phonon-polaritons because of magnons’ small
additional oscillator strength. It is also shown in Fig. 4ab, where green
lines were measured for a magnon-phonon matching condition ($T=361$ K), while
the violet curve shows pure phonon-polaritons with magnon mode around 0.96
THz. A small change of the gap $d$ between the slabs compensates for the drift
of dielectric constants of both materials with temperature. In Fig. 4cd, we
show pure magnon-polariton modes around 0.89 THz, and pure phonon mode around
0.92 THz.

Figure 4:  Reflection spectra showing the influence of NiO magnon on phonon-
polaritons. Magnon frequency sweep was achieved by changing NiO crystal
temperature. Magnitude (a) and phase (b) of reflection spectra. Green lines
guide the eye showing phonon-polaritons, and a violet line indicates the
magnon mode. Phasor plots (c,d) in the range of 985-940 GHz for different
temperatures, show a double loop indicating phonon polaritons (green arrows),
and a feature indicating magnon (violet arrow). Phasor plot (e) showing
phonon-magnon polaritons.

Phase data can also be analyzed with amplitude data, using phasor plots of the
imaginary part of the electric field as a function of its real part, as in
Fig. 4ce. Normally, the signal circles the origin with rising frequency, which
reflects phase progression. Magnon and phonon modes are very weakly
distinguishable as phase progression slowdowns and signal amplitude drops
related to absorption. In cavity modes, the circle progresses faster and
closer to the origin. In some cavity modes, the loop becomes so tight that it
no longer circles the origin, resulting in a reversed-phase progression, such
as the modes around 1.05 THz in Fig. 3b (more details in the Appendix C). The
phasor signature of the strong coupling is two intertwined loops [56], and
neither of them encircles the origin point. Thus, the phase progression
reverses twice during this signal evolution (green arrows indicate phonon-
polaritons in Fig. 4c). We also demonstrate in detail the formation of phonon-
polaritons phasor signature in Appendix C, when the magnon mode was spectrally
distant from the phonon. In Fig. 4c-e, we show the influence of magnon mode
(features indicated by violet arrows) on phonon polaritons (loops indicated by
green arrows) The loops indicating the strong coupling are much smaller for
phonon-polaritons (green arrows in Fig. 4cd) than for phonon-magnon-polaritons
(Fig. 4e). This is related to the modes narrowing and deepening, caused by the
hybridization with magnon. In contrast to Fig. 4cd, in Fig. 4e, the magnon
signature cannot be any more distinguished otherwise than the increase of the
loops, which indicates the hybridization of phonons and magnons.

Having demonstrated the strong coupling between magnons, phonons, and Fabry-
Pérot resonances experimentally, we now quantitatively investigate the mixing
between them. This mixing is best described within the framework of a quantum
model which treats the system as coupled quantum harmonic oscillators. This
model has been detailed in the Appendix A. Starting from a general quantum
electrodynamical description the total Hamiltonian of the tripartite system is
derived to be:

$\mathcal{H}=\hbar\omega_{m}(a^{\dagger}a+1/2)+\hbar\omega_{LO}d^{\dagger}d+\hbar\omega_{\alpha}c^{\dagger}c\\\
+i\hbar\Omega_{p}(a-a^{\dagger})(d+d^{\dagger})+i\hbar\Omega_{m}(c-c^{\dagger})(a+a^{\dagger}),$
(1)

This Hamiltonian was derived by considering the coupling of the electric field
of the FP resonance with the phonon polarization. In contrast, the magnetic
field is coupled to the magnon magnetization field through a Zeeman
interaction. Here $\omega_{m}$ is the frequency of the Fabry-Pérot resonance,
$\omega_{LO}$ is the LO optical phonon frequency of the CBO slab, and
$\omega_{\alpha}$ is the $\alpha-$magnon frequency. The operator $a$ is the
bosonic destruction operator for the photon mode and $d$ and $c$ are
collective boson operators describing the phonon and magnon excitations in the
two slabs. In the above equation, we have introduced the phonon-photon
coupling rate $\Omega_{p}$ and the magnon-photon coupling rate, $\Omega_{m}$.
As shown in Annex A they are provided by:

$\hbar\Omega_{p}=\frac{1}{2}\hbar\eta_{p}\sqrt{\frac{(\omega_{LO}^{2}-\omega_{TO}^{2})\omega_{m}}{\omega_{LO}}}$
(2)

and

$\hbar\Omega_{m}=\frac{1}{2}\left(\frac{H_{Ax}}{H_{E}}\right)^{\frac{1}{4}}\eta_{m}g\mu_{B}\sqrt{\frac{\hbar\omega_{m}\mu_{0}}{2}\left(\frac{N}{V_{NiO}}\right)}.$
(3)

We recall the transverse and longitudinal optical photon frequencies are
linked by the relation
$\omega_{TO}=\omega_{LO}\sqrt{\frac{\epsilon_{\infty}}{\epsilon_{st}}}$, where
the dielectric constants have been determined experimentally as described
further.

The collective nature of the magnon-light coupling is obvious from the
expression of the coupling constant, Eq.(3), which is proportional to the
square of the spin density $N/V_{NiO}$. In NiO, $\mu_{0}H_{Ax}\approx 0.5$ T
is the easyplane anisotropy field at room temperature, $\mu_{0}H_{E}=698$ T is
the exchange field [54], and spin density $N/V_{NiO}=5.489\times 10^{28}$ m-3.
[57, 58] By diagonalizing the Hamiltonian (1) through the Hopfield-Bogoliubov
procedure we obtain the dispersion relation for the polariton frequencies and
closed expressions for the mixing coefficients (Appendix A). To match the
coupled mode frequencies with the experimental ones, we have considered the
overlap factors $\eta_{p},\eta_{m}$ as the fitting parameters.

Figure 5:  Measured phase differential to frequency of reflection spectra
plotted as a function of NiO temperature $T$, for a matching mode frequency
(a). Predictions of the electromagnetic model (b). Predictions of interacting
mode frequencies, and magnon, phonon, and photon contents using the
microscopic model (Eq. 1) (c). Hopfield coefficients for photon, magnon, and
phonon (d).

In Fig. 5c, we use effective coefficients $\eta_{p}=.4$, and $\eta_{m}=0.6$,
which better reproduce the frequencies observed in the experimental spectra
(Fig. 5a). Calculations based on the microscopic model show that the middle
mode, under the frequency matching condition $f_{c}=f_{m}=f_{p}$, is composed
in about 80% of magnon and in 20% of phonon, while it has a less than 0. 3%
photon content. Unequal contents of the phonon and the magnon in this mode
result from their unequal coupling strengths to a cavity mode. This property
does not change much when tuning the magnon frequency within about 50% of the
phonon-polariton gap. However, in the experiment (Fig. 5a) this mode has a
small positive amplitude, which is not predicted using the microscopic model.
At the same condition of frequency matching, two other modes contain
approximately equally 45% of phonon and photon, while only about 10% of
magnon. Modes of larger amplitude and with more equal contents of phonon and
magnon can be observed when the magnon mode has a frequency close to either
upper or lower phonon-polaritons, that is, in Fig. 5, $|f_{m}-f_{p}|\approx
12$ GHz. Then, two phonon-magnon polariton modes with a smaller splitting are
formed, which both contain about 40% phonon and magnon and only 20% photons.
Due to a smaller splitting, they are strongly tunable with magnon frequency in
terms of magnon and phonon content.

Figure 6:  Measured phase differential to the frequency of reflection spectra
plotted as a function of NiO temperature $T$, for 7 different gaps $d$ between
the slabs (a-g). Predicted spectra of reflection phase differential to
frequency, as a function of the detuning of the NiO magnon from the frequency
of the CBO phonon, for 7 different gaps $d$ between the slabs (h-n).

The reflection spectra can be modeled using the expressions for the dielectric
and magnetic susceptibilities of the layers, which were determined using
transmission measurements of each of the layers (Appendix B). This approach
allows a direct comparison of the amplitude and phase-resolved measurement for
the electric field as obtained by the THz-TDs experiment. The compatibility of
this approach and the quantum description are commented in the Appendix A. We
analytically calculated reflection from a series of parallel-plane slabs using
a characteristic matrix method for isotropic media [59], details of which are
explained in Ref. 40. We assumed a Lorentzian shape of both magnon and phonon
modes, and otherwise constant background permittivity in the range of interest
(0.5-1.5 THz). The CBO sample is ceramics, so it is right to assume that it is
isotropic. In the case of the NiO crystal, $\mu$ does not contain any off-
diagonal terms, therefore, we can consider it isotropic in our model. Details
of assumed permittivity and permeability are shown in Appendix A.

The predicted reflection phase spectra are shown in Fig. 6(h-n) as a function
of magnon-phonon frequency detuning ($f_{p}-f_{m}$) for a few different gaps
between the modeled slabs. We used a simplified model where cavities are
independent of the detuning (temperature in the experiment). Qualitatively,
this model predicts a similar effect as measured as a function of temperature
(Fig. 6(a-g)). The model predicts that the gap, at which cavity mode matches
phonon frequency, is different from measured values by 29 $\mu$m, which might
stem from inaccuracies of thickness measurement, gap measurement, and
dielectric constant determination. The model correctly predicts a change in
the phase jump sign of polaritons when they interact strongly with either
phonons or magnons. No magnon-phonon coupling is predicted when none of the
cavity modes matches the phonon mode (Fig. 6gn). In Fig. 7, we present some
predicted spectra, which are chosen for parameters approximately reproducing
the experimental spectra presented in Fig. 4de and Fig. 2. Panels a and b of
Fig. 7 show the influence of the magnon on the phonon-polariton modes, with
spectra similar to those in Fig. 2. The influence of magnon on phonon-
polaritons is visible as a narrowing of the modes and a small increase in the
splitting. Panels c and d of Fig. 7 show magnon-polaritons at about 0.88 THz,
away from the phonon mode at 0.92 THz. Here phonon does not interact, because
there is no cavity mode at a close frequency.

Figure 7:  The predicted reflection spectra are calculated using model
parameters that roughly reproduce the experimental data presented in Fig. 4.
Magnitude (a) and phase (b) spectra show phonon-polaritons (violet) for the
magnon frequency away from the phonon and phonon-magnon-polaritons (green) for
zero detuning. Magnitude (c) and phase (d) spectra, showing magnon-polaritons
around 0.88 THz.

In summary, we presented systematic experimental studies on the interaction of
two distinct polariton systems as a function of temperature and separation
distance. The cooperative interaction of magnons in the slab of nickel oxide
(NiO) and phonons in the distant slab of CuB2O4 (CBO) is mediated by a cavity
mode, which leads to the formation of phonon-magnon polaritons states. The
experiments were performed at room temperatures, and the separation between
the crystals was in the millimeter range, i.e. showing the potential of this
coupling scheme for applications. We used the electromagnetic model, which
describes the observed features with realistic assumptions for the parameter
values. Cavity-mediated coupling enables control of antiferromagnetic
polaritons by modulating the resonator surroundings. It shows the possibility
of designing hybrid resonators operating in the THz range (few meV) that share
properties of their components, such as width, amplitude, and dependence on
external parameters.

More generally, our results present an even more intriguing concept of
composite states made of optically active quasiparticles of distinct nature,
which exchange the same cavity photon in the strong coupling regime. These
excitations could be plasmons in quantum confined electron gases [60, 61], or
intersubband cavity polaritons that can reach the ultra-strong coupling regime
[62]. A very appealing perspective is realizing such hybridization with THz
meta-atom resonators that can provide strong confinement for both the electric
and magnetic field [63]. Furthermore, the preparation of single magnon states
through coupling with superconductor qubits has been demonstrated in the GHz
range [25]. In principle, it is possible to realize a similar concept in the
THz range by employing high Tc superconductors [64]. Thus it could become
feasible to transfer quantum non-linearities into hybrid-magnon-quantum well
polariton systems in the ultra-strong coupling regime at 1 THz, overcoming the
intrinsic limitations owed to the bosonic character of these systems. We can
also benefit from technologies available at this frequency range, such as
semiconductor optoelectronics and THz-time domain to conceive original quantum
optical experiments. This is another venue for bringing quantum technologies
into the THz frequency range [65, 66], and paving the way towards innovative
tunable THz devices operating at room temperature.

###### Acknowledgements.

We acknowledge the funding from Pasific2 of the Polish Academy of Sciences
sponsored by the European Union’s Horizon 2020 research and innovation program
under the Marie Skłodowska-Curie grant agreement No. 847639 and from the
Ministry of Education and Science of Poland. This work was also partially
supported by the “International Research Agendas” program of the Foundation
for Polish Science, co-financed by the European Union under the European
Regional Development Fund (No. MAB/2018/9). We also acknowledge funding from
ERC-COG-863487 “UNIQUE”.

##  Appendix A Quantum model

Figure 8: Schematic of the tripartite system composed of phonon, magnon, and
mediating cavity photon.

This Appendix provides an effective quantum model for the interacting
tripartite system (Fig. 8). We have two objectives: the first one is to
provide a quantum Hamiltonian that is compatible with effective dielectric and
magnetic functions, where the polariton frequencies derived from the model can
be directly compared with the data. Our second objective is to provide
analytical expressions for the Hopfield mixing coefficients that allow
evaluating the mixing between magnons and phonons mediated by the
electromagnetic field.

###  A.1 Full Hamiltonian and photon field

We start with a very general Hamiltonian which describes photons coupled both
to electrically and magnetically active matter:

$\displaystyle\mathcal{\hat{H}}=\mathcal{\hat{H}}_{phonon}+\mathcal{\hat{H}}_{magnon}$
$\displaystyle+\int_{V}\Bigg{[}\frac{\hat{\mathbf{D}}^{2}(\mathbf{r})}{2\varepsilon_{\infty}(z)\varepsilon_{0}}+\frac{\mu_{0}\hat{\mathbf{H}}^{2}(\mathbf{r})}{2}\Bigg{]}d^{3}\mathbf{r}$
$\displaystyle-\int_{V}\frac{\hat{\mathbf{D}}(\mathbf{r})\hat{\mathbf{P}}(\mathbf{r})}{\varepsilon_{\infty}(z)\varepsilon_{0}}d^{3}\mathbf{r}+\int_{V}\frac{\hat{\mathbf{P}}^{2}(\mathbf{r})}{2\varepsilon(z)\varepsilon_{0}}d^{3}\mathbf{r}$
$\displaystyle-\int_{V}\mu_{0}\hat{\mathbf{H}}(\mathbf{r})\hat{\mathbf{M}}(\mathbf{r})d^{3}\mathbf{r}$
(4)

The light-matter interaction Hamiltonian is expressed here in the Power-
Zienau-Wooley representation of Quantum Electrodynamics [67]. The first line
contains the Hamiltonians of the phonon and magnon excitations that will be
detailed further. The second line of Eq. (4) describes the energy of the
electromagnetic field. The third line describes the light-matter interaction
with the dipolar excitations of interest, which in the present case are the
phonons of the CBO slab. Respectively, the field
$\hat{\mathbf{P}}(\mathbf{r})$ is the polarization field of the CBO phonons.
In that case the dielectric function $\varepsilon_{\infty}(z)$ is a piecewise
constant function that does not contain the phonon contribution and is
expressed as:

$\displaystyle\varepsilon_{\infty}(z)=\varepsilon^{CBO}_{\infty},\phantom{Q}z\in
CBO,$ $\displaystyle=\varepsilon^{NiO}_{\infty},\phantom{Q}z\in NiO$
$\displaystyle=1,\phantom{Q}z\in\mathrm{Elsewhere}$ (5)

The last line of Eq. (4) describes the light-matter interaction with the
$\alpha$ magnon of NiO [54], with $\hat{\mathbf{M}}(\mathbf{r})$ the
corresponding magnetization that is provided explicitly further.

We will consider the case where the matter systems interact with an
electromagnetic wave that is polarized along the dielectric layers, along the
unit vector $\mathbf{u}_{x}$, and propagates in the $z$-direction. We can then
write a general expression of the electromagnetic field:

$\displaystyle\hat{\mathbf{H}}(\mathbf{r})=\mathbf{u}_{y}\sqrt{\frac{\hbar\omega_{m}\varepsilon_{0}c^{2}}{2SI}}h(z)(a+a^{\dagger})$
(6)
$\displaystyle\hat{\mathbf{D}}(\mathbf{r})=-\mathbf{u}_{x}\sqrt{\frac{\hbar\varepsilon_{0}c^{2}}{2\omega_{m}SI}}\partial_{z}h(z)i(a-a^{\dagger})$
(7) $\displaystyle I=\int h(z)^{2}dz$ (8)

Here $\omega_{m}$ is the frequency of the photon, $c$ is the speed of light in
vacuum, $S$ is the surface of the system in the $xy$ plane, and $h(z)$ is a
function that describes the solution of Maxwell’s equation for a wave that
propagates in the multilayered medium described by the piecewise function from
Eq. (5). The quantity $I$ defined in Eq. (8) is a normalization coefficient
that allows the definition of the rising and lowering bosonic operators $a$
and $a^{\dagger}$. We assume that we can define a finite volume of space
containing the system and that contains the majority of the electromagnetic
energy. The corresponding length along the $z$ direction will be noted $L$.
This assumption is questionable for the type of waveguide modes that are used
in the experiment, however, we will see further how this formalism can still
be applied by introducing effective overlap factors and linking the polariton
dispersion to the effective dielectric and magnetic susceptibilities of the
layers. Finally, in accordance with Maxwell’s equation the function $h(z)$
must obey the identity:

$\int_{L}\frac{(\partial_{z}h(z))^{2}}{\varepsilon_{\infty}(z)}dz=\frac{\omega_{m}^{2}}{c^{2}}\int_{L}h(z)^{2}dz$
(9)

This identity states that the total magnetic energy is equal to the total
electric energy in the system, as expected for guided modes.

###  A.2 Light-phonon coupling in the CBO slab

A very general expression for the phonon Hamiltonian and polarization is
discussed in the Ref[68]. In particular, it can be shown for the optical
phonon the following results hold, in the long wavelength limit:

$\displaystyle\mathcal{\hat{H}}_{phonon}+\int_{V_{2}}\frac{\hat{\mathbf{P}}^{2}(\mathbf{r})}{2\varepsilon_{\infty}(z)\varepsilon_{0}}d^{3}\mathbf{r}=\sum_{\mathbf{Q},\mathbf{n}}\hbar\omega_{LO}d^{\dagger}_{\mathbf{Q}\mathbf{n}}d_{\mathbf{Q}\mathbf{n}}$
(10) $\displaystyle\hat{\mathbf{P}}(\mathbf{r})=$
$\displaystyle\sum_{\mathbf{Q},\mathbf{n}}\sqrt{\frac{\hbar\varepsilon^{CBO}_{\infty}\varepsilon_{0}(\omega_{LO}^{2}-\omega_{TO}^{2})}{2\omega_{LO}V_{2}}}\mathbf{n}e^{i\mathbf{Q}\mathbf{r}}(d_{\mathbf{Q}\mathbf{n}}+d^{\dagger}_{-\mathbf{Q}\mathbf{n}})$
(11)

The phonon excitations are expanded in a plane-wave basis of bosonic basis set
described by the operators $d_{\mathbf{Q}\mathbf{n}}$ and
$d^{\dagger}_{\mathbf{Q}\mathbf{n}}$ , where the $\mathbf{Q}$ is the
propagation wavevector and $\mathbf{n}$ is a unit vector defined as
$\mathbf{n.Q}=0$ which fixes the polarization of the phonon wave. Here $V_{2}$
is the volume of the CBO slab of thickness $\Delta_{2}$, and the two phonon
frequencies satisfy the Lyddane–Sachs–Teller relation:

$\varepsilon^{CBO}_{s}\omega_{TO}^{2}=\varepsilon^{CBO}_{\infty}\omega_{LO}^{2}$
(12)

The plane waves are normalized such as:

$\int_{V_{2}}e^{i\mathbf{Q}\mathbf{r}}d^{3}\mathbf{r}=V_{2}\delta_{\mathbf{Q}\mathbf{0}}$
(13)

We use Eq.(11) together with Eq.(7) to express the light-matter coupling term,
where the integral runs over the volume of the CBO slab alone:

$\displaystyle-\int_{V_{2}}\frac{\hat{\mathbf{D}}(\mathbf{r})\hat{\mathbf{P}}(\mathbf{r})}{\varepsilon^{CBO}_{\infty}\varepsilon_{0}}d^{3}\mathbf{r}=$
$\displaystyle i(a-a^{\dagger})\frac{\hbar
c}{2\sqrt{\varepsilon^{CBO}_{\infty}}}\sqrt{\frac{\omega_{LO}^{2}-\omega_{TO}^{2}}{\omega_{LO}\omega_{m}}}\sqrt{\frac{1}{I\Delta_{2}}}\times$
$\displaystyle\sum_{Q_{z}}\int_{\Delta_{2}}\partial_{z}h(z)e^{iQ_{z}z}dz(d_{Q_{z}\mathbf{u}_{x}}+d^{\dagger}_{-Q_{z}\mathbf{u}_{x}})$
(14)

Here we used the normalization condition Eq. (13). We then introduce new
bosonic operators for the phonons which integrate over all wavevectors:

$\displaystyle
d=\sum_{Q_{z}}\int_{\Delta_{2}}\partial_{z}h(z)e^{iQ_{z}z}dzd_{Q_{z}\mathbf{u}_{x}}\mathcal{N}$
(15) $\displaystyle
d^{\dagger}=\sum_{Q_{z}}\int_{\Delta_{2}}\partial_{z}h(z)e^{-iQ_{z}z}dzd_{Q_{z}\mathbf{u}_{x}}\mathcal{N}$
(16) $\displaystyle\mathcal{N}^{-1}=\sqrt{\int_{V2}(\partial h(z))^{2}dz}$
(17)

The normalization factor in the last equation ensures that the new bosonic
operators satisfy $[d,d^{\dagger}]=1$. Next, we use Eq.(9) to express the
integral over the function $(\partial h(z))^{2}$ that appears above. For that,
we make an approximation that the balance between electric and magnetic energy
is also satisfied within the CBO slab alone:

$\int_{\Delta_{2}}\frac{(\partial_{z}h(z))^{2}}{\varepsilon^{CBO}_{\infty}}dz\approx\frac{\omega_{m}^{2}}{c^{2}}\int_{\Delta_{2}}h(z)^{2}dz$
(18)

With that assumption, the light-phonon coupling Hamiltonian becomes:

$\displaystyle-\int_{V_{2}}\frac{\hat{\mathbf{D}}(\mathbf{r})\hat{\mathbf{P}}(\mathbf{r})}{\varepsilon^{CBO}_{\infty}\varepsilon_{0}}d^{3}\mathbf{r}=$
$\displaystyle i\hbar\Omega_{P}(a-a^{\dagger})(d+d^{\dagger})$ (19)
$\displaystyle\hbar\Omega_{P}=\eta_{p}\frac{1}{2}\hbar\sqrt{\omega_{LO}^{2}-\omega_{TO}^{2}}\sqrt{\frac{\omega_{m}}{\omega_{LO}}}$
(20)
$\displaystyle\eta_{p}=\sqrt{\int_{\Delta_{2}}h(z)^{2}dz\Big{/}\int_{L}h(z)^{2}dz}$
(21)

In the above expression, we introduced a dimensionless overlap factor
$0<\eta_{p}<1$ that quantifies the phonon coupling with the guided waves.
Following Ref. [68] we can further show how we can derive a dielectric
susceptibility for the CBO slab that includes the phonon contribution. We will
come back to this topic in the last section of this Appendix.

###  A.3 Light-magnon coupling in the NiO slab

In order to express the magnetization field in the NiO slab we use the
theoretical description of magnons described in the review paper Ref. [54]. In
this section, we will use the notations from that reference. Magnons of NiO
are treated in the framework of easy-plane antiferromagnets in the AF phase
(section D in Ref. [54]) that are described by the following Hamiltonian:

$\mathcal{\hat{H}}_{magnon}=\sum_{k}[\hbar\omega_{\alpha
k}\alpha_{k}^{\dagger}\alpha_{k}+\hbar\omega_{\beta
k}\beta_{k}^{\dagger}\beta_{k}]$ (22)

where $k$ is a tri-dimensional excitation vector of the system, and
$\gamma=g\mu_{B}/\hbar$ is the gyromagnetic ratio, and $\mu_{B}$ is the Bohr
magneton. In the long-wavelength limit, the frequencies of the two magnon
modes appearing in the above Hamiltonian are provided by:

$\displaystyle\omega_{\alpha,0}\approx\gamma\sqrt{2H_{E}H_{Ax}}$ (23)
$\displaystyle\omega_{\beta,0}\approx\gamma\sqrt{2H_{E}H_{Az}}$ (24)

Here we introduced the anisotropy fields of NiO which satisfy:

$H_{E}\gg H_{Ax}\gg H_{Az}$ (25)

Consequently, the highest frequency mode is the $\alpha$ mode which is the one
observed around 1 THz in the experiments described in the main text.

Using the inequalities (25) as well as the long-wavelength approximation we
can use the results from Ref. [54] in order to express the total spin along
the $y$-direction of the $\alpha$-mode, as well as the corresponding
magnetization vector as a function of the bosonic operators
$\alpha_{k}^{\dagger}$ and $\alpha_{k}$:

$M_{y}(\mathbf{r})=-\frac{g\mu_{B}}{2}\Big{(}\frac{H_{Ax}}{H_{E}}\Big{)}^{1/4}\frac{N}{V_{1}}\sum_{k}\frac{e^{i\mathbf{kr}}}{\sqrt{N}}i(\alpha_{k}-\alpha^{\dagger}_{-k})$
(26)

The $\beta$ mode does not contribute to this expression. Here we introduced
the volume $V_{1}=S\Delta_{1}$ of the NiO slab. The wavevectors are normalized
here with respect to the number of spins $N$ in the spin sublattice:

$\sum_{i}e^{i\mathbf{kr}_{i}}=N\delta_{\mathbf{k0}}$ (27)

Using Eq.(26) together with Eq. (6) in order to express the last line in Eq.
(4) we have:

$\displaystyle-\int_{V}\mu_{0}\hat{\mathbf{H}}(\mathbf{r})\hat{\mathbf{M}}(\mathbf{r})d^{3}\mathbf{r}=$
$\displaystyle i\hbar\Omega_{M}(c-c^{\dagger})(a+a^{\dagger})$ (28)
$\displaystyle\hbar\Omega_{M}=\eta_{m}\frac{g\mu_{B}}{2}\sqrt{\frac{\hbar\omega_{m}\mu_{0}N}{2V_{1}}}$
(29)
$\displaystyle\eta_{m}=\sqrt{\int_{\Delta_{1}}h(z)^{2}dz\Big{/}\int_{L}h(z)^{2}dz}$
(30)

Similar results have been obtained by Boventer et al. [69]. In the above
expression, we introduced the collective bosonic operators:

$\displaystyle
c=\frac{1}{\sqrt{\Delta_{1}\int_{\Delta_{1}}h(z)^{2}dz}}\sum_{k}h(z)e^{-ik_{z}z}dz\alpha_{k}$
(31)

that were derived in a very analogous manner to the collective phonon
operators from the previous section. As a consequence, we obtain the
expression of the overlap coefficient between photons and magnons $\eta_{m}$
(Eq. (30)) that is defined in a similar manner as the overlap coefficient
$\eta_{p}$ (Eq. (21)) between photons and phonon in the CBO slab.

###  A.4 Hopfield diagonalization

The full light-matter coupling Hamiltonian now becomes:

$\displaystyle\mathcal{H}=\hbar\omega_{m}a^{\dagger}a+\hbar\omega_{LO}d^{\dagger}d+\hbar\omega_{\alpha
0}c^{\dagger}c$ $\displaystyle+i\hbar\Omega_{P}(a-a^{\dagger})(d+d^{\dagger})$
$\displaystyle+i\hbar\Omega_{M}(a+a^{\dagger})(c-c^{\dagger})$ (32)

This is a bosonic Hamiltonian for three coupled harmonic oscillators, where
photons and magnons are both coupled to the photon of a frequency
$\omega_{m}$. In order to find the coupled modes of the system we will use the
Hopfield-Bogoliubov approach [Ref], where we introduce a mixed polariton
bosonic operator:

$\Pi=xa+ya^{\dagger}+md+td^{\dagger}+nc+sc^{\dagger}$ (33)

The coupled mode frequencies $\omega$ and the boson operators $\Pi$ must
satisfy the commutator relation $[\Pi,\mathcal{H}]=\hbar\omega\Pi$. The latter
can be cast in a matrix form:

$\mathbf{M}\mathbf{V}=\omega\mathbf{V}$ (34)

Here $\mathbf{V}^{T}=(x,y,m,t,n,s)$ is the vector of Hopfield coefficients,
and the matrix $\mathbf{M}$ is provided by:

$\mathbf{M}=\begin{bmatrix}\omega_{m}&0&i\Omega_{P}&-i\Omega_{P}&-i\Omega_{M}&-i\Omega_{M}\\\
0&-\omega_{m}&-i\Omega_{P}&i\Omega_{P}&-i\Omega_{M}&-i\Omega_{M}\\\
-i\Omega_{P}&-i\Omega_{P}&\omega_{LO}&0&0&0\\\
-i\Omega_{P}&-i\Omega_{P}&0&-\omega_{LO}&0&0\\\
i\Omega_{M}&-i\Omega_{M}&0&0&\omega_{\alpha 0}&0\\\
-i\Omega_{M}&i\Omega_{M}&0&0&0&-\omega_{\alpha 0}\end{bmatrix}$ (35)

Furthermore, the eigenvectors must satisfy the normalization condition:

$|x|^{2}-|y|^{2}+|m|^{2}-|t|^{2}+|n|^{2}-|s|^{2}=1$ (36)

which stems from the requirement $[\Pi,\Pi^{\dagger}]=1$. We thus can define
the coefficients that describe respectively the photon, phonon, and magnon
weight in the polariton states:

$\displaystyle W_{photon}=|x|^{2}-|y|^{2}$ (37) $\displaystyle
W_{phonon}=|m|^{2}-|t|^{2}$ (38) $\displaystyle W_{magnon}=|n|^{2}-|s|^{2}$
(39)

The matrix problem stated above can be solved numerically, thus providing
eigenfrequencies and the Hopfield coefficients. These results have been
illustrated in Fig. 5 of the main text. Actually, thanks to the particular
form of the matrix the equation providing the polariton frequencies can be
provided explicitly. This equation, which states
$\mathrm{det}|\mathbf{M}-\omega\mathbf{I}|=0$ is obtained to be:

$\displaystyle\Big{(}\omega_{m}^{2}-\omega^{2}-\frac{4\Omega_{M}^{2}\omega_{\alpha
0}\omega_{m}}{\omega_{\alpha
0}^{2}-\omega^{2}}\Big{)}\Big{(}\omega_{m}^{2}-\omega^{2}-\frac{4\Omega_{P}^{2}\omega_{LO}\omega_{m}}{\omega_{LO}^{2}-\omega^{2}}\Big{)}$
$\displaystyle=\frac{16\Omega_{M}^{2}\Omega_{P}^{2}\omega_{\alpha
0}\omega_{LO}\omega^{2}}{(\omega_{\alpha
0}^{2}-\omega^{2})(\omega_{LO}^{2}-\omega^{2})}$ (40)

This equation is used to fit the experimental data. It is interesting to
consider a case where $\omega_{m}\gg\omega_{\alpha 0},\omega_{LO}$. In that
case the first line of the above equation simplifies to $\omega^{4}$ and we
can write an effective coupled equation:

$(\omega_{\alpha
0}^{2}-\omega^{2})(\omega_{LO}^{2}-\omega^{2})\approx\frac{16\Omega_{M}^{2}\Omega_{P}^{2}\omega_{\alpha
0}\omega_{LO}\omega^{2}}{\omega_{m}^{4}}$ (41)

This equation illustrates that the photon mode introduces and effective
coupling between the phonon and magnon excitation.

###  A.5 Effective dielectric and magnetic susceptibilities

Let us now discuss how the coupled mode equation Eq.(40) leads to effective
susceptibilities. First, we consider the case where waves propagate only in
the CBO slab (we can extend the thickness of this slab to infinity), such as
$\eta_{p}=1$ and $\Omega_{m}=0$. The coupled wave equation thus becomes, using
Eqs. (20) and (40):

$\omega_{m}^{2}-\omega^{2}-\frac{(\omega_{LO}^{2}-\omega_{TO}^{2})\omega_{m}^{2}}{\omega_{LO}^{2}-\omega^{2}}=0$
(42)

This equation can be recast in the form:

$\varepsilon^{CBO}_{\infty}\omega_{m}^{2}=\varepsilon(\omega)\omega^{2}$ (43)

where we introduced the dielectric function $\varepsilon(\omega)$ of the
homogeneous CBO slab provided by:

$\frac{\varepsilon^{CBO}_{\infty}}{\varepsilon(\omega)}=1+\frac{\omega_{LO}^{2}-\omega_{TO}^{2}}{\omega^{LO}-\omega^{2}}$
(44)

This is exactly the dielectric susceptibility for a medium with an optical
phonon. Indeed, by posing
$\omega_{m}=c|\mathbf{Q}|/\sqrt{\varepsilon^{CBO}_{\infty}}$ Eq.(43) becomes
the Helmholtz dispersion relation for waves propagating in a CBO medium with a
wavector $\mathbf{Q}$ [70].

Likewise, we consider an infinite NiO slab, where the Helmholtz relation
becomes:

$\displaystyle\omega_{m}^{2}=\mu(\omega)\omega^{2}$ (45)

By setting $\Omega_{P}=0$ and $\eta_{m}=1$ in Eq.(29) from Eq. (40) we obtain:

$\displaystyle\frac{1}{\mu(\omega)}=1+\frac{W_{A}^{2}}{\omega_{\alpha
0}^{2}-\omega^{2}}$ (46) $\displaystyle W_{A}^{2}=\omega_{\alpha
0}^{2}\frac{n\mu_{0}\mu_{B}g}{2\sqrt{2}H_{E}}$ (47)

Here we introduced the density of spins $n=N/V_{2}$. Clearly, the above
expressions establish the magnetic susceptibility of the NiO material as a
function of microscopic quantities.

We know the material susceptibilities from Eqs. (44) and (46) together with
the full polariton equation Eq.(40) in order to define effective
susceptibilites. The latter is defined such as Eq.(40) can be rewritten in the
form:

$\varepsilon^{CBO}_{\infty}\omega_{m}^{2}=\varepsilon_{eff}(\omega)\mu_{eff}(\omega)\omega^{2}$
(48)

After some algebra the following functions are obtained:

$\displaystyle\frac{1}{\mu_{eff}(\omega)}=1+\eta_{m}^{2}\frac{W_{A}^{2}}{\omega_{\alpha
0}^{2}-\omega^{2}}$ (49)
$\displaystyle\frac{1}{\varepsilon_{eff}(\omega)}=\frac{1}{\varepsilon^{CBO}_{\infty}}\Big{[}1+\eta_{p}^{2}\frac{(\omega_{LO}^{2}-\omega_{TO}^{2})\omega_{m}^{2}}{\omega_{LO}^{2}-\omega^{2}}\Big{]}$
(50)

The definition of $\mu_{eff}(\omega)$ is acceptable, as the only magnetic
contribution arises from the NiO layer. However, the definition of
$\varepsilon_{eff}(\omega)$ is problematic, as even after the removal of the
CBO layer by setting $\eta_{p}=0$ we still have a medium with a dielectric
constant $\varepsilon^{CBO}_{\infty}$, which ignores the presence of NiO and
air dielectric layers. The following definition solves this problem:

$\displaystyle\frac{1}{\varepsilon_{eff}(\omega)}=\frac{1}{\bar{\varepsilon}_{\infty}}\Big{[}1+\bar{\eta}_{p}^{2}\frac{(\omega_{LO}^{2}-\omega_{TO}^{2})\omega_{m}^{2}}{\omega_{LO}^{2}-\omega^{2}}\Big{]}$
(51)
$\displaystyle\frac{1}{\bar{\varepsilon}_{\infty}}=1-\eta_{p}^{2}-\eta_{m}^{2}+\frac{\eta_{m}^{2}}{\varepsilon_{NiO}^{\infty}}+\frac{\eta_{p}^{2}}{\varepsilon^{CBO}_{\infty}}$
(52)
$\displaystyle\bar{\eta}_{p}^{2}=\eta_{p}^{2}\frac{\bar{\varepsilon}_{\infty}}{\varepsilon^{CBO}_{\infty}}$
(53)

It is thus sufficient to use the re-normalized value $\bar{\eta}_{p}$ for the
overlap factor in the light-phonon coupling constant defined in Eq.(20)
instead of $\eta_{p}$ in order to obtain a polariton dispersion relation that
takes into account all layers in the system.

##  Appendix B Dielectric characterization of samples

We measured transmission through each of our samples placed on a metal
aperture, normalized to transmission through an empty aperture (Fig. 9). We
fitted the results with a single-layer model and obtained background
dielectric constants of these materials and a resonant form of permeability in
the NiO crystal and phonon response in the CBO ceramics. In the case of the
CBO phonon, we assumed Lorentzian lineshape; however, this model (Fig. 9bc)
does not reproduce perfectly measured transmission around the resonance (Fig.
9a). This may be caused by other than the Lorentzian lineshape of the phonon
mode in CBO. The parameters obtained in this way were used to calculate the
predictions of electromagnetic wave reflection from structures composed of the
CBO, NiO, and tunable air gap.

Figure 9:  Transmission spectra of a standalone CBO with fits (a), under the
assumption of a dielectric function with real (b) and imaginary (c) parts,
permeability was assumed to be $1$. Transmission spectra of stand-alone NiO
with fitted transmission (d), under the assumption of permeability with real
(e) and imaginary (f) parts, permittivity was assumed constant.

Thus, for our calculations, we can assume the isotropic dielectric constants
of NiO and CBO, and $\epsilon_{air}=1=\mu_{air}$ for dry air. The details of
the calculation methods are presented in Ref. 40. Here, we take into account
the antiferromagnetic resonance in the permeability in NiO by writing

$\mu=1+\frac{\Delta\mu f_{m}^{2}}{f_{m}^{2}-f^{2}-ifg_{m}},$ (54)

where $f_{m}$ is magnon frequency, $g_{m}=8$ GHz is its width and its strength
$\Delta\mu=4.2\times 10^{-4}$ as determined in our transmission results (Fig.
9). The permittivity of NiO was $\epsilon_{Nio}=11.63+i0.02$. The resonance
strength $\Delta\mu$ is responsible for the strength of the coupling of magnon
with FP modes. We take the CBO layer as a dielectric resonator

$\epsilon=\epsilon_{\infty}+\frac{\Delta\epsilon
f_{p}^{2}}{f_{p}^{2}-f^{2}-ifg_{p}},$ (55)

where $f_{p}=0.92$ THz is the phonon frequency, $g_{p}=50$ GHz is its width,
its strength $\Delta\epsilon=0.021$, and $\epsilon_{\infty}=5.52+i0.05$, as
determined in our transmission results (Fig. 9). We assume $\mu_{CBO}=1$,
neglecting that CBO is a paramagnet at room temperature [71]. In this
electrodynamic model, we take into account the reflection from the metal
mirror as a 4th layer of 100 nm thickness characterized by very high real and
imaginary parts of its dielectric constant ($8\times 10^{8}+i10^{3}$).

##  Appendix C Phase

Figure 10: Evolution of a cavity mode near transition point from positive to
negative phase jump. (a) amplitude and (b) phase spectra for different gaps
between the slabs. Phasor plots of typical negative phase jumps encircling the
origin (c), near transition condition (d), and with positive phase jumps not
encircling the origin (e).

Here, we show in detail the evolution of a cavity mode near the transition
point between positive and negative phase jump at around $f=1.04$ THz and
around gap $d=400$ $\mu$m. We show in Fig. 10a standard reflection amplitude.
At larger distances of 407 and 398 $\mu$m, the cavity has a typical negative
phase jumps sign (Fig. 10b, in the direction of normal phase change with
rising frequency. In Fig. 10c, we present phasor plots for these two
distances, which show that in this case, the signal circles close to the
origin, reflecting the signal’s lower amplitude at the cavity mode. For the
case of $d=386$ $\mu$m, the transition is visible when cavity amplitude
reaches almost zero and phase jump becomes undefined. Phasor trajectory for
this case (Fig. 10e) passes very close to the origin. When we further
decreased the gap between the crystals to $d=371$ and 356 $\mu$m, the mode
amplitude dropped again, but its phase jumps flipped sign to positive. This
can be understood in phasor plot Fig. 10e as the signal no longer circling
around the origin, in a much smaller loop.

Since the data in Fig. 10 is an extract of the data presented in Fig. 3, we
notice that the condition for the change of the phase jump sign is related to
the frequency of the cavity mode, since they appear periodically for each
cavity modes, when they are at the same frequency range of about 1.05 THz.
Therefore this effect must be related to interference conditions, and as shown
in Fig. 10 it is not related to topological changes in the electric field
phasor.

Figure 11: Evolution of phonon-polariton phasor from weak to strong coupling.
(a) amplitude and (b) phase spectra for different gaps between the slabs. In
the weak coupling regime, the phasor shows a single loop (c), while in the
strong coupling regime, a topological change of the phasor into a double loop
is observed (d).

The change of mode phase jump also appears near the polariton formation
regions. However, in this case not only a phase sign change is observed, but
also a topological change in the phasor. In Fig. 11, we show in detail the
evolution of phonon-polariton formation at $T=300$ K, when the magnon is away
from the interaction region. Magnitude and phase spectra were collected for
different gaps between the slabs, which tuned the cavity mode frequency. The
strongest coupling was observed at $d=419$ $\mu$m when the cavity mode reached
the phonon frequency. In Fig. 11c, we show phasors of a weakly interacting
cavity mode, which forms a loop around the origin. The magnon can is
represented as a weak deformation of the loop for higher $d$. In Fig. 11d, we
show that a single loop is transformed into a double loop when entering the
strong coupling regime, which is a topological change of the phasor.

##  Appendix D Magnon and phonon contents using electromagnetism

Classical electrodynamics predicts some distant coupling between the
resonators, which can be observed in the reflection or transmission of
electromagnetic waves. Using our classical model, we determined the relative
phonon and magnon contents of polariton modes. We did that by calculating the
reflection amplitude derivatives to the small magnon strength change
$\Delta\mu$ (1. 0%) and the phonon strength change $\Delta\epsilon$ (0.1%).
Changes that give similar amplitude of reflection changes.

Figure 12:  Expected magnon and phonon content, for two different gaps, one
with cavity mode matching the phonon mode (a-c) and without matching (d-f).
Derivative to magnon strength (a,d) shows the magnon content of cavity modes,
and derivative to phonon strength (b,e) shows phonon content. The second
derivative shows modes coupled simultaneously to the magnon and phonon.

One sees that the magnon is visible in derivative to its strength (Fig. 12a)
around the diagonal, which is magnon frequency in our model. Signatures of the
magnon in reflection are either positive or negative, depending on the phase
of its interaction with cavity modes. Phonon is uncovered in around its
frequency in derivative to its strength (Fig. 12b). Here, the shape of the
phonon derivative spectrum is caused by interaction with a tuned cavity mode,
so that enhancing phonon strength increases the coupling rate, which results
in polariton modes more separated, which is visible as higher reflection (red-
yellow) in between the polaritons and lower amplitude outside (blue-black). In
this phonon derivative one already sees a trace of magnon interactions, in the
form of splittings of upper and lower phonon-polaritons. In particular, for
$f_{p}-f_{m}=0$, the splitting of the polariton modes is slightly larger, as
shown in Fig. 4 and 7, which is an indication of coupling of a cavity mode
coherently with magnon and phonon.

However, such a coherent coupling is not the only possibility, the model shows
(Fig. 12(d-f)) that even for modes not matching, there may be conditions in
which they interact with both the magnon and the phonon. This comes from the
fact that the phonon is such a strong oscillator that it couples even with
cavity modes detuning as much as about 50 GHz, or more [48]. In Fig. 12(d-f),
one sees magnon-polariton modes that are simultaneously coupled to the phonon;
however, due to a huge detuning, the coupling of these modes to the phonon is
weaker. Such modes give larger differential signals, which is caused by their
narrow line width, which in turn originates from their smaller phonon content
than in the case of modes presented in Fig. 12(a-c). We observed such modes
experimentally in Fig. 6, below 0.9 THz. However, in the experiment, we cannot
determine the phonon content of these modes, so we cannot show that they are
indeed coupled to the phonon, in another way than in the theoretical model.

##  Appendix E Phonon-magnon coupling at large distances

Figure 13:  Magnitude (a) and phase differential (b) are measured as a
function of the gap between the slabs with a step of 163 $\mu$m, which always
gives a cavity mode as the phonon frequency. Since the results were obtained
for $T=364$ K, there was also the magnon mode at the same frequency, leading
to the formation of phonon-magnon polaritons. The closing of the polariton gap
is visible with increasing distance $d$, up to about 4-5 mm.

Our mechanical systems allowed to control the gap between the sample slabs in
the range of tens of mm with $\mu$m precision. The focal point of the quasi-
optical systems was very elongated because we used a parabolic mirror of a
focal length of 101 mm. Therefore, we can widen the gap between the crystals
up to a few mm without any of the samples leaving the focal point and
remaining with high position accuracy. This allowed us to measure reflection
for the slabs separated by large gaps $d$, which were so chosen that for each
gap a cavity mode had a frequency at the phonon frequency $f_{p}$. The step
was 163 $\mu$m, which is half a wavelength of $f_{p}$. This means we were
stepping through subsequent cavity modes so that for each of these distances,
there was a cavity mode at $f_{p}$, which created conditions for cavity-
mediated phonon-magnon coupling. With rising distance the cavity volume
$V_{c}$ increases, but the number of oscillators (phonons and magnons) remains
constant. Since the coupling strength is proportional to $\sqrt{N/V_{c}}$,
this explains the drop in splitting with the rising gap between the slabs. The
rising gap finally leads to a collapse of the splitting, e.i. the transition
from the strong coupling regime (polaritons) to the weak coupling regime
(Purcell effect).

## References

  * Khitrova  _et al._ [2006] G. Khitrova, H. M. Gibbs, M. Kira, S. W. Koch, and A. Scherer, Vacuum rabi splitting in semiconductors, Nature Physics 2, 81 (2006).
  * Dovzhenko  _et al._ [2018] D. S. Dovzhenko, S. V. Ryabchuk, Y. P. Rakovich, and I. R. Nabiev, Light–matter interaction in the strong coupling regime: configurations, conditions, and applications, Nanoscale 10, 3589 (2018).
  * Rempe  _et al._ [1987] G. Rempe, H. Walther, and N. Klein, Observation of quantum collapse and revival in a one-atom maser, Phys. Rev. Lett. 58, 353 (1987).
  * Kasprzak  _et al._ [2006] J. Kasprzak, M. Richard, S. Kundermann, A. Baas, P. Jeambrun, J. M. J. Keeling, F. M. Marchetti, M. H. Szymańska, R. André, J. L. Staehli, V. Savona, P. B. Littlewood, B. Deveaud, and L. S. Dang, Bose–einstein condensation of exciton polaritons, Nature 443, 409 (2006).
  * Raizen  _et al._ [1989] M. G. Raizen, R. J. Thompson, R. J. Brecha, H. J. Kimble, and H. J. Carmichael, Normal-mode splitting and linewidth averaging for two-state atoms in an optical cavity, Phys. Rev. Lett. 63, 240 (1989).
  * Colombe  _et al._ [2007] Y. Colombe, T. Steinmetz, G. Dubois, F. Linke, D. Hunger, and J. Reichel, Strong atom–field coupling for bose–einstein condensates in an optical cavity on a chip, Nature 450, 272 (2007).
  * Delteil  _et al._ [2012] A. Delteil, A. Vasanelli, Y. Todorov, C. Feuillet Palma, M. Renaudat St-Jean, G. Beaudoin, I. Sagnes, and C. Sirtori, Charge-induced coherence between intersubband plasmons in a quantum structure, Phys. Rev. Lett. 109, 246808 (2012).
  * Törmä and Barnes [2014] P. Törmä and W. L. Barnes, Strong coupling between surface plasmon polaritons and emitters: a review, Reports on Progress in Physics 78, 013901 (2014).
  * Basov  _et al._ [2016] D. N. Basov, M. M. Fogler, and F. J. G. de Abajo, Polaritons in van der Waals materials, Science 354, aag1992 (2016).
  * Bayer  _et al._ [2017] A. Bayer, M. Pozimski, S. Schambeck, D. Schuh, R. Huber, D. Bougeard, and C. Lange, Terahertz light–matter interaction beyond unity coupling strength, Nano Letters 17, 6340 (2017).
  * Li  _et al._ [2018] X. Li, M. Bamba, N. Yuan, Q. Zhang, Y. Zhao, M. Xiang, K. Xu, Z. Jin, W. Ren, G. Ma, S. Cao, D. Turchinovich, and J. Kono, Observation of Dicke cooperativity in magnetic interactions, Science 361, 794 (2018).
  * Yahiaoui  _et al._ [2022] R. Yahiaoui, Z. A. Chase, C. Kyaw, F. Tay, A. Baydin, G. T. Noe, J. Song, J. Kono, A. Agrawal, M. Bamba, and T. A. Searles, Dicke-cooperativity-assisted ultrastrong coupling enhancement in terahertz metasurfaces, Nano Letters 10.1021/acs.nanolett.2c01892 (2022).
  * Roberts  _et al._ [1962] R. W. Roberts, B. A. Auld, and R. R. Schell, Magnetodynamic mode ferrite amplifier, Journal of Applied Physics 33, 1267 (1962).
  * Schuster  _et al._ [2010] D. I. Schuster, A. P. Sears, E. Ginossar, L. DiCarlo, L. Frunzio, J. J. L. Morton, H. Wu, G. A. D. Briggs, B. B. Buckley, D. D. Awschalom, and R. J. Schoelkopf, High-cooperativity coupling of electron-spin ensembles to superconducting cavities, Phys. Rev. Lett. 105, 140501 (2010).
  * Abe  _et al._ [2011] E. Abe, H. Wu, A. Ardavan, and J. J. L. Morton, Electron spin ensemble strongly coupled to a three-dimensional microwave cavity, Applied Physics Letters 98, 251108 (2011).
  * Huebl  _et al._ [2013] H. Huebl, C. W. Zollitsch, J. Lotze, F. Hocke, M. Greifenstein, A. Marx, R. Gross, and S. T. B. Goennenwein, High cooperativity in coupled microwave resonator ferrimagnetic insulator hybrids, Phys. Rev. Lett. 111, 127003 (2013).
  * Zhang  _et al._ [2014] X. Zhang, C.-L. Zou, L. Jiang, and H. X. Tang, Strongly coupled magnons and cavity microwave photons, Phys. Rev. Lett. 113, 156401 (2014).
  * Tabuchi  _et al._ [2014] Y. Tabuchi, S. Ishino, T. Ishikawa, R. Yamazaki, K. Usami, and Y. Nakamura, Hybridizing ferromagnetic magnons and microwave photons in the quantum limit, Phys. Rev. Lett. 113, 083603 (2014).
  * Tabuchi  _et al._ [2015] Y. Tabuchi, S. Ishino, A. Noguchi, T. Ishikawa, R. Yamazaki, K. Usami, and Y. Nakamura, Coherent coupling between a ferromagnetic magnon and a superconducting qubit, Science 349, 405 (2015).
  * Zhang  _et al._ [2015] X. Zhang, C.-L. Zou, N. Zhu, F. Marquardt, L. Jiang, and H. X. Tang, Magnon dark modes and gradient memory, Nature Communications 6, 8914 (2015).
  * Zhang  _et al._ [2016] X. Zhang, N. Zhu, C.-L. Zou, and H. X. Tang, Optomagnonic whispering gallery microresonators, Phys. Rev. Lett. 117, 123605 (2016).
  * Li  _et al._ [2019] Y. Li, T. Polakovic, Y.-L. Wang, J. Xu, S. Lendinez, Z. Zhang, J. Ding, T. Khaire, H. Saglam, R. Divan, J. Pearson, W.-K. Kwok, Z. Xiao, V. Novosad, A. Hoffmann, and W. Zhang, Strong coupling between magnons and microwave photons in on-chip ferromagnet-superconductor thin-film devices, Phys. Rev. Lett. 123, 107701 (2019).
  * Everts  _et al._ [2020] J. R. Everts, G. G. G. King, N. J. Lambert, S. Kocsis, S. Rogge, and J. J. Longdell, Ultrastrong coupling between a microwave resonator and antiferromagnetic resonances of rare-earth ion spins, Phys. Rev. B 101, 214414 (2020).
  * Potts and Davis [2020] C. A. Potts and J. P. Davis, Strong magnon–photon coupling within a tunable cryogenic microwave cavity, Applied Physics Letters 116, 263503 (2020).
  * Lachance-Quirion  _et al._ [2020] D. Lachance-Quirion, S. P. Wolski, Y. Tabuchi, S. Kono, K. Usami, and Y. Nakamura, Entanglement-based single-shot detection of a single magnon with a superconducting qubit, Science 367, 425 (2020).
  * Li  _et al._ [2020a] Y. Li, W. Zhang, V. Tyberkevych, W.-K. Kwok, A. Hoffmann, and V. Novosad, Hybrid magnonics: Physics, circuits, and applications for coherent information processing, Journal of Applied Physics 128, 130902 (2020a).
  * Bhoi  _et al._ [2021] B. Bhoi, S.-H. Jang, B. Kim, and S.-K. Kim, Broadband photon–magnon coupling using arrays of photon resonators, Journal of Applied Physics 129, 083904 (2021).
  * Awschalom and Flatté [2007] D. D. Awschalom and M. E. Flatté, Challenges for semiconductor spintronics, Nature Physics 3, 153 (2007).
  * Kockum  _et al._ [2019] A. F. Kockum, A. Miranowicz, S. De Liberato, S. Savasta, and F. Nori, Ultrastrong coupling between light and matter, Nature Reviews Physics 1, 19 (2019).
  * Roux  _et al._ [2020] K. Roux, H. Konishi, V. Helson, and J.-P. Brantut, Strongly correlated fermions strongly coupled to light, Nature Communications 11, 2974 (2020).
  * Yuan  _et al._ [2022] H. Yuan, Y. Cao, A. Kamra, R. A. Duine, and P. Yan, Quantum magnonics: When magnon spintronics meets quantum information science, Physics Reports 965, 1 (2022).
  * Spethmann  _et al._ [2016] N. Spethmann, J. Kohler, S. Schreppler, L. Buchmann, and D. M. Stamper-Kurn, Cavity-mediated coupling of mechanical oscillators limited by quantum back-action, Nature Physics 12, 27 (2016).
  * Ściesiek  _et al._ [2020] M. Ściesiek, K. Sawicki, W. Pacuski, K. Sobczak, T. Kazimierczuk, A. Golnik, and J. Suffczyński, Long-distance coupling and energy transfer between exciton states in magnetically controlled microcavities, Communications Materials 1, 78 (2020).
  * Xu  _et al._ [2019] P.-C. Xu, J. W. Rao, Y. S. Gui, X. Jin, and C.-M. Hu, Cavity-mediated dissipative coupling of distant magnetic moments: Theory and experiment, Phys. Rev. B 100, 094415 (2019).
  * Li  _et al._ [2022] Y. Li, V. G. Yefremenko, M. Lisovenko, C. Trevillian, T. Polakovic, T. W. Cecil, P. S. Barry, J. Pearson, R. Divan, V. Tyberkevych, C. L. Chang, U. Welp, W.-K. Kwok, and V. Novosad, Coherent coupling of two remote magnonic resonators mediated by superconducting circuits, Phys. Rev. Lett. 128, 047701 (2022).
  * Wang  _et al._ [2022] H. Wang, W. He, R. Yuan, Y. Wang, J. Wang, Y. Zhang, I. Medlej, J. Chen, G. Yu, X. Han, J.-P. Ansermet, and H. Yu, Hybridized propagating spin waves in a cofeb/irmn bilayer, Phys. Rev. B 106, 064410 (2022).
  * Harvey-Collard  _et al._ [2022] P. Harvey-Collard, J. Dijkema, G. Zheng, A. Sammak, G. Scappucci, and L. M. K. Vandersypen, Coherent spin-spin coupling mediated by virtual microwave photons, Phys. Rev. X 12, 021026 (2022).
  * Nair  _et al._ [2022] J. M. P. Nair, D. Mukhopadhyay, and G. S. Agarwal, Cavity-mediated level attraction and repulsion between magnons, Phys. Rev. B 105, 214418 (2022).
  * Yang  _et al._ [2022] Z.-B. Yang, H.-Y. Liu, and R.-C. Yang, Asymmetric quantum synchronization generation in antiferromagnet-cavity systems, The European Physical Journal Plus 137, 878 (2022).
  * Białek  _et al._ [2023] M. Białek, W. Knap, and J.-P. Ansermet, Cavity-mediated coupling of terahertz antiferromagnetic resonators, Phys. Rev. Appl. 19, 064007 (2023).
  * Berk  _et al._ [2019] C. Berk, M. Jaris, W. Yang, S. Dhuey, S. Cabrini, and H. Schmidt, Strongly coupled magnon–phonon dynamics in a single nanomagnet, Nature Communications 10, 2652 (2019).
  * Li  _et al._ [2021] Y. Li, C. Zhao, W. Zhang, A. Hoffmann, and V. Novosad, Advances in coherent coupling between magnons and acoustic phonons, APL Materials 9, 060902 (2021), https://pubs.aip.org/aip/apm/article-pdf/doi/10.1063/5.0047054/14565679/060902_1_online.pdf .
  * Khan  _et al._ [2020] P. Khan, M. Kanamaru, K. Matsumoto, T. Ito, and T. Satoh, Ultrafast light-driven simultaneous excitation of coherent terahertz magnons and phonons in multiferroic $\mathrm{BiFe}{\mathrm{o}}_{3}$, Phys. Rev. B 101, 134413 (2020).
  * Li  _et al._ [2020b] J. Li, H. T. Simensen, D. Reitz, Q. Sun, W. Yuan, C. Li, Y. Tserkovnyak, A. Brataas, and J. Shi, Observation of magnon polarons in a uniaxial antiferromagnetic insulator, Phys. Rev. Lett. 125, 217201 (2020b).
  * Liu  _et al._ [2021] S. Liu, A. Granados del Águila, D. Bhowmick, C. K. Gan, T. Thu Ha Do, M. A. Prosnikov, D. Sedmidubský, Z. Sofer, P. C. M. Christianen, P. Sengupta, and Q. Xiong, Direct observation of magnon-phonon strong coupling in two-dimensional antiferromagnet at high magnetic fields, Phys. Rev. Lett. 127, 097401 (2021).
  * Vaclavkova  _et al._ [2021] D. Vaclavkova, M. Palit, J. Wyzula, S. Ghosh, A. Delhomme, S. Maity, P. Kapuscinski, A. Ghosh, M. Veis, M. Grzeszczyk, C. Faugeras, M. Orlita, S. Datta, and M. Potemski, Magnon polarons in the van der waals antiferromagnet $\mathrm{Fe}{\mathrm{ps}}_{3}$, Phys. Rev. B 104, 134437 (2021).
  * Diederich  _et al._ [2022] G. M. Diederich, J. Cenker, Y. Ren, J. Fonseca, D. G. Chica, Y. J. Bae, X. Zhu, X. Roy, T. Cao, D. Xiao, and X. Xu, Tunable interaction between excitons and hybridized magnons in a layered semiconductor, Nature Nanotechnology 10.1038/s41565-022-01259-1 (2022).
  * Sivarajah  _et al._ [2019] P. Sivarajah, A. Steinbacher, B. Dastrup, J. Lu, M. Xiang, W. Ren, S. Kamba, S. Cao, and K. A. Nelson, THz-frequency magnon-phonon-polaritons in the collective strong-coupling regime, Journal of Applied Physics 125, 213103 (2019).
  * Marsh  _et al._ [2019] D. J. E. Marsh, K. C. Fong, E. W. Lentz, L. Šmejkal, and M. N. Ali, Proposal to detect dark matter using axionic topological antiferromagnets, Phys. Rev. Lett. 123, 121601 (2019).
  * Mitridate  _et al._ [2020] A. Mitridate, T. Trickle, Z. Zhang, and K. M. Zurek, Detectability of axion dark matter with phonon polaritons and magnons, Phys. Rev. D 102, 095005 (2020).
  * Wang  _et al._ [2011] J. Wang, R. Li, S.-C. Zhang, and X.-L. Qi, Topological magnetic insulators with corundum structure, Phys. Rev. Lett. 106, 126403 (2011).
  * Karaki  _et al._ [2022] M. J. Karaki, X. Yang, A. J. Williams, M. Nawwar, V. Doan-Nguyen, J. E. Goldberger, and Y.-M. Lu, An efficient material search for room temperature topological magnons (2022).
  * Marsh  _et al._ [2022] D. J. E. Marsh, J. I. McDonald, A. J. Millar, and J. Schütte-Engel, Axion detection with phonon-polaritons revisited (2022).
  * Rezende  _et al._ [2019] S. M. Rezende, A. Azevedo, and R. L. Rodríguez-Suárez, Introduction to antiferromagnetic magnons, Journal of Applied Physics 126, 151101 (2019), https://pubs.aip.org/aip/jap/article-pdf/doi/10.1063/1.5109132/13019547/151101_1_online.pdf .
  * Szwagierczak  _et al._ [2021] D. Szwagierczak, B. Synkiewicz-Musialska, J. Kulawik, and N. Pałka, Sintering, microstructure, and dielectric properties of copper borates for high frequency ltcc applications, Materials 14, 10.3390/ma14144017 (2021).
  * Thomas  _et al._ [2020] P. A. Thomas, W. J. Tan, H. A. Fernandez, and W. L. Barnes, A new signature for strong light–matter coupling using spectroscopic ellipsometry, Nano Letters 20, 6412 (2020).
  * Roth [1958] W. L. Roth, Magnetic structures of mno, feo, coo, and nio, Phys. Rev. 110, 1333 (1958).
  * Massarotti  _et al._ [1991] V. Massarotti, D. Capsoni, V. Berbenni, R. Riccardi, A. Marini, and E. Antolini, Structural characterization of nickel oxide, Zeitschrift für Naturforschung A 46, 503 (1991).
  * Born and Wolf [1980] M. Born and E. Wolf, Chapter 1 - basic properties of the electromagnetic field, in  _Principles of Optics (Sixth Edition)_, edited by M. Born and E. Wolf (Pergamon, 1980) sixth edition ed., pp. 1 – 70.
  * Todorov [2015] Y. Todorov, Dipolar quantum electrodynamics of the two-dimensional electron gas, Phys. Rev. B 91, 125409 (2015).
  * Ghosh  _et al._ [2021] A. Ghosh, M. Palit, S. Maity, V. Dwij, S. Rana, and S. Datta, Spin-phonon coupling and magnon scattering in few-layer antiferromagnetic ${\mathrm{FePS}}_{3}$, Phys. Rev. B 103, 064431 (2021).
  * Todorov  _et al._ [2010] Y. Todorov, A. M. Andrews, R. Colombelli, S. De Liberato, C. Ciuti, P. Klang, G. Strasser, and C. Sirtori, Ultrastrong light-matter coupling regime with polariton dots, Phys. Rev. Lett. 105, 196402 ((2010)).
  * Jeannin  _et al._ [2019] M. Jeannin, G. Mariotti Nesurini, S. Suffit, D. Gacemi, A. Vasanelli, L. Li, A. G. Davies, E. Linfield, C. Sirtori, and Y. Todorov, Ultrastrong light–matter coupling in deeply subwavelength thz lc resonators, ACS Photonics 6, 1207 (2019), https://doi.org/10.1021/acsphotonics.8b01778 .
  * Sanchez-Manzano  _et al._ [2022] D. Sanchez-Manzano, S. Mesoraca, F. A. Cuellar, M. Cabero, V. Rouco, G. Orfila, X. Palermo, A. Balan, L. Marcano, A. Sander, M. Rocci, J. Garcia-Barriocanal, F. Gallego, J. Tornos, A. Rivera, F. Mompean, M. Garcia-Hernandez, J. M. Gonzalez-Calbet, C. Leon, S. Valencia, C. Feuillet-Palma, N. Bergeal, A. I. Buzdin, J. Lesueur, J. E. Villegas, and J. Santamaria, Extremely long-range, high-temperature josephson coupling across a half-metallic ferromagnet, Nature Materials 21, 188 (2022).
  * Todorov  _et al._ [2024] Y. Todorov, S. Dhillon, and J. Mangeney, THz quantum gap: exploring potential approaches for generating and detecting non-classical states of THz light, Nanophotonics 13, 1681 (2024).
  * Frisk Kockum  _et al._ [2019] A. Frisk Kockum, A. Miranowicz, S. De Liberato, S. Savasta, and F. Nori, Ultrastrong coupling between light and matter, Nature Reviews Physics 1, 19–40 (2019).
  * Cohen-Tannoudji  _et al._ [2001] C. Cohen-Tannoudji, J. Dupont-Roc, and G. Grynberg,  _Photons et Atomes_ (EDP Sciences/CNRS Editions, 2001).
  * Rousseaux  _et al._ [2023] B. Rousseaux, Y. Todorov, A. Vasanelli, and C. Sirtori, Phonon-mediated dark to bright plasmon conversion, Phys. Rev. B 108, 125417 (2023).
  * Boventer  _et al._ [2023] I. Boventer, H. T. Simensen, B. Brekke, M. Weides, A. Anane, M. Kläui, A. Brataas, and R. Lebrun, Antiferromagnetic cavity magnon polaritons in collinear and canted phases of hematite, Phys. Rev. Appl. 19, 014071 (2023).
  * Todorov [2014] Y. Todorov, Dipolar quantum electrodynamics theory of the three-dimensional electron gas, Phys. Rev. B 89, 075115 (2014).
  * Boehm  _et al._ [2003] M. Boehm, B. Roessli, J. Schefer, A. S. Wills, B. Ouladdiaf, E. Lelièvre-Berna, U. Staub, and G. A. Petrakovskii, Complex magnetic ground state of ${\mathrm{cub}}_{2}{\mathrm{o}}_{4}$, Phys. Rev. B 68, 024405 (2003).