PRL 111, 053602 (2013)

week ending
2 AUGUST 2013

PHYSICAL REVIEW LETTERS

Nonlinear Interaction Effects in a Strongly Driven Optomechanical Cavity
Marc-Antoine Lemonde,1 Nicolas Didier,1,2 and Aashish A. Clerk1
1

2

Department of Physics, McGill University, Montréal, Quebec H3A 2T8, CanadaDépartement de Physique, Université de Sherbrooke, Sherbrooke, Quebec J1K 2R1, Canada(Received 19 April 2013; published 2 August 2013)

We consider how nonlinear interaction effects can manifest themselves and even be enhanced in astrongly driven optomechanical system. Using a Keldysh Green’s function approach, we calculate
modifications to the cavity density of states due to both linear and nonlinear optomechanical interactions,
showing that strong modifications can arise even for a weak nonlinear interaction. We show how thisquantity can be directly probed in an optomechanically induced transparency-type experiment. We also
show how the enhanced interaction can lead to nonclassical behavior, as evidenced by the behavior of g2
correlation functions.
DOI: 10.1103/PhysRevLett.111.053602

PACS numbers: 42.50.Wk, 07.10.Cm, 42.50.Lc, 42.65.k

Introduction.—The field of cavity optomechanics
involves understanding and exploiting the quantum interaction between a mechanical resonator and photons in adriven electromagnetic cavity. It holds immense promisefor both fundamental studies of large-scale quantum phenomena as well as applications to quantum information
processing and ultrasensitive detection, and has seen
remarkable progress in the past five years. Highlights
include the use of radiation pressure forces to cool amechanical resonator to close to its motional ground state[1,2] and experiments where the mechanical motion causessqueezing of the light leaving the cavity [3,4].
As remarkable as this progress has been, it has relied onstrongly driving the optomechanical cavity to enhance thebasic dispersive coupling between photons and mechanicalposition. While the resulting interaction can be made larger
than even the dissipative rates in the system [5–7], it ispurely bilinear in photon and phonon operators. As a result,it cannot convert Gaussian state inputs into nonclassicalstates or give rise to true photon-photon interactions.
Recent theoretical work has addressed the physics of thenonlinear interaction in weakly driven systems [8,9].
Unfortunately, one finds that the effects are suppressedby the small parameter g=!M .
In this Letter, we now consider nonlinear interaction
effects in an optomechanical system that (unlike
Refs. [8,9]) is also subject to a strong laser drive; we
consider the effects of this driving beyond simple linear
response. We find somewhat surprisingly that one can use
the strong drive to enhance the underlying single-photon
interaction. Using nonequilibrium many-body perturbation
theory [based on the Keldysh technique (see, e.g.,
Ref. [10])], we calculate how these effects modify thecavity density of states and hence the cavity’s response toan additional weak probe laser. This response is exactly
the quantity measured in so-called optomechanically induced transparency (OMIT) experiments [6,11–13]. We
find striking modifications of the OMIT spectrum, effects
0031-9007=13=111(5)=053602(5)which can be attributed to the nonlinear interaction causinga hybridization between one and two polariton states (withthe polaritons being joint mechanical-photonic excitations). We also find the possibility of enhanced polaritonpolariton interactions, which lead in turn to nonclassicalcorrelations (as measured by a g2 correlation function).
System.—The standard Hamiltonian of a driven optomechanical cavity is H^ ¼ H^ 0 þ H^ diss , with (@ ¼ 1)
^ a^ y a^H^ 0 ¼ !C a^ y a^ þ !M b^y b^ þ gðb^y þ bÞpﬃﬃﬃﬃ
þ ½ a in ðtÞa^ y þ H:c::(1)

Here, a^ is the cavity mode (frequency !C , damping rate ),b^ is the mechanical mode (frequency !M , damping rate ),and g is the optomechanical coupling. H^ diss describes thedissipation of photons and phonons by independent baths;
a in ðtÞ is the amplitude of the drive laser.
We consider the standard case of a continuous-wave
drive [i.e., a in ðtÞ / ei!L t ] and work in a rotating frame
at the laser frequency !L . We further make a displacement^ where a is thetransformation, writing a^ ¼ ei!L t ða þ dÞ,classical cavity amplitude induced by the laser drive.
Letting  ¼ !L  !C , the coherent Hamiltonian now
takes the form H^ 1 þ H^ 2 withH^ 1 ¼ d^y d^ þ !M b^y b^ þ Gðd^ þ d^y Þðb^ þ b^y Þ;(2)^ b^ þ b^y Þ:
H^ 2 ¼ gd^y dð(3)

G ¼ ga is the drive-enhanced many-photon optomechanical coupling; we set g, a > 0 without loss of generality.
The most studied regime of optomechanics is wherea  1 and g  , !M . It is then standard to neglect theeffects of H^ 2 . In the absence of any driving, a simple
perturbative estimate suggests that the effects of H^ 2 enteras g2 =!M , where the factor of !M corresponds to a virtualstate with one extra (or one less) phonon. This conclusion
can be made more precise by exactly solving the coherent,

053602-1

Ó 2013 American Physical Society

PRL 111, 053602 (2013)

week ending
2 AUGUST 2013

PHYSICAL REVIEW LETTERS

undriven system using a polaron transformation [8,9].
Thus, in this standard regime, one can ignore H^ 2 , leaving
only H^ 1 , which is easily diagonalized as H^ 1 ¼
P^ y c^  . Here, c^ þ; describe the two normal modes
¼ E cof the system. As these modes have both photon andphonon components, we refer to them as polaritons inwhat follows. Their energies are
qﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ1=2
1
E ¼ pﬃﬃﬃ 2 þ !2M  ð2  !2M Þ2  16G2 !M:
2(4)
For  ’ !M and G  , , the polariton energy splitting
can be resolved experimentally [5–7].
Polariton interactions.—Unlike previous work, we wish
to retain the effects of the nonlinear interaction H^ 2 but also
consider the effects of a large drive (and hence a large
many-photon coupling G). To proceed, we will treat theeffects of H^ 2 in perturbation theory. We use a Keldysh
Green’s function (GF) approach which is able to describe
the nonequilibrium nature of the system. The linear
Hamiltonian H^ 1 along with the dissipative terms in H^ dissdefine the free GFs of the system, which describe thepropagation of polaritons in the presence of dissipation.
Written in the polariton basis, the nonlinear interaction H^ 2
gives rise to number-nonconserving interactionsXðgA0 00 c^ y c^ y0 c^ y00 þ gB0 00 c^ y c^ y0 c^ 00 þ H:c:Þ;H^ 2 ¼
;0 ;00(5)where the coefficients gA=B
0 00 / g [14]. Note normal
ordering H^ 2 in terms of polariton operators introduces
small quadratic and linear terms which modify the diagonalized form of H^ 1 (see the Supplemental Material for
details [14]).
We start by considering how single-particle properties
are modified by the nonlinear interactions; such properties
can be directly probed by weakly driving the cavity witha second probe laser (i.e., an OMIT experiment [11–13])or by measuring the mechanical force susceptibility.
Understanding these properties amounts to calculating
the self-energy ½! of both the polaritons due to H^ 2 .
We have calculated all self-energy processes to second
order in g. Our approach captures both the modificationof spectral properties due to the interaction (i.e., the modification of the cavity and mechanical density of states), as
well as modifications of the occupancies of the mechanics
and cavity. While our approach is general, we will focus
on the most interesting case of a high mechanical quality
factor   !M , a cavity in the resolved sideband regime!M > , and a strong cavity drive G * .
Our full second-order self-energy calculation finds thatfor most choices of parameters, the polariton self-energies
scale as g2 =!M and thus have a negligible effect for the
typical case where g  !M . However, effects are much

