Quantum phase transitions in the driven dissipative Jaynes-Cummings oscillator: from the dispersive regime to resonance

We follow the passage from complex amplitude bistability to phase bistability in the driven dissipative Jaynes-Cummings oscillator. Quasidistribution functions in the steady state are employed, for varying qubit-cavity detuning and drive parameters, in order to track a first-order dissipative quantum phase transition up to the critical point marking a second-order transition and spontaneous symmetry breaking. We demonstrate the photon blockade breakdown in the dispersive regime, and find that the coexistence of cavity states in the regime of quantum bistability is accompanied by pronounced qubit-cavity entanglement. Focusing on the role of quantum fluctuations in the response of both coupled quantum degrees of freedom (cavity and qubit), we move from a region of minimal entanglement in the dispersive regime, where we derive analytical perturbative results, to the threshold behaviour of spontaneous dressed-state polarization at resonance.


I. INTRODUCTION
The Jaynes-Cummings (JC) oscillator is an archetypal source of intricate quantum nonlinear dynamics arising from the coupling of a quantized electromagnetic mode inside a resonator (cavity) to a two-level system (qubit) [1]. The behaviour of quantum nonlinear oscillators has been a subject of intense theoretical investigation (for an overview see Chapter 7 of [2]) providing at the same time the basis for numerous experiments in cavity and circuit quantum electrodynamics (see for example [3] where the extended JC oscillator is driven out of equilibrium in the presence of dissipation). In addition, controlled light-matter interaction has shifted the center of interest in phase transitions from condensed matter to quantum optics.
Amongst the most discussed light-matter quantum phase transitions in the literature are the Dicke phase transition [4,5], which is explicitly dissipative (with cooperative resonance fluorescence as its driven variant [6]), and the laser which exhibits a second-order phase transition out of equilibrium [7]. In comparison to those, however, the driven JC model is fundamentally different as it deals with the interaction of one field mode with one quantum object, the qubit, necessitating the reappraisal of the rôle of quantum fluctuations and a different definition of the thermodynamic limit [8][9][10][11]. The open coherently-driven dissipative qubit-cavity system yields a bistable response where quantum fluctuations are responsible for switching between two metastable states that are long-lived in relation to the characteristic cavity and qubit decay times [12]. The √ n splitting of the JC energy levels is a unique feature determining the nature of bistability both at resonance, where the cavity and qubit bare * t.mavrogordatos@ucl.ac.uk frequencies coincide, and in the dispersive regime, where the cavity and qubit are strongly detuned in relation to their dipole coupling strength [10,13,14].
At resonance, the mean-field nonlinearity diverges for zero photon number, while in the dispersive regime bistability builds up in a perturbative fashion with no associated threshold, unlike the laser. The perturbative approach becomes inadequate when the qubit participates actively in the bistable switching for stronger driving [15]. In that regime, the system response comprises an average over spontaneous switching between the metastable mean-field steady states where both the cavity and qubit are significantly excited. For complex amplitude bistability switching occurs between a dim (with lower n) and a bright (with higher n) state, while in phase bistability both states have the same magnitude and opposite phases following a transition from a discrete to a continuous spectrum in the system quasienergies [9,10,16]. In first-order dissipative phase transitions for the cavity field, weaker coupling implies a bigger photon number required for the nonlinearity to manifest itself, yielding a response which is a non-analytic function of the drive [2,10].
Motivated by the current experimental and theoretical interest in the response of quantum nonlinear oscillators, in this letter we track amplitude bistability, from its origin in the dispersive regime, up to a critical point at resonance, where phase bistability takes over. Mean-field results guide us to extract the relevant scaling parameters used to define the "thermodynamic limit" for this driven resonator in which the number of photons is not conserved. We present contour plots of quasidistribution functions for the cavity field, showing the passage from amplitude to phase bistability, and invoke the entanglement entropy, calculated via the reduced qubit density matrix, in order to demonstrate the active participation of both quantum degrees of freedom in the emerging bi- Average cavity photon number a † a as a function of the drive parameters. The dashed line indicates the driving frequency selected for (a). Parameters: g/δ = 0.14, 2κ/γ = 12, g/γ = 3347, n scale = 12.68. modality. We show further that, closer to resonance, enhanced multi-photon transitions appear for weak cavity excitation, followed by a breakdown of the photon blockade with stronger driving [10,11].

