Dissipative generation of significant amount of mechanical entanglement in a coupled optomechanical system

We propose an approach for generating steady-state mechanical entanglement in a coupled optomechanical system. By applying four-tone driving lasers with weighted amplitudes and specific frequencies, we obtain an effective Hamiltonian that couples the delocalized Bogoliubov modes of the two mechanical oscillators to the cavity modes via beam-splitter-like interactions. When the mechanical decay rate is small, the Bogoliubov modes can be effectively cooled by the dissipative dynamics of the cavity modes, generating steady-state entanglement of the mechanical modes. The mechanical entanglement obtained in the stationary regime is strongly dependent on the values of the ratio of the effective optomechanical coupling strengths. Numerical simulation with the full linearized Hamiltonian shows that significant amount of mechanical entanglement can indeed be obtained by balancing the opposing effects of varying the ratio and by carefully avoiding the system parameters that may lead to amplified oscillations of the mechanical mean values detrimental to the entanglement generation.

greatly enhancing the entanglement is to drive the coupled cavity modes with four-tone lasers of weighted amplitudes and specific frequencies so that we obtain an effective system Hamiltonian where two nonlocal Bogoliubov modes of the mechanical oscillators are coupled to the cavity modes via beam-splitter-like interactions. Notably, the Bogoliubov modes can be sufficiently cooled via swapping quanta with the cavity modes which interact with optical thermal baths with neglectful mean photons. In this way, after some time of dissipative dynamics, the mechanical modes are driven to close to a two-mode squeezed state which is, in fact, the joint vacuum of two Bogoliubov modes being cooled. The amount of entanglement is independent of initial states but is strongly dependent on the ratio of the effective optomechanical couplings rather than simply on their magnitudes. The change in the ratio will simultaneously have two confronting effects on creating mechanical entanglement. The entanglement can be maximized by balancing the confronting effects through choosing proper asymmetric driving amplitudes and as small mechanical decay rate as possible. Our numerical results with the full linearized Hamiltonian show that significant amount of mechanical entanglement can indeed be generated. In particular, we observe obvious amplified self-sustained oscillations of mechanical mean values for some system parameters, which may result from the intrinsic nonlinearity of the optomechanial interaction. The effects of these amplified oscillations in some parameter regimes, often unwanted in the generation of mechanical entanglement and largely unconsidered in many optomechanical schemes, have been numerically analyzed and been carefully avoided in our scheme.

Model and quantum Langevin equations
As illustrated in Fig. 1, we consider a coupled microtoroidal optomechanical system [39][40][41][42][43] where two phonon modes B 1 and B 2 respectively interact with two photon modes A 1 and A 2 which in turn are coupled via the photon tunneling. The arrangement is assumed to be symmetrical, i.e., the two mechanical oscillators have the same frequency ω m and damping rate γ m . The cavity modes with frequency ω c are driven by lasers with the frequency ω L and time-modulated amplitude E(t) through the tapered fibers. In the rotating frame with respect to laser frequency ω L , the Hamiltonian of our system reads (ħ = 1) where A j ( † A j ) and B j ( † B j ) are the annihilation (creation) operators of the jth photon mode and phonon mode, respectively. Δ 0 = ω c − ω L denotes the detuning between the cavity and the driving field. The parameters g and J represent the strengths of the single-photon optomechanical interaction and photon tunneling, respectively.
The dynamics of our system can be described by a set of quantum Langevin equations (QLEs) 44 : Here, κ is the cavity decay rate; a t ( ) j in stand for independent input vacuum noise operators with zero mean value and the following nonzero auto-correlation functions: where n b is the mean thermal occupancy of the mechanical baths.
In the presence of strong external pumping, we can write the system operators as A j = α j (t) + a j and B j = β j (t) + b j where a j and b j are quantum fluctuation operators with zero mean value around classical c-number amplitudes α j (t) and β j (t) of the system operators, respectively. Under the conditions |α , standard linearization techniques 4 can be applied by substituting A j = α j (t) + a j and B j = β j (t) + b j into Eq. (2). In this way, we obtain a set of nonlinear differential equations for the classical mean values α j (t) and β j (t) only (discarding the terms with quantum fluctuation operators and quantum noise operators) where we have assumed α 1 (t) = α 2 (t) = α(t) and β 1 (t) = β 2 (t) = β(t) considering the system symmetry. One can also get the following linearized QLEs for the quantum fluctuations by neglecting the terms containing classical mean values only and all nonlinear terms such as a 1 b 1 and † a b which correspond to a system Hamiltonian with linearized optomechanical interactions and G(t) = −gα(t) being the effective detuning and enhanced optomechanical coupling, respectively.

