Quantum noise properties of multiphoton transitions in driven nonlinear resonators

Coupling a driven quantum mechanical oscillator to environmental fluctuations allows the oscillator dynamics to reach a stationary state. In the stationary state, energy is coherently absorbed from the pump and leaks into the environment via random dissipative transitions, which inevitably induce noise in the resonator. This occurs even at zero temperature where only environmental zero-point quantum fluctuations (quantum noise) exist. The quantum noise properties of a nonlinear oscillator determine many fundamental nonequilibrium phenomena such as quantum heating [1–3] and quantum activation [4] (see below).

In recent years, considerable interest has been devoted to the nonequilibrium quantum dynamics of driven nonlinear quantum oscillators because of the tremendous progress in fabricating and thus controlling individual macroscopic quantum systems operating on the nanoscale. This includes superconducting Josephson junctions [5–8] in different variants and also nanomechanical devices which have been successfully realized in the deep quantum regime only recently [9–11]. In addition, quantum transport devices on the basis of molecular junctions have been realized where the interplay of charge transport and vibrational properties of the molecular bridge has been studied [12, 13]. An important aspect common to all these approaches is the fact that the nonlinear response and the noise properties of single macroscopic quantum systems can be addressed instead of measuring an ensemble of resonators where additional averaging is intrinsically involved.

In general, nonlinear oscillators are naturally used as basic elements for quantum state detection or amplification. Examples of those are the Josephson bifurcation amplifier [5, 14–18] and the cavity bifurcation amplifier [5, 19]. In this context, the noise properties of the resonator, which is used as a detector or amplifier, determine the backaction of the measurement or amplification on the system itself [20–22]. Clearly, it is desirable to keep the backaction as weak as possible, while on the other hand, a significant coupling of the amplification or measurement device to the system can be used to achieve a sufficiently strong detection or amplification efficiency. A fundamental lower limit for the disturbance introduced by the noisy detector into the qubit to be detected, however, will be set by the quantum noise acting in the detector. Hence, in order to design useful concepts for quantum state detection based on nonlinear resonators in the deep quantum regime, their quantum noise properties have to be addressed.

The Josephson bifurcation amplifier [5, 14–18] takes advantage of the dynamically induced bistability due to the intrinsic nonlinearity of the resonator. The eigenstates of the qubit whose state has to be measured are mapped onto the coexisting stable states of forced vibrations of the resonator, which have different amplitudes and phase relations relative to the phase set by the external drive. Hence, they allow for a large discriminating power. Up to now, these amplifying devices mostly operate in the semiclassical regime where many quanta in the resonator are excited. This implies that pure quantum fluctuations are typically small on average. Nevertheless, some experiments have been realized at low temperature where the relevant fluctuations are quantum mechanical in nature [5, 18]. The regime of weak fluctuations has been the subject of intense theoretical investigations [1–3, 23–27]. It has been shown [1–3] that a stationary distribution over the quasienergy states of the driven oscillator at zero temperature can be reached which has the form of an effective Boltzmann distribution with an associated nonzero effective temperature even when the statistical temperature is set to zero. Since only the zero-point fluctuations of the vacuum are responsible for this stationary distribution, this has been termed as quantum heating [3]. This stationary state is reached via activation-type transitions between discrete quasienergy states of forced vibrations which are induced by zero-temperature quantum noise [1–3] and are therefore called quantum activation transitions [4].

Signatures of the onset of quantum fluctuations can be seen in the relative intensities of the lines of the resonator noise spectrum [3, 23–26] and in the appearance of a fine structure in the spectral lines of resonators with comparatively large nonlinearities and large quality factors [3, 25]. Most importantly, it has been shown that the spectral fine structure of the noise power spectrum of a parametrically modulated oscillator yields detailed information on the population of the quasienergy states of the resonator in its stationary state [3, 25]. Since the noise power spectrum is, in principle, experimentally accessible, one can directly deduce the stationary nonequilibrium occupation distribution from this measurement signal. No other means is so far available to achieve this. Below, we also find a spectral fine structure in the noise power spectrum of the quantum Duffing oscillator which possesses a similar topology as the parametric oscillator.

Recently, it was proposed that nonlinear quantum detectors that operate in the regime of few quanta (deep quantum regime) would offer different advantages, such as a small backaction, a large discrimination power with an enhanced readout fidelity and a sufficiently large measurement efficiency [27]. In the deep quantum regime, the frequency-resolved nonlinear response of the oscillator to the external driving with frequency ω_ex shows a rich fine structure [1, 2, 28–30] which is mainly generated by few-photon transitions in the resonator. The splitting of the typical Lorentzian resonance of a harmonic oscillator into a series of non-Lorentzian resonances and antiresonances reflects the intrinsic nonequidistance of the energy levels E_n of a nonlinear oscillator. N-photon transitions with the resonance condition E_N − E₀ = Nℏω_ex, N = 1,2,... , and the subsequent drift down along the ladder of the few-photon Fock states generate a pronounced nonequilibrium quasienergy distribution which is strongly different from the Boltzmann-type [1, 2, 29]. Peaks or dips in the nonlinear response are a direct consequence of the nonequilibrium distribution over states with different oscillation amplitudes and phases [1, 2, 29]. The signatures of such a characteristic non-Lorentzian line shape of the response have been observed experimentally in a circuit quantum electrodynamics setup [31]. When such a nonlinear quantum resonator is used as a detector [27], it is important to understand and control its noise properties.

In this work, we investigate the quantum noise properties of modulated nonlinear oscillators in the deep quantum regime. As the simplest example, we consider a monostable anharmonic oscillator which has a quartic nonlinearity (quantum Duffing oscillator). Such a device can be parametrically coupled to a qubit [27] and thus its photon noise characteristics are of relevance when used as a parametric detector. Below, we analyse the power spectrum of the photon number fluctuations in the quantum Duffing oscillator. In the underdamped regime, we identify Lorentzian peaks in the power spectrum of the photon number noise which are associated with the multiphoton transitions in the quantum Duffing oscillator. Their intensities are determined by the stationary occupation probabilities of the quasienergy states. The latter is directly connected to the effective quantum temperature which can be identified in the stationary state and which is intimately connected to the quantum squeezing of the quasienergy states. Hence, measuring the power spectrum of the photon number fluctuations provides a direct and elegant way of determining the stationary occupation probabilities of the quasienergy states and thus the effective quantum temperature. As shown below, the latter is directly connected to the amount of quantum squeezing of the resobator state.

A weakly nonlinear Duffing oscillator has a remarkable symmetry: its energy levels E_n with n ⩽ N are pairwise resonant for the same driving frequency ω_ex, E_N−n − E_n = (N − 2n)ℏω_ex. An example of the energy spectrum for the case N = 3 is sketched in figure 1(a). After preparing the oscillator in its nth excited state n ⩽ N, it displays periodic quantum oscillations between the nth and the (N − n)th excited states. During these oscillations, |N − 2n| photons are being exchanged between the oscillator and the modulation field. The oscillations of the photon number $\hat {n}$ are usually referred to as multiphoton Rabi oscillations. Their characteristic frequency, the Rabi frequency Ω_n,N, depends on the intensity of the driving field and on the number of photons exchanged. The Rabi frequency Ω_0,N for the N-photon oscillations is the smallest Rabi frequency. The multiphoton Rabi oscillations with N − n photons involved are underdamped if their Rabi frequency Ω_n,N exceeds the dissipative rate of photon leaking into the environment. The latter is the oscillator relaxation rate γ. For γ ≪ Ω_0,N all the Rabi oscillations are, in general, underdamped. The periodically driven resonator reaches its stationary state on the time scale γ⁻¹.