more pronounced if one adjusts parameters so thatEþ ¼ 2E . This condition makes the term in H^ 2 which
scatters a þ polariton into two  polaritons (and viceversa) resonant. It can be achieved for any laser detuning in the range ( 2!M , !M =2) by tuning the amplitudea in of the driving laser so that the many-photon optomechanical coupling G ¼ Gres , where [15]qﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
pﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
Gres ½  172 !2M  4ð!4M þ 4 Þ=ð10 !M Þ: (6)In this regime, the dominant physics is well described by
the approximation H^ 0 ’ H^ eff withXE c^ y c^  þ g~ðc^ yþ c^  c^  þ H:c:Þ þ H^ NR ; (7)H^ eff ¼
¼

H^ NR ¼X c^ y c^  þ

X

U

0 ¼

¼

0

c^ y c^ y0 c^ 0 c^ 


:

(8)

The second term in H^ eff corresponds to making a rotatingwave approximation on the nonlinear interactionH^ 2 in Eq. (5), retaining only the resonant process;g~ ¼ gBþ / g is the corresponding interaction strength[Fig. 2(a) shows how g~ varies with ]. The terms in H^ NR
describe the small (i.e., / g2 =!M ) residual effects of
the nonresonant interaction terms in Eq. (5); we treat
them via straightforward second-order perturbation theory(i.e., a Schrieffer-Wolff transformation). They play no rolein the extreme good-cavity limit !M   [14].
Green’s functions for resonant nonlinear interactions.—
Focusing on the resonant interaction regime defined by
Eq. (6) and using the simplified Hamiltonian in Eq. (7),
we obtain simple expressions for the retarded GFs of thesystem. The retarded photon GF in the displaced, rotating
frame is defined as
Z1^ d^y ð0Þiei!t ;GRdd ½! ¼ idtðtÞh½dðtÞ;(9)
1with similar definitions for the polariton retarded GF
GR ½! ( ¼ ). As usual, d ½! ¼ ImGRdd ½!=
describes the cavity density of states; GRdd ½! also determines the reflection coefficient in an OMIT experiment(see Fig. 2). A standard linear response calculation [14]shows that the elastic OMIT reflection coefficient is given
by r½!pr  ¼ 1  icp GRdd ½!pr , where !pr is the frequencyof the weak probe beam and cp is the contribution to thetotal cavity  from the coupling to the drive port.
In the limit of interest where   E , there are no offdiagonal polariton GFs or self-energies P[14]. As a result,
GRdd ½! will be given as GRdd ½! ¼  ðC GR ½! þD fGR ½!g Þ, where the change-of-basis coefficientsC and D are given in Ref. [14]. The Dyson equations
for the polariton retarded GFs are

053602-2

GR ½! ¼ ð!  E þ i =2  R ½!Þ1 ;(10)

PRL 111, 053602 (2013)

PHYSICAL REVIEW LETTERS

where  is the effective damping rate of the  polariton[14]. Using the effective Hamiltonian in Eq. (7), a standard
Keldysh calculation yields that to second order in g, thepolariton self-energies take the simple forms2~
g2 ð1 þ 2n  Þ;!  2E þ i
4~g2 ðn   n þ Þ:R ½! ¼!  ðEþ  E Þ þ iðþ þ  Þ=2Rþþ ½! ¼

(11a)
(11b)

Here, n  is the effective thermal occupancy of the polariton [14]; for g~ ¼ 0, we have hc^ y c^  i ¼ n  . We
have taken the limit g=!M ! 0 and hence neglected theeffects of the nonresonant terms H^ NR in Eqs. (11);the explicit corrections due to these terms are given in
the Supplemental Material [14].
Equations (11) are central results of this work.
Equation (11a) describes the fact that a single þ polariton
can resonantly turn into two  polaritons and describes the
hybridization between these states that occurs for large
enough g. To see this explicitly, we consider the case ofexact resonance (i.e., Eþ ¼ 2E ) and write1
GRþ ½! ¼

X þ
1  i 24þ

2 ¼ !  Eþ þ i 24þþ þ þ

qﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
g2 ð1 þ 2n  Þ  ð2  þ Þ2 =16:þ ¼ 2~

;

(12)

(13)

For g~ * , we see that the þ polariton GF has two poles,
corresponding to the new hybridized eigenstates. We stressthat these eigenstates do not correspond to a fixed excitation number. Note that unlike the undriven system [8,9],
the effects of the nonlinear interaction can be significant
even if g  !M . Also note that the resonant coupling
between jþi and j i states is enhanced at finite
temperature by a standard stimulated emission factor(1 þ 2n  ). The form of this GF and self-energy are reminiscent to the photon GF for ordinary OMIT, where aphoton can resonantly turn into a phonon [11]; however,that effect does not involve any temperature-dependentenhancement.  in Eq. (11b) describes a process where
the propagating  polariton of interest interacts withan already present  polariton to turn into a þ. As thisprocess requires an existing density of polaritons, it isstrongly suppressed at low temperatures.
We note that it is possible to use resonance to enhance
the nonlinear optomechanical interaction without strong
driving, if one instead considers a system where two cavity
modes interact with a single mechanical resonator [16–18];
see also Ref. [19] for an alternate scheme based on two
cavity modes. Our approach has the benefit of only requiring a single cavity mode; further, for drive detunings near
 ¼ !M , it also has a natural resistance against
