Research ArticlePHYSICS

Critical slowing down in circuit quantum electrodynamics

See allHide authors and affiliations

Science Advances  19 May 2021:
Vol. 7, no. 21, eabe9492
DOI: 10.1126/sciadv.abe9492


Critical slowing down of the time it takes a system to reach equilibrium is a key signature of bistability in dissipative first-order phase transitions. Understanding and characterizing this process can shed light on the underlying many-body dynamics that occur close to such a transition. Here, we explore the rich quantum activation dynamics and the appearance of critical slowing down in an engineered superconducting quantum circuit. Specifically, we investigate the intermediate bistable regime of the generalized Jaynes-Cummings Hamiltonian (GJC), realized by a circuit quantum electrodynamics (cQED) system consisting of a transmon qubit coupled to a microwave cavity. We find a previously unidentified regime of quantum activation in which the critical slowing down reaches saturation and, by comparing our experimental results with a range of models, we shed light on the fundamental role played by the qubit in this regime.


The study of dissipative phase transitions has a long and interesting history not only due to their technological applications, such as in the construction of the laser (13), quantum limited amplifiers (4, 5), and optical switches (68) but also due to their theoretical interest since these phase transitions cannot be described by standard techniques such as mean-field theory (9). A key characteristic of first-order dissipative phase transitions is bistability (1013): close to the transition the two phases are metastable (14) and the dynamics of the system are highly sensitive to both its parameters and its initial state (1518). The steady state is reached via rare switching events during which the system transitions from one phase to the other (19, 20). This can be modeled using the theory of quantum activation in the case of dispersive optical bistability (21). Since the metastable states may be very long lived, this leads to critical slowing down in the equilibration time of the system. Critical slowing down has already been observed in a circuit quantum electrodynamics (cQED) lattice (22) and in an ensemble of nitrogen-vacancy centers coupled to a superconducting cavity (23) and has been modeled in the context of the Bose-Hubbard lattice (24).

Here, we observe critical slowing down in a cQED system with only two degrees of freedom: a transmon qubit (25) coupled to a three-dimensional (3D) microwave cavity (26). The nonlinearity introduced by the qubit causes the cavity to display bistability when a sufficiently strong microwave drive is applied. Within the bistable regime, the system divides its time between two metastable states, which are known as the bright and dim states according to the number of photons occupying the cavity. While the inherent nonlinearity of such a system has been exploited in (15) to achieve high-fidelity readout of the qubit state using high drive powers, here, we are interested in exploring and understanding the rich quantum dynamics happening at intermediate powers, close to the onset of bistability. We show that, at these powers, the cavity exhibits critical slowing down, reaching its steady state in a time much longer than the lifetimes of both the qubit and the cavity. We characterize the time scale of this slowdown as a function of driving frequency and power. We find a new regime of quantum activation in which the slowdown displays a saturation, which can only be explained by taking into account the full quantum description of the transmon. We demonstrate that even a simple superconducting circuit, consisting of only a qubit and a cavity, can be used to explore the rich physics of quantum phase transitions.

The device consists of a transmon qubit embedded in a superconducting aluminum 3D microwave cavity. Measurements of the transmitted signal through the cavity are performed using a standard cQED microwave setup described in Materials and Methods. This system can be described by the generalized Jaynes-Cummings (GJC) model, and its Hamiltonian can be written asH=ħn ωnnn+ħωcaa+ħm,n gm,nmn(a+a)+ħϵ(aedt+aedt)(1)a cavity mode of frequency ωc is coupled with strength gm,n to a transmon qubit whose unperturbed eigenstates can be written in terms of Mathieu functions (25). Here, we simply denote them by ∣n⟩ and their eigenenergies by ħωn. The cavity is represented using the annihilation(creation) operator a(a) and is driven by a monochromatic field of strength ϵ and frequency ωd. To model environmental noise, we use the Lindblad master equation (27)tρ=iħ[H,ρ]+(nc+1)κ D(a)ρ+ncκ D(a)ρ+γϕD(bb)ρ+(nt+1)γ D(b)ρ+ntγ D(b)ρ(2)where b is the ladder operator acting on the transmon and is defined by b=n=0n+1n n+1. The thermal occupations of the transmon and cavity baths are denoted by nt and nc, respectively, while γ and κ are the intrinsic transmon and cavity relaxation rates, and γϕ is the intrinsic transmon dephasing rate.

