Stroboscopic quantum optomechanics

We consider an optomechanical cavity that is driven stroboscopically by a train of short pulses. By suitably choosing the inter-pulse spacing we show that ground-state cooling and mechanical squeezing can be achieved, even in the presence of mechanical dissipation and for moderate radiation-pressure interaction. We provide a full quantum-mechanical treatment of stroboscopic backaction-evading measurements, for which we give a simple analytic insight, and discuss preparation and verification of squeezed mechanical states. We further consider stroboscopic driving of a pair of non-interacting mechanical resonators coupled to a common cavity field, and show that they can be simultaneously cooled and entangled. Stroboscopic quantum optomechanics extends measurement-based quantum control of mechanical systems beyond the good-cavity limit.


I. INTRODUCTION
Cavity optomechanics has proven extremely successful in controlling nanoscale and microscale mechanical motion at the quantum level [1]. Among the key achievements is the demonstration of ground state cooling [2,3], mechanical squeezing [4][5][6] and mechanical entanglement [7]. Most of these milestones have been obtained in sideband-resolved optomechanical systems operating in the continuous-wave or amplitude-modulated (two-tone) regime, where a notion of stationary regime can be defined, in some suitable rotating frame. Going beyond steady-state operation may be beneficial for several reasons, e.g. it allows to circumvent stability requirements. Sideband-resolved optomechanical systems driven by long pulses have been considered both for controlling mechanical motion [8][9][10][11] and as a model of quantum interface between flying quantum carriers; for instance, entanglement between microwave and mechanical degrees of freedom [12] and quantum state transfer [13] have been demonstrated in this regime.
Pulsed protocols can also lift the stringent requirement of sideband resolution. By employing pulses much shorter than the mechanical period, quantum state preparation and readout, e.g. of low-entropy and squeezed mechanical states [14][15][16][17][18][19], as well as opto-mechanical and all-mechanical entanglement [20] can in principle be achieved. However, in order to neglect non-unitary processes, coherent operations are restricted to very short times (also less than a single mechanical cycle). The conditional preparation of quantum states with few pulses also requires large interaction strengths, which has so far prevented pulsed optomechanics to attain the quantum regime [21]. Only very recently, pulsed operation close to the quantum regime has been demonstrated in a setup based on a photonic crystal nanobeam [22].
In this work we take a different approach and study the conditional dynamics of an optomechanical system driven by a train of pulses. We show that this new regime-stroboscopic quantum optomechanics-is effective to prepare and verify quantum states of mechanical motion beyond the sideband-resolved regime. In particular, by suitably choosing the spacing between the pulses, ground state cooling and squeezing of a single mechanical resonator can be achieved, as well as collective cooling and entanglement of two non-degenerate resonators (radiation-pressure coupled to a common cavity mode).
Compared to single-pulse protocols, our approach has the distinct advantage to allow for a cumulative effect of the measurements over many mechanical cycles, thus relaxing the requirement on the optomechanical coupling strength considerably. This however requires including mechanical dissipation in the description of the dynamics, as opposed to Refs. [15][16][17][18][19]. Due to the competition between radiation-pressure interaction and mechanical dissipation, the mechanical system eventually settles into a steady state, albeit a periodic one. For such a stroboscopic steady state, we provide simple analytic expressions for the conditional state. From this point of view, our work draws an interesting connection between the conditional dynamics of periodically measured systems and the Floquet theory of optomechanics [23,24].
Our study is inspired by early works in backaction-evading (BAE) measurements, where the stroboscopic dynamics of mechanical transducers was studied for the detection of weak classical signals [25]. We provide a full quantum-mechanical treatment of stroboscopic BAE measurements [26,27], which was so far missing. We discuss in details corrections to the ideal measurement regime stemming from thermal decoherence and the finite length of each pulse. Notably, we show that including the latter effect, usually considered detrimental, enables preparing pure mechanical squeezed states and opto-mechanical entanglement. In short, we show that stroboscopic quantum optomechanics bypasses the need for strong measurements and provides an effective and versatile tool for measurement-based quantum control of mechanical states.
The rest of the paper is organized as follows: in Sec. II we describe the system and derive an effective model of the dynamics based on stroboscopic measurements. The predictions of this model for stroboscopic squeezing and cooling of mechanical motion are presented in Sec. III and Sec. IV, respectively. In Sec. V we discuss engineering squeezed quantum states in connection with stroboscopic BAE measurements of mechanical motion. In Sec. VI we implement verification of the conditional state via retrodiction. In Sec. VII we extend stroboscopic quantum optomechanics to non-interacting mechanical resonators coupled to a common cavity field, and show that they can be simultaneously cooled and entangled. In Sec. VIII we discuss some experimentally relevant considerations for implementing our ideas. Finally, Sec. IX collects conclusive remarks and provides an outlook.

