Universal power law decay in the dynamic hysteresis of an optical cavity with non-instantaneous photon-photon interactions

We investigate, experimentally and theoretically, the dynamic optical hysteresis of a coherently driven cavity with non-instantaneous photon-photon interactions. By scanning the frequency detuning between the driving laser and the cavity resonance at different speeds across an optical bistability, we find a hysteresis area that is a non-monotonic function of the scanning speed. As the scanning speed increases and approaches the memory time of the photon-photon interactions, the hysteresis area decays following a power law with exponent -1. The exponent of this power law is independent of the system parameters. To reveal this universal scaling behavior theoretically, we introduce a memory kernel for the interaction term in the standard driven-dissipative Kerr model. Our results offer new perspectives for exploring non-Markovian dynamics of light using arrays of bistable cavities with low quality factors, driven by low laser powers, and at room temperature.

Descriptions of bistable optical cavities commonly assume instantaneous photon-photon interactions. In the mean-field equation of motion for the intracavity field α, this assumption manifests as a Kerr nonlinearity of the form |α| 2 α [25]. The same cubic nonlinearity is found in the Gross-Pitaevskii equation employed in atomic physics [26][27][28], in the Ginsburg-Landau theory of superconductivity [29], in the Lugiato-Lefever equation describing pattern formation in nonlinear optics [30], and in the force derived from Goldstone's Mexican hat potential V = −|φ| 2 + |φ| 4 for the scalar field φ at the heart of the Higgs mechanism [31]. In optics, some of the strongest Kerr nonlinearities arise in semiconductor cavities where exciton-exciton interactions are effectively instantaneous [32]. A drawback of those cavities is that optical bistability based on Kerr nonlinearities is typically only observed at cryogenic temperatures. In contrast, several optical resonators with relatively slow but strong thermal nonlinearities have routinely displayed bistabil-ity at room temperature [2,33,[35][36][37][38][39]. As bona fide bistable systems, thermo-optical resonators may open up new perspectives for classical Hamiltonian simulation and computation [21][22][23][24]40] at room-temperature. However, the influence of the thermal relaxation time on the dynamic hysteresis of bistable cavities remains to be addressed.
In this Letter, we demonstrate signatures of noninstantaneous photon-photon interactions in the dynamic hysteresis of a tunable micro-cavity. We investigate a laser-driven micro-cavity filled with oil as shown Fig. 1, operating at room-temperature. This cavity supports optical bistability at low driving powers P ∼ 70 µW. Scanning the cavity length under laser illumination, we observe an optical hysteresis that depends pronouncedly on the ratio of the scanning time to the memory time of the interactions. In contrast to previous reports of dynamic hysteresis in resonators with effectively instantaneous interactions [13,[41][42][43], we find a hysteresis area that is a non-monotonic function of the scanning speed. For fast scans, the hysteresis area decays following a power law with a universal exponent independent of the system parameters. Our results elucidate how the hysteretic behavior characterizing first-order phase transitions, and the boundary between phases, dynamically vanish when the nonlinearity has a finite memory time.
Figure 1 illustrates our system: a tunable Fabry-Pérot micro-cavity driven by a 532 nm continuous wave laser. The cavity is made by a concave and planar mirror, each comprising a distributed Bragg reflector (DBR) on a glass substrate. The mirrors have a peak reflectance of 99.9% at 530 nm, which is the center of the stop-band. The concave mirror was fabricated by milling a glass substrate with a focused-ion beam prior to the deposition of the DBR [44]. The Fig. 1 inset shows a chip containing concave mirrors suitable for making micro-cavities with different mode volumes. In this work, we use the mirror enclosed by the dashed circle in Fig. 1, which has a diameter of 7 µm and a radius of curvature of 12 µm. Thanks to the strong lateral confinement and high mirror reflectivity we can probe single optical modes across cavity length scans of several nanometers.
The chip containing the concave mirror is aligned parallel to the planar mirror using a hexapod nanopositioner. This nanopositioner controls all three translational (rotational) degrees of freedom of the concave mirror with nanometer (micro-degree) precision. The planar mirror is mounted on another actuator used to scan the cavity length. Optical excitation and collection are achieved through 10× microscope objectives with numerical aperture N A = 0.25. The cavity transmission is measured by a photodetector and an oscilloscope. Further details about our setup are included in supplemental information [45].
To endow the cavity with a nonlinear optical response, we placed a drop of olive oil inside. Oils are known for their thermo-optical nonlinearities [46][47][48]. Through zscan measurements we estimated the nonlinear refractive index n 2 of our olive oil to be ∼ −5 × 10 −8 cm 2 /W at 532 nm, consistent with Ref. [48]. Figure 2 shows the transmitted intensity through our oil-filled cavity averaged over 70 cycles and at three laser powers. Green and black data points correspond to opening and closing the cavity, respectively. For low powers P 20 µW, the cavity response is linear. The gray curve over the measurements for P = 20 µW is a Lorentzian fit, yielding a resonance linewidth of 0.104 ± 0.001 nm. For P = 70 µW, the transmission displays hysteresis (see arrows in Fig. 2) and a narrow bistability around a mirror position of 0.1 nm. The power needed for bistability in our cavity is similar to that in state-of-the-art monolithic semiconductor cavities [13,15,49], but at conveniently lower quality factors ( by a factor of ∼ 10) and operating at room-temperature instead of ∼ 5 K. For P = 150 µW, the bistability and hysteresis range enlarge as expected. An estimate of the maximum temperature rise in our bistable oil-filled cavity is provided in supplemental information [45].
All measurements in Fig. 2 correspond to linear ramps of the cavity length at 1.75 µm/s. Already for this relatively slow scan, a pronounced overshoot followed by a slow decay of the transmitted intensity arises when closing the cavity in the nonlinear regime. This overshoot is due to the finite thermo-optical response time of the cavity, which is not captured by the standard drivendissipative Kerr model for a single-mode cavity with instantaneous interactions [2]. The standard Kerr model for the intra-cavity meanfield α in a frame rotating at the driving frequency ω is is the laser-cavity detuning, with ω 0 the resonance frequency. U is the photon-photon interaction strength. F is the driving amplitude. The total loss rate Γ = κ 1 + κ 2 + γ is the sum of the input-output leakage rates through the two mirrors, κ 1,2 , and the intrinsic cavity loss rate γ due to absorption. The steady-state follows from settingα = 0 in Eq. 1. We attempted to fit the steady-state photon density |α| 2 calculated with Eq. 1 to the measurements for P = 150 µW in Fig. 2, with F as the only relevant adjustable parameter [45]. Γ is fixed by the resonance linewidth observed in the linear regime. Furthermore, since Eq. 1 is a mean-field model, the absolute number of photons |α| 2 or the interaction energy U alone are irrelevant; the spectral lineshape is determined by the ratio U |α| 2 /Γ. Therefore, any spectral lineshape can be obtained by varying F for any fixed U and Γ. Consequently, we adjusted F until obtaining the gray curve plotted over the measurements for P = 150 µW. Solid and dashed curves represent stable and unstable states, respectively, determined following Ref. 25. The fit is good far from resonance, but deviates from the data near the bistability. Next, we show that our data deviates more pronouncedly from predictions of the standard Kerr model as the scanning speed increases.
We performed hysteresis measurements for P = 150 µW and various scanning speeds. We selected a laser power far above the bistability threshold to limit the influence of noise on our measurements. Figure 3(a) shows average dynamic hysteresis measurements for three speeds. Top-to-bottom, the speed is ν, 7ν, and 49ν, with ν = 0.74 µm/s. The transmitted intensity is shown as a function of ∆/Γ, which we determined from the mirror position and the resonance linewidth in Fig. 2. Figure 3(a) shows how the hysteresis cycle qualitatively changes with the scanning speed. Increasing the speed from ν to 7ν makes the overshoot broader and the hysteresis wider. Interestingly, further increasing the speed to 49ν makes the overshoot broader but the hysteresis narrower. The measured lineshape for 49ν resembles a Lorentzian resonance for both scanning directions, although small deviations exist. This resemblance suggests that the response of the cavity is mostly linear for fast scans, regardless of the high power.
The behavior in Fig. 3(a) can be explained by consid-ering the finite heating and cooling time of our oil-filled micro-cavity; this makes photon-photon interactions noninstantaneous. Therefore, we modify Eq. 1 by letting (2) with the kernel function defined as K(t) = U τ e −t/τ . τ is the memory time of the nonlinearity, which corresponds to the thermal relaxation time of our microcavity.
Here we have followed the prescription of Mori [50] and Hänggi [51] for dealing with finite-time interactions. However, whereas Mori-type equations involve non-instantaneous dissipation, we introduced noninstantaneous photon-photon interactions.
Making the substitution 2 in Eq. 1 yields an integrodifferential equation, which can be conveniently written (for numerical simulation) as two coupled differential equations: Equations 2 and 3 imply that the state of the system depends on its entire past, weighted by the memory kernel K(t). Thus, photon-photon interactions are nonlocal in time. Note that when α(t) is constant and can be taken out of the integral in 2, we recover Eq. 1. Hence, steady-states are unchanged by K(t).
Figure 3(b) shows dynamic hysteresis calculations using Equations 3, with the same parameter values used for the steady-state calculations in Fig. 2. As for the experiments, we show scanning speeds a factor of 7 apart. The model faithfully reproduces all features observed in experiments. In the calculations, we set the memory time to τ = 10 4 Γ −1 and the slowest scanning speed to f = 5.97 × 10 −6 Γ 2 [45]. Relative to the experiments, the value of τ is smaller (details ahead) and the speed is larger. We rescaled the time scales to avoid unnecessarily long and memory-expensive calculations. Our mean-field calculations can be directly compared to experiments because we respect the hierarchy of time scales in experiments: Γ −1 τ T b , with T b the scanning time across the bistability. Moreover,the ratio T b /τ is similar for experiments and calculations.
Next, we analyze the hysteresis area across a range of scanning speeds. The hysteresis area is defined as A = T 0 |I ∆↓ − I ∆↑ |dt , with I ∆↓ and I ∆↑ the transmitted intensity when ∆ decreases and increases, respectively. T is the driving period, which exceeds T b . In Fig. 4(a) we plot the experimental average A. The corresponding calculations based on Equations 3 are presented in Fig. 4(b). In both measurements and calculations, A peaks at a certain scanning speed. This peak arises at the cross-over between two dynamical regimes. For slow scans, A increases with the speed because the cavity cannot adiabatically follow the driving force. This regime of dynamic hysteresis and the corresponding scaling laws for A have been previously explored experimentally and theoretically [13,41,42,52]. The second and new regime we investigate comprises speeds above the value for which A peaks. Therein, A decays with increasing speed because the nonlinearity does not have time to build up during the scan. Essentially, the second regime corresponds to a transition from nonlinear to linear dynamics. We interpret this transition as an effective reduction of the number of attractors in our system from two to one.
Our measurements are limited to scanning speeds between ∼ 0.5 µm/s and ∼ 40 µm/s. The upper speed limit is determined by the resonance frequency of our piezoelectric actuator. On the other end, we limited our measurements to speeds above 0.5 µm/s to avoid low-frequency mechanical noise in our setup. Despite these limitations, the good agreement between experiments and calculations in the available range encourages us to use our model to interpret the physics over an extended speed range.
In Fig. 5 we calculate A across a wide range of speeds for different F . At low speeds, the driving conditions determine the scaling of A [13]. At high speeds, we find that A decays following a power law with universal exponent -1. By 'universal' we mean that the exponent, i.e. the slope of the gray lines fitted to the data in Fig. 5, is independent of the system parameters. To assess whether our experiments show evidence of this behavior, in Fig. 4(a) we plot a power law with exponent -1 over our high-speed data points. The overlap between this power law and our experimental data for the highest speeds suggests that we reached the onset of the -1 power law regime. For comparison, we plot a -1 power law on top of the corresponding calculations in Fig. 4(b). In this case, the power law was fitted to the calculations in Fig. 5 over an extended range. As in experiments, we observe the onset of the -1 power law within the restricted speed range Fig. 4(b).
Recent calculations [42] and experiments [13] on dynamic hysteresis in cavities with instantaneous interactions observed a universal power law decay of A at low speeds. In that case, A decays due to the influence of quantum fluctuations. Coincidentally, the universal exponent discovered in Refs. [13,42] is also -1, as in the present work. However, the scaling behavior discovered herein has an entirely different origin (i.e., due to noninstantaneous interactions and unrelated to fluctuations) and arises in the opposite regime of fast scans.
Next, we estimate the experimental thermal relaxation time by comparing experimental and theoretical hysteresis cycles in Fig. 3. Since in theory we set τ and the speed at which ∆/Γ is scanned, τ can be converted to a range of ∆/Γ and viceversa. In Fig. 3(b) we indicate the range of ∆/Γ corresponding to τ by dashed gray lines for the three speeds considered. As expected, the range of ∆/Γ corresponding to τ increases with the speed. For the lowest speed f , the scanning time across the bistability range T b largely exceeds τ . In that regime, the overshoot observed when ∆ decreases is the most significant feature which is not captured by the standard Kerr model. The overshoot decays within a time ∼ τ , as expected. For the intermediate speed 7f , T b ∼ τ and the width of the overshoot approaches the bistability range. For the highest speed 49f , T b < τ and we have close-to-linear response. Experiments in Fig. 3(a) display the same behavior as the calculations. Hence, in a similar way we indicate the range of ∆/Γ corresponding to τ by two dashed gray lines in Fig. 3(a). Based on this range of ∆/Γ and our knowledge of the experimental scanning speed, all three measurements in Fig. 3(a) are consistent with a thermal relaxation time τ = 16 ± 1 µs.
In summary, we have shown how non-instantaneous photon-photon interactions influence the dynamic hysteresis of a coherently driven cavity. Non-instantaneous interactions maximize the hysteresis area at a finite scanning speed. At high speeds, the area decays following a universal power law with exponent -1. Through this scaling law, the hysteresis characterizing first-order phase transitions vanishes. Beyond single-cavity physics, our observation of optical bistability in oil-filled cavities paves the way for realizing bistable coupled cavities [53] and bistable cavity arrays at room-temperature. Such arrays could be used to probe Ising-type phase transitions [10] and solve combinatorial optimization problems [21][22][23], or to explore non-Markovian dynamics in complex optical networks. Unlike standard non-Markovian systems where finite-time dissipation makes the intracavity noise colored [50,51,[54][55][56][57], our system may support non-Markovian dynamics even for white intracavity noise. However, the cavity output noise spectrum is expected to be cut-off at high frequencies by the thermal relaxation time. This new regime of non-Markovian dynamics may be accessed by reducing the laser power so that existing noise in the cavity exerts a greater influence on the bistability, or by injecting tailored noise using modulators [49].  Figure S1 illustrates the setup we used to optically probe and characterize our oil-filled tunable microcavity. The cavity is driven by a single-mode continuous wave laser emitting at a wavelength λ = 532 nm. Excitation and collection are achieved through 10× microscope objectives with numerical aperture N A = 0.25. In all our experiments, we drive the fundamental transverse mode of the 9 th measurable longitudinal mode of the microcavity. Taking into account the electric field penetration into the distributed Bragg reflectors comprising our mirrors, the effective cavity length is ∼ 3 µm. To exclude multi-mode interference effects in our measurements, we optimized the in-coupling efficiency of the laser into the desired mode by finely adjusting the position of the concave mirror relative to the laser beam. Finally, the transmitted laser light was focused onto a photodetector by a f = 75 mm lens.

ACKNOWLEDGMENTS
The cavity length is modulated by displacing one of our mirrors with a piezoelectric actuator in closed-loop configuration. Through software, we specify the waveform, frequency, and travel range of the mirror. In all measurements we specified a linear ramp with a travel range of 350 nm. For modulation frequencies above ∼ 50 Hz, we suspected that the actual travel range of our actuator differed from the specified one. Therefore, we built a Michelson interferometer in the input arm of the setup [see Fig. S1] to characterize the mirror displacement. In particular, we measured the time-dependent intensity of a small section of the interferogram while the cavity length was being modulated. Next, we fitted the measured total intensity I t at the output of the interferometer with a two-beam interference equation of the form: I t = I 1 + I 2 + 2 √ I 1 I 2 cos(2kz + φ). I 1 and I 2 are the intensities in the two arms of the interferometer, k = 2π/λ is the angular wavenumber, z is the displacement of the actuator in one of the arms, and φ is a free parameter corresponding to the initial phase difference between the two arms. Through this analysis, we obtained the actual range traveled by our mirror. Finally, in combination with the scanning time, this travel range was used to calculate the scanning speeds along the horizontal axis of Fig. 4 in the main manuscript. Once the displacement of the mirror was properly calibrated, the beam splitter used for the Michelson interferometer was removed in order to avoid unwanted reflections which could disturb the hysteresis measurements.

Temperature rise in the oil-filled microcavity
In this section we estimate the temperature rise in our oil-filled microcavity in the nonlinear regime. To this end, let us first consider the mean-field equation of motion for the coherent field α in a driven-dissipative Kerr nonlinear cavity, i.e. Eq. 1 in the main manuscript. Calculating the steady-state solutions to that equation, one finds that the number of photons in the cavity N = |α| 2 satisfies: As in the main manuscript, Γ is the total loss rate, U is the photon-photon interaction strength, κ 1 is the inputoutput leakage rate, and F is the driving amplitude. In addition, we have defined the quantitỹ with ∆ = ω − ω 0 the frequency detuning between the driving laser and the cavity resonance, and U N the total interaction energy associated with a population of N photons. In Eq. 4, N appears as a Lorentzian function of the effective detuning∆. However, since N also enters into the right hand side of Eq. 4 via∆, N is a multi-valued function of the driving parameters F and ∆ . Indeed, Eq. 4 corresponds to a third order polynomial in N . This means that, in general, three steady-state solutions exist for a single driving condition. Of these three solutions, at most two are stable; this is known as bistability. The critical driving amplitude for observing bistability is F c = √ 3Γ 3 / (9κ 1 |U |). Next, we analyze the resonant response of our oil-filled microcavity. The cavity resonance frequency is with c the speed of light, q the longitudinal mode number, n the linear refractive index of the intra-cavity medium at ambient temperature, and L the cavity length. When the cavity length is scanned under laser illumination, the resonance frequency changes tõ ω 0 = qc 2 1 (n + δn) with δn the refractive index change due to laser-induced heating of the oil, and δL the change in cavity length due to the scan. Next, we give an approximate expression for Eq. 7 based on two observations: i) As Fig. 2 in the main article shows, the hysteresis range spans less than 0.2 nm in mirror displacement even for the highest laser power. Meanwhile, the initial cavity length is around 3 µm. Hence, δL L. ii) δn is on the order of 10 −4 [1]. Thus, we have δn n. Based on the above two inequalities,ω Next, we insert Eq. 8 into Eq. 5, and we take the value of L that satisfies qc/2nL = ω. Consequently, we obtaiñ The first term on the right hand side of Eq. 9 corresponds to the detuning ∆ in Eq. 5. The second term in Eq. 9 corresponds to the interaction energy, which is proportional to the intensity-induced refractive index change δn.
For our oil-filled cavity with thermo-optical nonlinearity, δn is given by with dn/dT a material constant and δT = T − T 0 the temperature change of the oil due to the light intensity. Next, we estimate δT based on energy conservation arguments, similar to Ref. 2. δT is related to the heat that goes inside and outside the cavity, i.e. q in and q out , via with C the heat capacity of the oil. Next, we assume thaṫ q in = AN , with N the number of photons in the cavity and A a constant describing the conversion of absorbed photons into heat. Furthermore, we assume thatq out = BδT , with B a constant describing heat dissipation into the environment. Hence, Eq. 11 becomes Therefore, in steady state (δ T = 0), the temperature rise δT is proportional to the photon number N : Combining Equations 9, 10, and 13, we find expressions for the rescaled detuning∆, the linear detuning ∆, and the thermo-optically induced interaction constant Based on the above analysis, we can estimate the temperature rise in the measurement of P = 150 µW shown in Fig. 2 of the main article. Equation 13 states that the highest temperature in the measurement corresponds to the largest N . From Eq. 4 the largest N corresponds tõ ∆ = 0. Using Eq. 9 and 10, we obtaiñ Finally, we solve the above expression for δT , and insert the parameter values corresponding to our experiment with a driving power of 150 µW. In particular, in Fig. 2 of the main article we observe δL = 0.177 nm. Furthermore, the cavity length is L ∼ 3 µm, and the linear refractive index of our oil at λ = 532 nm is n ≈ 1. 45. In addition, based on Ref. 1 we estimate dn/dT ∼ −4 × 10 −4 . Inserting these numbers in Eq. 15, we find the greatest temperature rise is δT = 0.2 • C.