mechanical heating, as mechanical contribution to thepolariton temperature scales as n th =, where n th is theweek ending
2 AUGUST 2013mechanical thermal occupancy. While a low temperature is
not essential for the density-of-states effects describedabove, it is essential for the correlation effects discussed
below. Finally, for superconducting microwave cavities,
the cavity linewidth  has a strong contribution fromtwo-level fluctuators and thus improves if one strongly
drives the cavity (as in our scheme).
Red-sideband drive.—For a detuning  ¼ !M , thepolariton resonance occurs when G ¼ 0:3!M . For thisdetuning, both polaritons are almost equal mixtures ofphoton and phonon operators. One finds  ¼ ð þ Þ=2
and that the resonant interaction strength g~ ’ 0:37g.
Because H^ 1 does not conserve the number of photons
and phonons, the polaritons are not eigenstates of d^y d^ þ
^ as a result, even at zero temperature, the effectiveb^y b;
thermal occupancies scale as n  / ðG=!M Þ2  1 [14].
The inset of Fig. 1 shows the evolution of the cavity densityof states for these parameters as g is increased from zero.
For g ¼ 0, one sees two symmetric peaks corresponding to
the two polaritons, i.e., the well-known normal-mode splitting [20,21]. As g increases, these peaks develop a marked
asymmetry. For g , a clear splitting of the þ peakoccurs, corresponding to the resonant hybridization of
one and two polariton states. Figure 1 also shows theresults of a numerical (but nonperturbative) masterequation simulation [14], showing our analytic approachis reliable even for moderately strong g.
Large-detuned drives.—The resonant-polariton interaction is also interesting for drives far from the red sideband,
where the value of Gres  !M . For a laser detuningnear the minimum possible value  ¼ 2!M at which

FIG. 1 (color online). Main: þ polariton resonance in thecavity density of states, for various values of the nonlinear
interaction strength g (as indicated), as obtained from Eq. (12)(with the inclusion of energy shifts from H^ NR [14]). For all plots,
the laser drive is at the red sideband  ¼ !M , and G ¼ 0:3!Mto ensure the resonance condition Eþ ¼ 2E ; we also take
!M = ¼ 50, T ¼ 0, and  ¼ 104 . The peak splitting signals
the hybridization of a þ polariton with two  polaritons. The
dashed curve is the result of a master-equation simulation for
g ¼  [14]. Inset: Full density of states, same parameters,
showing the asymmetry between þ and  polariton resonances.

053602-3

PRL 111, 053602 (2013)

PHYSICAL REVIEW LETTERS

resonance is possible (and setting G ¼ Gres ), the polaritons
are each either almost entirely phonon or photon,
implying a very small value of g~ / gG=!M . However, as
the  polariton is now almost purely phononic, its small
damping rate and potentially large thermal occupancy
enhance the self-energy in Eq. (11a). We can quantify theseeffects by considering the value of d ½! ¼ Eþ , which
will be suppressed by the hybridization physics described
here. One finds in general thatd ½Eþ  ¼

2= 1
4~
g2 ð1þ2n  Þ; Ceff ¼:þ 1þCeff
þ 

(14)

For a detuning  2!M ,  ’  þ ð8=9ÞG2 =!2M ,where the second term is the effective optical damping.
As a result, Ceff / ðg=Þ2 ð1 þ 2n  Þ. Even when Ceff  1,
the small density-of-states suppression is extremely sharpin frequency (with a width   ), leading to an observable feature. Further, the feature can be strongly enhancedby simply increasing the mechanical temperature.
Figure 2 shows the reflection coefficient jr½!pr j2 thatwould be measured in a standard OMIT experiment for
both the case of a control laser detuning near !M[Fig. 2(c)] and near 2!M [Fig. 2(d)]. As discussed, the!pr dependence of r reflects the structure in the cavity
density of states. The dashed curves in each panel are thepredictions of the usual linearized optomechanicalHamiltonian (i.e., g ¼ 0, G  0). For g , one sees aclear splitting of the reflection dip associated with the þpolariton resonance. In contrast, for g  , one canstill see effects of the nonlinear interaction by taking
(a)

(b)

(c)

(d)

FIG. 2 (color online). (a) Behavior of g~, Gres , and Ceff as afunction of detuning  of the main laser drive. (b) Schematic of
OMIT experiment, where a weak probe beam at a frequency !pr(defined in the lab frame) is reflected from the cavity. (c) OMIT
reflection coefficient, with , g given, G ¼ Gres ½, cp = ¼
0:25, and the remaining parameters identical to Fig. 1. The
dashed curve is the prediction of the linearized theory (g ¼ 0,
G  0). (d) Same as (c) but with a larger control laser detuning ’ 2!M . In this regime, the effect of even a small g   canbe resolved as it yields a very sharp feature.

week ending
2 AUGUST 2013 2!M , as one obtains an extremely narrow feature
in r½!pr . As an OMIT experiment involves measuring
average reflected power (as opposed to fluctuations),
such small, sharp features can be resolved. The valueof g= ¼ 0:01 used in Fig. 2(d) is comparable to thevalue 0.007 achieved in Ref. [22]; even larger values
of g= have been realized in cold-atom experiments(e.g., Ref. [23]).
Induced Kerr interaction.—The nonlinear interactionin the resonant regime defined by Eqs. (6) and (7) leads
to a strongly enhanced two-particle interaction between polaritons, mediated by the exchange of a þ polariton(Fig. 3). In a weakly driven optomechanical system,
Eq. (3) implies that phonons can mediate an effective
photon-photon interaction; however, as the virtual phonon
is off resonance, this interaction U / g2 =!M . In contrast,
the resonance condition Eþ ¼ 2E yields an induced
interaction Ures / g~2 =, an enhancement by a large
factor / !M =.
To assess the effects of the polariton-polariton interactions, we weakly drive our system with a second probe
tone and consider the g2 correlation functions g2u ¼^ u^ y ui^ 2 , where u ¼ b, d, cþ , c . g2u is a meahu^ y u^ y u^ ui=hsure of interaction induced correlations; g2 1 signifies
nonclassical correlation. Given the strong interaction experienced by  polaritons when the resonance condition
Eþ ¼ 2E is achieved, we expect that if the cavity isdriven at the E resonance, g2 will drop below 1. This
is indeed the result of a numerical, master-equation-based
calculation (see Fig. 3 and Ref. [14]). An analytic calculation based on a reduced state space (similar to that inRef. [18]) reproduces these results. For a weak probe driveat the E frequency, it yields [14]FIG. 3 (color online). Inset: Resonant interaction between polaritons. Main: Numerically calculated g2 correlation functionfor  polaritons (g2 ) and phonons (g2b ), in the presence of anadditional weak probe laser (frequency !pr ). Here, g ¼ ,  ¼
!M , and G ¼ 0:3!M ¼ Gres . We have taken !M = ! 1 tosuppress nonresonant interaction effects. The probe amplitude is¼ 0:2 (g2 ) and ¼ 0:3 (g2b ). Both phonon and polariton
g2 functions drop below 1 due to the interactions, indicating
nonclassical correlations despite the fact g=!M ’ 0. The dashedcurve is the result of an analytic theory (see Ref. [14]).

