# Cooling a Mechanical Resonator to Quantum Regime by heating it

###### Abstract

We consider a mechanical resonator made of diamond, which contains a nitrogen-vacancy center (NV center) locating at the end of the oscillator. A second order magnetic gradient is applied and inducing coupling between mechanical modes and the NV center. By applying proper external magnetic field, the energy difference between NV center electron spin levels can be tuned to match the energy difference between two mechanical modes and . A laser is used for continuously initializing the NV center electron spin. The mode with lower frequency is driven by a thermal bath. We find that the temperature of the mode is significantly cooled when the heating bath temperature is increased. We discuss the conditions that quantum regime cooling requires, and confirm the results by numerical simulation. Finally we give the intuitive physical explanation on this unusual effect.

###### pacs:

******## I introduction

Cooling and manipulating the motion of mechanical resonator is widely studied now ASP2014 . It has many applications in quantum information science Regal2011 ; Barzanjeh2012 ; Wang2012 ; Tian2012 ; Rips2013 ; Yin2009a ; Yin2015a ; Barzanjeh2015 , testing quantum effects in macroscopic systems PEN1996 ; Connell2010 ; Teufel2011 ; Poot2012 ; Pikovski2012 ; NIM2013 ; Bassi2013 ; Isart2011 ; Yin2013 ; Li2016 ; Barzanjeh2016 , ultra-sensitive sensing Caves1981 ; Geraci2010 ; Yin2011 ; Arvanitaki2013 ; Huang2013 ; Zhao2014 , etc. In order to cool mechanical resonator, we should couple it to other cold systems such as cavity mode (optomechanics) sideband1 ; sideband2 ; Liu2013 , electronic spin (Nitrogen-vacancy centers) RAB2009 ; MacQuarrie2016 ; Li2016b , etc. In optomechanics, laser cooling of mechanical resonators is a typical method. Intuitively, this is not hard to understand. A laser produces coherent light, which is well ordered (cold), and will make the mechanical system couple with a colder cavity bath sideband1 ; sideband2 . However, counter-intuitively, recent research showed that a thermal light, which is very hot , can also cool an opto-mechanical system MAR2012 . In this system, the mechanical resonator to be cooled is coupled with two optical modes Yin2009 . The low frequency mode is in contact with a hot thermal light while the high frequency mode is not, acting as an auxiliary part. After the system reaches the equilibrium state, the thermal phonon number in the mechanical mode will be smaller.

In a hybrid system that contains a NV center and a mechanical resonator, a magnetic gradient is applied to couple the NV center with the resonator. The NV center is continuously initialized by a laser, and could be used for cooling the mechanical resonator to the ground state RAB2009 . In previous literatures, the first order magnetic gradient was used for inducing coupling between the NV center and the mechanical mode RAB2009 ; YIN2015 ; Zhou2010 ; Xu2009 ; Chen2010 ; Arcizet2011 ; Kolkowitz2012 . It is quite naturally to ask whether the higher order magnetic gradient could be useful for coupling the NV center and the mechanical resonator and cooling it.

In this paper, we consider system that contains a mechanical cantilever and a NV center YIN2015 . A second order magnetic gradient is applied to induce the coupling between NV center and the two mechanical modes and . The electron spin transition frequency of the NV center is tuned to match the mechanical modes frequencies difference. The NV center is continuously resetted to the ground state by a laser. In this way, the NV center acts as an effective vacuum bath to cool the mode . An incoherent driving (thermal bath) is applied on the lower frequency mode , and resulting a higher temperature of it. The incoherent driving on mode also increases the effective coupling between mode and the NV center, and cools the mode . We find that the quantum regime cooling is possible under the current experimental conditions.

The paper is orgnized as follows. In section II, we introduce the scheme of the second order magnetic field gradient induced coupling between NV center and mechanical modes. In section III, we derived the analytical solution for the system. In Section IV, the analytical results are verified by numerical simulation. In section V, a discussion and conclusion are given.

## Ii The scheme

As shown in Fig. 1, we consider a diamond cantilever resonator that contains a NV center at the end of it. Two mechanical modes and are under consideration, with frequencies . Two spin states of the NV center electron are relevant here, and . By applying proper magnetic field, the energy split between the two states of the NV center, , is equal to the energy difference of the two mechanical modes, . Due to the second order magnetic gradient , the two mechanical modes and both couple with the NV center.

Initially, both mechanical modes and are in the thermal states with the average thermal phonon number and . The mode is driving by a thermal bath with average phonon number much higher than its environment. The NV center electron spin is initialized to the state at the beginning. The Hamiltonian of the whole system can be expressed as , where is the non-interaction part and describes the interaction between the mechanical resonator and the NV center electron spin.

