Entanglement dynamics of two nitrogen vacancy centers coupled by a nanomechanical resonator

TL;DR

This study investigates how entanglement between two remote nitrogen vacancy centers in diamond evolves over time when coupled via a nanomechanical resonator, revealing oscillations, sudden death, and revivals influenced by decoherence.

Contribution

It provides an analytical model of entanglement dynamics in NV centers coupled through a nanomechanical resonator, including effects of decoherence and controllable entanglement manipulation.

Findings

Entanglement oscillates and can be controlled by system parameters.

Sudden death and revival phenomena are observed.

Decoherence affects entanglement dynamics and revival patterns.

Abstract

In this paper we study the time evolution of the entanglement between two remote NV Centers (nitrogen vacancy in diamond) connected by a dual-mode nanomechanical resonator with magnetic tips on both sides. Calculating the negativity as a measure for the entanglement, we find that the entanglement between two spins oscillates with time and can be manipulated by varying the parameters of the system. We observed the phe- nomenon of a sudden death and the periodic revivals of entanglement in time. For the study of quantum deco- herence, we implement a Lindblad master equation. In spite of its complexity, the model is analytically solvable under fairly reasonable assumptions, and shows that the decoherence influences the entanglement, the sudden death, and the revivals in time.

Equations

$$m_1\frac{d^2{z}_1(t)}{d{t}^2}$$

$$m_2\frac{d^2{z}_2(t)}{d{t}^2}$$

$$m_1{\ddot{z}}_1=-(k_1+k_2)z_1+k_2{z}_2,$$

$$m_2{\ddot{z}}_2=-k_2{z}_2+k_2{z}_1,$$

$${\ddot{z}}_1=-{\omega }_1^2{z}_1+{\omega }_{21}^2{z}_2,$$

$${\ddot{z}}_2=-{\omega }_2^2(z_2+z_1),$$

$$z_1=A\cos (\Omega t+\varphi ),$$

$$z_2=B\cos (\Omega t+\varphi ).$$

$$({\omega }_1^2-{\Omega }^2)A-{\omega }_{21}^{2}B=0,$$

$$-{\omega }_2^{2}A+({\omega }_2^2-{\Omega }^2)B=0.$$

$${\Omega }_{1,2}^2=\frac{{\omega }_1^2+{\omega }_2^2\pm \sqrt{({\omega }_1^2-{\omega }_2^2{)}^2+4{\omega }_2^2{\omega }_{21}^2}}{2},$$

$$\displaystyle\biggl{(}\frac{B}{A}\biggr{)}_{\Omega_{1}}=\frac{\omega_{1}^{2}-\omega_{2}^{2}+\sqrt{(\omega_{1}^{2}-\omega_{2}^{2})^{2}+4\omega_{2}^{2}\omega_{21}^{2}}}{2\omega_{21}^{2}}\equiv\kappa_{1},$$

$$\displaystyle\biggl{(}\frac{B}{A}\biggr{)}_{\Omega_{2}}=\frac{\omega_{1}^{2}-\omega_{2}^{2}-\sqrt{(\omega_{1}^{2}-\omega_{2}^{2})^{2}+4\omega_{2}^{2}\omega_{21}^{2}}}{2\omega_{21}^{2}}\equiv\kappa_{2},$$

$$z_1$$

$$z_2$$

$$z_1=x+y,$$

$$z_2={\kappa }_{1}x+{\kappa }_{2}y,$$

$${\Omega }_{1,2}^2=\frac{k_1+k_2({\kappa }_{1,2}-1{)}^2}{m_1+m_2{\kappa }_{1,2}^2}.$$

$$\mathcal{H}$$

$$\mathcal{H}=\frac{1}{2M_1}{P}_x^2+\frac{1}{2}{M}_1{\Omega }_1^2{x}^2+\frac{1}{2M_2}{P}_y^2+\frac{1}{2}{M}_2{\Omega }_2^2{y}^2.$$

$$\mathcal{H}=\hslash {\Omega }_1({\hat{a}}_1^{†}{\hat{a}}_1+\frac{1}{2})+\hslash {\Omega }_2({\hat{a}}_2^{†}{\hat{a}}_2+\frac{1}{2}),$$