Effective Hamiltonian and the mechanism
In this paper, we focus on the weak optomechanical coupling regime, namely ω  g/ 1 m . In this case, approximate analytical solutions for Eq. (4) can be found by expanding the classical mean values α(t) and β(t) in powers of g as 29,30,38 Substituting these expressions into Eq. (4), one finds the equations for zero order of g When a four-tone driving laser One can follow similar procedures to derive higher-order corrections which are tedious and will not be presented here. In view of α(t) (1) In the asymptotic regime, the Hamiltonian in Eq. (6) then becomes To obtain the targeted Hamiltonian, we select the modulating frequencies ω k as specified in Fig. 2. We then rewrite the Eq. (11) in the interaction picture of and make the rotating-wave approximation by neglecting all fast oscillating terms under the conditions J > 2ω m and ω  G m k to arrive at the effective Assuming that we finally have the following kind of Hamiltonian with beam-splitter-like interactions The introduced Bogoliubov modes θ 1 and θ 2 are defined as unitary transformations of the mechanical modes b 1 and b 2 with a two-mode squeezed operator, respectively, Note that the joint ground state of θ 1 and θ 2 is the two-mode squeezed vacuum state of the mechanical modes b 1 and b 2 , which can be readily checked θ . For sufficiently small mechanical damping rate, the mechanical modes only weakly interact with the mechanical thermal baths with relatively large mean thermal occupancies. The dynamics of mechanical modes, i.e., the Bogoliubov modes, is dominated by the coupling to the cavity modes. In fact, the Bogoliubov modes θ 1 and θ 2 can be simultaneously cooled to near ground states after long enough time via the beam-splitter-like interactions with cavity modes a 1 and a 2 which in turn strongly couple to optical thermal baths with neglectable small thermal occupancies. In this way, the dissipation of the cavity modes is exploited to cool the Bogoliubov modes, in other words, to generate two-mode squeezed (thermal) states of the mechanical modes.

Evolution equation of the covariance matrix
The fact that the dynamics of our four-mode bosonic system is governed by a linearized Hamiltonian ensures that the evolved states are Gaussian states whose information-related properties are fully represented by the 8 × 8 covariance matrix (CM) σ with entries defined as [45][46][47] is a vector of dimensionless quadrature operators related to bosonic . By further introducing the vector of input noise quadrature operators we can transform the QLEs for the quantum fluctuations in Eq. (5) into a more compact form Here M is an 8 × 8 real coefficient matrix where G R (t) and G I (t) are respectively real and imaginary parts of the effective coupling G(t). Given that the quantum states of our system remain Gaussian throughout the evolution, the QLEs in Eq. (5) is equivalent to the equation of motion for the CM. From Eqs (3), (17) and (19), we can deduce a linear differential equation for the CM 1 T where D is a diffusion matrix whose components are associated with the input noise correlation functions in Eq. (3) In the following, we will utilize the Eq. (21) to study the time evolution of the mechanical entanglement. Note that the coefficient matrices in Eq. (20) correspond to the system Hamiltonian in Eq. (6) where the only approximation is the commonly used linearization techniques in optomechanics.

