Gauge-independent emission spectra and quantum correlations in the ultrastrong coupling regime of open system cavity-QED

1 Introduction

The intricate interactions between light and matter allow one to observe drastically different behavior depending on the relative strength of the light–matter coupling. In the weak-coupling regime, the losses in the system exceed the light–matter coupling strength, and energy in the system is primarily lost before it has the chance to coherently transfer between the matter and the light. Accessing this regime experimentally has allowed for breakthroughs in quantum technologies such as single-photon emitters [1], [2], [3], [4]. Beyond weak-coupling, in the strong-coupling regime the rate of decoherence is smaller than the rate of excitation exchange, allowing for the observation of vacuum Rabi oscillations: the coherent oscillatory exchange of energy between light and matter. The strong-coupling regime has helped initiate a second generation of quantum technologies [5, 6]. See Figure 1 for a simple schematic of a typical cavity-QED system with system-bath leakage.

Figure 1:

Schematic of a generic cavity-QED system. The optical cavity mode has quantized energy levels (in blue), with a decay rate κ. The matter system is a truncated TLS (in red), with a possible spontaneous emission decay rate γ. The two systems have a coherent coupling strength g. A coherent laser (in orange) drives the system with Rabi frequency Ω_d.

Around 2005, the “ultrastrong-coupling” (USC) regime was predicted for intersubband polaritons [7]. This regime is characterized not by still lower rates of decoherence, but by a coupling strength that is a comparable fraction of the bare energies of the system. The dimensionless parameter η = g/ω ₀ (i.e., the cavity-emitter coupling rate divided by the transition frequency) is used to quantify this coupling regime for cavity-QED. Typically, USC effects are expected when η ≳ 0.1, at which point the rotating wave approximation (RWA) used in the weak and strong regimes becomes invalid. Reported signs of USC emerged in 2009 with experiments involving quantum-well intersubband microcavities [8], achieving η ≈ 0.11. Terahertz-driven quantum wells have also demonstrated USC effects [9], and similar effects have been exploited to achieve carrier-wave Rabi flopping with strong optical pulses [10], [11], [12]. To date, many different systems have exhibited USC [13, 14]. Recently, using plasmonic nanoparticle crystals, η = 1.83 has been achieved, with potential to lead to η = 2.2 [15].

With experiments pushing the normalized coupling strength continuously higher, the interest in USC effects also continues to grow, helping to improve the underlying theories of light–matter interactions [16], even at arbitrarily high coupling strengths [17]. There have also been various predictions made about what novel technologies USC will bring about, including modifications to chemical or physical properties of various systems caused by their USC to light [7, 18], and the potential to create faster quantum gates and gain a high level of control over chemical reactions [13]. To push these advancements forward, it is essential to have a fundamental understanding of the physics involved with these systems and to accurately connect to experimental observables.

The cornerstone model in cavity-QED is a two-level system (TLS) interacting with a quantized cavity mode [19]. This model has been applied to atoms [20], [21], [22], [23], quantum dots [24], [25], [26], [27], and circuit QED [28], [29], [30], [31]. Outside the USC regime, this model is typically represented by the canonical Jaynes–Cummings (JC) Hamiltonian [32], which makes an RWA and can be easily diagonalized. In the USC regime, however, it is necessary to retain counter-rotating terms, giving rise to the quantum Rabi model (QRM) [13, 14, 33]. By detecting resonance fluorescence of light emitted from the cavity as quantified by the first-order degree of coherence correlation function (CF), the spectral content of these cavity-QED models can be explored, while the second-order intensity CF is fundamental to understanding the photon statistics as probed by intensity interferometry.