Zoom In Zoom Out Reset image size

In the stationary state, quantum noise induces—even at zero temperature—fluctuations in the photon number $\hat {n}$. The dynamics of these fluctuations is characterized by multiphoton oscillations which manifest themselves as peaks in the noise spectrum S(ω) of $\hat {n}$, located at plus/minus the Rabi frequencies Ω_n,N. In the underdamped regime, the dissipative dynamics of the driven oscillator is most appropriately described in terms of random transitions between the oscillator quasienergy states. When the driving is resonant, the pairs of oscillator Fock states with n- and N − n-photons are resonantly superposed. The corresponding oscillator quasienergy states are a symmetric and an antisymmetric superposition of the two Fock states. Their splitting in quasienergy is given by the Rabi frequency Ω_n,N. The corresponding peak in the noise spectrum at (−)Ω_n,N is due to random transitions from the state with the (highest) lowest to that with the (lowest) highest quasienergy of the doublet. The peak intensity is proportional to the stationary occupation probability of the initial quasienergy state. Therefore, the noise spectrum offers a convenient way of directly probing the stationary distribution over all the quasienergy states. Moreover, for weak driving and exactly zero detuning from the multiphoton resonance, the noise spectrum of the $\hat {n}$-photon transition is symmetric, i.e. S(ω) = S(−ω) and two inelastic peaks are signatures of an oscillatory decay of the fluctuations towards the stationary state. States belonging to a multiphoton doublet then have the same stationary occupation probabilities. For a weakly detuned modulation or a stronger driving, the spectrum becomes asymmetric. Besides, an additional quasielastic peak appears at zero frequency which represents incoherent relaxation of the fluctuations towards the stationary state. These features have some analogy to the spectral correlation function of a (static) quantum mechanical two-level system weakly coupled to a dissipative harmonic bath [32]. There, the spin correlation function is a sum of three Lorentzian peaks. The two inelastic peaks are symmetrically located at finite frequencies and their width determines the inverse of the dephasing time. In addition, the quasielastic peak at zero frequency represents incoherent relaxation with the inverse relaxation time given by its width. In the driven system, the appearance of a quasielastic peak depends on the intriguing interplay between the nonlinearity, the driving strength and the dissipation strength.

We consider a periodically modulated quantum oscillator with mass m, eigenfrequency ω₀ and a quartic (Kerr) nonlinearity described by the Hamiltonian

$$\begin{equation} H(t) = \frac{1}{2m}p^2 + \frac{1}{2}m \omega_0^2 x^2 + \frac{1}{4}\alpha x^4 +F x \cos(\omega_{\rm ex}t). \end{equation} \tag{ 1 }$$

The modulation amplitude F is assumed to be so small that it induces only weakly nonlinear vibrations. This is guaranteed by the condition αA² ≪ mω²₀, with A(F) being the typical amplitude of the nonlinear vibrations. The modulation frequency ω_ex is chosen to be close to the oscillator eigenfrequency ω₀ such that the detuning δω is small, i.e.

$$\begin{equation} \delta\omega\ll \omega_0,\qquad \delta\omega\equiv \omega_0 - \omega_{\rm ex}. \end{equation} \tag{ 2 }$$

Our theory applies to hard as well as to soft nonlinearities α≶0, but to be specific, we will focus in the rest of the paper on the case of a hard nonlinearity, α > 0.

The quantum dynamics of the weakly detuned and weakly nonlinear-driven oscillator is most conveniently described in terms of the oscillator ladder operators a and a^†, in a rotating frame determined by the unitary transformation

$$\begin{equation} R(t)=\exp[-\mathrm{i}\omega_{\rm ex}a^\dagger a t] . \end{equation} \tag{ 3 }$$

In the rotating frame, the typical time scale of the resonator dynamics is given by δω⁻¹, so that terms oscillating with frequencies ±2ω_ex and ±4ω_ex average out and can be neglected in the transformed Hamiltonian $R(t)H(t)R^\dagger (t) - \mathrm {i}\, \hbar R(t)\dot {R}^\dagger (t)$. With this rotating wave approximation (RWA), we obtain the Hamiltonian

$$\begin{equation} \skew3\tilde{ H} = \delta\omega \hat{n} + \nu \hat{n} (\hat{n} +1)/2 + f (a^\dagger + a)/2, \end{equation} \tag{ 4 }$$

where $\hat {n}\equiv a^\dagger a$ is the photon number operator, ν and f are the frequencies associated with the Kerr nonlinearity and the external field amplitude at the quantum scale $x_{\mathrm { ZPF}}=\sqrt {\hbar /m\omega _0}$, i.e. ν = 3αx⁴_ZPF/4ℏ and $f=Fx_{\mathrm { ZPF}}/\sqrt {2}\hbar$. In order to keep the notation compact we have set ℏ = 1 in equation (4) and in the remainder of the paper. The oscillator quasienergies ε_n and quasienergy states |ψ_n〉 are the eigenvalues and eigenvectors of the rotating wave Hamiltonian, $\skew3\tilde {H}|\psi _n\rangle =\varepsilon _n|\psi _n\rangle$. For vanishing driving, the quasienergy spectrum is given by

$$\begin{equation} \varepsilon_n=\delta\omega \, n + \nu n (n+1)/2\quad {\rm for }\ f \to 0 . \end{equation}\noindent \tag{ 5 }$$

We are primarily interested in studying the noise spectrum in the presence of multiple multiphoton resonances E_N−n − E_n = (N − 2n)ω_ex, wherein n ⩽ N and E_n are the energy levels of the undriven nonlinear oscillator, or equivalently ε_N−n = ε_n for f → 0. From equation (5) we find the resonant condition

$$\begin{equation} \delta\omega=-\nu (N+1)/2\equiv\delta\omega_N. \end{equation} \tag{ 6 }$$

Up to leading order in the driving, the quasienergy eigenstates |ψ_n〉 for n ⩽ N ≠ N/2 are a resonant superposition of the pair of oscillator Fock states |n〉 and |N − n〉, i.e.

$$\begin{eqnarray} &&|\psi_n\rangle\approx\left(|n\rangle\pm|N-n\rangle\right)/\sqrt{2} . \end{eqnarray} \tag{ 7 }$$

We choose the signs − and + for n < N/2 and N/2 < n ⩽ N, respectively. In the following, we refer to the resonant superposition of a pair of Fock states as resonant states or multiphoton states. The states |ψ_n〉 that are not involved in a multiphoton transition (n > N and n = N/2 for N even) can be approximated as the corresponding Fock states |ψ_n〉 ≈ |n〉. The Rabi frequency Ω_n,N−n of the multiphoton oscillations within the pair of Fock states |n〉 and |N − n〉 is given by the splitting of the corresponding levels ε_n and ε_N−n [29, 33] as