At sufficiently low drive powers, the transmon is confined to its ground state. The system can be treated as a single Duffing oscillator with a Kerr nonlinearity (28) whose Hamiltonian and master equation can be written as (note S1)H˜=ħωcaa+12ħKaaaa+ħϵ˜(aeiωdt+aeiωdt)(3)tρ=iħ[H,ρ]+(nc+1)κD(a)ρ+ncκD(a)ρ+κϕD(aa)ρ(4)

In this simplified model, the dispersive coupling with the transmon shifts the cavity frequency to ωc and introduces a Kerr nonlinearity K. The thermal occupation of the cavity bath is denoted by nc, and the cavity relaxation and dephasing rates are represented by κ and κϕ, respectively.

A further simplification can be made with the mean-field approximation, i.e. we assume the cavity to be in a coherent state ρ = ∣α⟩⟨α∣ and substitute this into ∂tα = Tr (a ∂tρ). We obtain the classical equation of motion in a frame rotating with the drivetα=(κ˜+i(ω˜cωd)+iKα2)αiϵ˜(5)and find the steady-state cavity amplitudeα=iϵ˜κ˜+i(Kα2+ω˜cωd)(6)

At weak drive powers, the occupation of the cavity is low, so the nonlinear term K∣α∣2 in Eq. 6 vanishes, and the standard Lorentzian response is obtained. However, when the number of photons in the cavity approaches the saturation number nsat = ∣ωc − ωd∣/K ∼ ∣α∣2, the nonlinearity becomes substantial, and the equation of motion may admit two stable steady-state solutions. The system enters the bistable regime.

In the mean-field approximation, the lifetimes of these states are infinitely long, but if fluctuations are taken into account, these states become metastable and the system may undergo rare escape events, switching from one state to the other. Both metastable states may coexist with each other over a range of drive amplitudes, and the time taken for the system to reach a steady state will be determined by their lifetimes. These lifetimes can be much greater than the lifetime of the cavity, and this gives rise to the phenomenon of critical slowing down. If an appropriate thermodynamic limit is taken, this time diverges and the model produces a first-order dissipative phase transition in which the two phases may coexist only at a single drive amplitude (13, 14, 29).


Cavity response in the bistable regime

We now show evidence of such a first-order phase transition in our system, containing only two degrees of freedom. We measure the signal transmitted through the cavity as a function of driving frequency (ωd) and power (Prf), as shown in Fig. 1A. We find that at low power, the cavity line is dispersively shifted to ωr/2π = 10.4960 GHz and has the Lorentzian shape that is typical of linear response. As the driving power increases, the lineshape shifts to lower frequencies and nonlinear features appear. The effective Kerr nonlinearity of the cavity is found to be K = −0.4221 MHz, and its relaxation rate is κ/2π=1.040 MHz. Above Prf = −29 dBm, a dip in the transmitted signal is observed. This indicates the presence of the bistable regime and is due to destructive interference between the two metastable states of the cavity (11, 12). The boundaries of this regime are modeled using the mean-field equations of motion derived from the Duffing approximation and are shown by the red lines. The bistable regime emerges just below the resonance frequency at a drive power of Prf = −35 dBm and opens up over a wider range of frequencies as the drive power increases.

Fig. 1 Phase transition in a system with two degrees of freedom.

