‘Mechano-optics’: an optomechanical quantum simulator

A widely-known paradigm in optomechanical systems involves coupling the square of the position of a mechanical oscillator to an electromagnetic field. We discuss how, in the so-called resolved sideband regime, this system allows to simulate dynamics similar to ordinary optomechanics, where the position of the oscillator is coupled to the field, but with the roles of the oscillator and the field interchanged. We show that realisation of this system is within reach, and that it opens the door to an otherwise inaccessible parameter regime.


I. INTRODUCTION
The field of optomechanics 1 has in recent years achieved some long sought-after milestones [1,2].Following the first observations of the cooling of a mechanical oscillator using radiation pressure [3] came that of strong coupling [4], and a host of architectures have appeared that incorporated mechanical elements in optical or microwave cavities.Cooling of a mechanical oscillator to its ground state first by cryogenic cooling [5] and then by means of radiation pressure [6] opened the door to experimenting with solid-state mechanical systems in the quantum regime, culminating in the observation of squeezed states of motion [7,8] and mechanical entanglement [9].
One of the long sought-after goals of optomechanics is to demonstrate manifestly quantum-mechanical behavior in the motion of a macroscopic mechanical oscillator.This has given rise to proposals discussing how, for example, one may observe jumps in the occupation number of the mechanical oscillator by monitoring the electromagnetic field leaking out of a cavity [10,11].The standard model within which this is explored is the so-called membrane-in-the-middle system [12], where a reflective membrane is placed at a node or antinode of a cavity field.All such models share a common interaction Hamiltonian that couples the photon number of the light field-â † â, with â being the corresponding annihilation operator-to the square of the position quadrature of the membrane, i.e., x2 .The interaction Hamiltonian therefore reads Ĥquad = g â † â x2 , wiere g is a, typically small, parameter that quantifies the strength of the interaction.This is the model that will form the basis of this paper.It stands in contrast to the more widely-studied linear optomechanics model [1], where the interaction Hamiltonian takes the form Ĥlin = g â † â x.The most significant limitation of this latter interaction is that g is typically very small compared to the other frequency scales of the problem.This requires that one must consider the case where the cavity field has a macroscopic coherent component α (assumed real and positive for simplicity), such that to lowest order Ĥlin ≈ G â + â † x, where G = αg g is an amplified coupling constant.The key drawback of operating under these conditions is that the resulting Ĥlin is, to a very good approximation, quadratic in the operators.As a result, initially Gaussian states (which are ubiquitous in nature and which tend to be quasi-classical) remain Gaussian at all times, making it exceedingly difficult to observe nonclassical behavior.This is the problem that we will tackle in this paper, by turning a quadratically-coupled optomechanical system into a quantum simulator (cf.also Ref. [13]).