$$\begin{equation} \Omega_{n,N-n}=|\varepsilon_n-\varepsilon_{N-n}|= f \left( \frac{f}{\nu} \right)^{N-2n-1}\frac{ {(N-n)!^{1/2}}}{{n!^{1/2}}(N-2n-1)!^2} . \end{equation} \tag{ 8 }$$

The resonance condition in equation (6) is not renormalized by a finite driving within the RWA. Only for a comparatively larger driving f ∼ ν ≪ ω₀, the multiphoton transitions have to be reinterpreted as tunneling transitions between semiclassical states [34, 35].

As we shall detail in section 4, the multiphoton Rabi oscillations induce peaks in the spectral densities of oscillator observables only when the Rabi frequency, Ω_0,N, for the multiphoton transition from the zero-photon ground state, is larger than the noise-induced level broadening of the relevant quasienergy levels ε₀ and ε_N. In the next section, we will pave the wave for the calculation of the noise spectrum in this regime, by formulating the master equation for a weakly nonlinear oscillator and by evaluating the stationary occupation populations over the quasienergy states.

In the presence of a weak bilinear coupling to the fluctuations of a bosonic bath, the assumptions of small detuning and weak nonlinearity that under the RWA naturally lead to a Liouville–von Neuman quantum master equation in the Lindblad form for the density matrix $\hat {\rho }$ of the weakly damped oscillator in the rotating frame [29, 30],

$$\begin{eqnarray} &&\begin{array}{l} \displaystyle \dot{\hat{\rho}} = {\cal L}\hat{\rho}\equiv -\mathrm{i}[\skew3\tilde{H},{\hat{\rho}}]+\gamma(\bar{n}+1){\cal D}[a]\hat{\rho}+ \gamma \bar{n}{\cal D}[a^\dagger]\hat{\rho},\\ \displaystyle {\cal D}[O]\hat{\rho}\equiv ([O\hat{\rho},O^\dagger]+[O,\hat{\rho} O^\dagger])/2. \end{array} \end{eqnarray} \tag{ 9 }$$

Here, ${\cal L}$ and ${\cal D}[O]$ are the Liouville and the Lindblad superoperators, respectively. Moreover, γ is the oscillator damping rate for which we assume that γ ≪ ω₀. It results from a standard Ohmic bath spectral density J(ω) = γω. In addition, $\bar {n}$ is the bosonic occupation number at frequency ω₀ and temperature T and is given by $\bar {n}=(\mathrm {e}^{\omega _0/k_{\mathrm {B}}T}-1)^{-1}$.

3.1. The stationary distribution

For long times, the density matrix in the rotating frame $\hat {\rho }$ relaxes to a stationary state $\hat {\rho }^{\infty }$, satisfying

$$\begin{equation} {\cal L} \hat{\rho}^{\infty} = 0. \end{equation} \tag{ 10 }$$

When the oscillator decay rate γ is larger than the driving, γ ≫ f, the width of the resonant quasienergy levels ε_n induced by the bath fluctuations is larger than the corresponding Rabi frequency Ω_n,N−n of the multiphoton transitions. Then, the multiphoton resonances are smeared out and the coherent effects associated with multiphoton oscillations are strongly suppressed. Hence, dissipation sets a lower limit for the driving strength, f ≫ γ, which has to be overcome in order to observe multiphoton transitions. On the other hand, for comparatively larger driving f ∼ ν, quantum fluctuations are significantly reduced and the oscillator is latched to a classical attractor at asymptotic times which are much larger than the typical relaxational time scale, γ⁻¹. In this regime, the (quasi)stationary distribution of the population over the quasienergy states assumes a Boltzmann form [1–3, 25, 26].

Here, we restrict our analysis to the deep quantum regime where the driving is larger than the damping but smaller than the nonlinearity, γ ≪ f ≪ ν. Thereby, we have implicitly assumed a comparatively large nonlinearity ν ≫ γ.

3.1.1. Underdamped regime

We start our discussion assuming that all Rabi oscillations are underdamped. Put differently, we assume that the smallest Rabi frequency Ω_0,N is larger than the relevant level broadening. We refer to this regime as the fully underdamped regime. Then, the off-diagonal matrix elements of $\hat {\rho }^\infty$ projected onto the quasienergy basis |ψ_l〉 are negligible and we can set them to zero, i.e. we perform a secular approximation,

$$\begin{equation} \rho^\infty_{lk}\equiv\langle\psi_l|\hat{\rho}^\infty|\psi_k\rangle= 0\quad {\rm for}\ l\neq k. \end{equation} \tag{ 11 }$$

Then, a balance equation for the stationary occupation probabilities ρ^∞_ll follows from equations (9) and (10) according to

$$\begin{equation} \gamma_l\rho^\infty_{ll}-\sum_{l\neq k}W_{lk}\rho^\infty_{kk}=0 . \end{equation} \tag{ 12 }$$

Here, W_lk is the transition rate from state |ψ_k〉 to state |ψ_l〉,

$$\begin{equation} \quad W_{lk}\equiv\gamma\left[(1+\bar{n}) |\langle \psi_l|a|\psi_k\rangle|^2+\bar{n}|\langle \psi_l|a^\dagger|\psi_k\rangle|^2\right], \end{equation} \tag{ 13 }$$

and γ_l is the width of quasienergy level ε_l given by $\gamma _l\equiv \sum _{k\neq l}W_{lk}$. We can now formulate more precisely the condition for underdamped Rabi oscillations to occur within the pair forming the narrowest resonance, which is Ω_0,N ≫ γ₀.

The solution for stationary occupation probabilities up to leading order in the small parameters f/ν and $\bar {n}$ is given in [29]: the pair of multiphoton states |ψ_n〉 and |ψ_N−n〉 in equation (7) have equal stationary population, i.e. ρ^∞_nn = ρ^∞_N−n,N−n. The pair with the narrowest resonance has the occupation probabilities ρ^∞₀₀ = ρ^∞_NN. The occupation probability grows algebraically with n < N/2 as

$$\begin{equation} \rho^\infty_{n+1,n+1}= \frac{N-n}{n+1}\rho^\infty_{n\,n}\quad {\rm for}\ n<N/2. \end{equation} \tag{ 14 }$$

The states |ψ_l〉 with l > N have vanishing occupation probability, ρ^∞_ll = 0. As follows from the discussion above, the degeneracy ρ^∞₀₀ = ρ^∞_NN is approximate and is lifted for higher order in f/ν.

3.1.2. Quasienergy distribution close to a multiphoton resonance

One can easily generalize the above expressions to the case when the detuning δω does not exactly match the resonant condition, δω ≠ δω_N. Since the Rabi frequencies for the different pairs of resonant transitions in equation (8) are exponentially different, we can choose |δω − δω_N| ≪ Ω_1,N−1, so that all the pairs of Fock states |n〉 and |N − n〉 with 1 < n < N/2 are still resonantly superposed, except for

$$\begin{eqnarray} &&|\psi_{0}\rangle = \cos \frac{\theta}{2} |0\rangle - \sin \frac{\theta}{2} |N\rangle \quad {\rm and} \quad |\psi_{N}\rangle = \sin \frac{\theta}{2}|0\rangle + \cos \frac{\theta}{2} |N\rangle , \end{eqnarray} \tag{ 15 }$$