(1) | ||||

where , , and , and is the zero field fluctuation for the mechanical mode , () is the effective mass for mode (). Experimentally, , can be realized (see Appendix).

Under the rotating wave approximation which shift the energy zero point by , we can change into

(2) |

Unavoidably, the system will go through some intrinsic loss. The mechanical decay and is in the order of (see Appendix). The intrinsic decay rate for NV center is much less than Hz at low temperature. By adding initializing laser, the decay rate for the NV center electron spin would be much larger than the decay rates for the mechanical modes and . The dephasing of NV center is neglected, as it could be eliminated by continuous dynamical decoupling Cai2012 . However, as we will use the same model to describe the decay of the spin and the vibration, here we should make sure that . As the evolution of the system contains decay and damping, it is convenient to describe it with quantum master equation OPS , which can be written out explicitly as

(3) | ||||

refers to the notation in Lindblad form, which is . The Eq.(3) will be solved in both analytical and numerical methods. The results of the both methods will be compared and discussed.

## Iii Analytical Solution

As mentioned above, the phonon bath for mode is very large, so the average phonon number of mode is almost always equal to the average phonon number in the bath, that is, . Moreover, as the phonon number of mode is much smaller than that of mode , it is safe to drop the term in . If we denote as the phonon occupation number fluctuation of mode , under the rotating wave picture, can be simplified to

(4) |

For the convenience of the notation, here we rename , .

Using standard method of quantum master equation QNO , and after adiabatically eliminating the mode Yin2007 ; Yang2015 , we can get the reduced master equation for mode and the electron spin mode (see Appendix)

(5) | ||||

where we have introduced two new Liouvillians for convenience

(6) | ||||

Eq. (6) can be written in a more compact form if we define the effective average excitation number of the NV center electron spin as

(7) |

In fact, is very close to zero because of the weak coupling requirement , in other words, due to the large energy split and the relative large damping rate . Then we can define an effective Liouvillian for the spin as well .

Now the equation of motion for the reduced density matrix becomes

(8) |

The corresponding adjoint equations for the phonon number operator and excitation number operator can be developed then (see Appendix).

(9) | ||||

Taking the average on both sides respect to the quantum state of the system, and under average field approximation which refers to , the stationary state for the average phonon number of mode can be solved from the quadratic equation below after we drop the nonphysical solution.

(10) |

where

(11) | ||||

In Fig. 2, the stationary phonon number of mode is plotted as a function of the average phonon number of mode . Unlike the previous Ref. MAR2012 , where there was an optimal heating bath temperature, here increasing heating bath temperature will always cooling the resonator, which is an outstanding property. As environment temperature increase, the quantum regime cooling can be realized by increasing both the mode temperature and NV center decay .

We can further look at how different values of can change the effect of cooling. The qualitative trend can be seen clearly from Fig. 3. For a given value of , when increases from a small value to a relative large value, there exists a certain for which the cooling effect is the largest. Before the optimal cooling point, is not large compared with the effective coupling (e.g., ), so increase will benefit the cooling process. However, after that, larger will result in smaller cooling rate. When increases, the optimal will increase, too. In the limit of , we can get a simple expression for under the weak coupling assumption

(12) |

As shown in Fig. 3, such an approximation is quite reasonable. In the limit of , the condition for quantum regime cooling () is

(13) |

which can be fulfilled by tuning with laser. These analytical results are based on rotating wave approximation (RWA). When thermal phonon number , the RWA is no longer valid. From the parameters used in Fig. 2, RWA requires , which is obviously fulfilled in Fig. 2 and the following figures in this paper.

## Iv Numerical Simulation

Now we discuss the numerical simulation of the cooling process, and compare the result with the analytical ones above. The idea is simple. We expand the quantum operators and quantum states in the uncoupled representation, that is, the Kronecker tensor product of the occupation number representations for the three modes involved. The matrix forms of Eq. (3) can then be calculated step by step with the classical four-order Runge-Kutta method. The average phonon number for mode is equal to , where is the phonon number operator of mode in the whole Hilbert space for the system, namely, . It should be pointed out, however, as the dimension of a matrix cannot be too large for a numerical simulation, the requirement is not able to be satisfied here. As a result, in Eq. (1) should be inserted in the equation of motion Eq. (3) instead of .