II. A SIMPLE MODEL OF STROBOSCOPIC CONDITIONAL DYNAMICS
We consider a standard optomechanical system where the positionx of a mechanical oscillator of frequency ω m modulates the frequency of a cavity modeâ of linewidth κ [1]. The cavity is illuminated with a train of short coherent pulses of length τ much smaller than the mechanical period, i.e. ω m τ 1. The number of photons N p in each pulse is large enough to make linearization of the optomechanical interaction an excellent approximation, such that we have ( = 1) where g is the pulsed coupling constant,X c = (â +â † )/ √ 2 and we expressed the mechanical position in terms of the slowly varying quadraturesX m ,P m . During the interaction time the coupling induces the unitary evolution For a very short pulse the harmonic motion can be neglected and the unitary evolution can be approximated aŝ where χ quantifies the strength of the interaction. An estimate of the latter for a fast cavity in the adiabatic limit κ τ −1 yields χ = 2g 0 Npτ κ , where g 0 is the single-photon optomechanical coupling. As we show in Appendix A, in this limit the interaction is de facto instantaneous, in that mixing between the two mechanical quadratures is fully neglected. Equation (3) realizes a quantum non-demolition (QND) gate between the optical and mechanical amplitudes [28]:X c ,X m are left untouched and information about them is acquired by the conjugate quadraturesÛ †P m(c)Û =P m(c) + χX c(m) . The probability of recording a value P c of the optical phase quadrature after such interaction is given by where the cavity starts off in the vacuum and the resonator in an arbitrary stateˆ m . In the second line of Eq. (4) we have rewritten the probability by introducing the family of Kraus operatorsΥ(P c ) = P c |Û |0 , elements of the positive operator-valued measure (POVM) {Υ †Υ } Pc , that satisfy Υ †Υ ≥ 0, dP cΥ †Υ = 1 m . An explicit expression forΥ is given byΥ which shows that the effect of the pulse on the mechanics is akin to a generalized position measurement. This expression has been first used to model momentum diffusion in continuous weak measurements [29]. Later, it was employed to model an optomechanical system driven by a single strong pulse, i.e., the regime of pulsed quantum optomechanics [14]. From Eq. (5) we notice that, when acting on a pure state, the measurement operator multiplies the wave function by a Gaussian function of width χ −2 and centered around the position P c /χ; by increasing the interaction strength χ, the wave function thus gets increasingly localized in position. Upon recording the outcome P c , the mechanical density matrix is transformed asˆ which is the conditional, or post-measurement state. When the pulse is off, environment-induced decoherence is affecting the otherwise free evolution of the mechanical resonator. The dynamics is governed by wheren and γ are the mean occupation and damping rate of the mechanical bath andb is the annihilation operator associated to mechanical quadratures. The evolution over a finite amount of time is given by the mapΦ th = e Lt . The pulsed interaction (6) and the free-evolution-plus-dissipation (7) form the 'unit cell' of our stroboscopic model, which can be thought of as a repetition of these two elementary steps, see Fig. 1(c). As we discuss below, when concatenating many such steps one is free to choose the spacing between two subsequent pulses. Over this time (i) the mechanical mode picks up a phase, which determines which quadrature is measured at the next interaction, and (ii) the mechanics exchanges phonons with the thermal environment. In particular, the presence of the latter contribution-neglected in previous studies [15][16][17][18][19]-competes with the measurement, eventually leading to a non-equilibrium steady state.