with θ = tan⁻¹[Ω_0,N/N(δω − δω_N)]. The corresponding solution for the stationary density matrix close to resonance is [29]

$$\begin{eqnarray} &&\begin{array}{l} \displaystyle \rho^\infty_{NN}=\rho^\infty_{00}\tan^4\frac{\theta}{2}, \qquad\rho^\infty_{11}=\rho^\infty_{00}N\tan^2\frac{\theta}{2},\\ \displaystyle \rho^\infty_{n+1,n+1}= \frac{N-n}{n+1}\rho^\infty_{n\,n}\quad {\rm for}\ 1\leqslant n<N/2. \end{array} \end{eqnarray} \tag{ 16 }$$

3.1.3. Partially underdamped regime

Next we consider a comparatively large relaxation rate γ, so that the narrowest Rabi resonance is overdamped but the remaining resonances are still underdamped, Ω_0,N ≪ Nγ ≪ Ω_1,N−1. We refer to this regime as the partially underdamped regime. Then, incoherent multiphoton transitions from the ground state |0〉 to state |N〉 with a small rate Ω²_0,N/(Nγ) and the subsequent emission of excitations into the bath determine a small but finite occupation of the resonant states ρ^∞_n n, n ⩾ 1. Formally, the stationary distribution $\hat {\rho }^\infty$ can be obtained by setting all the off-diagonal elements of ρ^∞_lk to zero except for ρ^∞_N0 and ρ^∞_0N and solving equation (10). Thereby, we find that

$$\begin{eqnarray} & &\langle0|\hat{\rho}^\infty|0\rangle\approx 1,\quad \langle N|\hat{\rho}^\infty|N\rangle\approx \Omega_{0,N}^2/(N^2\gamma^2), \nonumber\\ & &\rho^\infty_{11}=\Omega_{0,N}^2/(N\gamma^2)\quad {\rm for} \ \Omega_{0,N}^2/(N\gamma^2)\gg\exp[-\omega_0/(k_{\mathrm{B}} T)],\\ & &\rho^\infty_{n+1,n+1} = \frac{N-n}{n+1}\rho^\infty_{n\,n}\quad {\rm for}\ 1\leqslant n<N/2.\nonumber \end{eqnarray} \tag{ 17 }$$

The crossover between this solution and the fully underdamped solution (14) is given in [29]. Both stationary nonequilibrium distributions are determined by quantum fluctuations and are very different from the equilibrium Boltzmann-type distribution when a driven resonator is latched to a classical attractor.

3.2. The nonlinear response of the oscillator

In the steady-state regime, t ≫ γ⁻¹, the oscillator state is described by the time-independent density matrix $\hat {\rho }^\infty$ in the rotating frame and the oscillator dynamics is embedded in the time-dependent reference frame R(t). The mean value of an observable O is

$$\begin{equation} \langle O(t)\rangle_\infty\equiv\lim_{t\to\infty}\langle O(t)\rangle= \mathrm{Tr}\bigl\lbrace \hat{\rho}^\infty R^\dagger(t)OR(t)\bigr\rbrace. \end{equation} \tag{ 18 }$$

Therefore, the stationary oscillations of the position expectation value 〈x(t)〉_∞ are sinusoidal,

$$\begin{equation} \langle x(t)\rangle_{\rm \infty}= \sqrt{2}x_{\rm ZPF}\cos(\omega_{\rm ex}t+\varphi)|\langle a \rangle_\infty|,\qquad \langle a \rangle_\infty=\sum_{lk}\rho^\infty_{lk}\langle\psi_l|a|\psi_k\rangle . \end{equation} \tag{ 19 }$$

It has been shown that the nonlinear response 〈x(t)〉_∞ of the oscillator as a function of ω_ex shows resonances and antiresonances in the deep quantum regime [28–30]. The response is proportional to the transmitted amplitude in a heterodyne measurement scheme and has already been measured for a weakly nonlinear oscillator [31]. Clearly, such a measurement scheme, or more generally, any measurement scheme that probes stationary mean values as opposed to correlations does not allow one to resolve the different degenerate resonances separately. Nor do they allow us to access the stationary distribution ρ^∞_ll directly. This becomes possible only when correlations, e.g. via noise spectra, are measured. In the next section, we show that this can indeed be achieved by measuring the spectrum of the photon number noise.

4.1. Definition of the noise spectrum

The Lindblad master equation (9), in general, also allows one to investigate transient phenomena and correlation functions. Its formal solution for a given initial state $\hat {\rho }_0$ is given by $\hat {\rho }(t)=\exp [{\cal L} t]\,\hat {\rho }_0$. Moreover, a general correlator 〈O'(t')O(t)〉 can be evaluated as the mean value of the operator O' at time t' with the virtual operator R^†(t)OR(t)ρ(t) at time t. This view has been established several decades ago by the Lax formula [36, 37] according to

$$\begin{equation} \langle O'(t') O(t) \rangle = \mathrm{Tr} \lbrace R^\dagger(t')O'R(t') \mathrm{e}^{{\cal L}(t'-t)}R^\dagger(t)OR(t) \rho(t) \rbrace. \end{equation} \tag{ 20 }$$

For long times t ≫ γ⁻¹, we find that

$$\begin{eqnarray} &&\fl \langle O(t+\delta t) O'(t) \rangle_{\rm \infty}\equiv \lim_{t\to\infty} \langle O(t+\delta t) O'(t) \rangle=\mathrm{Tr} \lbrace R^\dagger(t+\delta t)OR(t+\delta t) \mathrm{e}^{{\cal L} \delta t} R^\dagger(t)O'R(t) \rho^{\rm \infty} \rbrace.\nonumber\\ \end{eqnarray} \noindent \tag{ 21 }$$

In general, such correlators are periodic functions of the preparation time t. The noise spectrum is defined as a double average over quantum fluctuations and the time t.

Here, we are specifically interested in the noise spectrum S(ω) of the autocorrelator of the photon number $\hat {n}$, $\langle \hat {n}( t+\delta t) \hat {n}(t) \rangle _{\mathrm { \infty }}$. From equations (3) and (21), we find that

$$\begin{eqnarray} & &\langle \hat{n}(t+\delta t) \hat{n}(t) \rangle_{\rm \infty}=\mathrm{Tr} \lbrace \hat{n} \mathrm{e}^{{\cal L} \delta t} \hat{n} \rho^{\rm \infty} \rbrace. \end{eqnarray} \tag{ 22 }$$

Since this correlator does not depend on the initial time t as a consequence of the RWA, we can define the noise spectrum in terms of a single average over quantum fluctuations according to

$$\begin{equation} S(\omega) = 2 \, \mathrm{Re} \int_0^{\infty} \mathrm{d} t \, \mathrm{e}^{\mathrm{i} \omega t}\langle \hat{n}( t) \hat{n}(0) \rangle_{\rm \infty}. \end{equation} \tag{ 23 }$$

It is useful to separate the contributions to S(ω) into those coming from the expectation value of $\hat {n}$, and those from its fluctuations, i.e.