II. THE DISPERSIVE JC MODEL
The JC Hamiltonian describes the interaction between a single resonant cavity mode and a qubit; however, it does not account directly for the coupling to the environment which is included only in the formulation of the Master Equation (ME) [17]. After adding dissipation in a frame rotating with the frequency ω d of the coherent driving field, the Lindblad ME for the reduced system density operator ρ (for a system Hamiltonian in the rotating wave approximation and setting = 1 for convenience) reads [9,18] In the right hand side of the above equation, ∆ω c,q = ω d − ω c,q are the detunings of the cavity resonance frequency ω c and the qubit bare frequency ω q from the frequency of the drive, coupled to a resonant cavity mode with photon annihilation and creation operators a and a † , respectively. The inversion operator σ z is related to the raising (lowering) operators σ + (σ − ) for the qubit with two states |g (ground), |e (excited) via σ z = 2σ + σ − − 1. The cavity mode is dipole-coupled to the qubit with strength g, while the classical coherent field (with very high photon occupancy) is coupled to the resonant cavity mode with strength ε d (also called drive amplitude). The cavity field is also coupled to a Markovian thermal bath at zero temperature, which induces a photon loss rate of 2κ. In addition to photon dissipation, there is also spontaneous emission to modes different than the resonant cavity mode, with rate γ. The strongly dispersive regime with weak spontaneous emission is defined through a qubit-cavity detuning such that δ ≡ |ω c − ω q | g 2κ γ. In our results we show complex amplitude bistability for the intracavity field in the following region of the drive phase space: 0 ≤ ∆ω c ≤ g 2 /δ and ε d < g 2 /δ < g (with γ/(2κ) ≈ 0.1). At resonance, where δ = 0, phase bistability is associated with a particular point in the phase space: (∆ω c = 0, ε d = g/2), and a threshold behaviour.  1) for the reduced cavity density matrix (Panel II). Parameters: g/δ = 0.14, 2κ/γ = 12, g/γ = 3347, n scale = 12.68.
Throughout our work, we solve numerically the Lindblad ME. We perform numerical simulations (employing the open-source Matlab software Quantum Optics Toolbox ) in a truncated Hilbert space for an initial product state with the qubit in the ground state and the cavity in a Fock state with zero photons [i.e. the initial ground state ρ(0) = (|n = 0 n = 0|) ⊗ (|g g|)]. The steadystate results are independent of the initial state, while convergence with respect to the number of photon levels has been ensured. The validity of the rotating wave approximation in the dispersive regime we are considering is checked against Fig. 1 of Ref. [19] and Fig. 1 of Ref. [20]. The former shows that the counter-propagating terms can be omitted without affecting the physical picture even for the maximum cavity-qubit detuning considered, since the coupling strength g remains sufficiently smaller than the bare cavity frequency (g/ω c ≈ 0.03), while the latter demonstrates good agreement between theory and experiment in a similar parameter regime.
In Fig. 1 [see both frames (a) and (b)] we depict the the intracavity photon field (in the coherent state space with |α ≡ |x + iy ) within a drive region where the quantum fluctuations are responsible for the deviation from the mean-field predictions. For low driving strengths, |α|, a † a and | a | coincide and the cavity is in the dim state, resembling a vacuum state with Gaussian distribution. With increasing ε d /κ, the bright state accrues probability resulting in the coherent cancellation we observe at the point A. As we follow the curve for | a |, probability transfers from the dim to the bright state crossing the boundary of a first-order dissipative quantum phase transition. The complex steady-state semi-classical intracavity amplitude α obeys the equation [16] withκ = κ − i∆ω c , which predicts two metastable states (dim and bright) and one unstable state that vanishes in the presence of fluctuations. In contrast to the meanfield prediction, the curve depicting | a | does not exhibit any bistability. Quantum fluctuations out of equilibrium manifest themselves through the absence of a Maxwell construction, since the line | a | does not cut the semiclassical curve in two equal areas. Furthermore, the curve for a † a does not exhibit the coherent cancellation dip, which is hence solely a quantum phase effect. These considerations hold also for the driven Duffing oscillator (the reader is referred to Fig. 1 of [18]) as a consequence of non-constant diffusion coefficients in the corresponding Fokker-Planck equation; here, however, we cannot formulate such an equation due to the active participation of the qubit [8,9].