We present several numerical results of the stationary phonon number with the results we get from Eq. (10), which is plotted in Fig. 4. As can be seen from the figure, when is satisfied, the analytical results match the numerical ones perfectly. But when and is comparable with , the two results no longer agree with each other quantitatively. However, it is quite reasonable, because in this case several approximations for the analytical results are not valid any more. When is not very large, the term does not dominate over in Eq. (3). As is also at the same scale with in this case, now neglecting the back action on mode is not a good approximation. When mode is being cooled, mode is being heated unavoidably, thus lowering the effect of cooling. The analytical results which neglect the back action thus exaggerate the cooling effect. This trend is shown clearly from the figure. But the analytical conclusions are still valid qualitatively. For instance, if we fix , the larger is, the larger the cooling effect will be. And if we fix , there exists a -dependent value of where the cooling effect is maximum.

## V Discussion and conclusion

From the above derivation and calculation, it is clear the higher frequency mechanical mode will be significantly cooled by heating the lower frequency mechanical mode. In fact, this trend can be seen directly from the Eq. (8), where acts as the cooling term. As the electron spin is much more likely to be in the spin state due to the continuous initialization, the effect of should be much larger than that of . So is the dominant operator for evolution, and the mode will be cooled. This can be understood that mode facilitate the energy transfer between mode and the spin. Quantized energy flows from the vibration mode to the NV center. Soon after that the excitation in the NV center electron spin decays into the environment through light. As this process is repeated, the heat of the mechanical mode is gradually removed.

In conclusion, we studied the second order magnetic gradient induced coupling between NV center and the mechanical modes. We proposed to drive one mechanical resonator mode with thermal bath to cool the other mode of the resonator. We solved the motion of the system analytically under the weak coupling assumption. It is found that under rotating wave approximation final thermal phonon number of the resonator approach to a minimum as the thermal bath temperature increases. The quantum regime cooling is feasible under the current experimental conditions. We numerically solved the equations, and found that the results fitted very well with the analytic results. The cooling by heating phenomena does not violate the laws of thermodynamics. This scheme can be understood as a thermal machine or a heat engine Zhang2014 ; Mari2015 ; Mitchison . The mode is an “engine”, and the NV center electron spin takes the role of a “condenser”. This work opens up the possibility that cooling the mechanical mode in solid systems by injecting the thermal phonon, which usually heats the system.

## Vi Acknowledgements

Z.Q.Y. is supported by National Natural Science Foundation of China Grant 61435007. W.L.Y. is supported by National Natural Science Foundation of China Grant 11574353 and 11274351. P.H and J.D. is supported by 973 Program (Grant No. 2013CB921800), the National Natural Science Foundation of China (Grant No. 11227901), and the CAS (Grant No. XDB01030400).

## Appendix A Experimentally realized paramters

In this section we discuss the values of parameters , , , , , , , , that can be realized in experiments.

As in former experiments (e.g. RAB2009 for a Si cantilever), the size of the resonator can be taken as around . Considering that the density of diamond is around 1.5 times the density of silicon Tao2014 , and the fact that the effective mass of a solid and normal shaped oscillator doing simple oscillation is of the same order as its mass, we can estimate that kg. The fundamental frequency of such a mechanical resonator is around 10 MHz RAB2009 . Although the zero-field splitting of the NV center is around 2.88 GHz Wrachtrup2006 , after an external magnetic field is exerted, the separation between and states can be reduced to the order of tens of MHz. As a result, it is possible to fulfill the requirement . Here we choose , for convenience. The quality factor of the resonator can reach up to Tao2014 . Second-order magnetic gradient can realize several times of in experiment Mamin2007 . By using finite element simulation, we got that the maximum could be around , if the distance between NV center and the magnetic particles is nm.

If is set to be , we will have , , . The decay rates of the mechanical resonator can be estimated by . We need to tune the decay rate of the NV center larger enough than the mechanical decay rate to reach ground state cooling, for instance, can equals to .

## Appendix B Derivation of the Reduced Master Equation

Here we provide a derivation of the reduced master equation Eq. (5) in main text. The reduced density matrix of mode and mode after mode is eliminated is

(14) |

We can adiabatically eliminate mode to get the evolution equation for the reduced density matrix QNO

(15) | ||||

where we denote

(16) | ||||

We expand , and to the second order and can get

(17) | ||||

Once we expand the commutation relations in Eq. (17), we can directly calculated it term by term.

(18) | ||||

What we need first is the term . It is easy to get from direct calculation.

(19) | ||||

Similarly, . From these we have , .

The following results are got using exactly the same method:

(20) | ||||

will then be calculated with the results above in mind.

(21) | ||||

It deserves notice that here we only keep the slow-time-varying term, so we neglect and above. Now we come back to Eq. (18), which can be simplified term by term. As means taking the average on , the trace of odd multiples of and equals to zero. Terms containing are neglected here because their effect is just changing the energy zero point.