053602-4

PRL 111, 053602 (2013)g2 ¼

PHYSICAL REVIEW LETTERS
1:
1 þ 4~
g2 =2

(15)

One also finds nonclassical correlations for photons andphonons. Shown in Fig. 3 is the phonon g2 function g2b(for the same parameters); it clearly drops below 1. The
double-peak structure of this curve is the result of the drive
inducing correlations between  and þ polaritons; it also
^ (see the Supplemental
occurs in the behavior of hb^y bi
Material for more details [14]).
Conclusions.—We have presented a systematic
approach for describing nonlinear interaction effects in adriven optomechanical system, identifying a regime wherea resonance enhances interactions between polaritons. We
have discussed how this would manifest itself in a OMITstyle experiment, as well as in g2 correlation functions. The
polariton interactions we describe could be extremelyinteresting when now considered in lattice systems orwhen considering the propagation of pulses.
We thank W. Chen and A. Nunnenkamp for useful discussions. This work was supported by CIFAR, NSERC, andthe DARPA ORCHID program under a grant from AFOSR.
Note added.—Recently, we became aware of a related
work by Børkje et al. [24].

[1] J. D. Teufel, T. Donner, D. Li, J. W. Harlow, M. S. Allman,
K. Cicak, A. J. Sirois, J. D. Whittaker, K. W. Lehnert, andR. W. Simmonds, Nature (London) 475, 359 (2011).
[2] J. Chan, T. P. M. Alegre, A. H. Safavi-Naeini, J. T. Hill,
A. Krause, S. Gröblacher, M. Aspelmeyer, and O. Painter,
Nature (London) 478, 89 (2011).
[3] D. W. C. Brooks, T. Botter, S. Schreppler, T. P. Purdy,
N. Brahms, and D. M. Stamper-Kurn, Nature (London)
488, 476 (2012).
[4] A. H. Safavi-Naeini, S. Groeblacher, J. T. Hill, J. Chan,
M. Aspelmeyer, and O. Painter, arXiv:1302.6179.
[5] S. Groeblacher, K. Hammerer, M. R. Vanner, and M.
Aspelmeyer, Nature (London) 460, 724 (2009).
[6] J. D. Teufel, D. Li, M. S. Allman, K. Cicak, A. J. Sirois,
J. D. Whittaker, and R. W. Simmonds, Nature (London)
471, 204 (2011).

week ending
2 AUGUST 2013[7] E. Verhagen, S. Deleglise, S. Weis, A. Schliesser, andT. Kippenberg, Nature (London) 482, 63 (2012).
[8] P. Rabl, Phys. Rev. Lett. 107, 063601 (2011).
[9] A. Nunnenkamp, K. Børkje, and S. M. Girvin, Phys. Rev.
Lett. 107, 063602 (2011).
[10] A. Kamenev and A. Levchenko, Adv. Phys. 58, 197(2009).
[11] G. S. Agarwal and S. Huang, Phys. Rev. A 81, 041803(R)(2010).
[12] S. Weis, R. Riviere, S. Deleglise, E. Gavartin, O. Arcizet,
A. Schliesser, and T. Kippenberg, Science 330, 1520
(2010).
[13] A. H. Safavi-Naeini, T. P. M. Alegre, J. Chan, M.
Eichenfield, M. Winger, Q. Lin, J. T. Hill, D. E. Chang,and O. Painter, Nature (London) 472, 69 (2011).
[14] See Supplemental Material at http://link.aps.org/
supplemental/10.1103/PhysRevLett.111.053602
foradditional
informationonthenormal-mode
transformation, the self-energy calculation, thenumerical master-equation simulation, and thecalculation of g2 correlation functions.
[15] Note that when the resonance condition Eþ ¼ 2E issatisfied (i.e., G ¼ Gres ½), the linear optomechanicalsystem is always stable, as can be confirmed by applying
the standard Routh-Hurwitz stability conditions.
[16] M. Ludwig, A. H. Safavi-Naeini, O. Painter, and F.
Marquardt, Phys. Rev. Lett. 109, 063601 (2012).
[17] K. Stannigel, P. Komar, S. J. M. Habraken, S. D. Bennett,
M. D. Lukin, P. Zoller, and P. Rabl, Phys. Rev. Lett. 109,013603 (2012).
[18] P. Komar, S. D. Bennett, K. Stannigel, S. J. M. Habraken,
P. Rabl, P. Zoller, and M. D. Lukin, arXiv:1210.4039v1.
[19] X.-Y. Lü, W.-M. Zhang, S. Ashhab, Y. Wu, and F. Nori,
arXiv:1210.8299v1.
[20] F. Marquardt, J. P. Chen, A. A. Clerk, and S. M. Girvin,
Phys. Rev. Lett. 99, 093902 (2007).
[21] J. M. Dobrindt, I. Wilson-Rae, and T. J. Kippenberg, Phys.
Rev. Lett. 101, 263602 (2008).
[22] J. Chan, A. H. Safavi-Naeini, J. T. Hill, S. Meenehan,and O. Painter, Appl. Phys. Lett. 101, 081115 (2012).
[23] K. W. Murch, K. L. Moore, S. Gupta, and D. M. StamperKurn, Nat. Phys. 4, 561 (2008).
[24] K. Børkje, A. Nunnenkamp, J. D. Teufel, and S. M. Girvin,following Letter, Phys. Rev. Lett. 111, 053603(2013).

053602-5

Supplemental information for “Nonlinear interaction effects in a strongly driven
optomechanical cavity”
Marc-Antoine Lemonde,1 Nicolas Didier,1, 2 and Aashish A. Clerk1
1

2

I.