The main contribution of this work is to present a self-consistent and unambiguous way to model observables in the USC regime of open system cavity-QED. Apart from addressing the subtle (and unknown) effects of dissipation, and excitation, and input-output, we show the striking influence of modelling experimentally relevant observables such as the emission spectra and quantum correlation functions. We also show how and why the form of the system-bath interactions matters, yet the form is gauge-invariant, if—and only if—treated properly (in contrast to the usual master equation approaches) using gauge-invariant master equations. We show equivalence between dipole gauge and Coulomb gauge master equations, if one applies gauge corrections in a consistent way, and we also demonstrate the drastic failure of currently adopted master equations in the USC regime. Our framework and formalism, to the best of our knowledge, constitutes a first way to do this, and can thus be applied to a wide range of measurements in the USC regime for open systems.

2 Gauge invariance and system-reservoir interactions

It was recently shown that extra care is needed when constructing gauge-independent theories [34], for computing experimental observables for suitably strong light–matter interactions. This development started with a series of papers dealing with so-called gauge ambiguities in the USC regime [35], [36], [37]. As a U(1) gauge theory, different gauges in QED manifest in different representations of the Hamiltonian of a given system, but these should be unitarily equivalent and give rise to equivalent physical observables. Without proper care, gauge invariance of cavity-QED theories can break down when considering USC [38]. This is due to the truncation of the matter system’s formally infinite Hilbert space to the two lowest eigenstates in forming the TLS—only keeping an infinite number of energy levels formally preserves gauge invariance [39]. Consequently, previous model predictions in the USC regime can be ambiguous since the predictions are impacted by the choice of gauge. While this issue has been known in general for several decades [40], only recently was this specific problem presented as rather insurmountable [38]. However, the issue has been resolved by using a self-consistent theory at the system Hamiltonian level [37, 41], restoring gauge invariance to the theory for systems with a finite Hilbert space.

Despite this, additional subtleties occur in the USC regime regarding the interaction of the cavity-QED system with its environment. To connect to experiments, one also requires an input–output model of dissipation from the cavity to external modes, requiring an open-system model of cavity-QED. In the USC regime, complications arise with this input–output formalism associated with approximations typically made outside of the USC regime. These complications originate from the hybridization of light and matter that occurs in USC, and as such the quanta of excitations inside the cavity-QED system have different quasiparticle representations than the photons actually emitted from the system. Moreover, the separation of operators into light and matter components becomes highly gauge-specific in the USC regime, and proper care must be taken to ensure self-consistency.

To fully synthesize these considerations with the restoration of gauge invariance, we present a dissipative and gauge-invariant master equation model, which is required to properly describe experimentally-observable quantities arising from output channels of the cavity. Key experiments to probe such observables include resonance fluorescence and two-photon detection schemes, and we make a direct connection to both of these. We also show how previous QRM master equations in the USC regime are ambiguous in general as they produce gauge-dependent results for observables, and we show how to fix such problems. Moreover, our theories can be used to explore the precise form of the system-bath interactions, which in fact yield different experimental signatures in the USC regime.

3 Model

In the dipole gauge (namely, the multipolar gauge in the dipole approximation), we can write the system Hamiltonian, using the QRM, as (ℏ = 1)

where ω ₀ (ω _c) is the TLS (cavity) transition frequency, σ ⁺ (σ ⁻) is the raising (lowering) operator for the TLS, and a ^† (a) is the cavity mode creation (annihilation) operator; g is the TLS-cavity coupling strength. We take ω _c = ω ₀ throughout. In contrast to the Coulomb gauge, straightforwardly truncating the dipole in the light–matter interaction to a TLS subspace does not break gauge invariance in the dipole gauge [37]. Making an RWA on Eq. (1) (i.e., neglecting counter-rotating terms a ^† σ ⁺ and aσ ⁻, which do not conserve excitation number), yields the simpler JC Hamiltonian.

Outside of the USC regime, the usual approach to include dissipation is with a Lindblad master equation [42],

where ρ is the reduced density matrix. The dissipation term, L bare ρ = κ 2 D [ a ] ρ , is the Lindbladian superoperator where D [ O ] ρ = 2 O ρ O † − ρ O † O − O † O ρ and κ is the cavity photon decay rate. Since dissipation is usually dominated by cavity decay, we neglect direct TLS relaxation and pure dephasing [43, 44]. However, the theory of how to include TLS dissipation is discussed in Appendix A.4