II. MODEL
We will consider a system consisting of one mode of the electromagnetic field coupled to one mechanical oscillator (Fig. 1).In the following, the operator â will denote the annihilation operator of the field and ω c its frequency.The operator x ( p) will denote the position (momentum) of the oscillator, the parameter m its mass, and ω m its frequency.The free Hamiltonian of the system can be written Ĥ Introducing the annihilation op-FIG.1.The system we consider in this work.(a) The simulator takes the form of a quadratically-coupled optomechanical system, such as a membrane-in-the-middle setup [12].(b) Suitable driving emulates a linearly-coupled optomechanical system with the roles of the light and mechanics reversed.(c) The full protocol we consider requires squeezing and displacement operations to simulate the "mechano-optical" Hamiltonian.Ĥint symbolises the respective interaction Hamiltonian.
erator for the mechanical field, b, through the relations We assume that the field is driven by a classical source of strength E at a frequency ω s .This is modelled by means of a Hamiltonian Ĥ dr = E * e ı ω s t â + E e −ı ω s t â † .It is convenient to transform to a frame rotating at the frequency of the source, and to neglect constant terms; we will drop the prime symbol to denote Hamiltonians in this rotating frame.This yields Ĥfree = − ∆ â † â + ω m b † b and Ĥdr = E * â + E â † , where ∆ = ω s − ω c , which could be negative or positive, is the detuning between the source and field frequencies.
We also assume that the interaction term between the field and oscillator is quadratic and reads Ĥquad where the coupling strength g is defined as g := g/(2 m ω m ).Next, define Ĥ := Ĥfree + Ĥdr + Ĥquad , which governs the unitary evolution of the system.The cavity field is assumed to couple to the external electromagnetic vacuum at a rate κ.This is conveniently modelled using an open systems formalism whereby the system is best described by means of its density matrix ρ, which obeys the master equation ρ A set of non-unitary terms similar to L[ρ] exists acting on the mechanical motion; however, as a simplifying assumption we will confine ourselves to evolution times significantly smaller than the decoherence time of the mechanical oscillator, which in turn allows us to ignore these terms in the master equation.It is useful to displace the field operators by a complex number α.For this purpose define D(α) := exp α â † − α * â , such that D † (α) â D(α) = â + α.For reasons that will become clearer in the next step, we will introduce a squeezing operator, Ŝ(z) := exp − 1 2 z * b2 − z b †2 , which acts on the mechanical state to yield, e.g., Ŝ † (z) b Ŝ(z) = cosh(r) b + e ıθ sinh(r) b † , where z = r e ıθ .
By means of a suitable choice for the phase reference for E, we can assume that α is real for convenience.We find that where g 0 := g α and Ω m := ω 2 m + 2 g ω m |α| 2 , and where we omitted constant terms.
We now argue that the terms in the second line of the above equation can be ignored.First, we will consider the situation where g is very small, such that a large |α| is used in order to effectively amplify the interaction; this is the situation most frequently encountered in present-day optomechanical systems [1].Terms of order |g| = |g 0 /α| |g 0 | can therefore be ignored safely.We also assume that |α| is, conversely, small enough to avoid any bistable or unstable dynamics (cf.Applications, below).Second, because we will also assume that Ω m κ ∼ Ω c it is possible to invoke the rotating-wave approximation (RWA) [14,15], also commonly used in optomechanics, to drop the second term.We note at this point that applying the rotating-wave approximation directly to Ĥquad neglects the mechanical squeezing produced by the mean field α, and will therefore give incorrect results when |α| is large.
In the next section we will show that the evolution operator derived from ĤDS is approximately equal to one derived from a greatly simplified model corresponding to the usual optomechanical Hamiltonian, but with the optical and mechanical fields interchanged.
From ĤDS to "mechano-optics."-Wewant to compute the time evolution operator ÛDS (t) := exp(−ı ĤDS t/ ) induced by the Hamiltonian (1).Operating under the assumption that g/Ω m 1, and that |α| 1, we can treat the term Ĥsmall := g â † â b † b + 1 2 b2 + b †2 as a small perturbation.Next, we define the "mechano-optical" Hamiltonian and the auxiliary term Ĥaux := 2 g 0 â + â † b2 + b †2 , such that ĤDS = ĤMO + Ĥaux + Ĥsmall .Operating under the assumption that Ĥsmall is a perturbation, it is possible to write a concise expression for ÛDS (t) in terms of the evolution operator ÛMO (t) := exp(−ı ĤMO t/ ); this calculation is detailed elsewhere (see Appendix).Our immediate aim is to quantify how similar the evolution of a state under the action of ÛMO (t) is to that under ÛDS (t).Let us introduce the fidelity , which measures the overlap between the states |ψ DS (t) = ÛDS (t)|ψ 0 and |ψ MO (t) = ÛMO (t)|ψ 0 for some arbitrary initial pure state |ψ 0 .We find, after some calculations (see Appendix), that to second order in the small parameter g 0 /Ω m , where we have introduced allowing g 0 to inherit an explicit time-dependence from α. Equation (3) may be simplified further for times t 1/ max t {|g 0 |}, yielding where ÊNO (t) is a normally-ordered sixth-order polynomial function of â, â † , b, and b † that is of order g 0 /Ω m 2 .The function F uni (t) ≥ 0 is a universal quantity independent of the initial state.For constant g 0 it can be written as  I.The numerical values used to illustrate the feasibility of the system presented here.Two sets of numbers are shown; the first ("Mechanics") refers to an optomechanical system with a macroscopic mechanical oscillator in an optical cavity, whereas the second ("cQED") is based on an electromechanical system that reproduces the quadratic optomechanical Hamiltonian, albeit with coupling strengths orders of magnitude larger than in the optical domain. with Starting from the vacuum state |ψ 0 = |0 we therefore obtain Our result (7) allows us, from now on, to consider situations where F (t) ≈ 1 during the relevant time period.Therefore, to a good approximation, we can claim that the Hamiltonian governing the system is ĤMO , which is identical in form to the usual linear optomechanical interaction Hamiltonian with g 0 playing the role of the single-photon coupling rate andcrucially-the roles of the optical field and mechanical oscillator reversed.This is the central result of this paper and yields what we will refer to as "mechano-optical" dynamics.Explicity, ĤDS ≈ ĤMO .Note also that the three frequencies Ω m , Ω c , and g 0 are all independent free parameters of the model.From this point on, we shall assume that g 0 is constant.Since Ŝ(r) and D(α) are unitary operations it follows that, as illustrated in Fig. 1, At this point we note two things about the effective Hamiltonian ĤMO .First, the sign and magnitude of the single-phonon coupling strength, g 0 , can now be set at will by means of an appropriate choice of the driving strength E. This was noted, but not discussed at length, in the context of a somewhat related model studied in Ref. [20].This freedom of choice allows us to implement experimental protocols that cannot be performed otherwise.For example, it allows us to simulate the non-equilibrium thermodynamics of optomechanical systems undergoing sudden quenches [21] by turning this interaction on or off as required.Second, the effective single-photon coupling strength in this model may be made very large.The system presented here therefore provides a much-needed shortcut towards simulating strong-coupling physics and makes it possible to effectively enter the single-photon strong-coupling regime of optomechanics; something which, despite significant theoretical [22][23][24] and experimental [25][26][27] progress in recent years, has thus far proven to be elusive, although atom-optomechanical systems have started approaching this regime [28].