Department of Physics, McGill University, Montréal, QC Canada H3A 2T8
Département de Physique, Université de Sherbrooke, Sherbrooke, QC Canada J1K 2R1

LANGEVIN EQUATIONS FOR LINEARIZED OPTOMECHANICAL SYSTEM IN THE POLARITONBASIS

We start by diagonalizing the linearized optomechanical Hamiltonian given in Eq. (2) of the main text, working as
always in a displaced interaction picture set by the laser drive on the cavity. Introducing
~ =
X



b̂ dˆ b̂† dˆ†

T

,



~ = ĉ− ĉ+ ĉ† ĉ† T ,
Y
−
+

(S1)

the diagonalization can be expressed in terms of a 4 × 4 real matrix U:~ = U · X.
~
Y

(S2)

U can be found by standard means, though its general form is both cumbersome and unenlightening; we define
V = U−1 . It is slightly less unwieldy in the special (but relevant) case of a drive at the red-detuned mechanicalsideband ∆ = −ωM . In this case, we have simply:hi1ĉ± = p(E± − ωM )(dˆ† ± b̂† ) + (E± + ωM )(dˆ ± b̂) ,(S3)
8ωM E±
while the inverse transformation is defined by:i
1 h − ω√
M +E+
M −E−
M −E+
M +E−
, ω√
, − ω√, ω√
V1j = √
E−
E+
E−
E+
8ωM
i
h
ωM +E−
1
M +E+
M −E−
M −E+
√
, ω√
, ω√
, ω√
V2j = √
E−
E+
E−
E+
8ωM

(S4)
(S5)

Next, we include the coupling of our system to the cavity and mechanical dissipative baths in the standard way,treating these baths as Markovian over frequencies of interest. Consider first the coupling to the cavity bath. Prior
to making displacement and interaction-picture transformations, the system-bath coupling has the form (see, e.g.,
Ref. 1):
r

X

κ Xˆ
(S6)
Ĥκ =
ωj fˆj† fˆj , Ĥκ,int = −i
fj − fˆj† â + â† ,2πρc jjwhere fˆj is the lowering operator for a bath mode, and ρc is the density of states of bath modes (which we treat
to be constant over frequencies of interest). In now moving to an interaction picture at the laser frequency, we
also transform the bathmodes,
i.e. fˆj → ifˆj e−iωL t . Formally, the interaction picture transformation involves ahPPunitary Û (t) = exp −iωL t â† â +fˆ† fˆjand transforms the bath Hamiltonian to Ĥκ =(ωj − ωL )fˆ† fˆj .
j

j

j

j

In this interaction picture, the counter-rotating terms in Ĥκ,int will explicitly oscillate at ±2ωL = ±2(ωc + ∆).
Even if we now write our photon operators â in terms of polariton operators, there is no possibility of having theseterms becoming resonant, as the cavity frequency is much larger than any other frequency scale in the problem (i.e.
ωc  |∆|, ωM , E+ , E− ). As such, one can safely make a rotating-wave approximation in the photon basis, resulting
in a standard system-bath interaction which is stationary in the interaction picture:rκ X  ˆ† ˆ ˆ† ˆ Ĥκ,int = i
fj d − d fj
2πρc jriκ X h ˆ† 
= −i
fj V21 ĉ− + V22 ĉ+ + V23 ĉ†− + V24 ĉ†+ − h.c.
(S7)2πρc j2
We have also made the displacement transformation â = dˆ+ ā as discussed in the main text. Note that in writing dˆ interms of polariton operators via Eq. (S2), we obtain terms of the form fˆj† ĉ†σ , which can cause polariton heating even
if the cavity bath is at zero temperature. Such terms are physical and must be retained. Formally, In the interaction
picture the bath now has negative frequency modes which can make such processes resonant. In more physical terms,
the combination of zero-point bath fluctuations with the cavity driving can excite the polaritons. This mechanism
has been discussed in other contexts under the name “quantum activation” by various authors [2, 3].
We turn now to the mechanical bath, where the basic interaction Hamiltonian can be written in an analogous way:
r


X
γ X
†
Ĥγ =
ωj ĝj ĝj , Ĥγ,int = −iĝj − ĝj† b̂ + b̂† .
(S8)2πρm jjHere, ĝj is a lowering operator for a mechanical bath oscillator, and ρm is the density of states of mechanical bathmodes (which we also treat to be constant over frequencies of interest). As there is no direct driving of the mechanicalresonator, the analysis here is somewhat simpler. We first re-write the phonon operator b̂ in the polariton basis, andthen make a rotating-wave approximation. The justification is that counter-rotating terms such as ĝj† ĉ†σ can never bemade resonant; as the mechanical resonator is not driven, there is no quantum activation mechanism involving themechanical bath. We thus obtainX r κm X  †σĝj ĉσ + h.c. ,(S9)Ĥγ,int '2πρm jσwhere (V = U−1 )2κM− = γ (V11 + V31 ) ,2κM+ = γ (V12 + V32 ) .

(S10)Having now established the correct system-bath coupling Hamiltonians in the polariton basis we wish to use, we
can derive the Heisenberg-Langevin equations for our system in the standard manner. As each bath couples to both+ and − polaritons, one obtains off-diagonal damping terms, (e.g. the mechanical bath produces a force on the +polariton that is proportional to the amplitude of the − polariton). Such terms will dynamically couple + and −polaritons, and can be included in our theory in a straightforward manner (i.e. by including Green functions that are
off-diagonal in the polariton index). However, in the weak-dissipation limit of interest (κ, γ  E+ , E− , |E+ − E− |),the mixing effects induced by such terms is strongly suppressed. As such, we will drop off-diagonal damping terms,
resulting in the form:

√
κσ 
(S11)
ĉσ (t) − κσ ξˆσ (t),∂t ĉσ (t) = − iEσ +2where the polariton damping rates κσ are given by:hi22κ− = κM
+
κ
(V
)
−
(V
)
,
21
23
−

h
i
2
2
κ+ = κM
+
κ
(V
)
−
(V
)
.
22
24
+

(S12)

The noise operators ξˆσ (t) are each linear combinations of the input noise emanating from the mechanical and cavity
baths. In the interaction picture we use, we will be sensitive to noise in the cavity bath at frequencies near ωc ,and noise in the mechanical bath at frequencies near ωM . In the limit of interest where the physical temperature
T  ~ωc /kB , there will be no thermal noise incident from the cavity bath at the frequencies of interest. Also, as
we focus on regimes where the polariton damping rates are much smaller than their energies, we can treat the noise
operators as being white noise (as is standard in input-output theory treatments). We thus haveDE
DEξˆσ† (t)ξˆσ (t0 ) = nσ δ (t − t0 ) ,ξˆσ (t)ξˆσ† (t0 ) = (1 + nσ ) δ (t − t0 ) ,(S13)
where the effective temperatures of the two noises are given by:


1  M
1  M
2
2
n− =
κ− nB [E− ] + (V23 ) κ , n+ =
κ+ nB [E+ ] + (V24 ) κ .
κ−
κ+

(S14)

nB [Eσ ] denotes the Bose-Einstein distribution function evaluated at energy Eσ and the mechanical bath temperature.
Finally, one also finds that ξˆ+ and ξˆ− are correlated with one another. Similar to the situation of off-diagonal
damping terms, such noise correlations could easily be included in our theory; however, as they play no role in theregime of interest where Eσ , |E+ − E− |  κ, γ, we drop them in what follows. One also finds that anomalous noisecorrelators can be non-zero (e.g. hξˆσ (t)ξˆσ0 (0)i). Again, while such terms can be retained in our theory, they play norole in the weak-damping regime of interest, and we hence drop them in what follows.

3
II.

UNPERTURBED POLARITONS GFS

The standard definitions of the three relevant GFs needed in the Keldysh technique are (σ stands for ±, but similarˆ
definitions hold for ĉσ = b̂ or d):
Z ∞
Dh
iE
GR
[ω]
≡
−i
dtθ(t) ĉσ (t), ĉ†σ0 (0) eiωt ,
(S15)
σσ 0
−∞
Z ∞
Dh
iE
dtθ(−t) ĉσ (t), ĉ†σ0 (0) eiωt ,
[ω]
≡
i
(S16)
GA
σσ 0
−∞
Z ∞ Dn
oE
GK
[ω]
≡
−i
dt ĉσ (t), ĉ†σ0 (0) eiωt .
(S17)
σσ 0
−∞

The retarded and advanced GFs keep track of spectral information, whereas the Keldysh Green function GK also
keeps track of the occupancy of states. For the linearized (non-interacting) theory, the GFs are easily obtained fromthe Langevin equations in Eq. (S11). These free GFs (which we denote by G) are diagonal (only non-zero for σ = σ 0 )and given by:
R
Gσσ
[ω] =

1
,
ω − Eσ + i κ2σ

(S18)

K
R
A
Gσσ
[ω] = (2nσ + 1) Gσσ[ω] − Gσσ
[ω] ,



(S19)

 R
∗
A
and Gσσ[ω] = Gσσ[ω] .
Finally, as there are no off-diagonal polariton Green functions, we can use Eq. (S2) to write the photon retarded
Green function GR
dd [ω] as:
X
 R
∗ 
GR
Cσ GR
,
(S20)
dd [ω] =σσ [ω] + Dσ Gσσ [−ω]
σ

where
2
2
2
2
C− = V21, C+ = V22, D− = V23, D+ = V24
.

(S21)

As the Green functions remain diagonal in the polariton index even with the nonlinear interaction (in the regime ofinterest, see below), the above relation also holds for the full Green functions (i.e. including the self-energy associatedwith g).

III.

Ĥ2 IN THE POLARITON BASIS

Using the change of basis matrix V = U−1 (c.f. Eq.(S2)), we can re-write the non-linear interaction Ĥ2 in Eq. (3)of the main text in the polariton basis via:
h
i h
i h
i
h
i 
~
~
~ + VY~Ĥ2 = g VY
VY
VY
.
(S22)
4

2

1

3

A,B
Expanding this equation allows one to obtain the interaction coefficients gσσ0 σ 00 in terms of g and matrix elements of† †V. In particular, the coefficient g̃ of the resonant interaction process ĉ− ĉ− ĉ+ will be given by:g̃ = g(V22 V21 (V11 + V13 ) + V23 V24 (V11 + V13 ) + V23 V21 (V14 + V12 )).

(S23)

As the normal-mode transformation described by U mixes raising and lowering operators, Ĥ2 will not be normalordered when written in terms of polariton operators (even though it is normal ordered when written in terms ofphoton and phonon operators). Normal-ordering Ĥ2 in the polariton basis yields the form:
Ĥ2 =: Ĥ2 : + (A− ĉ− + A+ ĉ+ + h.c.)(S24)

where the colons indicate normal-ordering in the polariton basis, and the constants Aσ ∝ g. The first term is thenormal-ordered polariton interaction written in Eq. (5) in the main text. We next make a unitary displacement

4

FIG. 1: (a) Feynman diagrams describing the dominant + and − polariton self-energies (to second order in g) in the casewhen the resonance condition E+ = 2E− is met. We have shown diagrams for the retarded self-energies; similar diagrams
also determine the advanced and Keldysh self-energies. The structure of propagators in Keldysh space is indicated by writing
explicit classical (cl) and quantum (q) indices at the ends of each Green function. (b) The non-resonant interaction terms induce
effective two-particle interactions between polaritons, as described by HN R in Eq. (8) of the main text. These two-particle
interactions modify the resonant self-energies depicted in Fig. (1); shown above are the diagrams describing the modificationof the + polariton self-energy by the effective interaction U−− . The summation of these ladder diagrams results in the simple
energy shift indicated in Eqs. (S27),(S28).

transformation of the form ĉσ → c̄σ + ĉσ , where the constants c̄σ are chosen to eliminate all linear-in-ĉσ terms in the
Hamiltonian; to leading order in g, c̄σ = −Aσ /Eσ . The resulting coherent Hamiltonian has the form:

 







ĉ
ĉ
ĉ− ĉ+ Z2 − + h.c. .
Ĥ2 = : Ĥ2 : + ĉ†− ĉ†+ Z1 − +
ĉ
ĉ
+

+

Here Z1 , Z2 are 4 × 4 matrices whose entries are all order g 2 /Eσ ∼ g 2 /ωM . The quadratic terms on the RHS describeP
corrects to the linear Hamiltonian Ĥ1 = σ Eσ ĉ†σ ĉσ arising from the displacement transformation. In principle, one
could combine these with the terms in Ĥ1 , re-diagonalize the resulting quadratic Hamiltonian, obtaining a new basisof non-interacting polaritons. However, to leading order in g 2 /ωM , all that is needed is to retain the diagonal elementsof Z1 , which simply shift the polariton energies but do not change their wavefunctions. These energy shifts (denotes
+ , − ) are just absorbed into the definition of our free Green functions. We are then left with a normal-orderedinteraction Hamiltonian Ĥ2 in the polariton basis which can be addressed perturbatively.
IV.