The Lindbladian can be derived by following the typical approach in which one neglects the TLS-cavity interaction when considering the coupling of these systems to the environment [29]. However, when moving into the USC regime, this approach fails, and the Lindbladian must be derived while self-consistently including the coupling between the subsystems. For sufficiently strong subsystem coupling, transitions occur between dressed eigenstates of the full Hamiltonian rather than between eigenstates of the individual free Hamiltonians [43].

In the USC regime, the system has transition operators |j⟩⟨k| which cause transitions between the dressed eigenstates of the system {|j⟩, |k⟩}. To obtain these transitions for the cavity mode operator, we use dressed operators [43],

and x − = ( x + ) † , where the sum is over states |j⟩ and |k⟩, with ω _k > ω _j , C _jk =⟨j|Π_C|k⟩, and we neglect thermal excitation effects; Π_C is an operator which couples linearly to dissipation channel modes which we assume proportional to the cavity electric field operator such that Π_C = i(a ^† − a). We then replace L bare in Eq. (2) with L dressed ρ = κ 2 D [ x + ] ρ , to arrive at the dressed state (DS) master equation. One can also use a generalized master equation to capture coupling to frequency-dependent reservoirs [43, 45]. See Appendix A for a derivation of the generalized master equation, and Section 5 for an example application using an Ohmic bath.

Beyond this dressing transformation, it has been shown that there exists a potential gauge ambiguity in the electric field operator which causes further problems when computing observables in the USC regime [37]; namely, Π_C corresponds to the Coulomb gauge electric field, but the QRM Hamiltonian is derived in the dipole gauge. The gauge transformation from the Coulomb gauge to the dipole gauge is generated by a unitary transformation, which for the restricted TLS subspace is given by the projected unitary operator [37] U = exp ( − i η ( a + a † ) σ x ) . The photon destruction operator transforms as a → U a U † = a + i η σ x [34]. Thus, to “gauge-correct” the master equation in the dipole gauge, we conduct the dressing operation as above, but with

where we take C j k ′ = 〈 j | U Π C U † | k 〉 = 〈 j | Π D | k 〉 = 〈 j | i ( a † − a ) + 2 η σ x | k 〉 ; see Appendix A for a derivation of the master equation in the dipole and Coulomb gauges and their equivalence.

To study the quantum dynamics and spectral resonances, we excite the system with an incoherent pump term, P inc D x GC − / 2 , or with a coherent laser drive, H drive ( t ) = ( Ω d / 2 ) ( x GC − e − i ω L t + x GC + e − i ω L t ) , added to H _QR, where Ω_d is the Rabi frequency and ω _L = ω _c is the laser frequency; thus, H _S = H _QR + H _drive. Note that the QRM with a coherent drive is time-dependent and oscillates around a pseudo-steady-state. In addition, because of the driving laser, the periodic nature of the system Hamiltonian means that in principle the QRM spectra, already quite rich, are modified further; however, we use Ω_d ≪ g, and neglect the influence of the coherent drive on the system eigenstates. The first few (lowest) energy eigenvalues are plotted for the QRM (dipole gauge) and JC model in Figure 2(a) for a range of normalized coupling strengths. Three transitions are shown, which we will refer to below.

Figure 2:

Example energy eigenvalues, as well as steady state excitation numbers and selected transitions rates (with and without gauge corrections). (a) The energy eigenvalues of the six lowest states of the QRM (blue, solid) and the JC model (red, dotted). Arrows mark transitions of interest, placed at arbitrary locations on the η-axis, (b) steady state excitation number for incoherent driving (cf. Figure 3), and (c) selected transition rates, with colors matching the arrows in (a). On the bottom two panels, solid (dashed) lines are with (without) the gauge correction in the dipole gauge. Note a sudden increase of N _cav near η ≈ 0.4 when states 2 and 3 cross.