(22) | ||||

(23) | ||||

(24) | ||||

(25) | ||||

## Appendix C Derivation of the Adjoint Equations for the Number Operators

In this part, the equations of motion for the number operators and , Eq. (9) in main text, is derived.

The adjoint equations for Eq. (9) in main text are

(26) | ||||

From simple calculations, we have

(27) | ||||

## Appendix D Detail of the Numerical Simulation

In this section, we discuss in detail the numerical simulation. For Fig. 4 in main text, we control the accuracy so that all numbers are accurate until the order of 0.001. The evolution time is , after which the state is at the stationary state up to the accuracy. As discussed in the Appendix before, we can renormalize the frequency parameters as , , , . As we want to cover a wide range of results, four values renormalized are chosen: 2, 15, 50 and 120. Small values for the phonon numbers and should be chosen due to the calculation cost.

Here we always choose as the initial condition. As to , we choose . The dimensions of the matrices in the Hilbert space are chosen to balance the accuracy and the cost of calculation resource. For the Hilbert space where the initial phonon number equals to 1, the dimension is chosen as . For the Hilbert space where the initial phonon number equals to 2, the dimension is chosen as . For the Hilbert space where the initial phonon number equals to 3, the dimension is chosen as . For the Hilbert space where the initial phonon number equals to 4, the dimension is chosen as . Obviously, the dimension of the Hilbert space for mode is always . The step length for the Runge-Kutta method is chosen to be .

Fig. 4 in the main text shows the results of only the case, even if for the example in the Appendix, should be renormalized as . In fact, we have tried that different will not influence the result of cooling, e.g., all give out same final phonon number for mode . So we are safe to use smaller here considering the stability of the simulation. This deduction also agrees with the analytical result because Eq. (10) in main text does not contain as long as is large enough for .

The time evolution of one set of parameters from numerical simulation is shown in Fig. 5. The phonon occupation number gradually goes down as the system evolves towards the stationary state. The cooling effect indeed happens. Also shown in the figure is the time evolution starting from the analytical result Eq. (9). The two lines are similar with the cooling effect of the analytical one a little bit more obvious, which also supports our findings.

## References