PERTURBATIVE TREATMENT

We have calculated the full Keldysh self-energy corresponding to Ĥ2 to order g 2 , without any further approximation.
We have also done this calculation in the case where G is not so large, such that one should work in the original basisof photons and phonons. Results of these full calculations will be presented elsewhere. Here, like in the main text,
we will focus on the resonant-interaction regime described in the main text, where the condition E+ = 2E− enhances
certain scattering processes. In this regime, the dominant self-energy processes for + and − polaritons are depicted
in Fig. 1. Using the standard rules of the Keldysh technique [4], these diagrams correspond to:Zn− − n+dω 0 KRA
K
2
G−− [ω 0 ] G++
g2, (S25)[ω 0 + ω] + G−−[ω 0 − ω] G++[ω 0 ] = 4e
ΣR
[ω]
=
2ie
g
−−
+
2π
ω − (E+ − E− ) + i κ− +κ2
Zdω 0 K(1 + 2n− )
2
ΣR
[ω]
=
2ie
g
G [ω 0 ] G R [ω − ω 0 ] = 2e
g2
.
(S26)
++
2π −−
ω − 2E− + iκ−Note that because of the resonance condition, these self-energies scale as g 2 /κ, whereas all other self-energy diagrams
are suppressed by an additional small parameter κ/ωM . They can thus be neglected in the limit κ/ωM → 0. For asmall but realistic value of κ/ωM , the small-energy shifts associated with the non-resonant interaction terms in Ĥ2
can shift the condition needed for resonance. To describe these small shifts (which can be important for realistic
parameters), it is sufficient to use standard second-order perturbation theory to treat the non-resonant terms. Thisis conveniently done via a Schrieffer-Wolff transformation, where a unitary transformation is performed to eliminate
the non-resonant terms to leading order in the Hamiltonian. The procedure generates new terms however at second
order in g. Keeping only such terms which do not change the total number of polaritons, we obtain the general form
given in the Hamiltonian ĤN R (Eq. (8) of the main text). Note that in this equation, the energy shifts δσ include
both the shifts arising from the non-resonant terms, as well as the energy shifts σ coming from the normal-ordering
procedure. Both such terms scale as g 2 /ωM .

5Including the energy shifts associated with the non-resonant processes (as described by ĤN R in Eq. (8) of the maintext), the self-energies above are modified to:n− − n+,Ẽ+− = E+ − E− + δ+ − δ− + U+− (n̄− − n̄+ ) ,+ω − Ẽ+− + i κ− +κ
2(1 + 2n− ),Ẽ−− = 2 (E− + δ− ) + 2U−− (1 + 2n̄− ) .
ΣRg2++ [ω] = 2e
ω − Ẽ−− + iκ−ΣRg2−− [ω] = 4e

(S27)
(S28)

The corrections due to the Kerr-type interaction constants U+− and U−− can be obtained by including ladder diagramsin the self-energies, as shown in Fig. 1b.
V.

OMIT REFLECTION COEFFICIENTIn an OMIT style-experiment, in addition to the main driving laser (which gives rise to the many-photon interactionG), a second weak drive tone (the “probe”) is applied at a frequency ωpr to the cavity. This driving is described bya term in the Hamiltonian:√Ĥpr = −i κcp dˆ† d¯in,pr e−iω̃pr t + h.c.
(S29)
where we work in the interaction picture determined by the main drive laser frequency, and hence ω̃pr = ωpr − ωL .
κcp parameterizes the coupling of the drive port to the cavity: for a one-port cavity, the total cavity damping rate
κ = κcp + κint , where κint describes internal cavity losses.
As the amplitude d¯in probe drive is weak, we can use standard linear response theory to describe its effects. The
Kubo formula thus tells us that to first order in d¯in,pr , the change in the cavity amplitude will be given by:
Z ∞
D
E


√
0 −iω˜pr t0
0 iω˜pr t0ˆ
dt0 d¯in,pr GR− d¯∗in,pr GR
(t
−
t
)e
(S30)
δ d(t) = −i κcp¯
dd (t − t )e
dd
−∞

√
iω̃pr t ¯∗
= −i κcp e−iω̃pr t d¯in,pr GR
din,pr GR(S31)
dd [ω̃pr ] − edd¯[−ω̃pr ] .
Here, GRdd [ω] is the retarded Green function of the cavity as defined above, calculated to zeroth order in Ĥpr , but
including the effects of the nonlinear interaction g. GR[ω] is the corresponding anomalous Green function defined as:
dd¯
Z
iE
Dh
ˆ d(0)
ˆGRdtθ(t) d(t),eiωt ,(S32)
dd¯ [ω] ≡ −i
Now, the standard input-output relation between input, output and cavity fields is [1]:√ˆdˆout (t) = dˆin (t) + κcp d(t)(S33)

Taking the average value of this equation, and defining the elastic amplitude reflection coefficient r[ωpr ] as theamplitude of hdˆout (t)i at the probe frequency divided by d¯in , we obtainr[ωpr ] = 1 − iκcp GR
dd [ωpr ],(S34)

as given in the main text.
VI.

MASTER EQUATION SIMULATION

Starting from the system-bath Hamiltonians written in the polariton basis (Eqs. (S7) and (S9)), one can trace overthe dissipative baths and derive a master equation for the reduced density matrix ρ̂ describing the polaritons using
standard Born-Markov approximations [5]. One obtains:hi∂t ρ̂ = −i Ĥ1 + Ĥ2 , ρ̂ + Lρ̂(S35)
where the coherent system Hamiltonian Ĥ1 + Ĥ2 is written without any approximation, and the super-operator L
describes the effects of the dissipative baths via standard Lindblad terms,
L = κ− (1 + n̄− )L[ĉ− ] + κ− n̄− L[ĉ†− ] + κ+ (1 + n̄+ )L[ĉ+ ] + κ+ n̄+ L[ĉ†+ ],

(S36)

6
with
1
L[ĉ]· = ĉ · ĉ† − {ĉ† ĉ, ·}.
2

(S37)