(A) Transmission spectroscopy of the cavity for a range of drive powers and frequencies. At low drive powers, we observe a resonance at ωd/2π = 10.4960 GHz. Nonlinearity causes the resonance to shift as power increases. The boundaries of the bistable regime are modeled using the mean-field equations of motion of the Duffing oscillator (red lines). (B) Transmitted signal as a function of drive power at ωd2π=10.4925GHz [vertical dashed line in (A)]. The measured transmission (black dashed line) is compared with the results of a Duffing model in green, for which we plot the cavity amplitude divided by the root of the photon saturation number (a/nsat). Insets I, II, and III show the Wigner function of the steady state of the cavity according to the Duffing model at the marked power, confirming that the transition is associated with a transfer of probability between two regions in phase space: the bright and dim states. a/nsat is also displayed for a range of saturation photon numbers: The transition becomes sharper as nsat increases toward the thermodynamic limit. (C) Asymptotic decay rate as a function of power for different values of nsat. The rate at which the system approaches the steady state drops orders of magnitude below its natural relaxation time κ˜, showing that the bistable regime is associated with critical slowing down. The inset shows that this slowdown increases exponentially as the thermodynamic limit is approached.

Figure 1B shows the signal transmitted through the cavity (black crosses) as a function of the drive power at ωd/2π = 10.4925 GHz. At this drive frequency, we calculate a saturation photon number of nsat = 8.4. We observe a sudden change in transmission at Prf = − 25 dBm in which the cavity switches from a low-amplitude (dim) state to a high-amplitude (bright) state. This transition is accurately modeled by the Duffing master equation (Eq. 4), the results of which are displayed by the green line.

We can make a connection between this behavior and the theory of phase transitions by defining the thermodynamic limit in which the saturation photon number nsat goes to infinity. This is achieved by rescaling the drive and nonlinearity according to ϵλϵ and KK/λ, which, in turn, gives nsat → λnsat (9, 13, 30). The simulated cavity amplitude is displayed for a range of values of nsat, and, as expected, we observe that the transition becomes sharper as the system moves toward the thermodynamic limit, which is typical of first-order phase transitions.

We are interested in exploring the system dynamics within this bistable regime, where the steady state consists of a mixture of bright and dim states. We can form an effective master equation by writing the state of the system asρ(t)=pb(t)ρb+pd(t)ρd(7)where ρb and ρd are the bright and dim states, and pb and pd are their occupation probabilities. We can then write a simple rate equation for evolution of the state(tpbtpd)=(ΓbdΓdbΓbdΓdb)(pbpd)(8)in which occasional fluctuations allow the system to switch between states at the rates Γb → d and Γd → b (12, 20). The system reaches a steady state when the occupation probabilities have reached equilibrium. The approach to this steady state is governed by(pbpd)=1Γad(ΓdbΓbd)+AeΓadt(11)(9)where the coefficient A is determined by the initial system conditions. Γad is referred to as the asymptotic decay rate and it is given byΓad=Γdb+Γbd(10)

Γad is calculated by extracting the gap in the Liouvillian superoperator (14, 22), derived from the Duffing master equation (Eq. 4). Figure 1C shows the asymptotic decay rate as a function of drive powers for different photon saturation numbers. We find that, within the bistable regime, the asymptotic decay rate drops far below the cavity decay rate (22). This effect is known as critical slowing down and is characteristic of the phase transition. Furthermore, the asymptotic decay rate is exponentially suppressed as nsat increases, indicating that ever larger fluctuations are required to cause switching between the states.

Critical slowing down

Measurements of critical slowing down are performed by recording the transient response of the cavity when a step function drive pulse is applied. Figure 2 shows the average cavity response outside (A) and inside (B) the bistable regime. Figure 2A shows the response at the cavity resonance at Prf = −40 dBm with the transmon initialized in either the ground state (blue line) or the first excited state (brown line). The time scale over which the cavity responds shows a clear dependence on the transmon state. When the transmon starts in the ground state, the cavity reaches equilibrium over a time scale 2 Tκ = 2/κ = 0.29 μs, set by the cavity relaxation rate κ, whereas when the transmon is initialized in the first excited state, the drive is initially off resonant with the cavity and the transmon must relax over a time T1 before the system can reach equilibrium.