$$\begin{eqnarray} &&\begin{array}{l} \displaystyle S(\omega) = \langle \hat{n}\rangle_{\rm \infty}^2\delta(\omega)+\delta S(\omega),\\ \displaystyle \delta S(\omega)\equiv 2 \, \mathrm{Re} \int_0^{\infty} \mathrm{d} t \, \mathrm{e}^{\mathrm{i} \omega t}\langle \delta\hat{n}( t) \delta\hat{n}(0) \rangle_{\rm \infty}. \end{array} \end{eqnarray} \tag{ 24 }$$

Here, $\delta \hat {n}$ is the operator for the photon number fluctuations, i.e. $\delta n=\hat {n}-\langle \hat {n}\rangle _{\mathrm { \infty }}$.

Our path to compute the noise spectrum consists of three steps. (i) We express the virtual preparation $\hat {n}\hat {\rho }^\infty$ in terms of right eigenvectors of the superoperator ${\cal L}$. (ii) We plug the resulting decomposition into equation (22). Then, each term decays exponentially with a different exponent which is given by the corresponding eigenvalue of ${\cal L}$. (iii) We compute the Fourier integral in equation (23), which thereby yields a sum over (overlapping) Lorentzian peaks.

The general expression, which is useful for a concrete numerical evaluation, for the noise spectrum given in terms of the eigenvectors and the eigenvalues of ${\cal L}$ is derived in the appendix. In the next section, we consider the special case of underdamped multiphoton Rabi oscillations.

4.2. Noise spectrum in the underdamped regime

When all the multiphoton Rabi oscillations are underdamped, Ω_N,0 ≫ Γ_N, the coherences |ψ_N−n〉〈ψ_n| and |ψ_n〉〈ψ_N−n| are approximate eigenvectors of the Liouvillian ${\cal L}$. Then,

$$\begin{eqnarray} &&\begin{array}{l} \displaystyle {\cal L}|\psi_n\rangle\langle\psi_{N-n}|= -(\Gamma_n-\mathrm{i}\Omega_{n,N-n})|\psi_n\rangle\langle\psi_{N-n}|\quad {\rm for} \ n<N/2,\\ \displaystyle {\cal L}|\psi_{N-n}\rangle\langle\psi_{n}|= -(\Gamma_n+\mathrm{i}\Omega_{n,N-n})|\psi_{N-n}\rangle\langle\psi_{n}|\quad {\rm for} \ n<N/2 \end{array} \end{eqnarray} \noindent \tag{ 25 }$$

with the level widths being given as $\Gamma _n=\gamma _n=\gamma (\bar {n}+1/2)N+\gamma \bar {n}$ for n < (N − 1)/2. For N odd, $\Gamma _{(N-1)/2}=\gamma (1+2\bar {n})(5N+1)/8+\gamma \bar {n}$. Up to leading order in f/ν, the decomposition of the virtual preparation $\hat {n}\hat {\rho }^\infty$ in terms of right eigenvectors of ${\cal L}$ has the simple expression

$$\begin{equation} \hat{n} \hat{\rho}^\infty\approx \left(N/2\right) \hat{\rho}^\infty -\sum_{n<N/2}(N/2-n)\rho^\infty_{n\,n} \left(|\psi_n\rangle\langle\psi_{N-n}|+|\psi_{N-n}\rangle\langle\psi_{n}|\right). \end{equation} \noindent \tag{ 26 }$$

Clearly, each term of the above decomposition yields a Lorentzian peak in the noise spectrum S(ω). The first term yields the contribution to S(ω) from the expectation value of $\hat {n}$, (N/2)²δ(ω). The remaining terms yield inelastic peaks associated with random transitions between quasienergy states belonging to the same multiphoton doublet. Since the populations ρ^∞_n n and ρ^∞_N−n,N−n are approximately equal, peaks at opposite frequency have approximately equal intensity. By putting together equations (22), (23), (25) and (26), we find that S(ω) = (N/2)²δ(ω) + δS(ω) with

$$\begin{equation} \delta S(\omega)\approx\sum_{n<N/2} S_n(\omega)+S_{N-n}(\omega), \end{equation} \tag{ 27 }$$

$$\begin{equation} S_n(\omega)=S_{N-n}(-\omega)= \frac{2\Gamma_n\rho_{n\,n}^{\infty}(N/2-n)^2}{(\omega-\Omega_{n,N-n})^2+ \Gamma_n^2}. \end{equation} \noindent \tag{ 28 }$$

Hence, the Lorentzian peaks are centred at the multiphoton Rabi frequencies Ω_n,N−n and have a resonance width of Γ_n. The factor (N − 2n)²/4 is the leading order expression for the squared matrix element $|\langle \psi _n|\hat {n}|\psi _{N-n}\rangle |^2$. Remarkably, the line intensities depend only weakly on the driving f and on the temperature through the stationary distribution ρ^∞_n n. Up to leading order, the driving f enters only in the splitting of the lines through the Rabi frequencies. Note that equation (27) is valid only in the vicinity of a multiphoton peak since terms of order γ are not taken into account. In order to evaluate the tails of the peaks more precisely, one has to take into account the contribution stemming from all eigenvectors of ${\cal L}$, see the appendix.

In the left and right panels of figure 2, we show the noise spectrum S(ω) for the cases N = 2 and 3, respectively. The noise spectrum for N = 2 shows a pair of symmetric peaks which correspond to the transitions $|\psi _0\rangle \leftrightarrow |\psi _{2}\rangle$. Likewise, the noise spectrum for N = 3 displays two pairs of symmetric peaks corresponding to the transitions $|\psi _0\rangle \leftrightarrow |\psi _{3}\rangle$ and $|\psi _1\rangle \leftrightarrow |\psi _{2}\rangle$. The green dashed lines mark the results from our approximate analytical formula in equation (27) while the yellow solid lines show the data obtained by numerically evaluating the expression in equation (A.2). Excellent agreement is found.

Zoom In Zoom Out Reset image size

In figure 2(a), additional smaller side peaks of the order of f/ν are also visible, see the grey lines representing a tenfold zoom. They are not associated with any resonant transition between multiphoton states and are thus not captured by the leading order expression given in equation (27). The particular subleading peaks in figure 2(a) belong to the transitions $|\psi _0\rangle \leftrightarrow |1\rangle$.

These features show a direct analogy to the spectral correlation function of a static quantum mechanical two-level system which is weakly coupled to a dissipative harmonic bath [32]. For a general biased two-state system with anticrossing energy levels, the pair correlation function is a sum of three Lorentzian peaks. The two inelastic peaks are symmetrically located at finite frequencies and their width determines the inverse of the dephasing time. For a biased static two-level system away from resonance, an additional quasielastic peak at zero frequency appears which represents incoherent relaxation with the inverse relaxation time given by its width. Since we consider here the case strictly at resonance (in the RWA), no zero-frequency peak is present.

4.3. Photon anti-bunching

In general, the photon emission characteristics of a quantum mechanical resonator can show peculiar nonclassical features. For instance, counterintuitive correlation phenomena such has photon antibunching can occur, where the photon number correlation function for short delay times is smaller than that for classical, uncorrelated photons. This implies that the probability for photons to arrive in pairs is suppressed [38]. Our approach provides a natural framework to investigate possible non-Poissonian statistics of the multiphoton events in the nonlinear resonator. Therefore, we consider the normalized photon number correlation function or second-order coherence function defined as