III. REALISABILITY
The key approximations made in the above derivation were two.First, we assumed that conditions are such that the terms in the second line of Eq. ( 1) can be ignored.Our fidelity calculations provide justification for the correctness of this statement.Second, we assumed that the entire protocol and evolution can be performed in a time much shorter than the mechanical decoherence time, such that no mechanical dissipative terms need to be included in the master equation.Ref. [8] demonstrates a mechanical oscillator with frequency ω m = 2π × 5.8 MHz and mechanical linewidth γ m = 2π × 8 Hz operating in a dilution refrigerator at a base temperature T = 10 mK, the number of phonons in steady-state is This yields a decoherence time of 1/(γ m n p ) ≈ 0.6 ms, which equates to n osc 20 000 oscillation periods.The greatest squeezing generated in Ref. [8] is reported as 4.7 ± 0.9 dB; the mean corresponds to a maximal squeezing parameter |r max | = −r max = − 1 2 ln 10 −4.7/10 ≈ 0.54.A mechanical oscillator closer in form to the one illustrated in Fig. 1 is explored in Ref. [16], which has a motional mass m = 1 ng, ω m = 2π × 140 kHz, γ m = 2π × 1.4 mHz and at T = 500 mK [18] (n p ≈ 74 000) has a decoherence time equating to n osc 1 300 oscillations.Temperatures down to 14 mK (n osc 47 000) are achievable [29], but uncommon, for this kind of system.Using the second derivative of the cavity mode frequency from Ref. [18] and the mechanical parameters from Ref. [16] yields g = 5.2 × 10 −4 /s for an effective quadratic coupling generated by exploiting avoided crossings between cavity resonances.The cavity is assumed to have a finesse of 60 000, consistent with the values reported in Ref. [17] for a similar setup, yielding a half-width at halfmaximum linewidth κ ≈ 2π × 70 kHz ω m .
Current mechanical realisations of quadratically-coupled systems suffer from a very small quadratic coupling strength normalised to the cavity linewidth, g/κ.A means for overcoming this problem was suggested in Ref. [19], using cavity QED ("cQED") techniques to realise an all-electronic analog of a quadratically-coupled optomechanical system.This allows to achieve relative coupling strengths many orders of magnitude larger than would otherwise be possible.Combining this technique with experimentally-achieved numbers  I) explicitly but the situation is very similar for the "Mechanics" system, with appropriately chosen parameters.The right-most figures should be compared with Figs.2(e) and 3(a) in Ref. [30], respectively; note that the optical and mechanical fields have had their roles switched.The parameters have been chosen such that g 0 = 0.5Ω c , and the initial state as product of coherent states |1 a ⊗ |1 b .from Ref. [8], it is possible to envisage a system with an effective mechanical frequency ω m = 2π × 300 MHz.Assuming a finesse of around 18 000 [8] for both resonators yields an optical linewidth κ = 2π × 330 kHz and a mechanical decay rate γ m = 2π × 17 kHz.Finally, using the data in Fig. 12 of Ref. [19], with m = 2, n = 0, 1, and Φ 0 ext /Φ 0 = 0.4, yields g ≈ 19×10 3 /s for this system.Operating at a base temperature of 10 mK yields a thermal population n p ≈ 0.3, which implies that the effective mechanical oscillator can be assumed to be in its ground state, and n osc 50 000, so that any decay and decoherence processes can be safely neglected.These values, summarised in Table I, were used for the examples that will be presented below.