4 Gauge-invariant observables

We first define the system excitation number,

and a quadrature operator matrix element squared,

which is proportional to the photodetection rate of cavity-emitted photons from the |j⟩ → |k⟩ transition [41]. In Figure 2(b), we show N _cav versus η, using incoherent driving (cf. Figure 3), where the solid curves show the effect of gauge corrections. Equivalent gauge-corrected results are obtained in the Coulomb gauge. With the correction, the population saturates, while the uncorrected population continues to increase superlinearly, and jumps when states 2 and 3 cross in energy, potentially related to the photon blockade [46]. With gauge corrections, we see a strong influence from the TLS operator physics. In Figure 2(c), we show | P j k ′ | 2 for the relevant transitions which are, for weak excitation, proportional to the transition linewidths; again, the solid lines show the gauge corrected results. Note that the corrected dipole gauge quadrature operator ( Π D ′ = i ( a † − a ) + 2 η σ x ) causes a major modification of the transitions, significantly impacting their behavior in the nonperturbative regime.

In Appendix E, we give analytical insight into these quadrature matrix elements using a Bloch–Siegert (BS) transformation, which analytically (to lowest order in η) predicts the following changes with gauge correction: | P I | 2 = 1 / 4 ( 1 + 3 η / 2 ) → 1 / 4 ( 1 − 5 η / 2 ) , and | P I I I | 2 = 1 / 4 ( 1 − 3 η / 2 ) → 1 / 4 ( 1 + 5 η / 2 ) , causing a reversed asymmetry with gauge corrections. Physically, this asymmetry arises from the BS shift of cavity and TLS resonances giving rise to photon-like and atom-like polariton branches; the composition of the Π operator (which is affected by the gauge correction) ultimately determines which state is more cavity-like, and thus has a greater decay rate (see Appendix E for details).

Next, we define the cavity-emitted spectrum,

where x GC , Δ ± = x GC ± − x GC ± and Ω = ω − ω _L. Beyond the spectrum, which uses a first-order quantum CF, we also compute the normalized second-order quantum CF,

which can quantify, for example, the likelihood of a photon being detected at (t + τ) if one was detected at t. We also introduce the time-averaged g ( 2 ) ( τ ) = ∫ t 1 t 1 + T g ( 2 ) ( t , τ ) d t / T , where t ₁ is an arbitrary time point at which the system has reached the pseudo-steady-state and T is the period of oscillation (see Appendix D). Note that without the gauge-correction, we use the uncorrected (corresponding to a Coulomb gauge representation) x ± , x Δ ± for computing the observables, and x ^± for incoherent or coherent driving (see Appendix A). All calculations use Python with the QuTiP package [47, 48].

For weak incoherent pumping, Figure 3 compares the computed spectra with and without the gauge correction (DGC: dipole-gauge-corrected and DG: dipole-gauge, respectively), for η ranging from 0.05 (strong coupling) to 0.5 (USC). For relatively small η = 0.05, the DGC (with gauge correction) spectra already begin to noticeably deviate from the DG spectra (usual QRM master equation solution).

Figure 3:

Cavity spectra outside the RWA (QRM) with DG model (orange dashed line), and DGC model (with gauge correction, blue line) for varying η and weak incoherent driving: P _inc = 0.01g. Spectra are normalized to have the same maxima. Other system parameters are κ = 0.25g, and ω _L = ω _c = ω ₀. Note a small change with the DG corrected model even below the USC regime (η = 0.05).