$$\begin{equation} g^{(2)}(\tau) = \frac{\langle a^\dagger(t)a^\dagger(t+\tau)a(t+\tau) a(t)\rangle_\infty}{\langle a^\dagger(t)a(t) \rangle_\infty \langle a^\dagger(t+\tau) a(t+\tau) \rangle_\infty} . \end{equation} \tag{ 29 }$$

For long delay times τ, the counts of two photons with a delay time τ are statistically independent events, g⁽²⁾(τ → ∞) = 1. For vanishing delay times, we have

$$\begin{equation} g^{(2)} (\tau \to 0) = 1+ \frac{\langle \hat{n}^2\rangle_\infty-\langle \hat{n}\rangle_\infty^2- \langle \hat{n}\rangle_\infty}{\langle n\rangle_\infty^2}. \end{equation} \tag{ 30 }$$

Photon antibunching corresponds to the case g⁽²⁾(τ = 0) < 1. For the fully underdamped case, we find the expression

$$\begin{equation} g^{(2)} (\tau = 0) = \frac{2N (N-1)+4\sum_{1}^{N-1}n(n-N)\rho_{n\,n}}{N^2}=1-\frac{1}{N}, \end{equation}\noindent \tag{ 31 }$$

which represents the known result of the second-order correlation function of the electromagnetic field [38]. Hence, the oscillator displays photon antibunching close to a multiphoton transition. The second-order coherence of the stationary state of the quantum Duffing oscillator at the Nth multiphoton resonance has the same value as the second-order coherence for an oscillator prepared in the single Fock state |N〉, in spite of its fluctuations over the quasienergy states.

4.4. Line shape of the noise spectrum close to a multiphoton resonance

In the presence of a small detuning from the multiphoton resonance, δω − δω_N ∼ Ω_0,N, the states |ψ₀〉 and |ψ_N〉 are no longer a resonant superposition of the Fock states |0〉 and |N〉. Hence, the corresponding stationary occupation probabilities ρ₀₀ and ρ_NN, given in equation (16), become significantly different. In turn, the pair of peaks S₀(ω) and S_N(ω), which are associated with the transitions $|\psi _0\rangle \leftrightarrow |\psi _{N}\rangle$, become asymmetric such that S₀(ω) ≠ S_N(−ω). This behaviour is shown in figure 3(a) for the case around the three-photon resonance. The peak line shapes can readily be evaluated and we find that

$$\begin{eqnarray} &&\fl S_{0}(\omega)= \frac{2\Gamma_0\rho_{00}^{\infty}N^2(\sin\theta\cos\theta)^2} {(\omega-\varepsilon_N+\varepsilon_0)^2 +\Gamma_0^2} =\rho_{00}^{\infty} \left( \frac{N^2\Omega_{0,N}(\delta\omega - \delta\omega_N)}{|\varepsilon_N - \varepsilon_0|^2} \right)^2 \frac{2\Gamma_0}{(\omega - \varepsilon_N + \varepsilon_0)^2 + \Gamma_0^2}, \end{eqnarray} \tag{ 32a }$$

$$\begin{eqnarray} &&\fl S_{N}(\omega) = \frac{2\Gamma_N\rho_{NN}^{\infty}N^2(\sin\theta\cos\theta)^2} {(\omega-\varepsilon_0+\varepsilon_N)^2 +\Gamma_N^2} =\rho_{NN}^{\infty} \left( \frac{N^2\Omega_{0,N}(\delta\omega - \delta\omega_N)}{|\varepsilon_N - \varepsilon_0|^2} \right)^2 \frac{2\Gamma_N}{(\omega - \varepsilon_0 + \varepsilon_N)^2 + \Gamma_N^2}. \end{eqnarray}\noindent \tag{ 32b }$$

Their distance increases with the quasienergy splitting, ε_N − ε₀ = sgn(δω − δω_N)(Ω²_0,N + N²|δω − δω_N|²)^1/2, whereas the peak width does not change close to the multiphoton resonance, δω − δω_N ∼ Ω_0,N. The asymmetry is determined by the stationary occupation probabilities ρ^∞₀₀ and ρ^∞_NN. From equation (16), we find that

$$\begin{equation} \frac{S_0(\omega)}{S_N(-\omega)}=\frac{\rho^\infty_{00}}{\rho^\infty_{NN}}= \cot^4\frac{\theta}{2}= \left[\frac{\Omega_{0,N}}{|\varepsilon_N-\varepsilon_0|-N(\delta\omega - \delta\omega_N)} \right]^4 . \end{equation}\noindent \tag{ 33 }$$

The above expression is valid for ω close to the centre of the largest peak, ω ∼ ε_N − ε₀, and |δω − δω_N| not too large such that S(±ω) ≫ γ.

Zoom In Zoom Out Reset image size

In addition to the peaks at finite frequencies (which induce decaying coherent multiphoton Rabi oscillations), a zero-frequency peak also appears. This quasielastic peak is associated with incoherent relaxational decay of the multiphoton Rabi oscillations and is also known for the noise correlation function of a static biased quantum two-level system [32].

In figure 3(b), we show the logarithm of the asymmetry ratio given in equation (33). The asymmetry shows a clear maximum at approximately ε₃ − ε₀.

To further illustrate the asymmetry in the peak heights, we show in figure 3(c) the peak maxima associated with the transitions |ψ₀〉 → |ψ₃〉 and |ψ₃〉 → |ψ₀〉. At the three-photon resonance (black dashed vertical line), both peaks are equal in height (symmetric noise spectrum). Away from the resonance, the low (high-)-frequency branch acquires more spectral weight for negative (positive) detuning.

4.5. Photon noise at zero frequency

Fluctuations of an oscillator (quasi)energy induce a broad (with width ∝γ) zero-frequency peak in the noise spectrum of an observable whose mean value depends on the (quasi)energy [39]. For weak driving f ≪ ν and at a resonance |δω − δω_N| ≪ Ω_0,N, the quasienergy states of the Duffing oscillator have large fluctuations as several quasienergy states have comparable occupation probabilities even at T = 0. However, the mean value of $\hat {n}$ becomes independent of the quasienergy, $\langle \psi _n|\hat {n}|\psi _n\rangle \approx N/2$ for n ⩽ N. As a consequence, the contribution to the noise spectrum of $\hat {n}$ coming from fluctuations δS(ω) does not have a peak at zero frequency since δS(0)∝γ. Close to resonance, when |δω − δω_N| ∼ Ω_0,N, two dynamical effects compete: on the one hand, the quasienergy fluctuations quickly decrease for increasing detuning, i.e. moving away from resonance as the occupation probability of the state |ψ₀〉 approaches one. On the other hand, the mean value of $\hat {n}$ becomes strongly dependent on the quasienergy. As a result of this competition, the intensity of the zero-frequency noise plotted as a function of δω has two maxima at the two opposite sides of the resonant value δω_N. In figure 4, we show the zero-frequency noise for the special case N = 2. The yellow solid line represents the intensity at zero frequency computed numerically, while the green dashed line is the leading order contribution (in f/ν)

$$\begin{equation} \delta S (\omega = 0 ) \approx \sin^2(2\theta)/4\gamma. \end{equation} \tag{ 34 }$$

Zoom In Zoom Out Reset image size

4.6. Noise spectrum towards the semiclassical regime

