Exceptional Point and Cross-Relaxation Effect in a Hybrid Quantum System

Exceptional points (EPs) are exotic degeneracies of non-Hermitian systems, where the eigenvalues and the corresponding eigenvectors simultaneously coalesce in parameter space, and these degeneracies are sensitive to tiny perturbations on the system. Here we report an experimental observation of the EP in a hybrid quantum system consisting of dense nitrogen (P1) centers in diamond coupled to a coplanar-waveguide resonator. These P1 centers can be divided into three subensembles of spins, and cross relaxation occurs among them. As a new method to demonstrate this EP, we pump a given spin subensemble with a drive field to tune the magnon-photon coupling in a wide range. We observe the EP in the middle spin subensemble coupled to the resonator mode, irrespective of which spin subensemble is actually driven. This robustness of the EP against pumping reveals the key role of the cross relaxation in P1 centers. It offers a novel way to convincingly prove the existence of the cross-relaxation effect via the EP.

For a small yttrium-iron-garnet (YIG) sphere in a threedimensional cavity [28], one can vary the magnon-photon coupling by moving the sphere in the cavity. However, this mechanical tunability is not precise, which limits its applications in quantum technologies. Also, Ref. [38] has proposed to control the interaction between magnons and photons via a tunable microwave cavity. Actually, there are various hybrid systems with a spin ensemble fixed on a coplanar-waveguide resonator [39][40][41][42][43][44][45]. For these on-chip systems, one also needs to explore a method to precisely tune the magnon-photon coupling. In this work, we experimentally investigate a hybrid quantum system consisting of dense nitrogen (P1) centers in diamond coupled to a coplanarwaveguide resonator. Different from the ferromagnetic spin ensemble in a YIG sphere [46][47][48][49][50], the P1 centers constitute a paramagnetic spin ensemble. We find that we can excite a large number of magnons with a drive field to precisely tune the magnon-photon coupling in a wide range and make it feasible to demonstrate the EP in this on-chip hybrid system. The term magnon is often used to describe collective excitations in a magnetically ordered system (e.g., a ferromagnetic system) with strong spinexchange interactions. In the long-wavelength limit, the effect of the exchange interactions can be ignored and the dipolar interactions become dominant [51]. This limiting case is analogous to paramagnetic systems where the dipolar interactions play the main role. In those paramagnetic systems, magnons can also be used to characterize the collective spin excitations (see, e.g., Refs. [52][53][54]). In the same manner, here we characterize the collective excitations of P1 centers by harnessing magnons.
The P1 centers in diamond can be divided into three spin subensembles and cross relaxation can occur among arXiv:2104.09811v1 [quant-ph] 20 Apr 2021 them [42,43,[55][56][57][58]. With the magnon frequency of the middle subensemble tuned in resonance with the resonator mode, we observe an EP by driving any of the three spin subensembles with a microwave field. Related to the middle subensemble, the other two subensembles are in the dispersive regime and each produces a frequency shift to the resonator mode. These two frequency shifts depend on the magnon occupations and can cancel each other when the magnon occupations in the other two subensembles become equal.
Here the EP is observed by pumping any of the three spin subensembles, indicating that the same magnon occupation is achieved for each subensemble, due to the cross-relaxation effect.
Also, it is known that an EP can enhance the sensitivity to a perturbation [10][11][12][13][14][15][16], while our observed transmission spectrum does not show any appreciable differences at the EP when pumping any of the three spin subensembles. This robustness of the EP against driving reveals the key role of the cross relaxation, which induces the same magnon occupation in each spin subensemble. To the best of our knowledge, this is the first experimental observation of the EP in this hybrid quantum system and it could potentially push the field of these hybrid quantum systems in a new direction. Moreover, it also offers a novel way of convincingly proving the existence of the cross-relaxation effect via the EP. In contrast to the cavity magnonics system consisting of a ferromagnetic spin ensemble coupled to a three-dimensional microwave cavity [21][22][23], this on-chip hybrid system has a more compact configuration and significantly reduces the system size.

II. EXCEPTIONAL POINT OF THE HYBRID QUANTUM SYSTEM
The hybrid system that we study is composed of P1 centers in diamond glued on a coplanar-waveguide resonator (see Fig. 1 and Appendix A). The [001] crystal axis of the diamond is perpendicular to the surface of the resonator and the static magnetic field B is applied along the [100] crystal axis. The P1 centers, which are the main defects in type-1b diamond synthesized under both high pressure and high temperature, are formed by substituting some carbon atoms with nitrogen atoms. These P1 centers can be divided into three subensembles of spins (s = 0, ±) with transition frequencies ω 0 = γ e B and ω ± = γ e B ± A , where γ e /2π = 28 GHz/T is the gyromagnetic ratio and A /2π = 94 MHz is due to the hyperfine interaction in P1 centers.
The collective spin excitations in each subensemble are magnons. As shown in the measured transmission spectrum [ Fig. 2(a)], there are three appreciable anticrossings, which reveal the formation of magnon polaritons by the strong coupling between the magnons of s = 0, ± subensembles and the photons in the superconducting coplanar-waveguide resonator.
For clarity, we first focus on the interaction between the where a and b are annihilation operators of the photons and magnons, with frequencies ω c and ω 0 , respectively. In ferromagnetic materials, spins are very polarized and less magnons are excited, so the effective magnon-photon coupling is often approximated as a constant [47][48][49][50] g eff ≈ g. Here, the P1 centers are paramagnetic and considerable numbers of magnons can be easily excited. Therefore, we need to consider the dependence of the magnon-photon coupling on the occupation number of magnons, where N is the number of P1 centers in the sample. Using the input-output theory [59], we can derive the transmission amplitude of the hybrid system, where γ is the damping rate of the magnon mode, κ i(o) is the decay rate of the resonator mode due to the input (output) port, and κ = κ i + κ o + κ int , with κ int being the intrinsic decay rate of the resonator. The damping of the magnon mode may result from the dipolar interactions of the P1 center with neighboring 13 C nuclei and other P1 centers [60] as well as the inhomogeneities of the external static magnetic field and the magnetic field of the resonator mode. For our coplanarwaveguide resonator, κ i ≈ κ o κ int . Away from the left or right anticrossing point, the polariton mode approaches the resonator mode and its line width gives κ/2π = 0.6 ± 0.05 MHz. When measuring the transmission spectrum in Fig. 2(a), we do not apply a drive tone to the system. Thus, b † b ≈ 0 at cryogenic temperatures. We obtain g eff /2π ≈ g/2π = 17.2 ± 0.5 MHz from the Rabi splitting at the middle anticrossing point. Finally, we use Eq. (1) to fit the measured transmission spectrum around the middle anticrossing point, which gives γ/2π = 11.9 ± 0.3 MHz.
When including both the decay rate of the resonator mode and the damping rate of the magnon mode [61,62], the Hamiltonian of the hybrid system can be effectively written as a non-Hermitian Hamiltonian: In matrix form, we can write the effective Hamiltonian as which has two complex eigenvalues, (4) When the magnon is tuned in resonance with the resonator mode (ω c = ω 0 ), the two eigenvalues of the magnon polaritons are reduced to For 2g eff > (γ − κ), Eq. (5) gives two separate polariton modes at the anticrossing point, each with line width 1 2 (κ + γ). However, these two polariton modes coalesce to one at 2g eff = (γ − κ), which is the EP of the hybrid system. When indicating that the two polariton modes have different line widths. Corresponding to the two eigenvalues in Eq. (5), the two eigenvectors also coalesce at the EP and each behaves differently on the two sides of the EP (cf. Appendix B).
To demonstrate the EP, we tune g eff via b † b by applying a drive tone of frequency ω d /2π = 3.093 GHz to the system, which is in resonance with the magnons of the s = 0 subensemble. To guarantee that the system is in its stationary state for each measurement, this driving tone lasts for 3000 s before the next measurement is performed (cf. Appendix C). Figure 2(b) shows the measured transmission spectrum of the system versus the power of the drive tone, where the static magnetic field is fixed at the middle anticrossing point in Fig. 2(a). This demonstrates the driving-power dependence of the Rabi splitting for the two polariton branches. It is clear that the EP occurs at P d ≈ −93.7 dBm. In the region of P d < −93.7 dBm, the two polariton branches are separated but when increasing P d , the separation (i.e., the Rabi splitting 2g eff ) decreases and, at the EP, the two peaks coalesce to one. We extract these data from the Rabi splitting and show the behavior of g eff versus the drive power P d in Fig. 2(c) (cf. the gray circles). There is a relation between the reduced magnon occupation χ ≡ b † b /N and the Rabi frequency Ω d (∝ √ P d ) of the drive tone (see Appendix D), where ξ = (1 − 2χ)/(1 − χ) and η = g 2 eff /κ 2 . We use this relation to fit the data of g eff [see the curve in Fig. 2(c)] and obtain Ω d = k √ P d , with the fitting parameter k = 1.54 × 10 14 . The theoretical value of k is k = κ/2hω c = 9.59 × 10 14 , which is larger than the above fitting value. This is due to the leakage of the drive power in the setup. Next, we use the fitted results [i.e., the curve in Fig. 2(c)] and Eq. (1) to simulate the transmission spectrum [ Fig. 2(d)]. Clearly, this simulated transmission spectrum agrees well with the experimental results in Fig. 2(b).
To make the observation of the EP more convincing, we compare the experimental data with the real and imaginary parts of ω 1,2 in Eq. (5) by varying the drive power P d . The two curves in Figs. 3(a) and 3(b) correspond to the real and imaginary parts of ω 1,2 in Eq. (5), where the dependence of the effective magnon-photon coupling g eff on the drive power is given by the fitting curve in Fig. 2(c). These two curves show the characteristics of the EP. In Fig. 3(a), the experimental data are extracted from the peak positions of the two polariton branches in Fig. 2(b); and the experimental data in Fig. 3(b) correspond to the line widths of the two magnon polaritons. At P d ≈ 93.7 dBm, the two polariton modes coalesce and only one peak is visible in the region of P d > 93.7 dBm. To obtain the line widths of the two modes in this region, we first fit the transmission spectrum and acquire the effective coupling strength g eff . Then, we use Eq. (5) and the fitted g eff to infer the line widths of the two modes, as in Ref. [17]. Indeed, the experimental data show a good agreement with the theoretical results, further confirming the observation of the EP in our hybrid system.

III. CROSS-RELAXATION EFFECT IN P1 CENTERS
Above, we only consider the coupling between the s = 0 spin subensemble and the resonator mode. In fact, as shown in Fig. 2(a), magnons in other two spin subensembles (s = ±) are also strongly coupled to the resonator mode. The results in Fig. 2(a) show that the magnon-photon coupling strengths are nearly equal for these three subensembles. Here, the static magnetic field is tuned to have the magnons of the s = 0 subensemble in resonance with the resonator mode, so the magnons in other two spin subensembles are very off resonant with the resonator mode. In such a dispersive regime, each of the s = ± subensembles yields a frequency shift to the resonator mode (Appendix E) and the frequency of the resonator mode is shifted from ω c tõ where g eff,± and δ ± = ω c − ω ± are the effective magnonphoton coupling strengths and frequency detunings between the resonator mode and the s = ± subensembles, respectively. In Fig. 2(b), the frequency of the drive tone is ω d /2π = 3.093 GHz, which is in resonance with the magnons of the s = 0 subensemble. Owing to the cross relaxation in P1 centers [56,57], identical occupations of magnons can be induced in other two spin ensembles when the drive tone is applied for a long time [43]. Thus, g eff,+ = g eff,− . If this is obeyed when varying the drive power, due to δ + = −δ − , the frequency shifts of the resonator mode induced by the s = ± subensembles then cancel each other andω c is reduced to ω c . Therefore, the frequency of the resonator mode is not modified by the presence of the s = ± spin subensembles, owing to the cross-relaxation effect that induces g eff,+ = g eff,− . Indeed, the observation of the EP in Fig. 2(b) reveals that the frequency of the resonator mode does not vary when tuning the drive power. Otherwise, the resonant conditioñ ω c = ω 0 for observing the EP cannot be satisfied. While the static magnetic field is tuned to have the magnons of the s = 0 subensemble in resonance with the resonator mode, we also manage to measure the transmission spectrum of the hybrid system by resonantly pumping the magnons in the s = + and s = − subensembles, respectively [Figs. 4(a) and 4(b)], instead of the magnons in the s = 0 subensemble. The results reveal that the transmission spectrum is also symmetric about the resonator frequency ω c and the EP occurs as well, as in Fig. 2(b). This means that g eff,+ = g eff,− , irrespective of which spin subensemble is pumped by the drive tone. This further proves the equal occupations of magnons in the s = ± subensembles and demonstrates the robustness of the EP against driving in this hybrid system. In Figs. 4(a) and 4(b), the simulated peak positions of the two polariton modes are shown, which also match the experimental results. Compared to the position of the EP in Fig. 2(b), the EP in Figs. 4(a) and 4(b) shifts to a higher drive power. This is because the drive tone now has a frequency detuning from the resonator mode and the relation between the magnon occupation and the drive power becomes different from Eq. (6) (cf. Appendix D).

IV. DISCUSSION AND CONCLUSIONS
We show how the EP can be removed by the varying frequency shift of the resonator mode when changing the drive power. The static magnetic field is tuned to have the magnons of the s = + spin subensemble nearly resonant with the resonator mode. Then, we resonantly pump the magnons in this subensemble with a drive tone. Figure 4(c) shows the measured transmission spectrum of the hybrid system versus the drive power. As in Fig. 2(b), this applied drive tone lasts for 3000 s before the next measurement is performed. The measured transmission spectrum becomes asymmetric and no EP is observed, in sharp contrast to Fig. 2(b). In the present case, the magnons of the s = 0, − subensembles are in the dispersive regime related to the resonator mode, yielding that the frequency of the resonator mode shifts from ω c tõ ω c = ω c + g 2 eff /δ 0 + g 2 eff,− /δ − , where the frequency detunings are δ 0 = ω c − ω 0 and δ − = ω c − ω − . The frequencies ω 1,2 of the two polariton modes are also given by Eq. (4) but with ω c , ω 0 and g eff replaced byω c , ω + and g eff,+ , respectively. Now, the terms (g 2 eff /δ 0 + g 2 eff,− /δ − ) inω c vary with the drive power, so the resonant conditionω c = ω + cannot always be obeyed. Therefore, we cannot reduce Eq. (4) to the simple form of Eq. (5) to exhibit the EP.
In Fig. 4(d), we also measure the transmission spectrum of the hybrid system versus the drive power but we tune the static magnetic field to have the magnons of the s = − subensemble nearly resonant with the resonator mode and then resonantly pump the magnons in this subensemble with a drive tone. As in Fig. 4(c), the transmission spectrum is asymmetric and no EP is observed. In this case, the magnons of the s = 0, + subensembles are in the dispersive regime related to the resonator mode and the frequency of the resonator mode is shifted from ω c toω c = ω c + g 2 eff /δ 0 + g 2 eff,+ /δ + , where δ 0 = ω c − ω 0 and δ + = ω c − ω + . Also, this modified frequency of the resonator mode varies with the drive power, so that the resonant conditionω c = ω − cannot always be satisfied. The frequencies ω 1,2 of the two polariton branches are given by Eq. (4) as well but with ω c , ω 0 , and g eff replaced byω c , ω − , and g eff,− , respectively. Compared to Fig. 2(b), the EP is removed in Figs. 4(c) and 4(d), due to the varying frequency shift of the resonator mode with the drive power.
Some experimental results on the drive-power dependence of the effective coupling have been presented in Refs. [43,63] but no EP has been observed there. In Ref. [43], the Rabi splitting versus the drive power has been shown only for the s = + subensemble. As discussed above, in our study the EP does not occur in this case [cf. Fig. 4(c)], because the resonant conditionω c = ω + cannot always be obeyed for the s = + subensemble. Also, as analyzed in Sec. III, only for the s = 0 subensemble can the EP occur in the hybrid quantum system considered. To confirm the observation of the EP, it is essential to both show the experimental results in Fig. 2(b) for the s = 0 subensemble and demonstrate the coalescence of the peak positions and line widths of the two polariton modes at the EP (i.e., Fig. 3). No such results have been shown in Refs. [43,63].
The results in Refs. [43,63] have been explained by a depolarization model involving the mixed states of the system. In our theory, the system is described by a non-Hermitian Hamiltonian H eff in Eq. (2). Then, the density operator ρ(t) of the system is governed by the following Liouvillian equation [61]: Thus, ρ(t) can be explicitly expressed as   4), but with ω c , ω 0 and g eff replaced byω c , ω ± and g eff,± respectively. In (c), the magnons of the s = + subensemble are tuned to have frequency ω + /2π = 3.106 GHz via the static magnetic field. The magnons of the s = 0 and − subensembles have frequencies ω 0 /2π = 3.012 GHz and ω − /2π = 2.918 GHz. The frequency of the resonator mode is found to be ω c /2π = 3.095 GHz, which is slightly blue shifted from ω c /2π = 3.093 GHz in Fig. 2 due to the static magnetic field (see Appendix C). In (d), the magnons of the s = − subensemble are tuned to have frequency ω − /2π = 3.080 GHz via the static magnetic field. The magnons of the s = 0 and + subensembles have frequencies ω 0 /2π = 3.175 GHz and ω + /2π = 3.269 GHz. The frequency of the resonator mode is found to be ω c /2π = 3.090 GHz, which is slightly red shifted from ω c /2π = 3.093 GHz in Fig. 2 due to the static magnetic field. The other parameters used are the same as in Fig. 2.
For pure states of the system, ρ 2 (t) = ρ(t), while ρ 2 (t) = ρ(t) in our case (see Appendix B). This reveals that the system is in the mixed state, even if the non-Hermitian Hamiltonian of the system has the simple form in Eq. (2). Therefore, the conclusion regarding the mixed states of the system is consistent with the depolarization model in Refs. [43,63]. Our model is based on the concepts and methods of non-Hermitian physics. Compared with the depolarization model, our model has the following three distinct merits. First, using the Holstein-Primakoff transformation and the meanfield approximation, we show that the effective coupling can be obtained using the quantum Langevin equation. By defining the spinpolarization factor in Refs. [43,63] as P eff = 1 − b † b /(N/2), the effective coupling strength in our model can be converted to g eff = g √ P eff in the depolarization model. Thus, our model provides a microscopic interpretation of the spinpolarization factor. Second, in our model we derive the effective non-Hermitian Hamiltonian of the hybrid system and convincingly explain the experimental results related to the EP by analyzing the real and imaginary parts of the eigenvalues of the effective non-Hermitian Hamiltonian, which is beyond the depolarization model. In contrast, the depolarization model in Refs. [43,63] only gives the relation between the effective coupling and the drive power but does not provide the effective non-Hermitian Hamiltonian and its eigenvalues (i.e., it cannot be used to describe the EP). Third, while our model uses the framework of the effective non-Hermitian Hamiltonian, which is different from the depolarization model, it may also be harnessed to study other phenomena of non-Hermitian physics [64]. Thus, our study sheds new light on power-induced depolarization phenomena in paramagnetic systems and will stimulate further work in this new direction.
In summary, we observe the EP in a hybrid quantum system consisting of P1 centers in diamond coupled to a coplanarwaveguide resonator and also convincingly show the crossrelaxation effect in P1 centers. Among various applications of EPs, enhancing the sensitivity of detection has been widely investigated [10][11][12][13][14][15][16]. In critical quantum metrology, it is important to precisely tune the coupling. For example, in Ref. [65], a metrological protocol is designed to probe one physical parameter (e.g., the spin frequency) of the Rabi model by slowly sweeping the coupling from zero to some desired value close to the critical point. Thus, the good tunability of the coupling in this hybrid quantum system may facilitate the application of EPs in metrology.
The cross relaxation is an important phenomenon for nitrogen impurities in diamond, because it can be used to produce population inversion in paramagnetic nitrogen donors to implement solid-state masers [56][57][58]. In the future, one can introduce an effective gain in the coplanar-waveguide resonator [66] to balance the loss and gain in the hybrid system. This can achieve a hybrid system with parity-time symmetry. Also, in the present work, only stationary-state properties of the hybrid system are studied. Thus, another future work can be focused on the time evolution of the non-Hermitian system. While coalescent eigenvectors at the EP are guaranteed by the non-Hermiticity of the Hamiltonian [1], it is interesting to probe the behavior of the eigenstates around the EP [6,8]. For example, with the good tunability of the coupling strength and the magnon-mode frequency, it is promising to explore the exotic topological phenomena related to the EP (e.g., nonreciprocal energy transfer [6] and asymmetric mode switching [8]) by dynamically encircling the EP via varying the coupling strength and the magnonmode frequency in the hybrid system. Moreover, monitoring the evolution of the non-Hermitian system can also reveal intriguing time-dependent responses of the system on the drive tone, including the time evolution of the magnon occupation in each spin subensemble. This can provide further information about the cross relaxation in P1 centers. Appendix A: EXPERIMENTAL SETUP The hybrid system consists of an ensemble of P1 centers in diamond coupled to a coplanar waveguide resonator [ Fig. 1(a)]. The P1 center is a substitutional nitrogen defect [ Fig. 1(b)]. Due to the hyperfine interaction with the host 1 4N nucleus, the P1 center has six energy levels and three allowed transitions [ Fig. 1(c)]. In the experiment, the sample is placed in a dilution refrigerator and cooled down to 20 mK [Fig. 5].
The coplanar waveguide resonator is fabricated by reactive ion etching of a 50-nm-thick d.c.-magnetron-sputtered niobium film on a thermally oxidized silicon substrate. The central conductor of the resonator is 20 µm wide and its gap to the ground plane is 11.6 µm, so that a 50 Ω characteristic impedance is obtained. The resonator has a length of 20 mm, defined by two nearly identical interdigital coupling capacitors with a capacitance of approximately 12 fF. The type-1b diamond used is synthesized under both high pressure and high temperature, in which the P1 centers are the main defects and constitute the spin ensembles harnessed in the experiment. The damping rate of the collective spin excitations (i.e., magnons) of the s = 0 spin subensemble in diamond is about γ/2π = 11.9 ± 0.3 MHz (half width at half maximum). This is determined via the line width Γ (half width at half maximum) of the polaritonic peaks under weak probe-field measurements [ Fig. 2(a)], with the relation Γ = (γ + κ)/2. The damping rates of the other two subensembles are nearly the same as that of the s = 0 subensemble. To perform the measurement, we first apply a drive tone of a given frequency on the coplanar waveguide resonator for a duration of time and then apply a fast and low probe-power signal (approximately 1 fW) to implement  Fig. 1(a) is placed in a cryogenic chamber.
the measurement of the transmission spectrum via a vectornetwork analyzer (VNA). The drive tone is generated by an analog signal generator and the duration time of the drive tone is set to be 3000 s to ensure that the system is in the stationary state. Then, we change the power of the drive tone and repeat the above process.