Calculation details
In this section we provide further details about the calculations in the main text, and we explain how parameter values were determined or selected. Let us first consider the steady-state calculations based on Eq. 1 and presented in Fig. 2 of the main text. In particular, we fitted the steady-state photon number |α| 2 to the experimental transmitted signal (∝ |α| 2 ) for a laser power P = 150 µW. The model parameters are the photon-photon interaction strength U , the total loss rate Γ, the driving amplitude F , and the input-output leakage rates κ 1,2 . At first sight, it may seem that these five parameters can be freely adjusted in order to fit the measured lineshape. However, as explained next, we do not have this freedom due to several considerations and constraints.
The starting point of our analysis is the realization that Eq. 1 in the main text is a mean-field model neglecting quantum fluctuations. Hence, the value of each parameter individually, or of |α| 2 , is irrelevant. The spectral lineshape is entirely determined by the ratio U |α| 2 /Γ. In particular, the linear regime is characterized by U |α| 2 Γ. Bistability emerges for U |α| 2 Γ.
Next, let us explain how the model parameters relevant to the fit in Fig. 2 were set. First, the value of Γ was determined by fitting a Lorentzian function to the measured lineshape in the linear regime (P = 20 µW). Second, note that κ 1 is just a multiplicative factor for F . In fact, we could have defined an effective driving amplitude F = √ κ 1 F in Eq. 1 and not introduced κ 1 at all.
We included κ 1,2 in our model for consistency with standard input-output theory, and to have the right units for F . Therefore, we set κ 1 = γ/2 without this choice having any impact on our analysis. Furthermore, the value of κ 2 and γ do not need to be specified in the calculation at all. Only Γ needs to be specified. Third, we set U = 0.005 Γ. This choice determines the number of photons |α| 2 involved in the bistability and the critical driving amplitude F c needed to reach the bistable regime. Note, however, that any spectral lineshape can be attained for any value of U by scaling F (which determines |α| 2 ) accordingly. Therefore, our choice of U doest not impact our analysis. Based on the above considerations and the choice of U , we are left with F as the only adjustable parameter used to fit the calculated lineshape to the measured lineshape. At this point, we would like to note an additional constraint related to the critical driving amplitude F c for which bistability emerges. In particular, for F = F c we have U |α| 2 ∼ Γ. Experimentally, we can estimate F c by performing dynamic hysteresis measurements at different powers. F c then corresponds to the minimum laser amplitude for which bistability is observed. Consequently, any laser amplitude can be referenced to F c . In this way, we find that the value F/F c = 1.52 (as obtained from the calculations) giving the best fit to the experimental data is fully consistent with our experimental estimate of F/F c ≈ 1.5. Finally, we note that the fact that F/F c is the relevant way to express the driving amplitude in our mean-field model, combined with the fact that κ 1 determines F c , further supports the statement that the value of κ 1 we selected is irrelevant to our analysis of the spectral lineshape.
For dynamic hysteresis calculations based on Equations 2 and 3 in the main text, we use all of the same parameter values mentioned above. In addition, we set the thermal relaxation time to τ = 10 4 Γ −1 . Our choice of the value of τ , which is somewhat shorter than the experimental one as explained in the main text, was based on two considerations. First, calculations using the experimental value of τ would be extremely long and memory-expensive. Second, as long as τ Γ −1 , the physics of dynamic hysteresis in our single-mode cavity remains qualitatively the same. A longer τ will simply shift the maximum hysteresis area to slower scans, but it will not change the shape of the hysteresis area curve as a function of the scanning speed. Therefore, we selected a sufficiently large value of τ that still allows us to perform the calculations within a reasonable (∼ days) time.
The scans we performed to calculate the hysteresis area consist of varying the detuning ∆ in a linear and symmetric way from −20 Γ to 20 Γ. Finally, since the equations of motion are deterministic, it is sufficient to calculate the response for only one period for each T . The calculation results are presented in Figures 3b, 4,