Next, we investigate the noise spectrum for larger driving strengths, f ≲ ν. In order to illustrate how the noise spectrum changes with increasing driving, we show the intensities of the brightest peaks as a function of the driving strength for the N = 5 photon resonance; see figure 5(a). In figure 5(b), we also show the quasienergy spectrum, and the noise spectrum for a comparatively large value of the driving amplitude f = ν is shown in figure 5(c). A peak in the noise spectrum at frequency ω = ε_l − ε_k is associated with a single transition |ψ_k〉 → |ψ_l〉 and is given by

$$\begin{equation} S(\omega)=\sum_{lk}\frac{2\rho_{ll}^{\infty}|\langle \psi_l |a^\dagger a|\psi_k\rangle|^2(\gamma(a_l-a_k)^2+\Gamma_l+\Gamma_k)}{ (\omega+\varepsilon_l-\varepsilon_k)^2+(\gamma(a_l-a_k)^2+\Gamma_l+\Gamma_k)^2} . \end{equation} \noindent \tag{ 35 }$$

Hence, the relative intensities of a pair of peaks at opposite frequencies are still related to the occupation probability of the corresponding initial states through S(ε_l − ε_k)/S(ε_k − ε_l) = ρ_kk/ρ_ll.

Zoom In Zoom Out Reset image size

For weak driving, we have three pairs of approximately symmetric peaks as described by equation (27). Each peak corresponds to a transition between two states belonging to a multiphoton doublet of quasidegenerate states: $|\psi _0\rangle \leftrightarrow |\psi _5\rangle$, $|\psi _1\rangle \leftrightarrow |\psi _4\rangle$ and $|\psi _2\rangle \leftrightarrow |\psi _3\rangle$. For increasing driving, the spectrum becomes increasingly asymmetric. For moderate values of the driving, the noise spectrum undergoes two major qualitative changes: (i) the peak at zero frequency becomes clearly visible; and (ii) a pair of peaks corresponding to the transitions $|\psi _1\rangle \leftrightarrow |\psi _3\rangle$ acquires a significant intensity. For f = ν, the peak associated with the transition |ψ₃〉 → |ψ₁〉 is even the second brightest peak.

These qualitative changes can be explained in terms of a semiclassical description valid beyond the weak driving limit. The RWA Hamiltonian in equation (4) can be rewritten in terms of rotating quadratures and interpreted as a quasienergy surface in phase space [34, 35]. It has the shape of a tilted Mexican hat and is sketched in figure 1(c) for two values of f. The larger the f, the stronger the induced tilt. The local maximum and the minimum of the quasienergy surface are the classical attractors. In the static frame, they describe stationary oscillations with a small and a large amplitude, respectively. In the vicinity of the attractors, the vibration amplitude and the slow part of the oscillation phase display slow vibrations with frequency ∝δω. In the absence of resonant transitions, each quasienergy state can be associated with a quantized quasiclassical orbit which lies on the internal surface around the local maximum, on the external surface or along the quasienergy well around the minimum. For very weak driving, $f\ll \nu /\sqrt {2(N+1)}$, the quantum mechanical Fock states |n〉 with n < N/2 are associated with quasiclassical trajectories on the internal surface around the local maximum, whereas the Fock states with photon number n larger than N/2 are associated with semiclassical orbits on the external surface. Within this representation, the multiphoton transitions can then be reinterpreted as tunneling transitions between the internal and the external part of the surface [34, 35]. For comparatively larger driving, the zero-point quasienergy associated with the slow vibrations around the minimum (∝δω) becomes smaller than the dynamical barrier height. Then, quasienergy states appear which are localized in the quasienergy well. In turn, the noise spectrum becomes qualitatively different from that for weak driving. The small quantum fluctuations around the minimum of the quasienergy surface can be described in terms of an effective auxiliary oscillator with ladder operators b and b^† and are given by

$$\begin{equation} a=a_h+b\cosh{r^*_h}-b^\dagger\sinh{r^*_h}. \end{equation} \tag{ 36 }$$

Here, a_h is the amplitude of the stationary oscillations rescaled by $\sqrt {2}x_{\mathrm { ZPF}}$ [25, 26]. They can be mimicked by a local effective quantum temperature T_e = (2k_B ln coth r*_h)⁻¹ which depends on the squeezing factor r*_h [1, 2, 25, 26]. For f = ν, the states |ψ₂〉, |ψ₃〉 and |ψ₁〉 can be identified with the ground state and the first two excited states of the auxiliary oscillator (but in the remainder of this discussion we keep the same labels for the states as in the weak driving limit). The level spacing ε₃ − ε₂ is of the order of the frequency of the slow classical oscillations of the amplitude and slow part of the phase.

Such oscillations appear in the noise spectral density of a classical oscillator as a pair of peaks. In a nonlinear quantum oscillator whose quasienergy levels are not equidistant and their distance exceeds the damping strength, the classical peaks have a 'quantum' fine structure [3]. In the present case of the Duffing oscillator, the classical noise peak is split into two peaks associated with nearest-neighbour transitions between the ground state and the first excited state, and the first and the second excited state, $|\psi _2\rangle \leftrightarrow |\psi _3\rangle$ and $|\psi _3\rangle \leftrightarrow |\psi _1\rangle$, respectively. Their peak height is proportional to the square of the rescaled vibration amplitude a_h and to the occupation of the initial state ρ^∞_n n. The latter, in particular, is governed by the quantum temperature T_e. For the ratio of the peak heights, we find that [3]

$$\begin{equation} \frac{ S(\varepsilon_3-\varepsilon_2)}{S(\varepsilon_2-\varepsilon_3 )}\approx \frac{\rho^\infty_{22}}{\rho^\infty_{33}}\approx\coth^2r^*_h\approx\frac{\rho^\infty_{33}}{\rho^\infty_{11}} \approx \frac{ S(\varepsilon_1-\varepsilon_3)}{S(\varepsilon_3-\varepsilon_1)}. \end{equation}\noindent \tag{ 37 }$$

Next-nearest-neighbour transitions can also yield peaks in the noise spectra of a Duffing oscillator [26]. In the present case, the transitions $|\psi _2\rangle \leftrightarrow |\psi _1\rangle$ yield a pair of dimmer peaks, but located at frequencies outside of the frequency range shown in figure 5.

In the weak damping, weak driving regime discussed so far, the quasienergy well around the minimum is still very shallow, and the oscillator can escape from the small amplitude attractor via tunneling. Therefore, the oscillator is not latched to any of the attractors and the noise spectral density also has peaks which are associated with intra-well transitions. In particular, the pair of peaks with the smallest splitting describes coherent tunneling oscillations between the internal and the external part of the quasienergy surface (coherent dynamical tunneling or multiphoton Rabi oscillations).

Before closing this section, we mention that for the stronger driving f = ν, a zero-frequency peak also appears in the noise spectrum, see figure 5(c), although the frequency detuning has been fixed at the five-photon resonance δω = δω₅. However, as discussed above, this resonance condition is only valid for small f ≪ ν, which is obviously not fulfilled. So the larger driving induces an effective small detuning away from the exact avoided quasienergy level crossing and generates an effective bias. Then, a relaxation pole appears in the relevant self-energy [32] which corresponds to a quasielastic relaxation peak at zero frequency.