Fig. 2 Averaged transient response of the cavity outside and inside the bistable regime.

The inset in (A) shows transmission spectroscopy of the cavity for a range of drive powers and frequencies during the first cooldown. The mean-field Duffing limits of the bistable regime are displayed in red, and the locations at which the data in (A) and (B) were taken are indicated by the white dots. (A) The cavity is driven at the low-power resonance ωd/2π = 10.4960 GHz and Prf = −40 dBm. The signal in blue (brown) is the transient response measured with the qubit initialized in its ground (first excited) state. The transient response is governed by the time scale 2Tκ = 0.29 μs if the transmon is in the ground state, whereas it is governed by T1 = 2.89 μs if the transmon is in the first excited state. (B) Transient responses for different initial qubit states in the bistable regime at Prf = −21 dBm. The inset shows spectroscopy of the cavity at this drive power with the drive frequency ωd/2π = 10.4898 GHz indicated by the dashed line. The transient response is divided into two parts. There is an initial fast response with a time scale ranging from Tκ to T1 depending on the initial transmon state, followed by a slow decay toward steady state over a time scale Ts = 73.2 μs, obtained from an exponential fit, which is much longer than both the transmon and cavity lifetimes. This critical slowing down allows us to distinguish the transients for different transmon states for more than 100 μs.

Figure 2B shows that the dynamics changes significantly if the system is driven at higher powers. The cavity is now driven at ωd/2π = 10.4898 GHz and Prf = −21 dBm. The inset displays the spectrum at this drive power and the dashed line indicates the drive frequency, which is chosen such that the system is in the bistable regime as signaled by the dip in transmission. The cavity response is now governed by multiple time scales. Initially, there is a fast rise in the cavity transmission over a time ranging from Tκ to T1 depending on the initial state of the qubit. However, after this fast response, we observe critical slowing down: a gradual decay toward equilibrium over a time much longer than both the cavity and qubit lifetimes. We label the time constant over which the system reaches equilibrium as Ts = 1/Γad. By initializing the transmon in a range of initial states, we show that the cavity retains a memory of the initial transmon state for over 100 μs, indicating our proximity to a phase transition.