Appendix B: EFFECTIVE HAMILTONIAN OF THE SYSTEM
The Hamiltonian of the P1 center includes the Zeeman energy and the hyperfine interaction between the electron (spin 1/2) and the host nucleus (spin 1) [67], where γ e /2π = 28 GHz/T is the gyromagnetic ratio, B is the static magnetic field, S ≡ (S x , S y , S z ) and I ≡ (I x , I y , I z ) are spin operators of the electron and the host nitrogen nucleus, respectively, and A is the hyperfine interaction tensor [68]. Note that in Eq. (B1), we leave out the Zeeman energy and the quadrupole interaction of the host nucleus, which only involve the operator I [69]. Due to the Jahn-Teller distortions, one of the four C − N bonds is elongated. When the Jahn-Teller axis (i.e., the elongated C − N bond) is along the z direction, the corresponding hyperfine tensor A is diagonal, i.e., A diag /2π = diag(81. 8, 81.8, 114.2) MHz [68]. In our experiment, the [100] crystal axis of the diamond sample is aligned along the external magnetic field B = B 0 e z . Because all C − N bonds have equivalent angles with B, the hyperfine interaction is the same for each of the P1 centers. In the coordinate frame with the z axis oriented along the [100] crystal axis, the hyperfine tensor A becomes off diagonal, which can be obtained via a transformation of the diagonal hyperfine tensor A diag [70]. The Hamiltonian of the P1 center can be rewritten as where A /2π = A ⊥ /2π ≈ 94 MHz are the diagonal elements of the hyperfine tensor A and the off-diagonal elements are included in the tensor A off . The hyperfine interaction is much smaller than the Zeeman energy and can be treated as a perturbation. From first-order perturbation theory, the Hamiltonian is reduced to where s = 0, ±1 are the three eigenvalues of the operator I z , with |s being the corresponding eigenvectors. According to the values of s, the ensemble of P1 centers can be divided into three subensembles of spins with transition frequencies ω + = γ e B 0 +A , ω 0 = γ e B 0 , and ω − = γ e B 0 −A , respectively. If we only take one subensemble into consideration (e.g., s = 0), the Hamiltonian of the hybrid system is where a (a † ) is the annihilation (creation) operator of the resonator mode with frequency ω c , S i ≡ (S i x , S i y , S i z ), and S ± i ≡ S i x ± iS i y are the operators of the ith spin in the s = 0 subensemble, and g s,i is the coupling strength between the ith spin and the resonator mode. For simplicity, the same coupling strength g s,i = g s is assumed for all spins in the subensemble. To describe the collective behavior of the spins, we define the macrospin operator J ≡ (J x , J y , J z ) = ∑ i S i and the collective coupling strength g ≡ √ Ng s , with N being the number of spins in the subensemble. Then, the Hamiltonian (B4) is reduced to with J ± ≡ J x ± iJ y .
To study the exceptional point (EP) of the hybrid system, we first use the Holstein-Primakoff transformation [71], to convert the Hamiltonian (B5) to is the annihilation (creation) operator of the magnons, which are the collective spin excitations in the s = 0 subensemble. Next, we linearize the above Hamiltonian under the mean-field approximation. Based on the Taylor's expansion in terms of b † b/N, the coupling term a † 1 − b † b/N b can be written as is the corresponding fluctuation. When keeping the terms up to first-order fluctuations, Neglecting the counter-rotating term a † b † under the rotatingwave approximation, we have Also, we have More generally, we have Therefore, Eq. (B8) can be approximately written as where When including the coefficient of the counter-rotating terms, the neglected counter-rotating part in the second term of the Taylor's expansion [cf. Eqs. (B7)-(B9)] is c 1 ga † b † , where c 1 = −β 2 /2N. For the third term of the Taylor's expansion, the neglected counter-rotating part is c 2 ga † b † , where c 2 = −|β | 2 β 2 /8N 2 . In fact, the neglected counter-rotating parts in all terms of the Taylor's expansion can be summed as Obviously, |C Σ | ≡ 1 − 1 − |β | 2 /N < 1. Because g ω c , ω 0 in our hybrid system, it is thus reasonable to perform the rotating-wave approximation in the above derivations, even under high pump powers. Since all high-order terms in the Taylor's expansion are considered in Eqs. (B8)-(B15), the mean-field approximation in Eq. (B14) is also valid for the high-excitation case. Substituting Eq. (B14) and its Hermitian conjugate into Eq. (B7), we obtain where g eff = g 1 − |β | 2 /(N/2). Obviously, the last term in the above equation can be absorbed into the drive Hamiltonian H d : Actually, the displacement term Ω b a † + Ω * b a in the Hamiltonian (B17) is a counter-rotating term. When the magnons are tuned to be nearly resonant with the resonator mode (i.e., ω 0 ∼ ω c ), this counter-rotating term can also be ignored and the Hamiltonian (B17) is reduced to H s = ω c a † a + ω 0 b † b + g eff (a † b + ab † ) in the near-resonance case. When the decay rates of the resonator mode and the spin ensemble, κ and γ, are included, the Hamiltonian of the system can be effectively written, in the non-Hermitian form [61,62], as With the obtained effective non-Hermitian Hamiltonian in Eq. (B19), the dynamics of the hybrid system is governed by the following Liouvillian equation [61]: where ρ(t) is the density operator of the system at time t.
Solving the Liouvillian equation, we can express the density operator as It is known in quantum mechanics that ρ 2 (t) = ρ(t) if the system is in a pure state. However, in the non-Hermitian case that we study, Due to H † eff = H eff , exp(iH † eff t) exp(−iH eff t) = 1, so there is ρ 2 (t) = ρ(t), even if the non-Hermitian Hamiltonian has a simple form in Eq. (B19). This means that the system is in a mixed state, instead of a pure state.
In matrix form, we can write the effective non-Hermitian Hamiltonian as Diagonalizing it, we obtain the two eigenvalues, this being Eq. (4) in the main text. In the resonant case (ω c = ω 0 ) that we study, Eq. (B24) reduces to Eq. (5) in the main text, i.e., Corresponding to these two eigenvalues, the two eigenvectors are where T denotes the matrix transpose and γ > κ in our hybrid system. The relative amplitudes A 1,2 can be written as and the relative phases φ 1,2 are Obviously, the two eigenvectors coalesce to |ψ 1 = |ψ 2 = (i, 1) T at the EP: g eff = (γ − κ)/2. On the left side of the EP, i.e., g eff < (γ − κ)/2, the relative amplitudes A 1,2 of the two eigenvectors are different, but the relative phases φ 1,2 are the same. On the contrary, on the right side of the EP, i.e., g eff > (γ − κ)/2, the relative amplitudes A 1,2 are the same, but the relative phases φ 1,2 are different (cf. Fig. 6). To measure the transmission spectra in Figs. 2(b), 4(a), and 4(b) (see the main text), we first tune the s = 0 subensemble in resonance with the resonator mode ω 0 /2π = ω c /2π = 3.093 GHz and then pump the system with a drive tone at ω d /2π = 3.093 GHz for a duration of time. Meanwhile, we monitor the transmission amplitude at 3.090 GHz and measure it once every 5 seconds. As shown in Fig. 7, the transmission amplitude becomes nearly stable when the duration time of the drive tone is increased to be approximately 2000 s. To ensure that the hybrid system is in the steady state, in both Fig. 2(b) and Fig. 4 the duration time of the drive tone is chosen to be 3000s before each measurement of the transmission amplitude is implemented (cf. the dashed vertical line in Fig. 7). To measure the transmission spectra in our experiment, the duration time of the drive tone is chosen to be 3000 s before each measurement of the transmission amplitude is implemented.
In the cases of Figs. 4(c) and 4(d), the bare frequencies of the resonator mode are fitted to be 3.095 GHz and 3.090 GHz, respectively. The different bare resonator-mode frequencies found in Figs. 4(c) and 4(d) are due to the effect of the static magnetic field on the superconducting waveguide resonator, because the applied static magnetic field can unavoidably affect the superconductivity of the resonator. To verify this, we measure the dependence of the transmission spectrum on the magnetic field in the same way as in Fig. 2(a), but the probe tone is chosen to be so intense that all three spin subensembles are almost decoupled from the resonator, i.e. the magnon occupation is achieved to be b † b = (N/2) for each subensemble to have g eff = g eff,± = 0. Indeed, as shown in Fig. 8(a), the bare frequency of the resonator mode is gradually shifted as the static magnetic field strengthens. The three dashed vertical lines (from the left to the right) correspond to the magnetic-field strengths used in Fig. 4(c), Fig. 2 In analyzing the experimental results, we assume that the decay rate of the resonator mode is independent of the drive power. To demonstrate this, in Fig. 8(b), we measure the decay rate of the resonator mode versus the drive power. It can be seen that the measured decay rate of the resonator fluctuates slightly around κ/2π = 0.6 MHz. This indicates that the decay rate of the resonator is nearly independent of the drive power. In the experiment, we apply a driving tone to achieve considerable magnon occupations in the considered spin subensemble. With the drive field included, the Hamiltonian of the system becomes in the rotating frame with respect to the frequency ω d of the drive field, where ∆ c = ω c − ω d (∆ s = ω s − ω d ) is the frequency detuning between the resonator mode (spin subensemble) and the drive field. When the decay rates κ and γ of the resonator mode and the spin subensemble are considered, the dynamics of the hybrid system can be described using the quantum Langevin equations [59]: From Eq. (D2), it follows that the mean values of a and J − satisfy where the mean-field approximation J z a = J z a is used for the two-operator term.
For the steady state of the system, ȧ = J − = 0, so The first equation in Eq. (D4) gives Substituting the above equation into the second equation in Eq. (D4), we have With the Holstein-Primakoff transformation and the meanfield approximation, J z and J − can be expressed as  8. (a) The dependence of the resonator-mode bare frequency on the static magnetic field. The transmission spectrum is obtained in the same way as in Fig. 2(a), but using a strong probe tone. (b) The measured line width of the resonator mode versus the drive power.
where ξ = (1 − 2χ)/(1 − χ) and η = g 2 (1 − 2χ)/(∆ 2 c + κ 2 ), with χ = b † b /N being the reduced magnon occupation. Since b 2 /N ( b 2 /N) * ≈ χ, multiplying the above equation with its complex conjugate, we obtain the steadystate solution for the magnon occupation: To measure the transmission spectrum in Fig. 2(b) of the main text, we use ω d = ω 0 = ω c , i.e. ∆ s = ∆ c = 0. In this resonant case, Eq. (D8) is reduced to Eq. (6) in the main text for fitting the effective magnon-photon coupling in Fig. 2(c), where η = g 2 (1 − 2χ)/κ 2 = g 2 eff /κ 2 . For the measured transmission spectra in Figs. 4(a) and 4(b) (see the main text), the drive tone is tuned to be resonant with the s = + and − subensembles, respectively, while the resonator mode is in resonance with the s = 0 subensemble (i.e., the resonator mode has a frequency detuning ∆ c from the drive tone). In the steady state, the drive-power dependence of the reduced magnon occupation χ = b † b /N as well as the corresponding effective magnon-photon coupling strength g eff,± = g √ 1 − 2χ can be described by where η = g 2 (1 − 2χ)/(∆ 2 c + κ 2 ). The experimental data (circles) in Figs. 9(a) and 9(b) are the effective magnonphoton couplings extracted from the Rabi splittings in Figs. 4(a) and 4(b), respectively. The solid curves are the corresponding effective magnon-photon couplings g eff,+ and g eff,− calculated using Eq. (D9) with ∆ c /2π = −94 MHz [ Fig. 9(a)] and ∆ c /2π = 94 MHz [ Fig. 9(b)], respectively. It is clear that the numerical results for both g eff,+ and g eff,− are in good agreement with the experimental data of g eff . This indicates that the cross relaxation occurs in the P1 centers, which yields g eff,± = g eff . Therefore, as explained in the main text,ω c = ω c + g 2 eff,+ /δ + + g 2 eff,− /δ − is reduced to ω c , because δ + = −δ − . This gives rise to the symmetric transmission spectra in Figs. 4(a) and 4(b) about the frequency of the resonator mode and the EP becomes observable in the hybrid system. Similar to Figs. 4(c) and 4(d) in the main text, in Fig. 10, we also show the transmission spectra without the EP, which are measured using an off-resonant drive tone. As in Fig. 4(c), we tune the s = + subensemble to be nearly resonant with the resonator mode. However, instead of driving the s = + subensemble, we use a drive tone to resonantly pump the s = 0 [ Fig. 10(a)] and s = − [ Fig. 10(b)] subensembles, respectively. Then, as in Fig. 4(d), we tune the s = − subensemble to be nearly resonant with the resonator mode, but use a drive tone to resonantly pump the s = 0 [ Fig. 10(c)] and s = + [ Fig. 10(d)] subensembles, respectively. The corresponding effective magnon-photon couplings are also calculated using Eq. (D9) and the two eigenvalues of the considered hybrid system can be then obtained. As in Figs. 4(c) and 4(d), the dashed curves are the simulated results for the real parts of the two eigenvalues, which also match the peak positions of the two polariton branches in Fig. 10.  Fig. 4(a) (see the main text), where the s = 0 subensemble is in resonance with the resonator mode and the s = + subensemble is resonantly pumped by a drive tone. The solid curve is calculated for g eff,+ using Eq. (D9) with ∆ c /2π = −94 MHz. (b) The data for g eff (circles) extracted from the Rabi splitting in Fig. 4(b), where the s = 0 subensemble is in resonance with the resonator mode and the s = − subensemble is resonantly pumped by a drive tone. The solid curve is calculated for g eff,− using Eq. (D9) with ∆ c /2π = 94 MHz. The other parameters are the same as in Fig. 2(c).

The effective interaction Hamiltonian H (I)
eff contains the second-order off-diagonal terms. When neglecting these second-order terms in the dispersive regime |ω c − ω ± | g eff,± , the effective Hamiltonian of the system becomes H eff ≈ H Furthermore, we remove the commutative diagonal terms eff and include the decay rates of the resonator mode and the s = 0 subensemble. Then, the reduced Hamiltonian of the hybrid system can be effectively written, in the non-Hermitian form [61,62], as  The Hamiltonian (E10) has the same form as Eq. (2) in the main text and its two eigenvalues are also given by Eq. (4) in the main text, i.e., ω 1, 2 in Eq. (B24), but ω c therein is replaced by ω c .