A. From measurement-induced evolution to deterministic CP maps
A great simplification comes from assuming that both the measurement and the dissipation act on a Gaussian state, in which case their output is a Gaussian state too [28,30,31]. For the case of a Gaussian measurement, such as the quadrature measurement in Eq. (4), the post-measurement state (6) depends on the measurement outcome only through the first moments or, equivalently, the measurement-induced evolution of the second moments is deterministic; this is a general feature of Gaussian measurements [32]. Therefore, the effect of the measurement can be cast in the form of a deterministic map E Υ for the second statistical moments. The action of this map σ = E Υ (σ) on the mechanical covariance matrix σ (with variance σ Xm , σ Pm and covariance σ XmPm ) is given by We explicitly see that the stochastic component of the measurement (P c ) is absent from the above expressions. The first and second expression describe the reduction of the variance alongX m , and the increased fluctuations of the conjugate quadrature due to the quantum backaction, respectively.
The (commutative) action of dissipation and free evolution (7) on the covariance matrix is described by the map E th,φ (σ) = e −γt R φ σR φ T + (1 − e −γt )σ th , where σ th = (n + 1 2 )1 2 is the covariance matrix of a thermal state and R φ is the rotation matrix due to harmonic evolution. Equivalently, under E th,φ , the input state gets rotated and mixed with a thermal state via a beam splitter of effective transmissivity η = e −γt . Measurement and dissipation compete over time. The former tries to reduce the uncertainty in one quadrature (at the expense of the other), while the latter tries to restore isotropy. Crucially, the spacing between two pulses determines the amount of mixing between the quadratures from one measurement to the next one. This consideration applies to any sequence of equally spaced pulses; for example, one can obtain a recursion relation σ (N ) = (E Υ • E th,φ )σ (N −1) to model a short train of pulses. This operation regime has recently become experimentally relevant for quantum applications [22]. Here we focus on a different regime: when the action of the measurement is undone by the dissipation there is no net effect over a 'unit cell' and the system reaches a stroboscopic steady state [see Fig. 1 We stress that the two operations do not commute, so that in general E Υ • E th,φ = E th,φ • E Υ , as we shall see below. This is a novel regime for cavity optomechanics, which has focused either on steady state properties of continuously driven systems or in the finite-time dynamics, as in pulsed optomechanics.

III. STROBOSCOPIC SQUEEZING OF MECHANICAL MOTION
A. Stroboscopic BAE measurement The first case we consider is that of a stroboscopic BAE measurement, for which a classical treatment is discussed in Refs. [25][26][27]. By choosing pulses interspaced by a multiple of half the mechanical period, we can in principle realize a QND measurement of position. Indeed, one has [x(t),x(t + T )] = i sin(ω m T ), so that at the stroboscopic times T = kπ/ω m a sequence of precise position measurement is possible with no fundamental limit imposed by quantum mechanics (in the following we always take the shortest interval k = 1). However, due to the presence of the environment, the covariance matrix does not come back to itself half a period later. We then look for solutions where the combined action of the measurement and the environment leaves the state invariant. Solving for the stroboscopic steady state σ ss = (E Υ • E th,π )σ ss , we get σ Xm = 2n + 1 and σ XmPm ≡ 0, where we set z = (2n + 1)χ 2 . We stress that the knowledge of such state is conditioned on the stream of measurement results. These expressions can be considerably simplified for large values of the mechanical quality factor Q = ω m /γ. The leading terms in the expansion are given by This simple result provides a quantum-mechanical treatment of stroboscopic BAE measurement and proves that mechanical decoherence does not preclude the occurrence of squeez-ing at long times. Indeed, uncertainty may fall below the zeropoint value, which implies a squeezed state of the resonator. Mechanical squeezing [expressed in −10 log 10 (2σ Xm ) Decibel (dB)] is plotted Fig. 2(a). We stress that different Q entail different characteristic times to approach the stroboscopic steady state. The solid lines are for the steady state relative to E th,π • E Υ , while the dashed for E Υ • E th,π . Physically, they correspond to the knowledge of the conditional state directly after or directly before the measurement. We can see discrepancies arising due to the non-commutative character of the two maps for low Q and large coupling values. In this parameter regime, if we start from a thermal state, the effects of measuring first are (partially) undone by the subsequent application of the thermal channel. On the other hand, by reversing the order (i.e., considering the map E Υ • E th,π ) E th,π acts as the identity, so the first measurement retains more conditioning power. The difference between the two cases thus boils down to an extra pulse, which has significative impact for large χ and explains the larger amount of squeezing. However, already for moderately large quality factors the two predictions coincide.
It is interesting to compare the condition for mechanical squeezing enforced by Eq. (13) with that required by pulsed optomechanics, i.e. by applying a single pulse (5). For a single pulse, values of the coupling χ > 1 are required to obtain squeezing (independently ofn), which has so far precluded reaching the quantum regime in pulsed optomechanics experiments. On the other hand, with stroboscopic driving approaching σ Xm < 1/2 only requires χ > 2π(n + 1/2)/Q, which can be considerably less demanding for large quality factors.
Finally, if we compare how the two variances in (13) scale with Q, it is clear that fluctuations increase faster inP m than they are reduced alongX m . This means that while getting squeezed, the resonator gets also heated up. This fact is highlighted in Fig. 2(b) where the mechanical purity µ = Tr ˆ 2 ss for the same cases of panel (a) is shown. In the large Q limit the purity takes the simple form Larger values of squeezing are accompanied by low purity. We will see in Sec. V that this conclusion gets drastically modified by considering the imperfect QND regime determined by mechanical evolution during the pulse.

B. Squeezed input pulses
Finally, we notice that the former results can be extended to the case where squeezed pulses, rather than coherent ones, are fed to the optomechanical cavity. In our simple model this observation amounts to replace the cavity vacuum seed state σ pulse = 1/2 [see Eq. (4)] with a squeezed state σ pulse = diag e r 2 , e −r 2 squeezed along the phase quadrature. This determines reduced fluctuations of the (measured) optical phase, which in turn enhances the conditioning effect of the measurement. One obtains results as in Eqs. (11), (12) with the substitution χ → e r 2 χ, namely a train of squeezed pulses magnifies the measurement strength by an exponential factor (in the degree of squeezing).

IV. STROBOSCOPIC GROUND STATE COOLING
Another interesting case is obtained by spacing the pulses by a quarter of a period. In this case the value of the variance alongX m andP m gets swapped by the free evolution, so that the measurement reduces both variances alternately. An exact expression for the stroboscopic steady state is available also in this case, although quite cumbersome. For convenience below we give the expansion for large Q The full expression of the functions F(χ,n), G(χ,n) is reported in Appendix B. Unlike Eq. (13), now fluctuations in both quadratures converge to a constant value for Q → ∞.
We also notice that there is a residual asymmetry between the two quadratures. The leading terms therefore describe a squeezed vacuum state, albeit one where the squeezing grows slowly with the coupling χ. For realistic values of the coupling the state is thus only weakly squeezed and has near-unit fidelity with the mechanical vacuum. The conditional purification of the mechanical state is also known as cooling-bymeasurement [21]. In the same spirit, we refer to this case as stroboscopic cooling. Of course for finite values of the quality factor the steady state will be mixed, but cooling close to the ground state is still possible. We show these features in Fig. 3.

V. IMPROVED DESCRIPTION AND NUMERICAL SIMULATIONS
In this Section we aim to provide a more accurate description of the stroboscopic conditional dynamics. We will focus on the case of stroboscopic BAE measurements but the analysis can be readily extended to the case of stroboscopic cooling. We expand along two directions: (i) we model the measurement as actually taking place outside the optical cavity and (ii) we evaluate the effects of the mechanical free evolution during the pulsed interaction. To this end, we consider the following extended Hamiltonian where, beside the term in Eq. (1), we also include an interaction with the continuum of electromagnetic modesâ in,t living outside the cavity. This stream of modes interacts with the system at time t and is otherwise uncorrelated [â in,t ,â † in,t ] = δ(t − t ). As customary, we assume they have Markovian cor- . For a short pulse of length τ (for now neglecting the free mechanical evolution) the corresponding propagator takes the formÛ whereX in ,P in are the proper (dimensionless) modes of the environment, i.e., [X in ,P in ] = i, which are being measured; homodyne detection of the phase quadrature corresponds to projection along |P in (see Appendix C for details). Formally, we can then proceed as in Sec. II to compute the conditional covariance matrix of the optomechanical system, include thermal decoherence, and enforce the stroboscopic steady-state condition. The full expression of the conditional state of the mechanical system is quite cumbersome, but in the large Q limit we get the following simple expressions with g = 2g 0 Np κτ . These expressions are to be seen as a refinement of Eq. (13); as we will show, they offer a useful comparison with numerical simulations.
Second, we include corrections to the ideal QND limit stemming from the mixing of the mechanical quadratures during a pulse of finite length. The ensuing unitary evolution contains two new terms (see Appendix A for the full expression): a squeezing term in the cavity amplitude, which however is O(g 2 /ω 2 m ), and a spurious term ∝X cPm , which spoils the QND nature of the interaction. The strength of this term is 2ω m /κ times the QND part, so that the QND limit is approximately recovered only for optomechanical systems deep in the bad-cavity regime. It is therefore important to address the corrections arising for finite values of the sideband parameter, which limit the amount of conditional squeezing attainable. Due to the presence of quadrature mixing, a closed expression of the conditional state can no longer be found. However, we can get a clear physical picture of the effects brought about by non-QND term in the following way. Consider the effective Hamiltonian generating the optomechanical evolution, first neglecting and then including the non-QND term (to faithfully model the measurement, we also include the interaction with the extra-cavity modes); the corresponding expressions are given by Eq. (1) and Eq. (A6), respectively.
The equations of motion for the two cases are schematized in Fig. 4(a) and (b), respectively, where an arrow connecting two terms means that the variable at the starting point drives the evolution of that at the ending point. Next, we incorporate the role of the measurement, which has a twofold effect: on the one hand, it enables to acquire information, i.e., to reduce the uncertainty about the mechanical quadratureX m ; this acquisition happens indirectly through the optmechanical coupling and requires that we keep track of the stochastic component. On the other hand, the measurement introduces disturbance, which directly affects the conjugate quadrature (X in ) and then, through the dynamics, reaches the mechanical system.
In the ideal QND case [see Fig. 4(a)] these two effects fully decouple. Fluctuations are reduced alongX m and increased inP m (backaction heating). Graphically, this corresponds to the fact that no arrow points towardX m , and hence no noise can drive it. Likewise, no arrow originates fromP m , which 'absorbs' all the backaction. Thus backaction confinement enables repeated measurements of the same quadrature with no added noise, which is the working principle of BAE measurements. When we take into account the finite mechanical evolution [cf. Fig. 4(b)], the non-QND terms open new paths (dashed arrows) for both backaction and conditioning to spread, with the following consequences: information is now acquired about both mechanical quadratures (and hence fluctuations of the conditional state are reduced in both directions) which entails that (i) the measurement purifies the state. Similarly, measurement backaction is no longer confined toP m but extends to both quadratures, i.e., (ii) the amount of squeezing is reduced with respect to the ideal case. Finally, information is simultaneously acquired about both the cavity and the mechanics (see multiple arrows incoming atP c ); such joint reduction of the uncertainty implies that (iii) correlations between cavity and mechanics are built. Depending on the occupancy of the mechanical bath, this may even lead to entanglement being established between the two resonators. We want to remark that, while the limitation (ii) posed by non-QND terms is known, their beneficial effects (i) and (iii) have not been previously appreciated. A similar situation is encountered in continuous BAE measurements, where RWA solution yields conditional squeezing with low purity, and the inclusion of counter-rotating terms lower the amount of squeezing but at the same time allows for larger purity and optomechanical entanglement [33].
To check the validity of these conclusions we numerically integrate the conditional evolution of the full optomechanical system subject to stroboscopic driving and continuous homodyne detection of the output phase quadrature (see Refs. [31,33] for details). Free mechanical evolution during each pulse is explicitly included in the simulation, i.e. we use the optomechanical interaction in Eq. (1). In Fig. 5(a) we show the numerical squeezing in the long-time limit (averaged over one period) and compare it with the prediction of Eq. (19). Our simple analytical formula shows excellent agreement except for large Q, where it does not capture the saturation of squeezing. Such saturation confirms our expec- tation (ii). Indeed, for a fixed duration of the pulse, the effects of the free mechanical evolution become more prominent for larger Q. From Fig. 5(b) we see that a realistic stroboscopic BAE measurement actually generates highly pure conditional squeezed states. Mixing of the two quadratures, present for any finite value of the sideband parameter, implies the simultaneous squeezing and cooling the mechanics by stroboscopic BAE measurement, as predicted in (i). Finally, in panel (c) the conditional optomechanical entanglement is displayed (quantified by the logarithmic negativity) which confirms (iii). Notice that entanglement is present in the high temperature regime.

VI. VERIFICATION OF THE MECHANICAL STATE: STROBOSCOPIC TOMOGRAPHY
Essential to any conditional protocol is a verification part. While in the previous sections we have focussed on state preparation through measurement, here we calculate how well a quadrature can be measured in a train of pulses. This is also known as retrodiction [34,35]. State verification via retrodiction has been recently employed to verify the quantum trajectory of a continuously driven optomechanical system [36]. The final result of a stroboscopic measurement is a measurement value with a given confidence interval. Repeatedly preparing and measuring a state allows for full tomography. Since it makes sense to keep measuring until the resonator is no longer correlated with its initial state, the resonator state at the end of the measurement is again a conditionally squeezed state as discussed above, independent of the initial state.
In a stroboscopic measurement of a harmonic oscillator, its position is measured at regular intervals. This results in a string of measurement results y = (y 0 , y 1 , y 2 , . . .), which are correlated with the actual position at that time y i = x i + m i . The measurement errors {m i } are normally distributed as the Kraus operator corresponding to the measurement (5) predicts. Specifically, if the oscillator is in a position eigenstate ρ = |x x|, the measurement probability distribution is . Thus, the {m i } are drawn from a Gaussian distribution of zero mean, with variance σ 2 m = 1/(2χ 2 ). In order to realize a QND measurement, the time t n between measurements has to be an integer multiple of half a period, t n = (n + 1)π/ω m . In this regime, the problem becomes classical, as the measurement backaction is evaded. Here we choose the time between measurements to be as short as possible, T = π/ω m , such that γT = π/Q. In between each measurement, the position of the oscillator decays, due to damping, and gets a random contribution from the thermal noise acting, x i = e −γT /2 x i−1 +d i . The random numbers d i also follow a normal distribution of mean zero and standard deviation σ 2 d = (n + 1/2)(1 − e −γT ) which follows from the fluctuation-dissipation theorem (or equivalently from the explicit discussion in Sec. II A).
In Appendix D we show that given a string of measurement results y, the conditional probability distribution inferred from Bayes' theorem is the normal distribution where the correlation matrix Q xx is given by Interestingly, as the number of measurements goes to infinity, the matrix Q xx can be inverted analytically (see Appendix D 2). The first element of the inverse, [Q −1 xx ] 11 , is the variance associated with the measured value, and as one would expect it coincides with the achieved squeezing [Eq. (13)]. While perhaps this could have been inferred from the results above, a very good approximation to the inverse and thus the variance can also be found for a finite number of measurements. Furthermore, this approach yields the weight each measurement value is associated with, although the formula is somewhat unenlightening (see Appendix D). Finally, we can also show that when taking into account many measurements before and after a certain point in time, i.e., using preparation and retrodiction, the associated variance is half of Eq. (13).

VII. COLLECTIVE ENTANGLEMENT AND COOLING OF TWO MECHANICAL RESONATORS
In this Section, we show how the previous results can be extended to the case of two non-degenerate mechanical resonators. In order to do that, we consider two mechanical resonators of frequency ω m,1 and ω m,2 coupled to a common cavity field. Collective BAE schemes have been proposed in this configuration for continuous and two-tone driving [37][38][39]. We introduce the mean and the relative mechanical frequency, respectively defined as ω = (ω m,1 + ω m,2 )/2 and Ω = (ω m,1 − ω m,2 )/2 (we assume ω m,1 > ω m,2 without loss of generality). We also define the collective mechanical variableŝ that satisfy [X ± ,P ± ] = i, [X ± ,P ∓ ] = 0. When the pulse is on, both mechanical resonators linearly couple to the common cavity amplitude, givinĝ j=1,2X c X m,j cos(ω m,j t) +P m,j sin(ω m,j t) For simplicity, we have considered the case of equal singlephoton optomechanical couplings and in the second line we have rewritten the interaction in terms of the rotated collective quadraturesX which still form a conjugate pair X ,Ŷ = i. Thanks to this change of variables we see that Eq. (25) has the same form as (1) and therefore we can rely on our previous analysis. In particular, for g(t) = g δ(t − kπ/ω) we recover the ideal case of stroboscopic QND interactionÛ ≈ e i √ 2χXmX (here k = 1). This corresponds to pulsing every half of the fundamental period 2T 1 T 2 /(T 1 + T 2 ), where T j are the single mechanical periods. Like in the single-mode case, we also include mechanical dissipation. For simplicity, in the following we consider equal mechanical damping rates and same occupancies for the two baths. For non-degenerate mechanical modes, these conditions may entail adjusting the local temperatures of the baths to achieve the same occupancy.
From the discussion of Sec. III A we conclude that the stroboscopic steady state is a squeezed thermal state in the collective variablesX andŶ , with the variance reduced alongX and heated up alongŶ by the backaction. We now want to express the state in terms of the original local variables. Fo this purpose, it is useful to parametrize the steady state σ ss of single mode BAE measurements [cf. Eqs. (11), (12)] as where an explicit expression of n eff , r eff can be obtained by inverting Eqs. (11) and (12) (one also needs to rescale g → g/ √ 2). Next we notice that the vector of original quadra-turesQ = (X m,1 ,P m,1 ,X m,2 ,P m,2 ) T and that of collective onesQ = (X,Ŷ ,Ŵ ,Ẑ) T are related via the following trans- 2b 2)Û BS ; hereẐ =X − cos Ωt +P + sin Ωt andŴ =P − cos Ωt − X + sin Ωt are the other two collective rotated quadratures and U BS is a beam splitter transformation. Note that the modes are rotated in opposite directions before getting mixed. By transforming state (27) accordingly we get where A = (n eff + 1 2 ) cosh 2r eff 1 2 and C = −(n eff + 1 2 ) sinh 2r eff σ z (with σ z the z-Pauli matrix), namely a twomode squeezed thermal state. This state is known to be entangled if and only if r eff > ln( √ 1 + 2n eff ). Therefore frequent measurements modeled by the pulses may induce entanglement between the two non-interacting resonators.
The same argument can be repeated for two-mode cooling, which corresponds to pulsing every quarter of the fundamental period 2T 1 T 2 /(T 1 +T 2 ). Indeed, Eqs. (15) and (16) are also in the form of a squeezed thermal state, although with very little squeezing. Therefore the two-mode stroboscopic steady state is still of the form (28), the difference being that now we have n eff ≈ 0 and smaller r eff (compared with the previous case). For both n eff , r eff ≈ 0 the steady state has a large overlap with the vacuum of the two modes. Notice however that for n eff = 0 the state is entangled for any value of r eff > 0.

