Preservation Macroscopic Entanglement of Optomechanical Systems in non-Markovian Environment

We investigate dynamics of an optomechanical system under the non-Markovian environment. In the weak optomechanical single-photon coupling regime, we provide an analytical approach fully taking into account the non-Markovian memory effects. When the cavity-bath coupling strength crosses a certain threshold, an oscillating memory state for the classical cavity field is formed. Due to the existence of the non-decay optical bound state, a nonequilibrium optomechanical thermal entanglement is preserved even without external driving laser. Our results provide a potential usage to generate and protect entanglement via non-Markovian environment.

The investigation of decoherence and dissipation process induced by environment is a fundamental issue in quantum physics [1][2][3][4] . Understanding the dynamics of such nonequilibrium open quantum system is a challenge topic which provides us the insight into the issue of quantum-classical transitions. Protecting the quantum property from decoherence is a key problem in quantum information science, therefore a lot of effort has been devoted to develope the methods for isolating systems from their destructive environment. Recently, people recognize that properly engineering quantum noise can counteract decoherence and can even be used in robust quantum state generation 5,6 . Meanwhile the features of the non-Markovian quantum process have sparked a great interest in both theoretical and experimental studies [7][8][9][10][11][12] . Numerous quantitative measures have been proposed to quantify non-Markovianity [13][14][15][16][17] .
As a promising candidate for the exploration of quantum mechanical features at mesoscopic and even macroscopic scales and for quantum information procession, cavity optomechanical systems come as a well-developed tool and have received a lot of attentions [18][19][20][21][22] . In the theoretical research of the cavity optomechanical system, the environment is often treated as a collective non-interacting harmonic oscillators, and the quantum Langevin equations 23 are developed to describe the radiation-pressure dynamic backaction phenomena. Significant progresses have been made in this framework 18,[24][25][26] . Almost all of these studies are focussing on the scenario of memoryless environment. However in many situations for optical microcavity system, the backaction of the environment and the memory effect of the bath play a significant role in the decoherence dynamics 27,28 . Quite recently, a nonorthodox decoherence phenomenon of the mechanical resonator is also observed in experiment 29 , which clearly reveals the non-Markovian nature of the dynamics. Therefore, it is necessary to investigate the non-Markovian dynamics for the nonlinear cavity optomechanical system so that we can use the memory effects to produce and protect coherence within it.
In the following, we investigate the cavity optomechanical dynamics under non-Markovian environment and put forward a method to solve the exact Heisenberg-Langevin equations where the non-local time-correlation of the environment is included. We find that when the cavity-bath coupling strength crosses a certain threshold, the optical bound state is formed, giving rise to the nonequilibrium dynamics of the entanglement. This remarkable result indicates the possibility of long-time protection of macroscopic entanglement via structured reservoirs.
Here ω c is the frequency of the cavity mode with bosonic operators â and ˆ † a satisfying =ˆ † a a [ , ] 1, while the quadratures q and p =q p i ([ , ] ) are associated to the mechanical mode with frequency ω m . The third term describes the optomechanical interaction at the single-photon level with coupling coefficient . The cavity is driven by an external laser with the center frequency ω 0 . The environment of such system can be described by a collection of independent harmonic oscillators 31 . The reservoir as well as the system-reservoir interaction is then given by The first term is the free energy of the cavity reservoir with the continuous frequency ω k as well as the hopping interaction between the cavity and the environment with the coupling strength g k . The second summation describes a mirror undergoing Brownian motion with the coupling through the reservoir momentum 23,31 . Here ω l is the reservoir energy of the mechanical mode, and γ l stands for the mirror-reservoir coupling. Dynamics of the system. To achieve a comprehensive understanding of the decoherence dynamics, one has to rely on precise model calculations. To this end, by making use of the reference frame rotating at the laser frequency, we can obtain the Heisenberg-Langevin equations  . When the intracavity photon number  a 1, we can apply the so-called linearization method 24,35 , which means the relevant quantum operators can be expanded about their respective mean values: a a q p ( , , , ) T . The superscript T represents the transpose operation. Then Eqs.(3) can be decomposed into two parts. The first is the classical part that describing the classical phase space orbits of the first moments of operators where, for simplicity, we have assumed In the single-photon weakly coupling regime, the coupling strength g 0 is the smallest parameter in Eqs. (4). We therefore perform the regular perturbation expansion in ascending powers of the rescaled dimensionless variable g 0 /ω m (for computational convenience one may set ω m = 1, and the other rescaled parameters are in units of ω m ). By substituting the expressions with rescaled g 0 (i.e., are the real and imaginary part of the Laplace transform of the self-energy correction respectively, and the Green's functions α 0 and q 0 obey the Dyson equations (18) with the initial conditions α . Base on Eq. (5), we can see that the non-vanishing intracavity field α 0 (∞ ) may induce an equilibrium position q 1 (∞ ) for the oscillator. This leads to the effective , which also alters the asymptotic dynamics of the cavity field. Accordingly, the interplay between the non-Markovian and nonlinear effects can be described more precisely in this way. Within the parameter space of our consideration, ∆ ≈ ∆ c , the validity of the power series assumption is guaranteed by the numerical simulations (see the Supplemental Material).
For general bosonic environments, the spectral density should be a Poisson-type distribution function 37 . We consider that the spectrum is of the where η a is a dimensionless coupling constant between the system and the environment, and ω  a is a high-frequency cutoff 2,37 . The parameter s a classifies the environment as sub-Ohmic (0 < s a < 1), Ohmic (s a = 1), and super-Ohmic (s a > 1). Using the modified Laplace transformation, one can give an analytical solution for the nonequilibrium Green's function . The first term survives only when ω ω and the pole ω r is located at the real z axis. This term corresponds to 'localized mode' 36,38 , which means that the cavity field oscillates with frequency ω r and does not decay. It is seem that the photons are 'trapped' in the cavity due to the backflow of the non-Markovian environment and do not diffuse. It is a term that determines the asymptotic dynamics of the optical field. Physically, this is equivalent to generate a bound state of the joint cavity-reservoir system 38 , which can be also determined by solving the energy eigenstates of the total Hamiltonian 39 , and such bound state is actually a stationary state with a vanishing decay rate during time evolution. The second term corresponds to nonexponential decays. In the long time limit, the bound state as well as the driving laser give rise to the non-vanishing intracavity photon numbers with . On the other hand, for the mechanical mode, however, due to the discontinuity of the self-energy correction at the real z axis, q 0 (∞ ) = 0. We find that the corresponding mechanical bound state can be formed only if some sharp cutoff appear in the spectral density (see section 'Methods' for details). Considering the feasibility, the single photon coupling g 0 (key parameter) is set close to that of recently performed optomechanical experiments 35 . Figure 1(a) is the density plot of the maximum value of |α 0 | in the long-time limit. It maps out the regions in parameter space where localized bound state occurs. As a result, a threshold characterizing the transition from weak to strong non-Markovian regions can be defined, and it is marked by the green-dashed line, which satisfying ω c + K c (− ω 0 ) = 0. In the weakly non-Markovian region 0 < η c < 0.01, as shown in the inset of Fig. 1(a), the red-detuned laser gives rise to a strong stationary amplitude values, which is determined by . Figure 1(b,c) show the dynamic evolution of |α 0 | and q 0 , where the first order solutions are shown in the insets. In Fig. 1(b) the optical bound state is formed when η c is above the threshold, while for Fig. 1(c), the coordinate average value of the oscillator is no longer zero, which reveals that radiation pressure push the oscillator to a new equilibrium position. The evolution in the phase space is shown in Fig. 1(d,e). The limit cycles that characterize the properties of the steady state depends on the initial states, which indicates the non-Markovian property and reflects the memory effect. The initial information of the system is maintained.
We now turn to the quantum fluctuation of operators that describes the actual quantum dynamics, which are deduced as t 0 By assuming the formal solution ( ), and substituting it into Eq. (7), we have subjected to the initial conditions is interpreted as a noise term 23,24 that depends on the initial states of the environments. It is easy to obtain the quadrature operator The first term is the projection of the quadratures on the system's Hilbert space. The second term characterizes the magnitudes of the input-noise that satisfies the non-local time correlation relations 23 . The last two terms describe the effect of initial system-reservoir correlations, which had been identified as an important factor in the decoherence dynamics 40,41 . For the sake of simplicity, here we assume as usual the system and the reservoirs are initially uncorrelated, and the reservoirs are in thermal states. Then the noise vector ξ t ( ) obeys the non-Markovian self- The thermal correlation functions are defined as

Preservation of entanglement.
To bring quantum effects to the macroscopic level, one important way is the creation of entanglement between the optical mode and the mechanical mode. If the initial state of the system is Gaussian, then Eq. (7) will preserve the Gaussian character. The entanglement can therefore be quantified via the logarithmic negativity 42  . In order to show clearly the phenomenon of entanglement preservation by non-Markovianity, we shut the classical pumping field all the time by setting E = 0. We see in Fig. 2 that the shape of the Ohmic spectrum characterized by s m , has negligible impact on the initial stage of evolution (see the insets). This is because the time scale of the mechanical oscillator as well as its bath is much larger than that of the bath of the cavity, which means the non-Markovian memory effects induced by the bath of the oscillator can be ignored completely when  t d 1/ , where d is the bandwidth of the oscillator's bath. For the short-time and long-time scale, s m slightly affect the amplitude and evolution period of the entanglement. The non-Markovian environments affect the entanglement through the bath spectrum as well as the specific type of system-reservoir interaction, which is clearly embodied in the correction functions f c (t) and f m (t), the form of which directly determines the entanglement dynamics. The Ohmic spectrum of the photon bath plays a decisive role for the time evolution of the entanglement, as it allows the existence of optical bound state. If the environmental spectral density is not a finite band, then the optical bound state is unable to form 36 , and the memory effect is not strong enough to compensate the loss of photons. Consequently the entanglement will be destroyed completely. However for the phonon bath, the Ohmic spectrum is unable to form a bound state between the mechanical mode and the environment unless some band gaps appear in the spectrum (see section 'Methods' for details). Therefore the resonator may only present weak non-Markovian effect, and it plays a indecisive role for the entanglement (it has limited influence on the evolution period and the value of entanglement). By changing the cutoff frequency of the phonon spectrum, we verified that the cutoff frequency do not play a major role in the entanglement dynamics, which also manifests the weak non-Markovian effects of the resonator. In addition, although the sudden death and rebirth of entanglement is also observed, it differ with the case of Markovian environment 44 where the photons would rapidly dissipate to the memoryless environment. Here, due to the present of the bound state, the entanglement can be produced and can be preserved in non-Markovian environment after long time even without external drive. This provides a way to decoherence control of optomechanical systems, in which, a controllable quantum environment indeed have the ability to protect the quantum correlation of the internal system.
We finally discuss how to detect the generated optomechanical entanglement. Although we show the evolution of the dynamics, it still might be difficult for detection of the temporal entanglement. Fortunately, we can obtain a larger entanglement for long time, for example, ω m t around 1040 than short time ω m t < 50, which means that we can detect the entanglement at long time evolution. For measurement the logarithmic negativity, one can detect all independent entries of the covariance matrix which may be achievable by utilizing the Q-switching technology 45 . If the Q-switch is off, the composite system (i.e., the cavity and it's environment) is closed, and the system has no output. When the Q-switch is on, the measurement can be performed by homodyning the cavity output. The duration of the Q-switched pulse should be short enough, so that it has a negligible effect on the system. The mechanical mode can be detected by employing the method put forward by Vitali et al. 18 , where an assistant fixed mirror and the mechanical oscillator form an additional 'probe' cavity mode. If this additional cavity is driven by a much weaker intracavity field so that its back-action on the mechanical mode can be neglected. In addition, the interaction between the additional cavity and it's bath is also weak enough, so it can be well treated in the Markovian region. By adjusting parameters, the probe mode adiabatically follows the dynamics of the mechanical mode, therefore the output field gives a direct measurement of the mechanical mode. Finally the covariance matrix can be determined by changing the phases of the corresponding local oscillator and measuring the correlations between the two cavity outputs. Then one can numerically extract the logarithmic negativity.

Conclusion
In conclusion, we have put forward a scheme to preserve the entanglement of optomechanical system in non-Markovian environment. An analytical approach for describing non-Markovian memory effects that impact on the decoherence dynamics of an optomechanical system is presented. The exact Heisenberg-Langevin equations are derived, and the perturbation solution is given in the weak single-photon coupling regime. Employing the analytical solution, we have shown that, the system dynamics change dramatically when the cavity-bath coupling strength crosses a certain threshold, which corresponds to dissipationless non-Markovian dynamics. The interplay between non-Markovian and nonlinear effects can be also explained though the perturbative method. As a quantum device which may subjected to dissipative and decoherence effects, however, our results show that the surroundings of such physical setting can protect the quantum entanglement, rather than destroy it even in the long-time scales. Our research provides a new approach to explore non-Markovian dynamics for the cavity optomechanical systems. The classical phase space orbits. In order to understand the classical dynamics that describes an optomechanical system embedded in memory environments, we apply a perturbative solution for Eq. (4). In realistic physical settings, the single-photon coupling strength is extremely weak, thus the perturbation expansion approach is a good approximation. The corresponding variables in the zeroth-order is given by is given in Eq. (6).
The dissipationless non-Markovian dynamics usually requires ∞ ≠ f ( ) 0 (in absence of any external driving). In other word, it is referred to a process of nonthermal stabilization 46 , in which the initial information of the system partially maintains. In this case, there should be poles exist on the real z axis for  f z ( ). Due to the discontinuity of the Laplace transform of the self-energy correction, however, the real poles exist only in the frequency regions that all spectral density vanishes. In this region, . From the physical point of view, it can be explained by the bound state generated between the system and its environment [47][48][49] , which is actually a stationary state with a vanishing decay rate. It usually occurs when the environments has band gaps or a finite band.
The reservoir and its coupling to the system are fully characterized by the spectral density J(ω). For thermal bosonic (photon and phonon) baths, the most general spectral density as introduced in 37 should be a Poisson-type distribution function with some frequency cutoff. From the experimental point of view, it is often not possible to model the environment in an accurate way, because their density of states is unknown. However one can in principle extract the information of the spectral density of the heat baths by ingenious experiment. Recently, this work is accomplished by Gröblacher et al. 29 . The main experimental results show that the specific geometry of the slab cause strongly sub-Ohmic spectral densities. To theoretically explore the decoherence dynamics, we therefore consider the Ohmic spectrum.
On the other hand, as discussed above, the imaginary part of ω ∑ ( ) with ω ω ≤ − r 0 . For the mechanical mode, the imaginary part of ω ∑ ( ) m vanishes only at the zero point, although an Ohmic spectrum with a finite band is considered. The main reason lies in their different internal memory kernels, the structure of which is entirely determined by the specific type of interactions between the system components and the corresponding environment. For the zero point, note that . Thus it is unable to form a bound state between the mechanical mode and the environment unless some band gaps appear in the spectrum 38 .