III. PERTURBATIVE EXPANSION FOR THE CAVITY BISTABILITY
At first we will examine the birth of dispersive amplitude bistability for a driving frequency in the region ∆ω c g 2 /δ and weak drive strength. In the strongly dispersive limit the presence of the small term g/δ precludes the divergence of nonlinearity at low intracavity field amplitudes. When the length of the Bloch vector is conserved, in the absence of spontaneous emission (γ = 0), The driving field has a phase difference of −π/2 with respect to the drive in Fig. 1(a), leading to a rotation of the distribution by that angle in the x-y plane, as expected from Eq. 2.
We apply the dispersive transformation to diagonalize the JC Hamiltonian, generating the term δ 1 + n s /n scale , where n s = a † a + σ + σ − is the operator of system excitations (see [14] and [20] for more details). Expanding Eq. (3) to the lowest order in n s /n scale (with |α| 2 the semiclassical analogue of n s ) yields in agreement with Eq. (32) of [10]. In the dressedcavity Duffing approximation we retain only terms up to the second-order in n s /n scale , and the reduced Hamiltonian acquires the quartic correction (g 4 /δ 3 )σ z a † 2 a 2 (in agreement with the semi-classical prediction of Eq. 4 for σ z = −1). The series expansion also renormalizes progressively the driving phase space such that the effective drive strength and frequency are functions of the system operators [14]. We can then derive the Wigner function for the effective dressed Duffing oscillator [18], calculated via the generalized P -representation [21][22][23] where 0 F 1 (λ; x) and 0 F 2 (λ, µ; x) are generalized hypergeometric functions of the variable x, with parameters λ, µ.
The effective detuning ∆ω c = ∆ω c + (g 2 /δ)σ z − (g 4 /δ 3 )(2σ z + 1) accounts for the correction by the dispersive shift and higher-order terms. We note that the perturbative expression (5) reproduces the Gaussian form of the distribution function corresponding to a vacuum state, W = (2/π)e −2|α| 2 , for very low driving strengths ε d /κ, and allows us to track the progressive participation of the various nonlinear terms arising from the hypergeometric function 0 F 1 . It is therefore more instructive to write a perturbation series expansion for the numerator: c · (c + 1) · . . . · (c + m − 1), showing explicitly the development of nonlinearity for increasing drive strength. In the regime where n s /n scale ≈ 1, the Duffing approximation breaks down as the qubit vector becomes increasingly entangled to the cavity mode, moving towards the equatorial plane in the Bloch sphere representation [16]. In Fig. 2(a) we depict the cavity field distribution in the absence of entanglement with the qubit. The excitation pathways 'flow' around the nodes of the Wigner function in a spiral-like fashion, as the departure from the Gaussian form becomes more apparent. These perturbative distributions approximate very well the exact ME results, in which the qubit is included as an independent degree of freedom [see Fig. 2(b)]. Hence, the agreement verifies the fact that the qubit participates only in dressing the cavity with a Kerr term depending on σ z = σ z = −1. Treating σ z as a constant of motion when solving Hamilton's equations (for arbitrary excitation but for a short time scale in comparison to γ −1 ) underlies the method followed by the authors of Ref. [20], who provide a semiclassical expression for |α| 2 in the dispersive regime. The cavity response is described therein by a skewed Lorentzian curve with an amplitude dependent frequency shift, approaching the value σ z g 2 /δ for |α| → 0, and restoring the linear regime limit [24]. Conversely, for large intracavity amplitudes and δ = 0, their expression yields the drive detuning corresponding to the two resonant paths, as revealed by the nonlinear equation that explains the split Lorentzian response at resonance. This distinct split has a direct relation to the origin of phase bistability, namely the formation of two quasiindependent excitation ladders with a vanishing connection between them [10]. The region of high powers in the strongly dispersive regime can be accessed for nondemolition qubit readout with ∆ω c ∼ κ [20,25,26]. Dynamical mapping of the qubit to the photon states has been also proposed in [27] as an optimized protocol, exploiting the first-order phase transition by means of which the photon blockade breaks down (see also the following section).