Discussion
The entanglement of two mechanical oscillators can be calculated from the two-mode CM σ m which is the first four rows and columns of the CM σ for the whole system. When σ m is arranged in the following block form with each V j being a 2 × 2 matrix, the entanglement of the mechanical modes b 1 and b 2 called logarithmic negativity can then be calculated 48,49 3 . To demonstrate the mechanism of generating mechanical entanglement via cavity dissipation discussed priviously, we plot in Fig. 3 the time evolution of the entanglement E N between two mechanical modes b 1 and b 2 with all mechanical and cavity modes initially in thermal equilibrium with their baths. The results are numerically evaluated with the full linearized Hamiltonian in Eq. (6), using a set of experimentally achievable parameters 39,42,50,51 . Obviously, there is no entanglement between b 1 and b 2 until the Bogoliubov modes θ 1 and θ 2 have been sufficiently cooled after some time. Then, following a dramatic increase, E N eventually tends to be saturated with small vibrations which derive from the effects of non-resonant terms. The steady-state mechanical entanglement achieved (E N ~ 1.6) in our scheme is much larger than the optomechanical entanglement (E N ~ 0.63) generated in many previous schemes 1,4,21,36 which are based on the coherent parametric interactions and subjected to the stability constraint. For larger mechanical decay rate, one has stronger interactions of the mechanical modes with the mechanical thermal baths, which raises the final effective temperature of Bogoliubov modes, accordingly reducing the steady-state mechanical entanglement as illustrated in Fig. 3.
As shown in related three-mode cases 30,33 , the amount of stationary entanglement is a nonmonotonic function of the ratio of the effective couplings G + /G − . The increase in G + (holding G − constant) has two competing effects. On the one hand, it can increase the squeezing parameter = − + − r G G tanh ( / ) 1 of the two-mode squeezed thermal states for the mechanical modes in the stationary regime, i.e., enhance the stationary entanglement. On the other hand, it will weaken the cooling effects of the Bogoliubov modes due to the declining coupling strength between the cavity mode a 1(2) and the Bogoliubov mode θ 1 (2) . The achievable mechanical entanglement is determined by balancing these opposing effects.
Another effect of varying G + /G − that needs particular attention is its influence on the dynamics of the mechanical mean values. All previous analyses are based on the assumptions that we have a stationary and well-behaved system dynamics which does not enter the strongly oscillating, unstable, or chaotic regime [52][53][54][55][56] . From Eqs (10), (12) and (14) we have Assuming that all system parameters except the driving amplitudes are kept fixed, we have G + /G − proportional to the ratio of the driving amplitudes E 2 /E 1 . Figure 4 shows the time evolution of the classical mechanical mean values β(t) for two selected values of G + /G − , which are numerically calculated via Eq. (4). After some transient time, β(t) reaches self-sustained oscillations. The oscillation amplitudes can be apparently different depending on the values of the ratio G + /G − . Generally, for those driving amplitudes E 2 /E 1 corresponding to G + /G − that are not  very close to 1 (such as G + /G − = 0.8), as demonstrated in Fig. 4(a) we have small mecahincal oscillations in the asymptotic regime so that β β ω and ∆ ∆  t ( ) 0 , ensuring that the effective Hamiltonian derived in Eq. (15) is valid. However, for some larger E 2 /E 1 associated with G + /G − only slightly less than 1 (such as G + /G − = 0.98), β(t) oscillates with amplified amplitudes 55,56 as shown in Fig. 4(b). In this case, the condition ∆ ∆  t ( ) 0 and accordingly the derived Hamiltonian in Eq. (15) for generating mechanical entanglement are no longer effective, thus no significant entanglement can be obtained.
Finally we would like to briefly discuss the experimental feasibility of our scheme. The microtoroidal cavities in our system shown in Fig. 1 can be fabricated at the edges of two separate chips. The whispering-gallery modes of the two microtoroidal cavities couple to mechanical modes via the radiation pressure, while the direct coupling between the two cavities can be achieved by placing them on nanopositioning systems capable of precisely controlling the distance 39,40 . An electro-optic modulator is exploited to generate the required time-modulated amplitudes of lasers pumping the whispering-gallery modes of the resonators by means of taper couplings. A set of optomechanical parameters comparable to the adopted ones in our scheme is reported in recent experiments with mechanical resonance frequency ω m = 2π × 78 MHz, single-photon optomechanical coupling rate g = 4.35 × 10 −5 ω m , cavity decay rate κ = 0.09 ω m , mechanical decay rate γ m = 1.5 × 10 −4 ω m , and thermal occupancy ~2 using cooling technique 42,50 . Although we use coupled microtoroidal optomechanical system in our scheme, the mechanism for generating distant mechanical entanglement can be extended to other coupled optomechanical systems, such as microwave-circuit optomechanical systems where similar parameters are achieved in recent state-of-the-art experiments 57,58 .

Conclusion
In summary, we have proposed a scheme to generate steady-state mechanical entanglement in a coupled optomechanical system. By applying four-tone driving lasers with weighted amplitudes and specific frequencies, we can get beam-splitter-like interactions between the cavity modes and the delocalized Bogoliubov modes of the two mechanical oscillators, which enables the cooling of the Bogoliubov modes by the cavity decay as long as the mechanical decay rate is small. The mechanical oscillators are then driven to close to two-mode squeezed states when the Bogoliubov modes have been sufficiently cooled to near vacuum after some time of dissipative dynamics. By balancing the opposing effects of varying the ratio of the effective optomechanical couplings and carefully avoiding the system parameters that may lead to the unwanted amplified oscillations, we obtain steady-state mechanical entanglement that is significantly larger than the generated entanglement based on the coherent parametric interactions in many previous schemes.