- (1) M. Aspelmeyer, T. J. Kippenberg, and F. Marquardt, Reviews of Modern Physics 86, 1391 (2014).
- (2) C. Regal and K. Lehnert, in Journal of Physics: Conference Series (IOP Publishing, 2011), p. 012025.
- (3) S. Barzanjeh, M. Abdi, G. J. Milburn, P. Tombesi, and D. Vitali, Phys. Rev. Lett. 109, 130503 (2012).
- (4) Y. Wang and A. A. Clerk, Phys. Rev. Lett. 108, 153603 (2012).
- (5) L. Tian, Phys. Rev. Lett. 108, 153604 (2012).
- (6) S. Rips and M. J. Hartmann, Phys. Rev. Lett. 110, 120503 (2013).
- (7) Z. Yin and Y. Han, Physical Review A 79, 024301 (2009).
- (8) Z. Yin, W. Yang, L. Sun, and L. Duan, Physical Review A 91, 012333 (2015).
- (9) S. Barzanjeh, S. Guha, C. Weedbrook, D. Vitali, J. H. Shapiro, and S. Pirandola, Phys. Rev. Lett. 114, 080503 (2015).
- (10) R. Penrose, General relativity and gravitation 28, 581 (1996).
- (11) A. D. O’Connell, M. Hofheinz, M. Ansmann, R. C. Bialczak, M. Lenander, E. Lucero, M. Neeley, D. Sank, H. Wang, and M. Weides, Nature 464, 697 (2010).
- (12) J. Teufel, T. Donner, D. Li, J. Harlow, M. Allman, K. Cicak, A. Sirois, J. D. Whittaker, K. Lehnert, and R. W. Simmonds, Nature 475, 359 (2011).
- (13) M. Poot and van der Zant, Herre SJ, Physics Reports 511, 273 (2012).
- (14) I. Pikovski, M. R. Vanner, M. Aspelmeyer, M. Kim, and C̆. Brukner, Nature Physics 8, 393 (2012).
- (15) S. Nimmrichter and K. Hornberger, Phys. Rev. Lett. 110, 160403 (2013).
- (16) A. Bassi, K. Lochan, S. Satin, T. P. Singh, and H. Ulbricht, Reviews of Modern Physics 85, 471 (2013).
- (17) O. Romero-Isart, A. C. Pflanzer, F. Blaser, R. Kaltenbaek, N. Kiesel, M. Aspelmeyer, and J. I. Cirac, Phys. Rev. Lett. 107, 020405 (2011).
- (18) Z. Yin, T. Li, X. Zhang, and L. Duan, Physical Review A 88, 033614 (2013).
- (19) T. Li and Z. Yin, Science Bulletin 61, 163 (2016).
- (20) S. Barzanjeh and D. Vitali, Physical Review A 93, 033846 (2016).
- (21) C. M. Caves, Physical Review D 23, 1693 (1981).
- (22) A. A. Geraci, S. B. Papp, and J. Kitching, Phys. Rev. Lett. 105, 101101 (2010).
- (23) Z. Yin, T. Li, and M. Feng, Physical Review A 83, 013816 (2011).
- (24) A. Arvanitaki and A. A. Geraci, Phys. Rev. Lett. 110, 071105 (2013).
- (25) P. Huang, P. Wang, J. Zhou, Z. Wang, C. Ju, Z. Wang, Y. Shen, C. Duan, and J. Du, Phys. Rev. Lett. 110, 227202 (2013).
- (26) N. Zhao and Z. Yin, Physical Review A 90, 042118 (2014).
- (27) I. Wilson-Rae, N. Nooshi, W. Zwerger, and T. J. Kippenberg, Phys. Rev. Lett. 99, 093901 (2007).
- (28) F. Marquardt, J. P. Chen, A. Clerk, and S. Girvin, Phys. Rev. Lett. 99, 093902 (2007).
- (29) Liu Yong-Chun, Hu Yu-Wen, Wong Chee We, and Xiao Yun-Feng, Chinese Physics B 22, 114213 (2013).
- (30) P. Rabl, P. Cappellaro, M. G. Dutt, L. Jiang, J. Maze, and M. D. Lukin, Physical Review B 79, 041302 (2009).
- (31) E. MacQuarrie, M. Otten, S. Gray, and G. Fuchs, arXiv preprint arXiv:1605.07131 (2016).
- (32) P. Li, Z. Xiang, P. Rabl, and F. Nori, Phys. Rev. Lett. 117, 015502 (2016).
- (33) A. Mari and J. Eisert, Phys. Rev. Lett. 108, 120602 (2012).
- (34) Z. Yin, Physical Review A 80, 033821 (2009).
- (35) L. Zhou, L. Wei, M. Gao, and X. Wang, Physical Review A 81, 042323 (2010).
- (36) Z. Xu, Y. Hu, W. Yang, M. Feng, and J. Du, Physical Review A 80, 022335 (2009).
- (37) Q. Chen, Z. Xu, and M. Feng, Physical Review A 82, 014302 (2010).
- (38) O. Arcizet, V. Jacques, A. Siria, P. Poncharal, P. Vincent, and S. Seidelin, Nature Physics 7, 879 (2011).
- (39) S. Kolkowitz, A. C. Jayich, Q. P. Unterreithmeier, S. D. Bennett, P. Rabl, J. G. Harris, and M. D. Lukin, Science 335, 1603 (2012).
- (40) Z. Yin, N. Zhao, and T. Li, Science China Physics, Mechanics & Astronomy 58, 1 (2015).
- (41) J. Cai, B. Naydenov, R. Pfeiffer, L. P. McGuinness, K. D. Jahnke, F. Jelezko, M. B. Plenio, and A. Retzker, New Journal of Physics 14, 113023 (2012).
- (42) H. Breuer and F. Petruccione, The theory of open quantum systems (Oxford University Press on Demand, 2002).
- (43) C. Gardiner and P. Zoller, Quantum noise: a handbook of Markovian and non-Markovian quantum stochastic methods with applications to quantum optics (Springer Science & Business Media, 2004), 56.
- (44) Z. Yin, F. Li, and P. Peng, Physical Review A 76, 062311 (2007).
- (45) C. Yang, J. An, W. Yang, and Y. Li, Physical Review A 92, 062311 (2015).
- (46) J. Wrachtrup and F. Jelezko, Journal of Physics: Condensed Matter 18, S807 (2006).
- (47) Y. Tao, J. Boss, B. Moores, and C. Degen, Nature communications 5 (2014).
- (48) H. Mamin, M. Poggio, C. Degen, and D. Rugar, Nature nanotechnology 2, 301 (2007).
- (49) K. Zhang, F. Bariani, and P. Meystre, Phys. Rev. Lett. 112, 150602 (2014).
- (50) A. Mari, A. Farace, and V. Giovannetti, Journal of Physics B: Atomic, Molecular and Optical Physics 48, 175501 (2015).
- (51) M. T. Mitchison, M. Huber, J. Prior, M. P. Woods, and M. B. Plenio, arXiv preprint arXiv:1603.02082 (2016).