VIII. EXPERIMENTAL CONSIDERATIONS
The principal considerations for implementing stroboscopic optomechanics concern the appropriate hierarchy of time scales κ −1 τ ω −1 m and sufficient measurement strength χ. Low bath occupationn and high quality factors Q facilitate access to the quantum regime. However, experimental non-idealities have to be taken into account as well. These can include sub-unity detection efficiency, optical absorption heating, mechanical frequency drift, and spurious mechanical modes, among others.
Given the recent progress with measurement-based quantum state preparation [36,40] with membrane-in-the-middle optomechanical systems, we discuss this platform first.
With MHz resonance frequencies, sufficiently short (submicrosecond) pulses are readily implemented using standard modulation techniques. Such pulses could be accommodated in short (L ∼ 1 mm), medium-finesse resonators with κ/2π 15 MHz, and the cavity output detected with high efficiency, as previously demonstrated [36,40]. Optical power levels tolerated in continuous-wave operation [40] suggest χ = gτ ∼ 0.1 can be achieved without significant device heating at a temperature of ∼ 10 K, orn ∼ 10 5 . Since stroboscopic operation would lower the thermal load from optical absorption by a factor of order ω m τ < 1, even higher χ may be possible, provided instabilities are avoided and the cavity lock maintained. If the experiment were implemented with soft-clamped membrane resonators [36,40,41], very high quality factors Q ∼ 10 9 are available (for comparison with some of the results presented so far, gentle laser pre-cooling can be assumed to trade equivalent bath occupancy with quality factor, leaving the productnQ constant). With this set of parameters, significant levels of squeezing can be achieved, see Fig. 5.
In practice, however, other mechanical modes at harmonics of the stroboscopic sampling frequency contribute to the measured signal, in an effect known as aliasing in the context of periodically sampled data. This would lead to spurious noise and interactions, and a degradation of the prepared state. In contrast to pulsed optomechanics [14], which offers virtually no spectral discrimination of mechanical modes at all, the stroboscopic protocols are only sensitive to spurious modes that coincide with a harmonic of the sampling frequency. Yet soft-clampled membrane resonators with their high density of states outside the bandgap, in which the high-Q modes lie, would be strongly affected by this effect. Membrane-in-themiddle setups with mechanical systems that feature a sparser mode spectrum, such as trampoline resonators [42,43], may therefore be preferable, provided sufficiently high Q-factors and/or low temperatures can be achieved.
Nanophotonic structures could be an alternative platform of interest. Their already sparse mechanical mode spectrum could conceivably be engineered to be sufficiently anharmonic. Measurement strengths as high as χ ≈ 0.1 have already been demonstrated with single optical pulses much shorter (τ = 20 ns) than the mechanical period (2π/ω m ≈ 300 ns) [22]. In combination with high-efficiency readout and efficient heat removal [44], stroboscopic optomechanics might also allow generation of squeezed and entangled states in such systems.