IV. PHASE TRANSITION CROSS-OVER
We will now delineate the defining features of complex amplitude bistability with the aim of approaching phase bistability, the occurrence of which signals the ultimate difference between the effective Duffing and the full JC nonlinearity. For that purpose we plot the Q function in the steady state, Q(x + iy) = (1/π) x + iy|ρ c |x + iy , for the reduced cavity matrix ρ c and the coherent state |x + iy . The region of coexisting states with probabilities of the same order of magnitude marks the boundary in the phase space of the drive where quantum fluctuations induce equiprobable transitions between the metastable states [12,28]. The resulting critical slowing down is a direct consequence of nonlinear dynamics (see Chapter 5 of [8]) and the departure from a Gaussian probability distribution. As the authors of Ref. [12] note, "bistability is a macroscopic phenomenon reached in the limit n scale → ∞". In our case, the bimodal distributions identify distinct states that are long-lived on the time scale γ −1 (and consequently on the scale (2κ) −1 ) even for n scale = 12.68. Bimodality is depicted in Fig. 3 [Panel I, frames (a)-(c)], associated with maximal qubitcavity entanglement, as we will see later on. The bright state is quadrature-squeezed along the mean-field direction (in a similar way to resonance fluorescence [8,29]), another display of the JC nonlinearity in this regime.
In the Panel I of Fig. 3 we plot quasidistribution functions showing coexistent metastable states in a region where the qubit is significantly excited and approaches progressively the equator in the Bloch sphere. For increasing drive amplitude ε d /κ we observe a growing separation of the two metastable state distributions followed by a change in their orientation [frames (a)-(c)]. The quasidistribution function in Fig. 3 (d) of Panel I illustrates a precursor of phase bistability for ε d = g/2, lacking nevertheless complete symmetry with respect to the horizontal axis (and hence having peaks of unequal height) because δ = 0. We build upon this theme in the Panel II of Fig. 3, where we track the emergence of phase bistability for decreasing values of δ/g, and ∆ω c = 0, ε d = g/2. As δ/g → 0, nonlinearity is triggered by lower photon numbers and the two peaks ap-proach each other (for the same values of ε d /κ), while complete symmetry is restored only when δ = 0.
We proceed now to the study of entanglement as a measure of the joint participation of both quantum degrees of freedom, employing the von Neumann entropy for the reduced qubit density matrix ρ q = Tr c ρ in the steady-state (where Tr c denotes the partial trace over the cavity field states), defined as S q = −Tr[ρ q ln ρ q ] = − i=1,2 λ i ln λ i . The eigenvalues λ i of the reduced qubit matrix ρ q = (ρ gg , ρ ge ; ρ * ge , ρ ee ) are given by the expression [30,31]: The entropy S q quantifies the entanglement between the two quantum oscillators, assessing the purity of the reduced quantum state for the qubit in the steady state (following the evolution of the open system from an initial pure state). In the dispersive regime there is still appreciable entanglement between the cavity and qubit despite their strong detuning. It has recently been shown that entanglement is also present in the linear region [32], which we have neglected when setting σ z = −1 in our analytical mapping to the Duffing oscillator. The entanglement entropy tracing a first-order phase transition in the drive phase space is shown in Fig. 4(a). From the linear region, where entanglement is very weak [light blue region in Fig. 4(a) appearing at ∆ω c = g 2 /δ], we move to the nonlinear regime where the maximum shifts to the left with a very steep drop, in a similar manner to the average photon number n = a † a , due to the presence of growing amplitude bistability [compare to Figs. 1 and 3(a) of [10]]. In Fig. 4(b) we plot the second-order correlation function for zero time delay, g (2) (τ = 0), defined via the relation g (2) (0) = n(n − 1) /( n 2 ), in order to reveal the effect of quantum fluctuations. The peak of quantum correlations is shifted relatively to the entanglement entropy maximum, with the two curves (blue and red, respectively) intersecting closer to the position of the coherent cancellation dip in the cavity amplitude | a | [similar to the point A in Fig. 1(a)] and the pseudospin projection | σ − |. The aforementioned dip has a purely quantum origin at zero temperature, which explains the amplification of quantum fluctuations in that region [18,33]. On the other hand, the maximum of the von Neumann entropy occurs at the frequency where the dim and bright state distributions attain peaks of equal height, as we can observe in Fig. 5(a). When δ/g → 0 the system response becomes highly nonlinear for low drive strengths, as n scale decreases. We observe enhanced resonant multi-photon transitions [inset of Fig.  5(b)] gradually disappearing in the region of high drive strengths [main panel of Fig. 5(b)]. This phenomenon is referred to as breakdown of the photon blockade [see Figs. 2(a) and 5(a) of [11], and [10] for an extensive discussion at resonance -δ = 0] accompanied by the appearance of amplitude bistability [see the Q function plot in the bottom inset of Fig. 5(b)]. (a) Growing entanglement entropy Sq as a function of the driving frequency for g/δ = 0.14 (with n scale = 12.68) and five equispaced drive amplitudes in the range ε d /γ = [20,100]. The red curve corresponds to ε d /γ = 100 at the peak of which we plot the function Q(x + iy) for the intracavity amplitude distribution (inset). The green curve depicts the relative difference of the Q function peak values h1, h2 (with h1 > h2), defined as r = (h1 − h2)/h1. Changing excitation paths at resonance. Entanglement entropy Sq in the region of phase bistability (δ = 0) as a function of the driving frequency past the threshold ε d = g/2 (compare with Fig. 1 of Ref. [10] and see Eq. (7), in the absence of spontaneous emission). The insets depict the joint quasidistribution function Q(x + iy) for the marked points A, B and C. Parameters: 2ε d /g = 1.06, 2κ/γ = 500, g/γ = 3347.