IV. APPLICATIONS
In this section we will outline two applications of our techniques, which have far-reaching consequences.First, we can apply the system we discussed to simulating the dynamics of optomechanical systems in the single-photon strong-coupling regime.For concreteness, we discuss explicitly the "cQED" system in Table I; the "Mechanics" system can be treated similarly but requires stronger squeezing to overcome the decay of the optical field.Setting Ω c = 2g 0 to approach the singlephoton strong-coupling regime and choosing our parameters self-consistently (α 80 752) we can reproduce dynamics reminiscent of the study in Ref. [30]; see Fig. 2 for further details.Of interest are two facts.First, despite the rather small bare coupling coefficient g we have obtained effective strongcoupling dynamics.Second, the mechanical state acquires a strongly non-classical character, as shown by the strongly negative Wigner function.This is different from the standard single-photon strong-coupling regime in optomechanics, where under simple driving it is the optical field that acquires a non-classical character.We neglect decay processes in the figure because the last time-step shown is at 2π/Ω c , which is smaller than the decay time 1/κ.In the example illustrated in Fig. 2 the fidelity between the dynamics generated by ĤMO and ĤDS is ≥ 99.58%.
The second application we propose is aimed at probing the instabilities that naturally arise in this system.It is not difficult to show that the eigenvalues λ n,l of the mechano-optical Hamiltonian ĤMO have the form where n, l ∈ N, and the corresponding eigenstates |λ n,l have the form |λ n,l := D † l g 0 Ω c |n a ⊗ |l b .For all l ≥ l max := Ω m Ω c /g 2 0 , these eigenvalues are negative for some values of n, which introduces instabilities if the state of the oscillator has significant overlap with |l ≥ l max b .Our scheme allows to access these regimes, which are otherwise inaccessible-using the parameters in Ref. [8], for example, the instability is expected to become noticeable at temperatures around 0.6 GK, since l max ∼ 10 12 .In our case, however, one can bring l max to a low value and trigger instabilities in the system since it no longer has a well-defined Gibbs state.We note in passing that for large enough displacements any realisable mechanical potential becomes anharmonic; this instability is therefore not easily observable otherwise.Indeed, adding a non-linearity of the form χ( b † b) 2 to ĤMO , with χ ≥ g 2 0 /Ω c , restores the positivity of all eigenvalues.