With increasing η, notably, the DGC and DG spectra are substantially different above η = 0.1: the DGC spectra still show a reversed asymmetry, with a significant narrowing of the lower polariton resonance (I) and a broadening of the upper polariton resonance (III); the ratio of higher-lower polariton peak areas under weak excitation changes from 1 − 3 η + O ( η 2 ) to 1 + 5 η + O ( η 2 ) with gauge correction—a dramatic change even for η < 0.1 (see Appendix E). These peaks can be identified as resulting from the |1⟩ → |0⟩ (olive arrow on Figure 2(a)) and |2⟩ → |0⟩ (brown arrow) transitions, respectively. Since | P j , k | 2 contributes to photon emission directly through the κ decay channel [41], the narrowing (broadening) of peak I (III) with increasing η can be explained with Figure 2(c). Without the correction, the opposite trend is observed, which is again consistent with Figure 2(c) (dashed lines). At η = 0.5, there is also a noticeable resonance (II) around ω = 0.8g, showing a deep mixing of the TLS and cavity dynamics in the USC regime. We can identify this energy difference with the |3⟩ → |1⟩ transition, pink arrow on Figure 2(a), which also has reduced broadening with η, cf. Figure 2(c).

We have shown how the gauge correction manifests in modified linewidths and drastically different spectral weights in comparison to the usual QRM—even so far as to result in a complete reversal of the asymmetry predicted from a non-gauge-corrected model [49] (Figure 9, η = 0.5). We now demonstrate how this gauge correction manifests in the Coulomb gauge. To do this, we display results for the cavity-emitted spectrum and CFs with coherent and incoherent pumping, using the discussed dipole gauge and the Coulomb gauge master equation.

In the Coulomb gauge, the standard system Hamiltonian for the QRM is [37]

where g _C = g _D ω ₀/ω _c and D is the strength of the diamagnetic term. Using the Thomas–Reiche–Kuhn sum rule [47], then D ≥ g C 2 / ω 0 , and for our simulations we take D = g C 2 / ω 0 . Thus, with ω ₀ = ω _c and g _D ≡ g, we have D = η ² ω ₀. Unfortunately, this form does not satisfy the gauge principle, and produces the wrong eigenenergies and eigenstates in the USC regime [35, 37, 41]. Instead, the corrected Coulomb gauge uses a different system Hamiltonian [37],

which contains field operators to all orders, and the C′ superscript indicates we are using the corrected form for the system Hamiltonian. In the Coulomb gauge, the gauge-invariant dissipator term is (see Appendix A)

where x C + = ∑ j , k > j C j k C | j 〉 〈 k | with C j k C = ⟨ j | Π C | k ⟩ , and we now compute the dressed states in the Coulomb gauge, using both uncorrected and corrected forms of the system Hamiltonian.

Figure 4 (top) shows the coherent and incoherent spectra at η = 0.5, showing that the gauge correction results in a profound effect in either case. For coherent driving, using Ω_d = 0.1g, there is a significant sharpening of the resonances. The Coulomb gauge result without the gauge correction corresponds to a minimal coupling Hamiltonian naively truncated to a TLS, which results in incorrect energy levels for the dressed-state master equation [37]. This effect of having the incorrect energy levels and eigenstates is clearly shown in the uncorrected Coulomb gauge results in Figure 4, which is especially wrong with coherent pumping, since the system is effectively being pumped off resonance, because of the diamagnetic term. For coherent pumping, additional Rabi field strengths are shown in Appendix C, where we also show simulations with and without an RWA for the pump field.

Figure 4:

Direct comparison between master equation results using the dipole and Coulomb gauges at η = 0.5, for both coherent and incoherent excitation, showing the profound effect of the gauge correction and how this manifests in identical spectra (top) and g ⁽²⁾(τ) correlation functions (bottom). Solid and dashed curves are with and without the gauge correction, respectively. For the coherent drive (left), we use Ω_d = 0.1g, and the incoherent pumping (right) is the same as in Figure 3 (P _inc = 0.01g).

Next, in Figure 4 (bottom), we examine the second-order coherence, which is important for characterising the generation of non-classical light. In all cases shown, we observe photon bunching at short time-delays. With the gauge correction, there is a significant reduction in the level of bunching, and the usual USC master equations significantly overestimate the bunching characteristics. Moreover, the dynamics are qualitatively different, and thus the non-GC master equations results clearly fail in the USC regime. In all cases, we confirm full agreement between the corrected dipole gauge and corrected Coulomb gauge results, since these are the correct gauge-invariant solutions, and thus produce identical results.