Next, we characterize the dependence of the critical slowing down Ts on both driving frequency and power within the bistable regime, and we model our findings using a range of approaches. We can either exploit the analytical theory of quantum activation for a Duffing oscillator, as provided in (21), or we can use the master equations for the Duffing and GJC models. In this latter approach, we rewrite the GJC master equation as ∂tρ = ℒρ where ℒ is the Liouvillian superoperator. The eigenbasis of ℒ can be used to express the state of the system asρ(t)=ncne(Γn+iωn)tρn(11)where ℒρn = −(Γn + iωnn. At long times, the state will be dominated by the steady-state ρss and the asymptotically decaying state ρad. In this limit, we write the state asρ(t)=ρss+cadeΓadtρad(12)

While for the Duffing model, it is possible to extract Γad by diagonalizing the Liouvillian, due to the larger Hilbert space size of the GJC model, Γad is extracted by integrating the master equation (Eq. 2) for a sufficiently long time. We are hence able to compare the experimentally attained values of critical slowing down time with simulated values. Figure 3A shows the measured valued of Ts (red circles) as a function of the drive frequency ωd along the dashed line in the spectroscopy inset, which is located at a drive power of Prf = −17 dBm. We find that Ts reaches a maximum close to the dip in the cavity spectrum. This is expected, since the dip is a key signature of the coexistence of the two phases in our transition.

Fig. 3 Measuring and modeling the critical slowing down time.

(A) Critical slowing down time Ts in the bistable regime as a function of driving frequency at Prf = − 17 dBm during the second cooldown, during which the qubit frequency had shifted to 8.7965 GHz and the low-power resonance of the cavity had shifted to 10.4761 GHz. The red points represent the experimental data, which we compare with the results of master equation calculations applied to the Duffing oscillator (blue line) and the GJC model with transmon dephasing (green line) and without (purple line). We also display the results of previous analytical theory of switching rates for the Duffing oscillator (orange line) (21). At this power, both the master equation and the analytical calculations qualitatively reproduce the experimental values of Ts. The horizontal dashed line in the inset shows the location of our measurements within the overall cavity spectrum. (B) Maximum value of Ts for different drive amplitudes (red points). These data were collected along the diagonal dashed line in the inset of (A). As the drive power increases beyond −17 dBm, Ts reaches a saturation at a value of ≈100 μs, which is consistent with the simulations based on the GJC model with transmon dephasing (green line). Removing the dephasing by setting γϕ = 0 (purple line) does not change the power at which saturation occurs but it does raise the upper limit on Ts. Meanwhile, analytical (orange line) and master equation (blue line) calculations with the Duffing approximation predict that Ts rises exponentially with drive amplitude, as can be seen using the logarithmic scale of the inset.

We compare our measurements with master equation simulations applied to a Duffing oscillator model (blue line) and to the GJC model with (green line) and without transmon dephasing (purple line). We also compare our data to results attained from the analytical theory of quantum activation for a Duffing oscillator (orange line) (21). At this drive power, we find that all of our models give at least qualitative agreement with the measured dependence of critical slowing down on frequency, but only the GJC model is able to model this effect quantitatively.

However, if we plot how the maximum value of Ts varies with the amplitude of the drive, as shown in Fig. 3B, we observe a significant divergence between our data and the values attained using the Duffing model. Whereas the theory of the Duffing oscillator predicts that Ts should increase exponentially with the drive, we instead observe that, at sufficiently strong drive amplitudes, Ts saturates. To account for this difference, we require the full GJC model. We find that the master equation predicts the same saturation in Ts as found in experiment when we explicitly include the transmon in the simulation.


The system of two strongly coupled oscillators is governed by essentially different activation dynamics. To shed light on the dissimilarity between the Duffing model and the GJC model, we examine the switching rates Γb → d and Γd → b more closely. To obtain these rates, it is necessary to first find the steady state and extract the occupation probabilities pb(t → ∞ ) and pd(t → ∞ ). We can then calculate the switching rates Γd → b and Γb → d according toΓd(b)b(d)=pd(b)(t)Γad(13)

The occupation probabilities can be extracted using the asymptotically decaying and steady states (note S2). The resulting rates are displayed in Fig. 4B, which shows Γd → b and Γb → d as a function of frequency at a driving power of Prf = −14 dBm for both the Duffing and the GJC model. Whereas Γd → b is in close agreement between the two models, Γb → d is significantly greater in the GJC model. This limits the critical slowing down time according to Eq. 10 and leads to the saturation seen in Fig. 3B. It also indicates that the bright state has a shorter lifetime in the GJC model. Prior work suggests that the instability of the bright state increases with the nonlinearity of the ladder of states in the vicinity of the wave packet (17). The instability of the bright state in the GJC model may be due to extra nonlinearity that is present when the transmon becomes excited.

Fig. 4 Relating the steady-state Wigner to the switching rates.

(A) Steady-state cavity Wigner functions produced using the GJC model at Prf = −14 dBm. At ωd/2π = 10.4709 GHz the steady state consists mainly of the dim state, which corresponds to the peak near the origin. However, as we increase the drive frequency, the occupation of the bright state increases as well. At ωd/2π = 10.4720 GHz, the occupations of the bistable states are approximately equal, and, at ωd/2π = 10.4723 GHz, the bright state is dominant. (B) Switching rates between metastable states as a function of driving frequency. Whereas Γd → b is similar in both the Duffing and the GJC models, Γb → d is significantly different. In the GJC model, Γb → d is much greater compared to the Duffing model. This explains the saturation we observe in Ts.

The Wigner function for three cavity steady states using the GJC model is shown in Fig. 4A. We observe that, at low drive frequencies, the dim state is the main contributor to the overall state of the system; this is consistent with Γb → d dominating over Γd → b. At higher drive frequencies, the reverse is true, while at some intermediate drive, frequency the two bistable states are equally occupied and the switching rates are balanced.

In summary, we have explored the rich quantum activation dynamics happening at an intermediate driving regime in cQED. We have observed a phase transition in our system, which contains only two degrees of freedom: a transmon qubit coupled to a microwave cavity. A key signature of this transition is critical slowing down, in which the time taken for the system to reach a steady state can extend far beyond the natural lifetime of the qubit or cavity. We have measured the slowdown time for a wide range of powers and frequencies, and we have compared our results with simulations. We found that the transition and its associated critical slowing down are well modeled by the Duffing approximation at low drive powers. However, at higher drive powers, we observed a saturation in the critical slowing down time, which can only be captured by the full GJC model.

It is known that in this regime the transmon becomes highly excited and starts to participate in the dynamics (31) so it is no longer valid to apply the Duffing approximation. An accurate model must include the quantum fluctuations of the qubit and the resulting destabilization of the bright state that this causes. Currently, there exists no analytical theory for the switching rates in the bistable regime of a cavity coupled to spins or multilevel systems, and this suggests that one avenue of future work could focus on extending the existing theory for the Duffing oscillator to these models. Moreover, we link the critical slowing down to the switching rates between the two metastable states and show how these differ for the GJC and Duffing models. This experiment is a powerful demonstration of the versatility of superconducting circuits, showing that even with few degrees of freedom, it is possible to explore rich nonlinear physics and phenomena such as dissipative phase transitions.


A two-port Al microwave cavity holds a lithographically patterned Al transmon qubit, which is fabricated on a sapphire substrate. The transmon qubit consists of two Al pads of dimensions 350 μm by 450 μm connected by an Al/AlOx/Al Josephson junction, which is patterned using standard e-beam lithography and double-shadow evaporation techniques. The cavity is thermally anchored to the 10 mK plate of a dilution refrigerator. Input signals are heavily cryogenically attenuated to reduce thermal noise, and measurements of the signal transmitted through the cavity are made via cryogenic circulators and a low-noise high-electron-mobility transistor amplifier, with the signal lastly being recorded as a voltage with an analog-to-digital converter.

Measurements were collected during two successive cooldowns. During the first cooldown (C1), spectroscopy of the transmon reveals its lowest transition to be ω01/2π = 9.1932 GHz with anharmonicity α/2π = −203.6 MHz. The bare resonance frequency of the cavity is ωc/2π = 10.4263 GHz and its quality factor is found to be Q = 7900. The relaxation and dephasing times of the transmon are T1 = 2.89 μs and T2 = 2.37 μs, respectively. During the second cooldown (C2), the system is described by a GJC model with the following parameters: ωc/2π = 10.423 GHz, g0/2π = 295 MHz, κ = 1.432 MHz, γ = 33 kHz, γϕ = 1 kHz, nc = 0.01, and nt = 0.02. The eigenstates of the transmon were produced using a Josephson energy of EJ/2π = 46.7 GHz and a charging energy of EC/2π = 221 MHz (25). Applying the Duffing approximation to this system, we find ωc/2π=10.4761 GHz, K/2π = − 0.152 MHz, κ/2π=1.387 MHz, κϕ/2π=1.02 Hz, and nc=0.0100.


Supplementary material for this article is available at

This is an open-access article distributed under the terms of the Creative Commons Attribution license, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


Acknowledgments: Funding: This work has received funding from the EPSRC under grant nos. EP/J001821/1, EP/J01350/1, EP/M013243/1, and EP/K003623/2. Author contributions: G.T. performed the experiments; P.B. provided modeling and support on the theory; J.R. fabricated the sample; A.D.P. optimized the FPGA measurement setup; M.E. provided support for analyzing the results; T.K.M. contributed to understanding the theory; E.G., P.J.L., and M.H.S. supervised the project. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. The data underlying the paper can be found at:

Stay Connected to Science Advances

Navigate This Article