V. CONCLUSION
We have introduced a mechanical quantum simulator based on a quadratically-coupled optomechanical system.This system can effectively reproduce the dynamics of a standard optomechanical system, but where the roles of the optical and mechanical fields are switched, and where the single-photon coupling strength can be chosen by driving the system appropriately.As an example, we have shown how to apply our system to simulate the single-photon strong-coupling regime, and have discussed its role in exploring instabilities in quan-tum systems.Our work opens the door to an entirely new use case for optomechanical systems.

VI. ACKNOWLEDGMENT
to zeroth order in 1/α and first order in g 0 /Ω m .The full formula from which this expression was derived reads From this formula we can see that the entire exponential in the time-ordered part is multiplied by g 0 Ω m which is assumed to be significantly smaller than 1.Therefore, this exponential contributes to at least first order.More importantly, from the reasoning below, we argue that the genuine second order corrections, those denoted by ĥ2 , will not contribute.Therefore, we can safely approximate all functions to their zeroth order, which means F ± = F (2)  b ≈ 0 and F b ≈ η.All our results are valid for Ω m t 1.
Appendix C: Error from neglecting correction in ÛDS (t) First, introduce |ψ(t) := Û(t)|ψ(0) and |χ(t) := Û0 |ψ(0) for some arbitrary initial state |ψ(0) .We want to compute an approximate expression for the fidelity F (t) := | χ(t)|ψ(t) | 2 .We consider the rather generic case where where the parameter is such that 1, and where ĥi (i = 1, 2, 3) are generic hermitian operators.Notice, now, that there are two contributions to second order in in the above expression.The first comes from the square of the first-order contribution.The second is a genuine second-order contribution.
For convenience of presentation we rewrite Eq. (C1) as Proceeding, it is straightforward to see that, to second order in , Let us recall, however, that we want to study the particular scenario where |χ(t) = Û0 (t)|ψ(0) .Then, Eq. (C4) simplifies dramatically to which also guarantees formally that the fidelity F (t) satisfies F (t) ≤ 1.
To simplify this equation, notice that . The proof of this statement is rather straightforward.Either it is clear from the unitarity definition of the time evolution and the fact that it must satisfy U(t) U(t) † at all orders, or d dt meaning that these two quantities are identical up to a constant.Since all these operators are equal to zero for t = 0, then this constant must also be equal to zero.This proves the statement.Therefore we have which is non-negative, always smaller than or equal to 1, and manifestly real-as expected and required.
To proceed we define a dimensionless time variable η = Ω m t (and, similarly, η = Ω m t ).Given the form of our time evolution operator Û(t) → ÛDS (t), and that = g 0 /Ω m 1, we have which is correct only for (g 0 /Ω m )η 1.
Let us simplify Eq. (C11).We start by computing Ê1 (η).We have Ê1 When g 0 is constant these can be calculated explicitly to yield , and (C16) The second term in Eq. (C11) is more complicated and requires the evaluation of sixteen terms; its explicit expression is not illuminating.However, we note that, after some calculations, it is possible to see that it has the general form where ÊNO (η) is an operator made of normal-ordered sestic, quartic, and quadratic combinations of the operators â † , â, b † , and b, while F uni (η) is a function of η but independent of the initial state that can be obtained after some algebra.In particular, if we start with the vacuum for both modes â and b, i.e., |ψ(0) = |0 , we obtain the vacuum fidelity F vac (η) that reads Furthermore, it can be shown that This expression is used, in terms of the time variable t, in the main text.

FIG. 2 .
FIG. 2. (Color online.)Simulating the single-photon strong-coupling regime of optomechanics.The first (second) row represents the Wigner function of the optical (mechanical) field; from left to right the images show t = 0, π/Ω c , and 2π/Ω c .Here we show the "cQED" system (cf.TableI) explicitly but the situation is very similar for the "Mechanics" system, with appropriately chosen parameters.The right-most figures should be compared with Figs.2(e) and 3(a) in Ref.[30], respectively; note that the optical and mechanical fields have had their roles switched.The parameters have been chosen such that g 0 = 0.5Ω c , and the initial state as product of coherent states |1 a ⊗ |1 b .