V. PHASE BISTABILITY
Let us finally link the increasing entanglement entropy to the appearance of phase bistability past the threshold set by the critical point of the second-order quantum dissipative phase transition: (∆ω c = δ = 0, ε d = g/2). At resonance, the nonlinearity can be triggered by low photon numbers with a different scaling parameter, associated with a strong-coupling limit [10], as opposed to the strongly dispersive regime. Fig. 6 shows the development of a phase-bimodal distribution as we cross the line ∆ω c = 0, where the entanglement entropy has a local maximum. For growing drive strength, the entropy at point B increases and the two peaks of the Q function move further apart compared to their threshold position, always remaining symmetrical with respect to the horizontal axis.
At this stage, it is instructive to invoke for a final time the solution above threshold of the so-called neoclassical equations, i.e. the semiclassical equations that conserve the length of the Bloch vector [10] (in the absence of spontaneous emission) which are also combined to derive the steady-state expression of Eq. 3 in the dispersive regime. Neoclassical theory predicts a parity-breaking transition at resonance, according to the equation as well as a bistable qubit vector lying on the equatorial plane (ζ ≡ σ z = 0) with ν ≡ σ − = ±α/(2|α|). In that regard, phase bistability corresponds to maximally entangled states of the two coupled quantum degrees of freedom, in which the qubit polarization and the cavity field are not enslaved to the external drive, as already predicted by the mean-field analysis of Ref. [17].

VI. CONCLUSION
In this letter we have examined the interplay of qubit-cavity entanglement and cavity bimodality when connecting the dispersive and the resonance regimes in the driven dissipative Jaynes-Cummings model for varying qubit-cavity detuning. For the assessment of the cavity nonlinearity we have employed both the mean-field and the Master Equation treatment including quantum fluctuations. We have followed the change of the intracavity field quasidistribution functions from the strongly dispersive regime to the gates of a critical point related to a second-order quantum phase transition at resonance. We have also included in our discussion the complex amplitude bistability encountered in the driven dissipative Duffing oscillator, adopting a perturbative approach for weak driving fields. This is a region of minimal entanglement and very weak qubit involvement in the formation of the system nonlinearity, for the quantum description of which we have employed an analytical form of the Wigner quasidistribution function. The growing participation of both coupled quantum degrees of freedom marks the passage from a first-order to a second-order dissipative quantum phase transition.