Fixing the rotating-wave approximation for strongly-detuned quantum oscillators

Periodically-driven oscillators are commonly described in a frame co-rotating with the drive and using the rotating-wave approximation (RWA). This description, however, is known to induce errors for off-resonant driving. Here we show that the standard quantum description, using creation and annihilation of particles with the oscillator's natural frequency, necessarily leads to incorrect results when combined with the RWA. We demonstrate this on the simple harmonic oscillator and present an alternative operator basis which reconciles the RWA with off-resonant driving. The approach is also applicable to more complex models, where it accounts for known discrepancies. As an example, we demonstrate the advantage of our scheme on the driven quantum Duffing oscillator.

Nature displays many degrees of freedom that are coupled to one another. For simplicity, we describe them as linear resonators that drive each other as they oscillate with a bare harmonic around their respective equilibrium points. Hence, finding the response of a system to an applied drive is an essential task in all areas of science, ranging from fundamental effects in classical and quantum physics to applications in spectroscopy and diagnostics of our daily lives.
The bare harmonic translates to a drive that is periodic in time with a given frequency. To describe the system's response to the drive, a common approach involves moving to a rotating frame, where the system's behavior appears stationary on short timescales. In the rotating frame, the so-called rotating-wave approximation (RWA) is used, whereby any remaining rapidly oscillating terms are dropped. This procedure becomes less precise when the system and drive frequencies are far detuned from one another, i.e., when the neglected terms are less rapid. Nevertheless, the RWA is a standard approach to find the response of driven systems [1,2], such as nanomechanical resonators [3][4][5], optical cavities [6][7][8][9][10][11][12], phononic and magnonic modes [13][14][15][16][17][18][19] and superconducting junctions [20][21][22]. Similarly, it is a common starting point for Floquet engineering and perturbative expansions of nonlinear driven systems [23][24][25][26][27].
Despite its ubiquity, the RWA is known to induce errors at large detuning [28]. There the neglected oscillating terms become significant and must be included as perturbations that affect the system's response [29][30][31][32][33]. Yet, such a perturbative treatment becomes very tasking in many cases that rely on detuned driving, such as optomechanical cooling [34,35] and the formation of frequency combs [36][37][38][39][40]. These cases are therefore commonly treated using the RWA. The issue is particularly salient in the case of the simple quantum harmonic oscillator, where the solution obtained with the RWA disagrees with the exactly soluble classical limit.
In this work, we pinpoint the source for the quantumto-classical discrepancy in the driven harmonic oscillator and correct for it. Specifically, we show the description of the oscillator in terms of its bare excitations is illsuited for the RWA. Instead, we propose a non-standard operator basis (frame) in which the RWA yields the correct solution for the driven case. We furthermore demonstrate that this basis choice is also beneficial for describing driven nonlinear systems, as it significantly improves the fidelity of the RWA in the driven Duffing oscillator at large detuning. We expect that our approach will lead to a significant improvement in the description of numerous scenarios where the response of an oscillator to a coherent drive is calculated.
The Hamiltonian of a periodically driven harmonic oscillator reads with position x, its conjugate momentum p, mass m, and resonance frequency ω 0 . The driving force has amplitude F 0 and oscillates at frequency ω. The time evolution of the system is given by Hamilton's equations of motion (EOMs),ṗ = −mω 2 0 x + F 0 cos(ωt) ,ẋ = p/m .
These EOMs can be solved exactly by transforming to the Fourier domain, where the driving term appears as a Dirac delta function. As a result, the stationary response of the oscillator to the drive reads [41] x = X cos(ωt) , p = −mωX sin(ωt) , . We can parameterize the oscillator's phase space by x and p/mω 0 . The trajectory described by Eq. (3) is then an ellipse with vertices at X and Xω/ω 0 , see Fig. 1(a).
Within the framework of quantum mechanics, the Hamiltonian (1) is an operator, written in terms of a pair of non-commuting operatorsx andp. To find the oscillator's response in the quantum realm, one commonly introduces raising and lowering operators a and a † that diagonalize the time-independent part of H, with F a = F 0 /(2 √ 2mω 0 ), see Fig. 1(b). The driving renders the Hamiltonian (5) time-dependent and several methods exist to find the resulting response. For example, one approach involves time-dependent perturbation theory [42], where we move to the interaction picture using the transformation V (t) ≡ e −ia † a ω0t , and diagrammatically expand the remaining interaction Hamiltonian H int . While this approach is possible for the harmonic oscillator, it rapidly becomes infeasible when additional terms are included in H. Hence, in many cases only linear response is considered, e.g., using the input-output method [43].
Another common approach involves the rotating-wave approximation (RWA). Here, we aim at removing the oscillatory time-dependence altogether. To do so, we first use the unitary transformation U a (t) = e −iωta † a . In the resulting rotating frame, we obtain the (still timedependent) effective Hamiltoniañ whereã ≡ U a (t) † aU a (t) denotes the operator a in the rotating frame, and ∆ = ω − ω 0 is the so-called detuning away from resonance. In the rotating frame, the remaining time-dependent terms are far-detuned since 2ω |∆| and therefore elicit negligible response. Hence, we can apply the RWA and remove all of these from the Hamiltonian (6) to obtaiñ H RWA . We may then find the response by solving Heisenberg's EOM, i d dt ã ≡ [ã,H RWA ] / = −∆ ã − F a , and obtain the stationary solution ã = −F a /∆. Then, by inverting the transformation U a (t) we have a = ã e −iωt . An analogous EOM is similarly solved for ã † . Combining the two solutions [cf. Eq. (4)], we obtain in the classical limit Crucially, the response amplitude found with the RWA does not match the analytical solution in Eq. (3), see Fig. 1(a). Indeed, the ratio only approaches unity for ∆ → 0, which demonstrates the limitation of the RWA to near-resonant drives. A comparison between the exact and RWA results in phase space is also illustrative. To this end, we take a stationary amplitude ã ≡ | ã |e iθ , and using Eq. (4), we find the equivalent trajectory in  sin(ωt − θ). The resulting path is always circular, unlike the elliptical path of the exact solution, cf. Eq. (3) and Fig. 1(a). Indeed, the dropped timedependent oscillating terms, dubbed micromotion, add on to the circular path to produce the full elliptical motion, cf. small circle in Fig. 1(a). This implies that to describe the elliptical path, ã must be time-dependent, which cannot be obtained only by time-independent corrections from a high-frequency expansion [23][24][25][26][27]. In other words, the assertion of a stationary ã violates Hamilton's EOMs, specifically the relation p = mẋ.
As the driven quantum harmonic oscillator is exactly solvable, the RWA was not necessary in the aforementioned analysis. Without the approximation, we would still have obtained the exact result. However, when the system becomes nonlinear, its EOMs are usually not exactly solvable. The RWA, being effectively the lowest order of a perturbative expansion, is then a common method to find the stationary response. As we have seen, however, it neglects crucial corrections due to the detuned drive frequency. We identify that this systemic mismatch stems from the usual starting point (4), which relies on operators that diagonalize the bare oscillator Hamiltonian and describes excitations with energies ω 0 , see Fig. 1(b).
We now turn to our scheme for an exact treatment of the driven harmonic oscillator. We seek operators akin to a and a † which recover the ellipsoidal phase-space path corresponding to drive frequency ω. A natural choice is to replace ω 0 by ω in Eq. (4), i.e., to define new operators b and b † viâ which satisfy [b, b † ] = 1, and describe the system in terms of excitations with the drive's frequency ω. The driven harmonic oscillator Hamiltonian (1) now reads with Fig. 1(c). Note that we effectively rotated from the a to the b operators using a unitary (Bogoliubov) transformation [44]. With the new operators, Eq. (10) contains squeezing terms that do not conserve the particle number, since b and b † do not diagonalize the bare harmonic oscillator unless ω = ω 0 .
We repeat the same procedure in the new basis. We first apply the unitary transformation U b (t) = e −iωtb † b to move to a rotating frame, and then obtain the corresponding Heisenberg's EOMs (11) Crucially, when we now search for a stationary solution for b , we obtain which upon inverting U b (t) matches the exact result in Eq. (3). Importantly, the correct result is obtained regardless of whether we take the RWA or not. In other words, when the operators b, b † are used, the stationary solution in the rotating frame is exact and satisfies Eq. (11). This is the main result of this work. Example: Duffing oscillator.-The RWA is a common starting point for dealing with time-dependent systems that are not exactly solvable. Our result (12) implies that a suitable choice of operators made prior to applying the RWA can significantly reduce errors at large detuning. To demonstrate this, we compare the results obtained using the RWA with the two operator definitions, a, a † and b, b † . We consider a Duffing oscillator, described by the Hamiltonian where α is the Duffing nonlinearity. We plug both operator definitions [Eqs. (4) and (9)] into Eq. (13), move to a rotating frame [cf. Eq. (6)], and apply the RWA in both procedures. Note that the RWA drops multiple oscillating terms that describe frequency conversion processes due to the Duffing nonlinearity. We furthermore apply a semiclassical mean-field ansatz [45] to obtain Heisenberg's EOMs We can now search for stationary solutions in either of the rotating frames [ d dt ã = 0 or d dt b = 0]. Both Eqs. (14) and (15) generate a cubic polynomial condition which has up to three solutions [4,46,47], see Fig. 2(a). Both approaches produce the expected tail-shaped response, i.e., a single solution in the ω < ω 0 regime that bends up towards a high-amplitude solution and a coexistence region in the ω > ω 0 regime, where both lowand high-amplitude solutions appear. The coexistence region manifests mathematically as a bifurcation point for the roots of the polynomial condition. Despite qualitative agreement, the solutions in the two frames differ quantitatively, i.e., in their amplitudes as a function of detuning and in the positions of their bifurcation points.
To find out which of the two approximate procedures describes the driven Duffing oscillator (13) more accurately, we compare their predictions with a "numerical experiment" in the classical picture. Specifically, we numerically evolve Hamilton's equations for the Hamiltonian (13) until we reach stationary motion for different detunings and initial conditions and then adiabatically sweep the detuning to explore the solution landscape. Note that we add an infinitesimal dissipation term to enforce convergence in the simulation, which negligibly shifts the stationary outcome. Additionally, the numerical time-trace describes the oscillatory response at the drive frequency ω including high harmonic generation. Hence, the relevant benchmark for the RWA solutions is the Fourier component x ω at frequency ω of the time trace. In Figs. 2(b) and (c), we plot the relative discrepancy δ x = |(x ω − x ω,RWA ) / x ω | between the numerically obtained x ω and the two RWA results in the high-and low-amplitude solution branches, respectively. We observe that our new RWA approach agrees much better for all values of ω. In the low-amplitude regimes, the system responds quasi-linearly and our scheme performs better as discussed above in the harmonic oscillator case. Slight deviations are seen in the high-amplitude regime, where nonlinear effects dominate the physics.
In driven-dissipative systems, phase diagrams depict the number or type of stationary solutions as a function of system parameters [48,49]. Here, we focus on the location of the bifurcation point as a function of detuning and drive amplitude. The bifurcation point is marked in Figs. 2(a) and (c), where the better performance of our RWA approach manifests as a clear change in its position. In Fig. 2(d), we draw phase diagrams obtained by the two RWA approaches, and observe a distinct difference in the phase boundaries. The exact numerical solution matches our scheme much better, see Fig. 2(e).
Conclusion.-In summary, we have presented an operator basis for periodically-driven systems which anticipates a response at the driving frequency and is thus better suited for using the RWA. For both the harmonic and the Duffing oscillators, this basis choice significantly improves the fidelity of the RWA while retaining its simplicity. Our approach is applicable to systems in many areas of physics beyond mechanical oscillators described by x and p [3]. In quantum optics [43] and optomechanics [50], a driven cavity mode is described by Eq. (1) with the magnetic and electric fields playing the roles of x and p. For models with inherent detuning such as optomechanical cooling and frequency-comb generation, our approach may result in significant numerical correc-tions. A specific case of interest are models coupling a cavity to one or more two-level systems [51,52]; here the so-called general RWA is used, where a basis diagonalizing the spin-cavity interaction reduces the RWA error [53]. As general RWA uses the standard operators a and a † , we expect that it can be improved further by using our scheme. Detuned driven systems also appear in circuit QED where the magnetic flux and charge appear in place of x and p [20]. Overall, a sizable body of work exists elaborating on the errors due to the RWA at large detuning [54][55][56][57]. Our work suggests that in many cases, these errors can be dramatically reduced by choosing an appropriate operator basis.