Note that this master equation corresponds to each polariton seeing independent thermal baths, in direct analogy to
the form of the quantum Langevin equations in Eq. (S11).
While Eq. (S35) is not a convenient starting point for deriving analytic results, it does allow us to numerically study
the system without having to assume a small value of g. Using Eq. (S35) and the quantum regression theorem [5],we have numerically calculated the cavity density of states ρd [ω] = −Im GR
dd [ω]/π, finding good agreement with ouranalytic perturbative results even for g as large as κ (see Fig. 1 in main text). Note that to find agreement with thesenumerical results, it was crucial to include in the analytic theory the corrections associated with the non-resonantinteraction processes, ĤN R (c.f. Eqs.(S27) and (S28)).
Finally, we have also numerically studied our system using a more conventional master equation, in which thedissipation is described by Lindblad superoperators which act in the photon and phonon basis, i.e. replace L in
Eq. (S35) with L0 , where:ˆ + γ(1 + n̄M )L[b̂] + γ n̄M L[b̂† ],L0 = κL[d](S38)

and n̄M is a Bose-Einstein distribution evaluated at the mechanical frequency ωM and mechanical bath temperature.
For the parameters studied in the paper, this conventional master equation yields results very similar to those obtained
from Eq. (S35).
VII.

g2 CORRELATION FUNCTIONSAs discussed in the main text, the resonant two-particle interaction between − polaritons can lead to non-classicalvalues of the g2 correlation functions. In Fig. 3 of the main text, we have used a numerical solution of our masterequation and the quantum regression theorem to calculate g2 functions of both photons and − polaritons, taking
g = κ, ∆ = −ωM and tuning G to Gres = 0.3ωM to ensure the resonance condition E+ = 2E− . For simplicity, we also
take the large ωM limit (i.e. ωM /κ → ∞, g/ωM → 0), so that corrections due to non-resonant interaction terms in Ĥ2
can be neglected. We have also included a weak probe drive on the cavity at a frequency ωpr close to E− by including
a term of the form Ĥpr (c.f. Eq.(S29)). When writing this driving term in terms of polarities, we will have terms that
effectively drive both − and + polaritons. However, in the large ωM limit and for ωpr ' E− , the direct driving termon the + polaritons will be strongly off-resonance and can be neglected. The probe field driving Hamiltonian thus
reduces to:Ĥpr = e−iωpr t ĉ†− + h.c.

(S39)

Eq. (15) of the main text and the dashed curve in Fig. 3 are the results of a simple analytic theory which accurately
describes the − polariton g2 function in the limit of an extremely weak probe drive amplitude . For a weak drive(and only keeping the resonant interaction process in Ĥ2 ), we can restrict attention to the four lowest eigenstates of
the coherent Hamiltonian. These are:|ψ0 i = |0, 0i,|ψ1 i = |1, 0i,

E0 = 0,E 1 = E− ,|ψ2 i = − sin θ|0, 1i + cos θ|2, 0i,

E2 = 12 E+ + E− −q( 12 E+ − E− )2 + 2g̃ 2 ,(S42)

E3 = 12 E+ + E− +q( 12 E+ − E− )2 + 2g̃ 2 ,(S43)|ψ3 i = cos θ|0, 1i + sin θ|2, 0i,

(S40)
(S41)

√
where |n− , n+ i denotes a Fock state having n− − polaritons and n+ + polaritons, and tan 2θ = 2 2g̃/(E+ − 2E− ).
The system is driven by adding the probe drive Hamiltonian Ĥpr given in Eq. (S39). Including the drive, the wave
function is written in terms of the above states as:|ψi = A00 |0, 0i + A10 |1, 0i + A20 |2, 0i + A01 |0, 1i.

(S44)

We take into account the dissipation with the anti-Hermitian Hamiltonian Ĥdamping = − 12 i(κ− ĉ†− ĉ− + κ+ ĉ†+ ĉ+ ). In††

the rotating frame corresponding to the unitary operator ei(ĉ− ĉ− +2ĉ+ ĉ+ )ωpr t , the Hamiltonian reads
HRWA = (∆− − 2i κ− )ĉ†− ĉ− + (∆+ − 2i κ+ )ĉ†+ ĉ+ + g̃(ĉ†+ ĉ− ĉ− + ĉ+ ĉ†− ĉ†− ) + (ĉ†− + ĉ− ),

(S45)

7

FIG. 2: Plot of average polariton, photon and phonon numbers versus probe frequency ωpr , for parameters identical to thoseused in Fig. 3 of the main text. We have taken a drive strength of  = 0.2κ, where  is defined in Eq. (S39). The double-peakstructure in the phonon population is due to higher-order effects of the probe drive which induce correlations between + and− polarities.

with the detuning ∆− = E− − ωpr and ∆+ = E+ − 2ωpr = 2∆− at the resonance E+ = 2E− . The dynamics isgiven by ∂t |ψi = −iĤ|ψi. Ignoring other states than those involved in |ψi, the steady state number and correlation
g2− = hĉ†− ĉ†− ĉ− ĉ− i/hĉ†− ĉ− i2 arehĉ†− ĉ− i =2,∆2− + (κ− /2)2[∆2− + (κ− /2)2 ][4∆2− + (κ− /2)2 ]g2− =,[2∆2− − g̃ 2 − (κ− /2)2 ]2 + 9∆2− (κ− /2)2

(S46)
(S47)

where we approximate A00 = 1, take κ+ = κ− and keep the leading terms in . The minimum value of g2− is foundat ωpr = E− with g2− = κ4− /(4g̃ 2 + κ2− )2 and the correlation tends toward unity for infinite detuning.
We comment briefly on the unusual double-peak resonance exhibited by g2b [ωpr ] in Fig. 3 of the main text. This
structure is also present in the behaviour of the average phonon number versus probe drive but not in the behaviour of
the average photon number, see Fig. 2. This difference is a higher-order effect of the probe drive-field that would not
be captured in linear-response (i.e. contributions to hb̂† b̂i beyond order 2 ). At higher orders, the probe field causes
non-zero correlations between + and − polaritons, e.g. hĉ†+ ĉ− i =6 0. It follows from the normal-mode transformationof Eq. (S2) that such averages contribute to the photon and phonon populations with opposite sign; this ultimately
explains the strong difference in the behaviour of photon and phonon number versus ωpr .

[1] A. A. Clerk, M. H. Devoret, S. M. Girvin, F. Marquardt, and R. J. Schoelkopf, Rev. Mod. Phys. 82, 1155 (2010).
[2] M. Marthaler and M. Dykman, Phys. Rev. A 73, 042108 (2006).
[3] I. Serban and F. K. Wilhelm, Phys. Rev. Lett. 99, 137001 (2007).
[4] A. Kamenev and A. Levchenko, Advances in Physics 58, 197 (2009).
[5] C. W. Gardiner and P. Zoller, Quantum Noise (Springer, Berlin, 2000).