5 Influence of the spectral bath function on the gauge correction and gauge-invariant spectra with an Ohmic bath

In the simulations above, for simplicity, we used a flat density of states (DOS) for the spectral bath function; namely, the DOS was assumed to be constant relative to the energy scale of the resonances. This helps to better identify intrinsic spectral asymmetries related to gauge correcting.

For completeness, here we explicitly show an example numerical solution without invoking the approximation that κ(ω) is frequency independent. Specifically, we compute the emitted spectra when κ(ω) = κ as well as κ(ω) = κω/ω _c (Ohmic bath). We use the same example as in Figure 3 with incoherent driving at η = 0.5. These numerical solutions are obtained from the generalized master Eq. (A9), described in Appendix A.

As can be seen in Figure 5, clearly the form of the spectral bath function does not affect any of our general conclusions, as the gauge correction is, in both cases, dramatic, and of course produces exactly the same result for both the dipole gauge and the Coulomb gauge. To be clear, if we plot these together, then they are indistinguishable, which also confirms that our numerical results are well converged in terms of basis size and time steps.

Figure 5:

The computed cavity spectra using the dipole gauge (left) and Coulomb gauge (right), with a flat DOS [κ(ω) = κ, panels (a), (c)] and an Ohmic DOS [κ(ω) = κω/ω _c , panels (b), (d)] using the generalized master equation [see Eqs. (A9) and (A10) in Appendix A]. In both cases, the effect of the gauge correction (solid lines versus dashed lines) is dramatic. We use the same parameters as in Figure 2 of the main text, with incoherent driving, and parameters η = 0.5 and κ = 0.25g. Notably, in all cases, regardless of the spectral function, the corrected dipole gauge and corrected Coulomb gauge results are identical.

6 Conclusions

We have presented a gauge-invariant master equation approach and calculations for the cavity emission spectra in the USC regime, and shown how the usual QRM in the dipole gauge fails, yielding effects that are just as pronounced (or even more pronounced) as counter-rotating wave effects in this regime. We have demonstrated how the gauge correction significantly affects the intensity CF and cavity excitation number. We have also shown how the gauge correction modifies results in the Coulomb gauge compared to typically used models. Apart from yielding new insights into the nature of cavity-QED system-bath interactions and presenting gauge-invariant master equations that can be used to explore a wide range of light–matter interaction in the USC regime, our results show that currently adopted master equations in the USC regime produce ambiguous results since they do not satisfy gauge invariance.

While we have shown explicit results for the cavity spectrum and intensity CF, the gauge correction causes profound effects on any observable that is computed from the master equations in the same coupling regimes. The nature of the system-bath coupling is also very important, which must also be related to the quadrature coupling to the external fields and the observables to ensure a gauge-invariant master equation. For example, it may be more appropriate to use Π_C = a + a ^† (vector potential coupling) rather than Π_C = i(a ^† − a) (electric field coupling) for the interaction (in the Coulomb gauge), or some linear combination of the two; this change affects the dissipators, incoherent pumping, and coherent excitation in a way that still yields gauge-independent results, but the observables are different. By unitary equivalence, the form of the quadrature coupling used in the system Hamiltonian is thus also not arbitrary, which is in stark contrast to the JC model, where both these coupling forms yield identical results. These two coupling forms are widely used in the USC literature and are assumed to lead to the same result; however, they differ significantly, which reinforces the need, highlighted recently [50, 51], to go beyond the usual phenomenological formulation of system-environment coupling Hamiltonians in the USC regime of cavity QED in favor of a general fundamental microscopic derivation. Solutions to such problems can likely be rigorously addressed using quantized quasinormal modes [52], [53], [54], [55], [56], which even apply to cavities and media in the presence of gain [57, 58].