IX. CONCLUSIONS
In this work we provided a description of the conditional dynamics of an optomechanical system driven by a train of pulses. We showed that the resulting framework-dubbed stroboscopic quantum optomechanics-provides a versatile toolbox for measurement-based quantum control of optomechanical systems in the bad-cavity regime, ranging from ground state cooling to mechanical squeezing, and applicable to single as well as multimode optomechanics. Crucially, it enables the generation and characterization of measurementbased squeezing and entanglement. Stroboscopic driving alleviates the requirements of pulsed protocols based on a single (or a few) pulse(s). In the following we derive the complete expression of the unitary evolution Eq. (2) induced by the pulsed interaction and discuss when the simple QND expression Eq. (3) is recovered. At a formal level, the evolutionÛ (t, t 0 ) can be equivalently described by the Magnus expansion which comes in the form of a non-ordered exponential. Compared to Eq. (2), the complexity of the expression has just been shifted to the argument of the exponential. The first two terms of the Magnus expansion are given bŷ whereĤ I (t) is the linearized interaction in the rotating frame, as given in Eq. (1). For concreteness, let us consider a rectangular pulse of length τ centered at the origin, described by the normalized profile ε(t) When the pulse drives the optomechanical cavity on resonance, the evolution of the field amplitude inside the cavity is well approximated byα = − κ 2 α + κN p ε(t), where we neglected the mechanical response during the short interaction time. The solution reads where we set f ± (t) = 1−e − κ 2 (t± τ 2 ) . The expression captures the build-up and the decay of the coherent field inside the cavity. The time-dependent optomechanical coupling in Eq. (1) is g(t) = g 0 α(t). Notice that in the fast cavity limit the expression reduced to a rectangular pulse of height g ad = 2g 0 Np κτ . We can now compute the expressions (A2a), (A2b) for this profile. For late times t τ we get where the prefactors read In the above expression we set g ad = 2g 0 Np κτ . In the adiabatic limit χ reduces to χ ad = g ad τ ≡ 2g 0 Npτ κ , which is the expression given in the main text. Also notice that, as expected for the problem at hand, all the nested commutators corresponding to higher orderΩ k≥3 identically vanish. The evolution thus takes the following exact expression Two extra terms have appeared compared to Eq. (3). A singlemode operator that is responsible for squeezing of the cavity field and an interaction term that spoils the QND character of the quadratureX m . Note that both terms are present for any length of the pulse, even though the spurious term is suppressed by a factor 2ω m /κ. The QND limit can then be recovered only for a vanishing sideband parameter ω m /κ.
We would like to model a series of BAE measurements on a damped harmonic oscillator in a thermal environment. In order to do so, we need to find the probability distribution of x given some measurement record y, P ( x| y), assuming that the parameters σ m , σ d , γ are known.
From the description in the main text, we can determine the conditional probability distribution for the measurement outcomes as well as our prior (D2) The initial variance for x 0 could be from a thermal state. Optionally, we could already have performed measurements at that point, in which case the initial variance and mean are given by the resulting state. We can thus write down the joint probability distribution P ( x, y) = P ( y| x)P ( x).
For now, we are only interested in the mean and variance of the first entry. We thus need to determine and Note that Σ xx commutes with the matrix in round brackets. We calculate the inverse Q −1 xx in Appendix D 2.

Explicit calculation of the inverse
Consider the matrix 1 a 0 · · · 0 0 0 a b a · · · 0 0 0 0 a b · · · 0 0 0 . . . . . . . . . . . . . . . . . . . . . 0 0 0 · · · b a 0 0 0 0 · · · a b a 0 0 0 · · · 0 a 1 The entries of the inverted matrix are [45] where θ and φ fulfil certain recurrence relations. In our case, they are in fact the same, and we have The solution can be found from the characteristic polynomial x 2 − bx + a 2 = 0, which has roots Fitting the general solution to the boundary conditions, we find but with θ n determined by Eq. (D12) above. This allows us to write the inverse of the above matrix in an exact analytical, albeit unwieldy manner This exact formula for the inverse represents the central result of this appendix. We consider simplifications that arise in certain cases below.

a. Physically relevant case
The relevant case is b 2 /4 > a 2 . In this case ξ + > ξ − , such that in the limit of a large number of measurements n → ∞, the formula for the inverse turns into (D16)