4.7. Dependence of the noise spectrum on damping and temperature

So far, we have analysed the case of zero temperature and small damping, $\bar {n}\ll 1$ and γ ≪ Ω_0,N. In this section, we briefly address how the noise spectrum is modified for larger damping and finite temperature by presenting the numerical results of the spectrum in a broad parameter range.

In figure 6(a), we show S(ω) for different values of the damping for the three-photon resonance where δω = δω₃. As expected, the peaks in the noise spectrum get broader for increased damping. Outside the fully underdamped regime, the two peaks of the pair associated with the transitions $|\psi _0\rangle \leftrightarrow |\psi _{3}\rangle$ start to overlap and eventually merge into a single peak at zero frequency. Thereby, the zero-frequency noise is no longer suppressed S(ω ≈ 0)∝γ⁻¹, since incoherent relaxation prevails over coherent decay for large damping. The peaks associated with the underdamped transitions $|\psi _1\rangle \leftrightarrow |\psi _{2}\rangle$ are still described by equation (27) even when the spectrum has a peak at zero frequency. The decrease in peak intensities reflects the decrease of the populations ρ^∞₁₁ and ρ^∞₂₂ in the partially underdamped regime.

Zoom In Zoom Out Reset image size

The dependence of the noise spectrum on temperature is shown in figure 6(b) and behaves qualitatively similarly. For small temperatures $\bar {n}\ll 1$, the spectrum is described by equation (27). The temperature dependence enters the line widths of the quasienergy levels as well as in the stationary distribution ρ^∞_n n. For higher temperatures, the two low-frequency peaks merge into a single peak at zero frequency and the side peaks become increasingly broader as expected.

In recent years, the rich phenomenology of driven and damped nonlinear quantum oscillators has been impressively consolidated, including their nonlinear response behaviour in the form of resonant and antiresonant amplification, quantum coherent multiphoton Rabi oscillations, quantum activation and quantum heating. Gradually, the nontrivial effects visible in noise correlation functions have also moved to the focus of interest. These become relevant whenever a nonlinear quantum oscillator is used as a central element in an amplifier or quantum measurement device. In this work, we have analysed the quantum noise properties of the quantum Duffing oscillator in the regime where only a few quanta are excited. Then, the nonlinear response shows pronounced multiphoton peaks which are associated with resonant multiphoton Rabi oscillations. The noise properties of these multiphoton transitions show a rich phenomenology with a fine structure of the noise power spectrum. To obtain the noise spectrum by analytic means, we invoke the Lax formula for the autocorrelation function of the photon number at different times and calculate its Fourier transform. Exactly at a multiphoton resonance, the noise spectrum consists of a collection of pairs of related resonances which are located at opposite frequencies and which are equal in height. Each pair is associated with a multiphoton doublet. In spite of large fluctuations over the oscillator quasienergy, no quasielastic peak occurs at zero frequency. This is a consequence of a special symmetry of the quantum Duffing oscillator: all quasienergy states that are associated with a multiphoton doublet have the same mean value of the photon number $\hat {n}$.

Slightly away from a multiphoton resonance, the noise spectrum becomes asymmetric and the two resonances are no longer equal in height. In addition, as the mean values of $\hat {n}$ become different for quasienergy states with comparable occupations, the quasielastic peak emerges. Since the quasienergy fluctuations are suppressed away from a multiphoton resonance, the intensity of the quasielastic peak as a function of the detuning displays a maximum at the two opposite sides of the resonant value δω_N.

We have shown that the amplitudes of the resonances of the photon noise power spectrum can be used to directly determine the asymptotic occupation probabilities of the quasienergy states of the quantum Duffing resonator. Similar results have been obtained for the parametric oscillator [3]. No other scheme has been proposed for this task up to now. Once these populations are measured via the noise power spectrum, the quantum temperature as well as the amount of quantum squeezing of the resonator state (including its nonclassical photon statistics) follow immediately in the semiclassical regime by increasing the photon number by a larger driving amplitude. Then, a quasiclassical quasipotential landscape in phase space is a convenient tool for understanding the stationary nonequilibrium dynamics. This view directly leads to the quantum mechanical squeezed states that exist close to the local minimum of the quasienergy landscape. A harmonic expansion allows us to characterize the quantum fluctuations via an effective quantum temperature. At larger (real) temperature and damping strengths, all these quantum coherent features are washed out.

Although the time-resolved detection of noise properties of quantum observables of driven resonators requires considerably more experimental effort, we are confident that future experiments will soon elucidate the importance of quantum noise in these systems.

This work was supported by the DAAD (German Academic Exchange Service) through research grant no. A/08/73659 and by the DFG through SFB 925 'Light induced dynamics and control of correlated quantum systems' (project C8). VP was supported by the NSF (grant no. EMT/QIS 082985). We thank P Nalbach, M Marthaler and M Dykman for valuable discussions.

In order to numerically compute the noise spectrum S(ω), we diagonalize the Liouvillian superoperator ${\cal L}$,

$$\begin{equation} {\cal L}\hat{k}_{\mathrm{R}}=\lambda_k\hat{k}_{\mathrm{R}},\qquad \hat{k}_{\mathrm{L}} {\cal L}=\lambda_k\hat{k}_{\mathrm{L}}. \end{equation} \tag{ A.1 }$$

Here, $\{\hat {k}_{\mathrm {L}}\}$ and $\{\hat {k}_{\mathrm {R}}\}$ are the sets of the left and right eigenvectors of ${\cal L}$, respectively. They constitute a biorthogonal system implying that ${\mathrm { Tr}} \lbrace \hat {k}^\dagger _{\mathrm {L}}\hat {k}'_{\mathrm {R}}\rbrace = \delta _{kk'}$. The corresponding eigenvalues λ_k are, in general, complex. The stationary density matrix $\hat {\rho }^\infty$ is the only eigenvector with zero eigenvalue. All other eigenvectors represent transient dissipative processes. The corresponding eigenvalues have a negative real part. Note that we have implicitly assumed that ${\cal L}$ admits a spectral decomposition. Put differently, we have ruled out that any Jordan block has a dimension larger than one.

With this at hand, the noise spectrum in equation (23) can be rewritten in terms of $\hat {k}_{\mathrm {L}}$, $\hat {k}_{\mathrm {R}}$ and λ_k as a sum of partial spectra characterized by Lorentzian peaks according to

$$\begin{equation} S(\omega)=\frac{2\sum_k(\mathrm{Tr}\lbrace \hat{n}\hat{k}_{\mathrm{R}}\rbrace \mathrm{Tr}\lbrace \hat{k}_{\mathrm{L}}^\dagger\hat{n}\hat{\rho}^\infty\rbrace){\rm Re}\,\lambda_k}{(\omega -{\rm Im}\, \lambda_k)^2+({\rm Re}\, \lambda_k)^2}. \end{equation} \tag{ A.2 }$$

The sum extends over those eigenvectors $\hat {k}_{\mathrm {R}}$ that belong to nonzero eigenvalues λ_k. Thereby, we have not included the elastic Rayleigh peak $\langle \hat {n}\rangle _{\infty }^2\delta (\omega )$ which trivially comes from the stationary state $\hat {k}_{\mathrm {R}}=\hat {\rho }^\infty$.