$${\mathcal{H}}_{nr}=\hslash {\Omega }_1{\hat{a}}_1^{†}{\hat{a}}_1+\hslash {\Omega }_2{\hat{a}}_2^{†}{\hat{a}}_2.$$

$${\mathcal{H}}_{sr}=\hslash ({\lambda }_1({\hat{a}}_1+{\hat{a}}_1^{†})+{\lambda }_2({\hat{a}}_2+{\hat{a}}_2^{†}))S_z,$$

$$H_{NV}=\sum_{i=\pm}\biggl{(}-\hbar\delta_{i}|i\rangle\langle i|+\frac{\hbar\Omega_{i}}{2}\bigl{(}|0\rangle\langle i|+(|i\rangle\langle 0|\bigr{)}\biggr{)},$$

$${\mathcal{H}}_{\mathrm{NV}}=\hslash {\omega }_g|g⟩⟨g|+\hslash {\omega }_e|e⟩⟨e|+\hslash {\omega }_d|d_1⟩⟨d_1|.$$

$${\mathcal{H}}_s={\mathcal{H}}_{nr}+{\mathcal{H}}_{NV}+{\mathcal{H}}_{sr},$$

$${\mathcal{H}}_{NV}=\hslash {\omega }_g{R}_{11}+\hslash {\omega }_d{R}_{22}+\hslash {\omega }_e{R}_{33},$$

$${\mathcal{H}}_{sr}=\hslash ({\lambda }_g({\hat{a}}_1{R}_{21}+{\hat{a}}_1^{†}{R}_{12})+{\lambda }_e({\hat{a}}_2{R}_{32}+{\hat{a}}_2^{†}{R}_{23})).$$

$$\mathcal{H}={\mathcal{H}}_0+V_0,$$

$${\mathcal{H}}_0$$

$${\mathcal{H}}_{NV}^i$$

$$V_0$$

$${\mathcal{H}}_{sr}^i$$

$${\mathcal{H}}_{\mathrm{eff}}=\frac{i}{2}{\int }_{-\infty }^{0}d{t}^{\prime }[V_0(t^{\prime }),V_0(0)]$$

$$V_0(t)=e^{-i{\mathcal{H}}_{0}t}{V}_0(0)e^{i{\mathcal{H}}_{0}t}.$$

$${\hat{\mathcal{H}}}_{\mathrm{eff}}={\hat{\mathcal{H}}}_0+\hat{V},$$

Full text

I Introduction

Coherent steering of the dynamics of quantum systems has always been a subject of intense research due to its wide spread applications in devices and setups where quantum coherence is relevant. Notable examples are nano-electromechanical devices r1 ; r2 , quantum resonatorsr3 ; r4 ; r5 , two- and three-level quantum systems r6 ; r7 ; r8 ; r9 ; r10 , and confined cold atoms r11 ; r12 ; r13 ; r14 ; r15 . In recent years, much attention has been paid to nitrogen-vacancy (NV) impurities in diamonds and to the possibilities of their usage in quantum computing and entanglement r16 ; r17 ; r18 ; r19 ; r20 ; r21 ; r22 ; r23 ; r24 ; r25 ; r26 ; r26-1 . These defects can be considered as a three-level system or spin triplets ($S=1$). The advantage of NV-based systems is that they have a long coherence time even at room temperatures and they are easily manipulated by light, electric or magnetic fields.

Recent developments in fabrication technologies of nanomechanical resonators made it possible to implement the coherent connection between the motion of the magnetized resonator tip and the individual point defect centers in diamond r19 . The resonators also allow to link two NV centers of diamonds. r20 ; r26 Previous works studied mainly the entanglement between two NV centers, which are connected by a resonator in vibration. The resonator is assumed to have only one fundamental frequency of vibration and only two states of the NV center spin triplets are entangled with the resonator. Modern technologies allow producing nanomechanical resonators with two or several fundamental frequencies of vibrations r27 ; r28 ; r30 . The dual-mode resonator link allows the participation of all the three spin states of NV centers. In this paper, we study the quantum correlation between two remote NV Centers through a dual-mode nanomechanical resonator. Two magnetic tips are attached to both sides of the resonator, which is placed symmetrically between two NV centers.

The physical set-up consists of two cantilevers in series, both excited at their effective resonant frequencies to produce a tip response as a superposition of a low- and a high-frequency oscillation state. Since the resonator is attached to the top of two magnetic tips, the vibration causes a time varying magnetic field. The vibrating resonator tip can be considered as the superposition of two oscillations, so that the time-varying magnetic field is a superposition of two different frequencies. Such a magnetic field connects the two remote NV centers and causes an indirect interaction between them through a nanomechanical resonator. This interaction involve all the three levels of the NV centers. In previous studies, the "dark" superposition state was decoupled. Using the dual-mode resonator allow considering the case when both "bright" and "dark" superposition states are involved in the dynamics.

We will demonstrate that in spite of the complexity, with realistic assumptions, this model is analytically solvable even in the presence of quantum decoherence. The paper is organized as follows: In section II, we study the dynamics of the dual-mode cantilever, in section III we present an exact analytical solution for the case of one and two NV centers. We derive the Hamiltonian of the non-direct interaction between NV centers mediated by the dual-mode cantilever. In section IV, we study the time dependence of the system’s entanglement, and in the last section we study the problem of quantum decoherence using a Lindblad master equation.

II Dynamics of the dual-mode resonator

A prototype model of a nanomechanical resonator connected to two NV centers of diamond is shown in Fig. 1. A small high-frequency cantilever is attached at the end of a larger low-frequency cantilever, which in turn is attached to a base. Two magnetic tips are attached to the small cantilever and two spin-one systems (NV centers) denoted by ‘1’ and ‘2’ are placed at distances $h_{1}$ and $h_{2}$ from the respective magnetic tips.

The model represents a union of two consecutive resonators. A mechanical analogue in terms of springs and masses is shown in Fig.2. The equation of motion of two rods in the mechanical model can be represented by the following equations r28 :

[TABLE]

where $m_{1}$ and $m_{2}$ are the effective masses of the first and the second rod, $k_{1}$ and $k_{2}$ are spring constants (vibrating rods obey Hooke’s law). $z_{1}(t)$ and $z_{2}(t)$ are the coordinates of the first and second rod tips, $Q_{1}$ and $Q_{2}$ are quality factors, $\nu_{1-0}$ and $\nu_{2-0}$ are free-resonant frequencies of the first and the second rod. $F_{ts}[z_{2}(t)]$ is the external force acting on the second rod tip. The system can be excited by flapping the stem base, $z_{0}$ is the location of the first stem base at time $t$.

We will study the dynamical equation Eq. 1 by employing several assumptions: First, let us assume that the base is stopped at $z_{0}=0$ and there is no magnetic interaction between magnetic tips and spins: $F_{ts}[z_{2}(t)]=0$. Before proceeding further, it is useful to estimate the range of realistic parameters for the cantilever r26 . A typical cantilever fabricated out of Si(100) has a Young’s modulus $Y=130$ GPa and a density $\rho=2.33\times 10^{3}$ kg/m3. In the dual-mode hybrid cantilever (see Fig. 1), the length of the two different parts is $L_{1}\approx 15\cdot 10^{3}$nm, $L_{2}\approx 9\cdot 10^{3}$ nm, width $w_{1}=$300 nm, $w_{2}=$200 nm and thickness $d_{1}=30$ nm, $d_{2}=20$nm, masses $m_{1}=3.5\times 10^{-16}$ kg and $m_{2}=10^{-16}$ kg. Using equations $2\pi\nu_{1-0}=1.8^{2}\sqrt{\frac{Yd_{1}^{2}}{3\rho}}\frac{1}{2L_{1}^{2}}$ and $2\pi\nu_{2-0}=1.8^{2}\sqrt{\frac{Yd_{2}^{2}}{\rho}}\frac{1}{2L_{2}^{2}}$, we therefore obtain the frequencies $\pi\nu_{1-0}=0.54$ MHz, $\pi\nu_{2-0}=1.7$ MHz. The estimated spring constants are: $k_{1}=\frac{Y\times w_{1}\times d_{1}^{3}}{L_{1}^{3}}=3\cdot 10^{-4}$ kg/s2, $k_{2}=\frac{Y\times w_{2}\times d_{2}^{3}}{L_{2}^{3}}=2\cdot 10^{-4}$ kg/s2. The quality factor is large r19 , $Q_{1,2}\approx 10^{5}$, and the frequencies are $\omega_{1}=\sqrt{\frac{k_{1}+k_{2}}{m_{1}}}=1.2$MHz, $\omega_{2}=\sqrt{\frac{k_{2}}{m_{2}}}=1.4$MHz. Thus, by designing geometrical characteristics of the hybrid dual cantilever, one can easily achieve the following conditions: $\frac{2\pi\nu_{1-0}}{Q_{1}}\ll\omega_{1}$, $\frac{2\pi\nu_{2-0}}{Q_{2}}\ll\omega_{2}$. Given these estimates we can rewrite Eq. (1) as

[TABLE]

which can further be simplified to

[TABLE]

where $\omega_{1}=\sqrt{(k_{1}+k_{2})/m_{1}}$, $\omega_{21}=\sqrt{k_{2}/m_{1}}$ and $\omega_{2}=\sqrt{k_{2}/m_{2}}$.

The general solution of the system of equations Eq. 3 has the following form:

[TABLE]

Here $A$, $B$ and $\Omega$ are constants and can be calculated by substituting Eq. 4 into Eq. 3 as

[TABLE]

The non-trivial solution for $\Omega$, obtained by solving the above equations, is

[TABLE]

where $\Omega_{1}$ and $\Omega_{2}$ are the system’s own frequencies. They differ from the resonance frequencies $\omega_{1}$ and $\omega_{2}$ of the system because of the coupling term $\omega_{2}\omega_{21}$ in Eq. 6.

The ratio of amplitudes $A$ and $B$ at $\Omega=\Omega_{1}$ can be written as

[TABLE]

where $\kappa_{1}$ is completely determined by the parameters of the system and does not depend on the initial conditions. It is called the distribution coefficient amplitude for $\Omega_{1}$. Similarly we find the distribution coefficient for $\Omega_{2}$ as

[TABLE]

and the general solution of the coupled equations Eq. 4 can be written as:

[TABLE]

This means that the oscillations of each stem tip are the superposition of two harmonic oscillations with frequencies of vibration given by Eq. 6.

Let us define two normal coordinates $x$ and $y$ such that $x=A_{1}\cos(\Omega_{1}t+\phi_{1})$ and $y=A_{2}\cos(\Omega_{2}t+\phi_{2})$. In terms of these normal coordinates, the solution takes on a simpler form,

[TABLE]

and $\Omega_{1,2}$ can be written as

[TABLE]

The Hamiltonian of the system can be written as

[TABLE]

Introducing the masses $M_{1,2}=m_{1}+m_{2}\kappa_{1,2}^{2}$ and momenta $P_{x}=M_{1}\dot{x}$, $P_{y}=M_{2}\dot{y}$, the Hamiltonian can be expressed as

[TABLE]

In this way, the Hamiltonian of the nanomechanical resonator can be considered as the sum of two noninteracting harmonic oscillators. The quantum-mechanical Hamiltonian operator can be written by replacing $x$, $y$, $P_{x}$ and $P_{y}$ by quantum-mechanical operators $\hat{x}$, $\hat{y}$, $\hat{P}_{x}$ and $\hat{P}_{y}$, respectively. A further transformation into creation and annihilation operatorsr30 leads to

[TABLE]

where $\hat{a}_{1}$ $\hat{a}_{2}$ have the usual meaning. If we neglect the zero point vibrations, the nanomechanical resonator may be described by the Hamiltonian:

[TABLE]

The magnetic tip movement generates a magnetic field $|\vec{B}_{\rm tip}|\simeq G_{m}\hat{z}$, where $G_{m}$ is the magnetic field gradient and $\hat{z}$ is the tip location operator. Since the tip peak vibration is a superposition of two harmonic waves, the $\hat{z}$ operator can be expressed as $\hat{z}=a_{10}(\hat{a}_{1}+\hat{a}_{1}^{\dagger})+a_{20}(\hat{a}_{2}+\hat{a}_{2}^{\dagger})$ and the interaction between the magnetic tip and the spin is expressed as

[TABLE]

where $\lambda_{i}=g_{s}\mu_{B}G_{m}a_{i0},~{}~{}i=1,2$, the gyromagnetic ratio $g_{s}\approx 2$, and $a_{i0}=\sqrt{\hbar/2M_{i}\Omega_{i}}$ are the amplitude zero-point fluctuations. The system is placed in an external magnetic field $B_{0}$ along the $z$ direction. If the conditions $\frac{2\pi\nu_{1-0}}{Q_{1}}\ll\omega_{1}$, $\frac{2\pi\nu_{2-0}}{Q_{2}}\ll\omega_{2}$ do not hold true, instead of Eq. (2) one needs to solve directly Eq. (1). The difference between Eq. (1) and Eq. (2) is the weak damping term and this damping of the cantilever’s oscillation leads to quantum decoherence. Usually, NV centers are characterized by a low decoherence rate. However, in order to take into account environmental effects, in section V we utilize a Lindblad master equation for NV centres and address the problem of quantum decoherence more completely.

III NV centers coupled to a dual-mode cantilever

In this section, we present an analytical solution for the single and two NV centers coupled to the dual-mode cantilever. We derive the Hamiltonian of non-direct interaction between NV centers mediated via the dual-mode cantilever.

III.1 One-spin case

The total spin of a nitrogen vacancy center in diamond is $S=1$, with the three spin sub-states, $m_{s}=-1,\ 0$ and $+1$ being separated from each other by the frequency $\omega_{0}/2\pi\simeq 2.88$ GHzr19 . An added external magnetic field $B_{0}$ shifts $|-1\rangle$ and $|+1\rangle$ states (Zeeman shift) proportionally to $B_{0}S_{z}$.

The spin part of the Hamiltonian $\mathcal{H}_{NV}$ reads

[TABLE]

where $\delta_{\pm}$ and $\Omega_{\pm}$ denote the detunings and the Rabi frequences of the two transitions. We consider here the case when the Rabi frequencies $\Omega_{\pm}$ are not equal, $(\Omega_{\pm}=\Omega_{0}\pm\Delta\Omega(B_{0}))$, where $\Omega_{0}$ is the Rabi transition frequency in zero magnetic field and $\Delta\Omega(B_{0})=\mu_{B}B_{0}$. In addition, we consider the case when the detunings are equalr19 $\delta_{+}=\delta_{-}=\delta$. The schematics of the transition is shown in Fig. 3(a). With this assumption, we can calculate the eigenfunctions of the Hamiltonian as $|d_{1}\rangle=\frac{1}{\Omega\sqrt{2}}(\Omega_{+}|-1\rangle-\Omega_{-}|+1\rangle)$, $|e\rangle=\cos\theta|b_{1}\rangle+\sin\theta|0\rangle$ and $|g\rangle=\cos\theta|0\rangle-\sin\theta|b_{1}\rangle$, where $|b_{1}\rangle=\frac{1}{\Omega\sqrt{2}}(\Omega_{-}|-1\rangle+\Omega_{+}|+1\rangle)$ and $\tan(2\theta)=-\sqrt{2}\Omega/\delta$, $2\Omega^{2}=\Omega_{+}^{2}+\Omega_{-}^{2}$. It should be noted that states $|b_{1}\rangle$ and $|d_{1}\rangle$ are superpositions of "bright" $|b\rangle=(|-1\rangle+|+1\rangle)/\sqrt{2}$ and "dark" $|d\rangle=(|-1\rangle-|+1\rangle)/\sqrt{2}$ statesr19 : $|b_{1}\rangle=\frac{\Omega_{0}|b\rangle-\Delta\Omega(B)|d\rangle}{\Omega}$, $|d_{1}\rangle=\frac{\Omega_{0}|d\rangle+\Delta\Omega(B)|b\rangle}{\Omega}$. The corresponding eigenfrequencies are $\omega_{d}=-\delta$, $\omega_{e}=\frac{-\delta+\sqrt{\delta^{2}+2\Omega^{2}}}{2}$ and $\omega_{g}=\frac{-\delta-\sqrt{\delta^{2}+2\Omega^{2}}}{2}$. Note that $\omega_{g}<\omega_{d}<\omega_{e}$.

In the diagonal basis, the Hamiltonian assumes the following form:

[TABLE]

In this diagonal basis, the z-component of the spin is S_{z}=-\frac{\Omega_{+}\Omega_{-}}{\Omega^{2}}\big{(}\sin\theta(|d_{1}\rangle\langle g|+|g\rangle\langle d_{1}|)-\cos\theta(|d_{1}\rangle\langle e|+|e\rangle\langle d_{1}|)\big{)}. Since $\omega_{g}<\omega_{d}<\omega_{e}$ we can relabel the energy levels $|g\rangle$, $|d_{1}\rangle$ and $|e\rangle$ as $|1\rangle$, $|2\rangle$, and $|3\rangle$, respectively. We also introduce the operators $R_{ij}=|i\rangle\langle j|$, obeying $R_{ij}R_{kl}=R_{il}\delta_{kj}$ where $i,j=1,2,3$. After implementing the operators $R_{ij}$ and using Eqs. (15),(16),(17), and (18), we rewrite the total Hamiltonian of the system in a more convenient form

[TABLE]

Here $\mathcal{H}_{nr}$ is the Hamiltonian of the nanomechanical resonator (see Eq.(15), $\mathcal{H}_{NV}$ is the Hamiltonian of the NV-center and $\mathcal{H}_{sr}$ is the Hamiltonian of interaction between them, obtained in the rotating wave approximation: $\lambda_{g}=-\frac{\lambda_{1}\Omega_{-}\Omega_{+}\sin\theta}{\Omega^{2}}$ and $\lambda_{e}=\frac{\lambda_{2}\Omega_{-}\Omega_{+}\cos\theta}{\Omega^{2}}$. The above Hamiltonian represents a generalized Jaynes–Cummings modelr7 ; r34 . The levels and transitions are shown in Fig. 3(b).

III.2 Two-spin case

Consider the case of two identical spins (NV centers) placed on the different sides of a resonator and coupled by it. Let us assume that the two spins have the same Rabi transition frequency, $\Omega_{\pm}^{i}=\Omega_{\pm}$, detuning $\delta_{\pm}^{i}=\delta$, and interaction constants with the resonator $\lambda_{g,e}^{i}=\lambda_{g,e}$, where the superscripts $i=1,2$ represent the spins. The Hamiltonian for the two spin case can be written as

[TABLE]

where

[TABLE]

and

[TABLE]

The superscripts $i=1,2$ stand for the first and the second spin, respectively. The interaction of the spins with the magnetic tips of the resonator results in an indirect coupling between the spins.

The Hamiltonian of the indirect interaction between the NV spins can be evaluated using the Fröhlich methodr35 (set $\hbar=1$):

[TABLE]

and

[TABLE]

In the rotating-wave approximation, we write Eqs. (III.2)- (24) as

[TABLE]

where

[TABLE]

and

[TABLE]

In the above $\hat{n}_{i}=\hat{a}^{\dagger}_{i}\hat{a}_{i}$ is the mean photon number operator, $\alpha=\lambda_{g}^{2}/\Delta_{1}$, $\beta=\lambda_{e}^{2}/\Delta_{2}$, $\Delta_{1}=\Omega_{1}-\omega_{e}$, $\Delta_{2}=\Omega_{2}+\omega_{g}$. In the next section we demonstrate the phenomenon of the early-stage disentanglement in the system (the entanglement sudden death). Importantly for a particular ratio between the parameters $\alpha=\beta$, meaning that \big{(}\Omega_{1}+\big{(}\sqrt{2\Omega^{2}+\delta^{2}}-\delta\big{)}/2=\tan^{2}\big{(}\theta\big{)}\big{(}\Omega_{2}-\big{(}\sqrt{2\Omega^{2}+\delta^{2}}+\delta\big{)}/2 the entanglement in the system acquires a persistent value. Here $\Omega_{1,2}$ are the frequencies of the dual mode cantilever and $\Omega$ is the Rabi frequency. Together with the detuning $\delta$ and the angle $\tan(2\theta)=-\sqrt{2}\Omega/\delta$ all these parameters can be controlled by the size of the cantilever and the applied magnetic field.

IV Measure of the entanglement: negativity

The amount of entanglement shared between the two spins can be measured by the "negativity" which is defined as

[TABLE]

where $\chi_{i}$’s are the eigenvalues of the partial transposed density matrix with respect to the spin1. If $\chi_{i}>0$ then $|\chi_{i}|-\chi_{i}=0$, however, if $\chi_{i}<0$, then $|\chi_{i}|-\chi_{i}=-2\chi_{i}$ and the negativity

[TABLE]

The time dependence of the negativity can be calculated using the Schrödinger equation (with $\hbar=1$):

[TABLE]

where the initial state can be considered in its most general form

[TABLE]

Here $|i\rangle|j\rangle\equiv|i\rangle\otimes|j\rangle$, the kets corresponding to the first and the second spin, and $a_{1}$, $a_{2}$ $\cdots$ $a_{9}$ are coefficients to be determined.

Using the Schrödinger equation, we can write:

[TABLE]

The equations for $a_{2}$ and $a_{4}$ as well as for $a_{6}$ and $a_{8}$ show the same oscillations, which is obvious as the pair of coefficients $a_{2}$ and $a_{4}$ relate states $|1\rangle|2\rangle$ and $|2\rangle|1\rangle$ such that $a_{2}|1\rangle|2\rangle\Leftrightarrow a_{4}|2\rangle|1\rangle$ and the pair of coefficients $a_{6}$ and $a_{8}$ relate states $|2\rangle|3\rangle$ and $|3\rangle|2\rangle$ such that $a_{6}|2\rangle|3\rangle\Leftrightarrow a_{8}|3\rangle|2\rangle$. The solutions of equations (32) are

[TABLE]

Using these coefficients, we can write the density matrix $\rho$ with elements $\rho_{nm}=a_{n}a_{m}^{*}$. For this density matrix, we can easily calculate the partial transpose matrix with respect to the spin1. The elements of the partial transpose matrix are connected to the elements $\rho_{nm}$ of the density matrix as follows:

[TABLE]

The eigenvalues of the matrix $\rho^{T_{\rm spin_{1}}}$ can be calculated and the negative eigenvalues used in Eq.(29) to calculate the time dependence of negativity. The value of the Negativity is highly dependent on the choice of initial conditions, i.e., the choice of $C_{i},i=1\cdots 9$.

Let us consider the initial state of the system to be

[TABLE]

with $a_{2}(0)=a_{6}(0)=1/\sqrt{2}$ and $a_{1}(0)=a_{3}(0)=a_{4}(0)=a_{5}(0)=a_{7}(0)=a_{8}(0)=a_{9}(0)=0$. For this choice of initial state, we can calculate the non-zero elements of the density and the eigenvalues of its partial transpose with respect to spin1. The partial transpose matrix has only one negative eigenvalue given as $\chi=-\sqrt{\rho_{24}\rho_{42}+\rho_{26}\rho_{62}+\rho_{48}\rho_{84}+\rho_{68}\rho_{86}}$ and the negativity from Eq. (29) comes out to be

[TABLE]

It should be noted that if $\alpha=\beta$ then $N(\rho)=0.5$. We obtain the same result if we take the initial condition $a_{4}(0)=a_{8}(0)=1/\sqrt{2}$ and $a_{1}(0)=a_{2}(0)=a_{3}(0)=a_{5}(0)=a_{6}(0)=a_{7}(0)=a_{9}(0)=0$. In Fig. 4, the time dependence of the negativity is shown for a choice of the initial conditions $a_{2}(0)=a_{6}(0)=1/\sqrt{2}$ or $a_{4}(0)=a_{8}(0)=1/\sqrt{2}$ and the rest of the coefficients being zero. We see that the negativity oscillates in time and the oscillations are more rapid as the difference between $\alpha$ and $\beta$ increases. In this case, the value of the negativity is bounded between [math] and $0.5$.

It is interesting to find the instances where the negativity vanishes resulting in zero entanglement between the spins. From Eq.(36) we infer

[TABLE]

It is easy to see that the above equation holds true when

[TABLE]

The solution of these equations delivers the condition for the choice of $\alpha$ and $\beta$ for a zero negativity as

[TABLE]

where $n$ and $m$ are integers. $n+m+1$ and $n-m$ have different parity. If the first is even, then the second is odd and vice versa. Thus, negativity becomes zero when the ratio $\alpha/\beta$ is equal to a fraction whose denominator is even and the numerator is odd or vice versa. In Fig.4, we see that $N(\rho)=0$ when $\frac{\alpha}{\beta}$ is equal to $\frac{1}{2}$ or $\frac{1}{4}$ and $t=\frac{\pi}{2}+\pi k,k=0,1,2...$

It should be noted here that if we take the initial conditions as $a_{2}(0)=a_{4}(0)=a_{6}(0)=a_{8}(0)=1/2$ and the rest of the $a_{i}$ to be zero, the negativity will be $N(\rho)=1/2$. For an initial state with all the $a_{i}$ being the same and equal to $1/3$, the eigenvalues of the partial transpose of the density matrix cannot be calculated analytically except for the case $\alpha=\beta$. For this case, the negativity is

[TABLE]

which can be expressed in the form

[TABLE]

where $k$ is an integer. The maximum of the negativity is $N_{\rm max}(\rho)=\frac{8}{9}$ at $t=-\frac{\pi}{2\alpha}+\frac{2\pi k}{\alpha}$ We find that the entanglement is zero during the time interval $\frac{2\pi k}{\alpha}\leq t\leq\frac{(2k+1)\pi}{\alpha}$, however, there is a partial entanglement sharing between the spins during the interval $\frac{(2k+1)\pi}{\alpha}<t<\frac{2(k+1)\pi}{\alpha}$.

V Lindblad equation

We consider two three-level atoms (NV centers), with the energy levels of each NV center as in Fig. (3). We have discussed above these two NV centers interacting with each other indirectly via a nanomechanical resonator, a scheme described in ref. 27. We now introduce damping by assuming that the excited levels $|2\rangle$ and $|3\rangle$ decay to the ground state $|1\rangle$, and a direct transition between the excited levels is not allowed. The time evolution of such a system is given by the master equation (Lindblad Equation) r38 ; r39 :

[TABLE]

where $\hat{V}$ is the (27), and $L\rho$ is the damping termr38 ; r39 :

[TABLE]

The spontaneous emission of atoms 1 and 2 from their excited states $|3\rangle$ to the ground states $|1\rangle$ is described by the spontaneous decay rate $\gamma_{e}$, similarly $\gamma_{d}$ is the spontaneous decay rate of excited $|2\rangle$ to the ground $|1\rangle$. For the derivation of the master equation, we choose the basis:

[TABLE]

and the $R_{ij}=|i\rangle\langle j|$ operators in this basis have the following forms:

[TABLE]

and $R_{ij}^{1}=R_{ij}\otimes I_{3}$, $R_{ij}^{2}=I_{3}\otimes R_{ij}$, with $I_{3}$ the three-dimensional unit matrix. The most general solution of the master equation (41) for an arbitrary initial state is given in the appendix. In fig. 5 we have shown the effect of damping on the entanglement between the two spins by varying all the parameters of interest. Fig. 5(a) shows the dynamics of entanglement for different values of $\gamma_{d}$ with the other parameters fixed as $\alpha=\beta=\gamma_{e}=1$. We see that the entanglement decays more sharply due to the increase in decoherence parameter $\gamma_{d}$. The same can be seen by varying the other decoherence parameter $\gamma_{e}$ in fig. 5(b). As we see, due to the decoherence, entanglement’s sudden death is smeared out.

VI Conclusions

One of the main challenges for the NV spin-based nanomechanical resonator is to achieve a high degree of entanglement between NV spins. The success of this proposal naturally depends on the strength of the coupling between NV spins. On the other hand, single NV spins $S=1$ are a triplet state with two characteristic transition frequencies between the three triplet states. The realization of the controlled transitions between all the three levels requires two-frequency nanomechanical resonator with a special type of dual frequency cantilever. The frequency characteristics of the cantilever depends on the particular choice of the cantilever’s geometry. A proper setup of the dual cantilever helps tuning the oscillation frequencies to the resonance frequencies. The dual cantilever is supplemented by the magnetic tips. Thus, the oscillation of the cantilever leads to the indirect interaction between the NV spins, while the direct interaction is small. Hence, the entanglement between the NV spins is fueled by the dynamics of the cantilever. In spite of the complexity of the system (two coupled NV centers and a dual frequency cantilever), the model is analytically solvable. In order to study decoherence, we utilized the Lindblad master equation. An exact analytical solution of the Lindblad equation shows the influence of the decoherence. A prominent effect of the entanglement’s sudden death is smeared out by decoherence.

Acknowledgements.

We acknowledge financial support from by Deutsche Forschungsgemeinschaft SFB 762.

VII Appendix

The most general solution of the Lindblad equation (41) for an arbitrary initial state is given as ($\rho_{ij}$ on the right-hand side are the initial values):

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]