Theory of nuclear spin dephasing and relaxation by optically illuminated nitrogen-vancy center

Dephasing and relaxation of the nuclear spins coupled to the nitrogen-vacancy (NV) center during optical initialization and readout is an important issue for various applications of this hybrid quantum register. Here we present both an analytical description and a numerical simulation for this process, which agree reasonably with the experimental measurements. For the NV center under cyclic optical transition, our analytical formula not only provide a clear physics picture, but also allows controlling the nuclear spin dissipation by tuning an external magnetic field. For more general optical pumping, our analytical formula reveals significant contribution to the nuclear spin dissipation due to electron random hopping into/out of the $m=0$ (or $m=\pm1$) subspace. This contribution is not suppressed even under saturated optical pumping and/or vanishing magnetic field, thus providing a possible solution to the puzzling observation of nuclear spin dephasing in zero perpendicular magnetic field [M. V. G. Dutt \textit{et al}., Science \textbf{316}, 1312 (2007)]. It also implies that enhancing the degree of spin polarization of the nitrogen-vacancy center can reduce the effect of optical induced nuclear spin dissipation.


I. INTRODUCTION
Diamond nitrogen-vacancy (NV) center 1 is a leading platform for various quantum technologies such as quantum communication, quantum computation, and nanoscale sensing [2][3][4][5][6][7][8] . The electronic spin of the NV center and a few surrounding nuclear spins form a hybrid quantum register [9][10][11] . Important advantages of this solid-state quantum register include the long electron and nuclear spin coherence time 12 , the capability of high-fidelity initialization, coherent manipulation, and projective readout of the electronic/nuclear spins 5,13 and even the entire quantum register 11,14,15 by optical and microwave (or radio frequency) illumination. However, during the optical illumination for initialization and readout 2,8,14,[16][17][18] , the dissipative spontaneous emission and non-radiative decay of the NV electron generates substantial noise on the nuclear spin qubits through the hyperfine interaction (HFI), which may significantly degrade the control precesion. This motivates widespread interest in using the NV center electron to engineer the nuclear spin dissipation, including pure dephasing and relaxation 16,17 .
In the past few years, the optically induced nuclear spin dissipation has been investigated in many works 16,17,[19][20][21] . Generally, the nuclear spin dissipation originates from the random fluctuation of the NV electron under optical illumination, which falls into two categories: one involving the flip of the NV electron spin and the other does not. The former is usually strongly suppressed by the large energy splitting of the NV electron unless the NV electron is tuned to the ground state or excited state level anticrossing 19,21 . The latter is energetically more favorable and dominates the nuclear spin dissipation in many situations, as confirmed by a series of experiments 2,12,17 . The theoretical investigation of this latter mechanism has been carried out in the framework of a phenomenological spinfluctuator model 16 . This work gives an intuitive understanding for the optically induced nuclear spin dissipation: the generation of a rapidly fluctuating effective magnetic field on the nuclear spins by the optically induced random hopping of the electron between different states. When the hopping is sufficiently fast and hence the noise correlation time is sufficiently short, the nuclear spin dissipation could be suppressed 16 in a way similar to the motional narrowing effect in NMR spectroscopy in liquids. This effect has been successfully used to significantly increase the nuclear spin coherence time 12 . Despite these remarkable success, this spin-fluctuator model still suffers from two drawbacks. First, its analytical form is qualitative, while obtaining quantitative results require numerical simulations. This not only complicates the calculation, but also smears the underlying physics picture. Second, the various parameters in this model are phenomenological, i.e., they are not directly related to the physical parameters of the NV center, but instead are obtained from fitting the experimental data. This precludes a straightforward guidance on controlling the nuclear spin dissipation by tuning various experimental parameters.
To bridge this gap between experimental observation and theoretical understanding, we present a microscopic and analytical theory on the nuclear spin dephasing and relaxation by an optically illuminated NV center at room temperature. In addition to performing numerical simulation of the coupled NV-nuclear spin evolution, we further derive analytically a closed Lindblad master equation for the nuclear spin by adiabatically eliminating the fast electron spin dynamics in the Born-Markovian approximation. We begin with the simplest case in which a single cyclic transitions (e.g., between the ground and excited m = 0 states) of the NV center is optically driven. Our analytical expressions for the nuclear spin dephasing and relaxation provide a quantitative description and a physically transparent interpretation that substantiates the previous analytical (but qualitative) and numerical results 16 .
They also demonstrate the possibility to control the nuclear spin dissipation by tuning the magnetic field 16 . Next we consider general optical illumination of the NV center incorporating finite inter-system crossing between m = 0 and m = ±1 subspaces. Our numerical results agree well with the experimental measurements 17 . Our analytical results shows that the random hopping between the m = 0 (or m = ±1) triplet states and the metastable singlet of the NV center could significantly contribute to nuclear spin dissipation. This contribution is not suppressed under saturated optical pumping and is nearly independent of the magnetic field. This provides a possible solution to the puzzling observation of nuclear spin dephasing in zero magnetic field 2 . An analytical formula for the nuclear spin dissipation in terms of the HFI tensors also allows us to measure the HFI tensor for the excited electron state, which is usually smeared by the short electron spontaneous emission lifetime.

A. Model
To begin with, we present a microscopic theory for the decoherence of an arbitrary nuclear spinÎ (e.g., 13 C, 15 N, or 14 N) by the electron of the NV center undergoing optically induced cyclic transition |g ↔ |e , e.g., |g = |0 and |e = |E y in the widely used setup for single-shot readout 14,22 . In the rotating frame, the electron dynamics is governed by the Liouville superoperator L e (·) ≡ −i[Ĥ e , (·)] + α γ α D(L α )(·), wherê H e = ∆σ e,e + Ω R 2 (σ e,g + h.c) is the electron Hamiltonian,σ i, j ≡ |i j|, ∆ is the detuning of the optical pumping, and γ α is the rate of the αth dissipation processL α in the Lindblad form D(L α )(·) ≡L α (·)L † α − {L † αLα , (·)}/2. Here we include the spontaneous emission L =σ g,e from |e to |g with rate γ 1 ≈ 1/(12 ns) and the pure dephasingL =σ e,e of the excited state |e with rate γ ϕ , which has a strong temperature dependence, from a few tens of MHz at low temperature up to 10 7 MHz at room temperature 23,24 . Including the electron-nuclear HFI (Ŝ g ·A g +Ŝ e ·A e )·Î ≡F·Î and the nuclear spin Zeeman term γ N B ·Î (γ N = −10.705 kHz/mT is the 13 C nuclear gyromagnetic ratio) under a magnetic field B, the electron-nuclear coupled system obeyṡ in the rotating frame of the pumping laser. There are two contributions to the nuclear spin dissipation. One involves the flip of the electron spin and hence is strongly suppressed by the large electron-nuclear energy mismatch away from the NV center ground state and excited state anticrossings. The other does not flip the electron spin and hence is energetically favorable in most situations. In our analytical derivation, we neglect the former contribution by dropping the off-diagonal electron spin flip terms inF and only keep the diagonal part:F ≈σ g,g ω g +σ e,e ω e , where ω g = g|Ŝ g |g · A g and ω e = e|Ŝ e |e · A e . The second term of Eq. (1) describes the precession of the nuclear spin with angular frequency γ N B + ω g and γ N B + ω e , respectively, conditioned on the electron state being |g and |e . When ω g ω e , the optically induced random hopping of the electron between |g and |e gives rise to random fluctuation of the nuclear spin precession frequency and hence nuclear spin dissipation: the fluctuation of the precession frequency orientation (magnitude) leads to nuclear spin relaxation (pure dephasing) 16 . Below we derive analytical a closed equation of motion of the nuclear spin to describe these effects.

B. Lindblad master equation for nuclear spin
The time scale for the optically pumped two-level NV center to reach its steady state is ∼ τ NV ≡ 1/(2R + γ 1 ) < 12 ns, where R = 2π(Ω R /2) 2 δ ((γ 1 +γ ϕ )/2) (∆) is the optical pumping rate from |g to |e and δ (γ) (x) = (γ/π)/(x 2 + γ 2 ) is the broadened δ-function. When τ NV is much shorter than the time scale of the nuclear spin dissipation, we can regard the NV center as always in its steady stateP as determined by L eP = 0, e.g., the steady state population on |e and |g are P e = R/(2R + γ 1 ) and P g = 1 − P e , respectively. Then we treat the dissipative NV center as a Markovian bath 25 and use Born-Markovian approximation to derive a Lindblad master equation for the reduced density matrix of the nuclear spinp(t) ≡ Tr eρ (t) (see appendix A for details): (2) whereω ≡ γ N B + P g ω g + P e ω e is the average precession frequency that defines the nuclear spin quantization axis e Z ≡ ω/|ω|. The last three terms describe the nuclear spin dissipation in the tilted cartesian coordinate e X = e x sin ϕ − e y cos ϕ, (3a) e Y = cos ϕ cos θe x + sin ϕ cos θe y − sin θe z , (3b) e Z =ω/|ω| = sin θ cos ϕe x + sin θ sin ϕe y + cos θe z , (3c) where θ (ϕ) is the polar (azimuth) anlge ofω in the conventional coordinate (e x , e y , e z ) with e z along the N-V symmetry axis. The nuclear spin dissipation include pure dephasing [the second term of Eq. (2)] due to the fluctuation ofF Z and relaxation [the last two terms of Eq. (2), withÎ ± ≡Î X ± iÎ Y ] due to the fluctuation ofF ± ≡F X ± iF Y . Typically the nuclear spin level splitting |ω| ≪ γ 1 , γ ϕ , so we obtain where O ⊥ ≡ O X e X + O Y e Y is the component perpendicular to the nuclear spin quantization axis, T = 1/R + 1/(γ 1 + R) is the duration of one electron hopping cycle (excitation time 1/R and de-excitation time 1/(γ 1 + R)), and is the uncertainty of the electron dwell time in the excited state in each hopping cycle. Here the last step of Eq. (5) holds at room temperature, where γ ϕ ∼ 10 7 MHz is much larger than typical γ 1 , R, and ∆. Equation (4) shows that nuclear spin dissipation vanishes when ω g = ω e , simply because in this case the nuclear spin precession frequency is not randomized by the optically induced electron hopping.

C. Physical picture
Equations (2)-(5) not only provide an quantitative and analytical description for the dissipative nuclear spin dynamics due to an optically pumped NV center, but also have a physically transparent interpretation that substantiates the previous analytical (but qualitative) and numerical results 16 . For example, the pure dephasing rate in Eq. (4a) is directly connected to the nuclear spin phase diffusion process by the optically induced random hopping of the electron between the ground state |g and the excited state |e 16 . To clearly see this, lets consider the phase accumulation of the nuclear spin during an interval [0, t]. Suppose that during this interval, the electron undergoes N hopping cycles, and that during the kth cycle, the electron stays in |g for an interval τ k , so the total dwell times in |g and |e are τ = N k=1 τ k and t − τ, respectively, and the nuclear spin accumulates a phase factor e −i(a g +γ N B) Z τ−i(a e +γ N B) Z (t−τ) . For t ≫ T , the number of hopping cycle N ≈ t/T ≫ 1, i.e., τ is the sum of many random variables {τ k }, so τ obeys Gaussian distribution centered at P g t with a standard deviation √ Nτ e , where τ e is the rms fluctuation of each τ k . Averaging the phase factor over this Gaussian distribution gives e −i|ω|t e −Γ ϕ t , where Γ ϕ coincides with Eq. (4a) as long as τ e is given in Eq. (5), e.g., at room temperature, for weak pumping R ≪ γ 1 , the uncertainty τ e ≈ √ 2/γ 1 of the dwell time in |e is dominated by the uncertainty in the spontaneous emission; while for saturated pumping, τ e ≈ 1/( √ 2R) is strongly suppressed by the rapid optically induced transition between |e and |g . The relaxation rate Γ ± in Eq. (4b) can be understood in a similar way.
analytical results Eqs. (4) provide a microscopic basis for the previous model 16 and experimental observations 2,12,26 , e.g., it clearly shows the initial increase of the dissipation rates Γ ϕ , Γ ± ∝ R under weak pumping R ≪ γ 1 and the motional narrowing Γ ϕ , Γ ± ∝ 1/R under saturated pumping R ≫ γ 1 . The former arises from the increase of T with decreasing R under weak pumping, while the latter comes from both the decrease of τ e ∼ 1/R and T ∼ 1/R under saturated pumping. Our analytical formula also demonstrate the possibility 16 to control Γ ϕ and Γ ± by using the magnetic field to tune the nuclear quantization axis e Z ∝ω, e.g., if we tune e Z to be perpendicular (parallel) to ω g − ω e , then we can eliminate nuclear spin pure dephasing (relaxation). Interestingly, the sum rule |ω e − ω g | 2 (6) suggests that reducing Γ ϕ (Γ ± ) inevitably increases Γ ± (Γ ϕ ) and it is impossible to suppress Γ ϕ and Γ ± simultaneously, unless the NV states are tuned such that ω g = ω e .

D. Connection to experimental observations
Equation (2) describes the dissipative evolution of the nuclear spin in the tilted coordinate (e X , e Y , e Z ) with e Z ∝ω. From Eq. (2), we obtain the Bloch equations for the average nuclear spin Î (t) ≡ TrÎp(t), where T 1 = 1/(Γ + + Γ − ) and T 2 = 1/(Γ ϕ + (Γ + + Γ − )/2). Then the sum rule in Eq. (6) implies 1/T 2 + 1/(2T 1 ) ∝ |ω e − ω g | 2 , i.e., tuning the magnetic field can prolong T 1 time (T 2 time) at the cost of reducing T 2 time (T 1 time). The above Bloch equations have simple solutions Î Z (t) = Î Z (0) e −t/T 1 and Î + (t) = e i|ω|t e −t/T 2 Î + (0) . However, nuclear spin initialization and measurement are usually performed in the conventional coordinate (e x , e y , e z ) with e z along the N-V axis, so T 1 and T 2 will be mixed in the observed signals. For example, Dutt et al. 2 initialize a strongly coupled 13 C nuclear spin-1/2 (hereafter referred to as 13 C b , according to the notation of Gali 27 ) into the eigenstate (| ↑ + | ↓ )/ √ 2 ofÎ x , let it evolve freely for an interval τ, and then measureÎ x through a π/2 pulse e −iπÎ y /2 followed by a fluorescence readout ofÎ z via the NV center. According to Eq. (3), the measured signal Î x (t) = α=X,Y,Z (e x · e α ) Î α (t) consists of a non-oscillatory term e −t/T 1 sin 2 θ cos 2 ϕ/2 that decays with a time scale T 1 and an oscillatory term e −t/T 2 (1 − cos 2 ϕ sin 2 θ) cos(|ω|t)/2 that decays with a time scale T 2 . The oscillating feature has been observed experimentally 2 . When the nuclear spin quantization axis e Z is parallel to the initial state polarization direction e x , i.e., θ = π/2 and ϕ = 0, the oscillatory feature disappears.
At room temperature, when the magnetic field is along the z axis, the optical transition is spin conserving. The cyclic transition between the m = 0 ground state |g = |0 g and excited state |e ≡ |0 e does not contribute to nuclear spin dissipation since ω g = ω e = 0. When B deviates from the z axis, its transverse component B T ≡ B x e x + B y e y mixes the m = 0 sublevels and the m = ±1 sublevels, so that ω g = −(2γ e /D gs )B T · A g and ω e = −(2γ e /D es )B T · A e , where D gs (D es ) is the zero-field splitting in the NV ground (excited) state and γ e = 28.025 MHz/mT is the gyromagnetic ratio of the NV electron. This can be understood as a hyperfine enhancement of the nuclear spin g-factor 2 (see the next section for more detailed discussion). As a result, the nuclear spin dissipation rates Γ ϕ , Γ ± are proportional to |B T | 2 , as observed experimentally 2,16 . For the 13 C b nucleus studied by Dutt et al. 2 , the HFI tensor has been obtained by ab initio calculations 27 as From this HFI tensor, we estimate |ω g − ω e | ∼ 0.3 MHz when |ω g | = 1 MHz. Under optical pumping rate R = γ 1 (the experimental condition 2 ), the two-level fluctuator model [Eq. (4a)] gives a nuclear spin dephasing rate 1/T 2 ≈ Γ ϕ ≈ (250 µs) −1 [≪ 1/τ NV ∼ (4 ns) −1 , so the NV center is a good Markovian bath], which is two orders of magnitudes smaller than the experimentally observed value ∼ (1 µs) −1 . Equivalently, to be consistent with the experiment 2 , the difference |ω g − ω e | must be assumed to be 10 times larger 16 . This large discrepancy suggests that the leakage from m = 0 subspace to m = ±1 subspace may plays an important role in determining the nuclear spin dissipation. In the next section, we shall demonstrate that the small leakage from |0 subspace to the | ± 1 subspace could introduce additional contributions that may dominates the nuclear spin dissipation.

A. Model
Now we consider general optical pumping of the NV center at room temperature, incorporating the finite intersystem crossing between m = 0 and m = ±1 subspaces. In this case, the electron-nuclear coupled system still obeys Eq. (1). The only difference is that now the Liouville superoperator L e for the optically pumped NV center includes seven energy levels: the ground triplet |m |g ≡ |m g (i.e., |0 g , | ± 1 g ), the excited triplet |m |e ≡ |m e (i.e., |0 e , | ± 1 e ), and the metastable singlet |S , where |g (|e ) denote the ground (excited) orbital. The unitary part of L e is a seven-level NV Hamiltonian (However, we haven't consider the effect of ionization of the NV center which would happen when the laser intensity is very strong 12 )Ĥ e = ∆σ e,e +Ĥ gs +Ĥ es + where ∆ is the optical detuning between the zero-phonon line and the laser frequency,Ĥ gs = D gsŜ 2 g,z + γ e B ·Ŝ g and H es = D esŜ 2 e,z + γ e B ·Ŝ e describe, respectively, the ground state triplet with zero-field splitting D gs = 2.87 GHz and the excited state triplet with zero-field splitting D es = 1.41 GHz. The dissipative part of L e includes various dissipation processes between the seven levels of the NV center as sketched in Fig. 1: the spontaneous emission from the excited orbital |e to the ground orbital |g with rate γ 1 = 1/(12 ns) 28 , the non-radiative decay from | ± 1 e to the metastable singlet |S with rate γ s1 ≈ γ 1 followed by the non-radiative decay from |S to |0 g with rate γ s = 1/(143 ns) 29 , the leakage from |0 e to | ± 1 g with equal rates γ s2 ≪ γ s1 , and the orbital dephasing of the excited stateL =σ e,e with rate γ ϕ ∼ 10 7 MHz 23,24 .
As discussed in the previous section, there are two processes contributing to nuclear spin dissipation. The one involving the electron spin flip is strongly suppressed away from the NV center ground state and excited state level anticrossing. Thus, in our analytical derivation below, we consider the other process that does not change the electron spin state, i.e., we drop the off-diagonal electron spin flip terms inF and only keep the diagonal part. The magnetic field component B T ≡ B x e x + B y e y perpendicular to the N-V axis (z axis) slightly shifts the electron levels and mixes the electron states from |m g and |m e to |m g and |m e (m = 0, ±1). For γ e |B T | ≪ D gs , D es , the level shift can be safely neglected, but the state mixing has a nontrivial influence on the diagonal part ofF, i.e., we need to keep the terms diagonal in the mixed basis |m g and |m e : F ≈ m (σm g ,m g m g |Ŝ g |m g · A g +σm e ,m e m e |Ŝ e |m e · A e ).
Up to the first order of the small quantities |γ e B T |/D gs and |γ e B T |/D es , we obtain (hereafter |m g/e stands for |m g/e ): The terms proportional to B T lead to hyperfine enhancement of the nuclear spin g-factor 2 , e.g., the termσ 0 g ,0 g 0 g |Ŝ g |0 g ·A g inF can be written asσ 0 is the correction to the nuclear gyromagnetic ratio by the HFI conditioned on the electron being in |0 g . Since γ e |B T | ≪ D gs , D es , the hyperfine enhancement terms in ±1 g |Ŝ g | ± 1 g and ±1 e |Ŝ e | ± 1 e can be safely dropped, so that where a g = −(2γ e /D gs )B T · A g , a e = −(2γ e /D es )B T · A e and b g/e = e z · A g/e . Now the second term of Eq. (1) describes the nuclear spin precession conditioned on the electron state: the precession frequency is γ N B ± b g (γ N B ± b e ) when the electron state is | ± 1 g (| ± 1 e ), or γ N B ± a g (γ N B ± a e ) when the electron state is |0 g (|0 e ), or γ N B when the electron is in the metastable singlet |S . The optically induced hopping of the electron between different states randomizes the precession frequency and leads to nuclear spin dissipation 16 . Below we derive analytically a closed equation of motion of the nuclear spin to describe these effects.

B. Lindblad master equation for nuclear spin
The time scale τ NV = 1/ min{γ s2, R} for the seven-level NV center to reach its steady state is determined by the time scale of the slowest process: the intersystem crossing from m = 0 subspace to m = ±1 subspaces if the optical pumping is strong, or the optical pumping rate R = Ω 2 R /γ ϕ from the ground orbital to the excited orbital if the optical pumping is weak. When τ NV is much shorter than the time scale T 1 , T 2 of the nuclear spin dissipation, we can follow exactly the same procedures as used in the previous section to derive the Lindblad master equation for the nuclear spin density matrix and the Bloch equation for the average nuclear spin angular momentum. The former (latter) has exactly the same form as Eq.
(2) [Eq. (7)] and describes the nuclear spin dissipation in the tilted cartesian coordinate (e X , e Y , e Z ) with e Z ≡ω/|ω| [see Eq. (3)] andω = γ N B + P 0 g a g + P 0 e a e , where P 0 g and P 0 e are steady state populations of the NV center on |0 g and |0 e . The detail expression of steady populations is given in appendix A. Now we discuss the analytical expressions for the nuclear spin pure dephasing rate Γ ϕ and relaxation rate Γ ± . The former comes from the fluctuation ofF Z , while the latter comes from the fluctuation ofF ± ≡F X ± iF Y . SinceF is a linear combination ofŜ g,z ,Ŝ e,z ,σ 0 g ,0 g , andσ 0 e ,0 e , the fluctuation ofF Z andF ± involve various cross-correlations among these four operators. Fortunately, due to the large orbital dephasing at room temperature, the optical pumping rate from the ground orbital to the excited orbital is nearly independent of the spin state and the coherence between electron states can be neglected. This allows us to neglect the cross correlation between {Ŝ g,z ,Ŝ e,z } and {σ 0 g ,0 g ,σ 0 e ,0 e } (see Appendix B for details). So Γ ϕ and Γ ± are the sum of the contributions Γ (1) ϕ , Γ (1) ± from the fluctuation ofŜ g,z ,Ŝ e,z associated with the m = ±1 subspace and the contributions Γ (0) ϕ , Γ (0) ± from the fluctuation ofσ 0 g ,0 g ,σ 0 e ,0 e associated with the m = 0 subspace. Unless explicitly specified, hereafter we consider a typical situation |ω| ≪ R, γ s1 .
The contribution from m = ±1 subspace is where P i is the steady-state population of the electron state |i . Formally Γ (1) ϕ and Γ (1) ± are independent of the magnetic field, but actually the components b g,Z , b g,⊥ ≡ b g,X e X + b g,Y e Y , etc. are defined in the tilted coordinate e X , e Y , e Z [see Eqs. (3)], which in turn depends on the magnetic field. Importantly, Γ (1) ϕ and Γ (1) ± do not vanish even in zero magnetic field. This provides a possible solution to the puzzling observation of nuclear spin dephasing in zero magnetic field 2 , which has been speculated to be due to the orbital fluctuation of the NV center in the excited state 16 .
Equations (9) and (10) 4)]. This is because in the two-level fluctuator model, the electron only hops between |g (with nuclear spin precession frequency γ N B + a g ) and |e (with nuclear spin precession frequency γ N B + a e ): when a g = a e , the electron hopping does not randomize the nuclear spin precession, so there is no nuclear spin dissipation. By contrast, in the seven-level fluctuator model, the electron can hop between seven energy levels, each of which corresponds to a different nuclear spin precession frequency (see the discussion at the end of the previous subsection). Therefore, even if the hyperfine of the excited state is the same as that of the ground state, the dissipation process also exists. This conclusion is different from the expectation 14,18 that the decoherence comes from the hyperfine difference between the ground state and excited state. The nuclear spin dissipation vanishes only when all these precession frequencies are equal, i.e., when a g = a e = b g = b e = 0. Second,Γ (1) ± and Γ (1) ϕ are proportional to the electron population P −1 e = P +1 e ∝ γ s2 in the | ± 1 e level, which vanishes when the leakage rate γ s2 from m = 0 subspace to m = ±1 subspace vanishes. Third, under weak pumping R ≪ γ 1 , γ s1 , we have Γ (1) ϕ , Γ (1) ± ∝ 1/R increasing with decreasing pumping strength, until the pumping is too weak for the Markovian assumption τ NV ≪ T 1 , T 2 , based on which our analytical formula are derived, to remain valid. Upon further decrease of the pumping strength, the NV center becomes a non-Markovian bath and the nuclear spin dissipation rates would show a maximum and then decrease (see the next subsection for more discussions). Finally, under saturated optical pumping, Γ (1) ϕ and Γ (1) ± are saturated instead of being suppressed: whereT ≡ 2/γ s1 +1/γ s2 +1/γ s ≈ 1/γ s2 is the average duration of one electron hopping cycle,τ 1 ≡ 2/γ s1 is the uncertainty of the dwell time in the | + 1 (| − 1 ) level, and the prefactor 2 accounts for the contribution from the m = +1 and m = −1 subspaces. Equations (11) can be understood as follows. First, under strong pumping, the hopping time between the ground orbital and the excited orbital is negligibly small, so the m = +1 (or m = −1) subspace effectively becomes a single energy level with nuclear spin precession frequency γ N B + (b g + b e )/2(or γ N B − (b g + b e )/2). Second, the durationT ≈ 1/γ s2 of one hopping cycle is ultimately limited by the slowest process: the intersystem crossing from m = 0 to m = ±1 subspace. Therefore, Eqs. (11) correspond to an effective two-level fluctuator model [cf. Eqs. (4)]: one state is the m = +1 (or m = −1) subspace with nuclear spin precession frequency γ N B+(b g +b e )/2 (or γ N B − (b g + b e )/2) and the other state is the subspace excluding | ± 1 subspace, which produce a nuclear spin precession frequency γ N B. Although strong optical pumping suppresses the randomization of the nuclear spin precession due to spin-conserving electron hopping between the ground orbital and excited orbital inside the m = +1 (or m = −1) subspace, there is extra contribution due to the random electron hopping between the m = ±1 subspace and the subspace excluding | ± 1 subspace.
The contributions from m = 0 subspace involve a g and a e in a quadratic form, so Γ (0) ϕ , Γ (0) ± ∝ |B T | 2 increases significantly with the magnetic field components perpendicular to the N-V axis. Due to the finite leakage from m = 0 into m = ±1 subspace, the analytical expressions for Γ (0) ϕ and Γ (0) ± are very tedious (see Appendix B), so here we discuss the limits of weak pumping and strong pumping. Under weak pumping, Γ (0) ϕ , Γ (0) ± ∝ 1/R decrease with increasing pumping strength (this behavior does not persists down to R ≪ 1/T 1 or 1/T 2 , where our Markovian assumption does not hold). Under strong optical pumping, they are saturated: is the uncertainty of the time for the electron dwelling at the m = 0 subspace. Similar to the contributions from the m = ±1 subspace, under strong optical pumping, the contributions from the m = 0 subspace correspond to an effective two-level fluctuator model: one is the m = 0 subspace with nuclear spin precession frequency γ N B + (a g + a e )/2, the other state is the subspace excluding the m = 0 subspace, which produce a nuclear spin precession frequency γ N B.
When the leakage from m = 0 subspace to m = ±1 subspace is neglected (i.e., γ s2 = 0), the steady-state populations in the m = ±1 subspace vanish, corresponding to perfect optical initialization of the NV center into the state |0 g . In this case, we have Γ (1) ϕ = Γ (1) ± = 0 and where P 0 e = 1 − P 0 g = R/(2R + γ 1 ). This recover the roomtemperature two-level fluctuator model [Eqs. (4) and (5)]. This can be easily understood: since the population is trapped in the m = 0 subspace, the fluctuation of the nuclear spin precession frequency could only come from the difference between a g and a e . Below we discuss two situations: (i) the magnetic field is along the N-V axis (z axis); (ii) the magnetic field is perpendicular to the N-V axis (z axis).

C. Magnetic field along N-V axis (z axis)
When the magnetic field is along the N-V symmetric axis (z axis), we have a g = a e = 0, so the average precession frequencyω = γ N B is along the −z axis, and the tilted axis (e X , e Y , e Z ) can be chosen as (e x , −e y , −e z ). Since a g = a e = 0, only m = ±1 subspace contribute to nuclear spin dissipation: Γ ϕ = Γ (1) ϕ and Γ ± = Γ (1) ± [see Eqs. (9) and (10)]. To begin with, we demonstrate the validity of our analytical formula Eqs. (9) and (10) by comparing them with the exact numerical results from directly solving the electron-nuclear coupled equations of motion [Eq. (1)]. we estimate the typical nuclear spin dissipation time ∼ T 1 , T 2 ≪ τ NV for 13 C b . In this case, the NV center is a highly non-Markovian bath beyond the description of our analytical formula. To see how our analytical fomula becomes progressively applicable when going from the non-Markovian regime to the Markovian regime, we manually scale down A g and A e by a factor η = 1, 10, and 100 to decrease the nuclear spin dissipation. Figures 2(a) and 2(b) show three features: (i) Both the exact results (dashed lines) and our analytical results (solid lines) tend to saturate at large R, even for very strong HFI (η = 1), where the NV center is highly non-Markovian. (ii) With increasing η and/or optical pumping rate R, the nuclear spin dissipation rates 1/T 1,2 decrease and/or the electron dissipation rate 1/τ NV increases, thus our analytical results begin to agree with the exact numerical results. (iii) For successively small R, the analytical dissipation rates 1/T 1,2 (solid lines) tend to diverge, while the numerical results (dashed lines) exhibit a maximum value ∼ 1/τ NV . This is because at sufficiently small R, the time scale of the NV dissipation τ NV ∼ 1/R is longer than the nuclear spin dissipation and the NV center becomes a non-Markovian bath. In this case, the electron-induced nuclear spin dissipation rates 1/T 1 and 1/T 2 are upper limited by the electron dissipation rate ∼ 1/τ NV .
In Figs. 2(c) and 2(d), both the exact numerical results (dashed lines) and our analytical formula (solid lines) show that the nuclear spin dissipation rates 1/T 1,2 increase rapidly with increasing leakage rate γ s2 from m = 0 to m = ±1 subspaces, due to the rapid increase of the population P −1 e [see Eqs. (9) and (10)]. When γ s2 = 0, the population P −1 e = 0, so our analytical formula gives vanishing nuclear spin dissipation rates, while the exact numerical results give a extremely small dissipation rates. This residue dissipation comes from the process involving the electron spin flip, which have been neglected in our theory since it is strongly suppressed by the large electron-nuclear energy mismatch away from the ground state and excited state anticrossing. Nevertheness, for extremely small γ s2 (= γ s1 /1000), it is responsible for the small difference between the analytical results (blue solid line) and the exact numerical results (blue dashed line) at large optical pumping rate in Fig. 2(d).
Next we study the magnetic field dependence of T 1 and T 2 and compare them with the experimental measurements 17 . For the nuclear spin dephasing time T 2 in Fig. 3(a), away from the ground state and excited anticrossing of the NV center [indicated by arrows in Fig. 3(a)], our analytical formula agree well with the exact numerical results, whether or not the electron spin flip terms inF is included. This indicates that the contribution involving the electron spin flip is negligibly small compared with the contribution not involving the electron spin flip.
In deriving Eq. (10) for the nuclear spin relaxation rates, we have neglected the small nuclear spin level splitting |γ N B|. When this effect is included, the analytical expressions for Γ ± are given in Eq. (A2), which shows a Lorentzian dependence on the magnetic field Γ ± ∝ 1/(|B| 2 + δ 2 B ) with a characteristic width Under saturated pumping, as is usually used for optical readout, this width ∼ γ 1 /|γ N | ∼ 500 mT. By contrast, the contribution involving the electron spin flip also has a Lorentzian dependence on the magnetic field, but with a much smaller characteristic width ∼ γ 1 /γ e ∼ 1 mT. The magnetic dependence of the relaxation time of 13 C nucleus has been measured by Dreau et al. 17 . They found that the anisotropic components A g,zx = A g,xz and A g,zy = A g,yz of the ground state HFI significantly contribute to the nuclear spin relaxation. The compoment A g,z,z has been measured to be 0.25 MHz, while the other components are not clear. Here we assume A e = A g with an isotropic diagonal components A g,x,x = A g,y,y = A g,z,z = 0.25 MHz and a small anisotropic component A g,z,x = A g,x,z = 1.5 kHz and A g,zy = A g,yz = 0. Figure 3(b) shows that the exact numerical results obtained by directly solving Eq. (1) agree well with the experimentally measured T 1 time 17 . As discussed previously, the exact numerical results contain two contributions: the one not involving the electron spin flip (which is treated by our analytical formula) and the one involving the electron spin flip (which is not treated by our analytical formula). Figure 3(b) shows that our analytical formula provides an accurate description to the former contribution, although in the present case the latter contribution dominates because of the much larger isotropic HFI ∼ 0.25 MHz compared with the anisotropic HFI ∼ 1 kHz. Finally, for the nuclear spin at lattice O as reported in the supplement of Ref. 17, it has a much shorter relaxation time ∼ 40 ms at 200 mT. Such short relaxation time is obviously dominated by the mechanism of Eq.(10) ,from which, we can estimate the anisotropic HFI component of this nuclear spin to be ∼ 20 kHz.

D. Magnetic field perpendicular to N-V axis
Without loosing generality, we consider the magnetic field B = B y e y along the y axis of the conventional coordinate. In this case, the precession frequencies a g = −(2γ e B y /D gs )e y · A g and a e = −(2γ e B y /D es )e y · A e are proportional to the magntic field. The nuclear spin precession frequencyω ≡ γ N B + P 0 g a g + P 0 e a e deviates from the z axis. In this case, both Γ (1) ϕ , Γ (1) ± [see Eqs. (9) and (10) coordinate e X , e Y , e Z ≡ω/|ω| that differs from the conventional coordinate (e x , e y , e z ).
First, we compare our analytical formula for the nuclear spin 1/T 1 and 1/T 2 to the exact numerical results from directly solving the electron-nuclear coupled equations of motion [Eq.
(1)]. To see how our analytical fomula becomes progressively applicable when going from the non-Markovian regime to the Markovian regime, we start from the strongly coupled nuclear spin 13 C b and downscale its HFI tensors A g ( 13 C b ) and A e ( 13 C b ) [see Eqs. (8)] by a factor η = 1, 10, and 100 to decrease the nuclear spin dissipation. The nuclear spin 1/T 2 and 1/T 1 shown in Fig. 4 show very similar behaviors to the case when the magnetic field is along the N-V axis [cf. Fig. (2)], including the saturation at large optical pumping rate R and the improved agreement between the analytical results and the numerical results with increasing η and/or R. In particular, Figs. 4(c) and 4(d) show that 1/T 1,2 increase rapidly with the leakage rate γ s2 , indicating that in addition to the contributions Γ (1) ϕ and Γ (1) ± from m = ±1 subspace, the contributions Γ (0) ϕ and Γ (0) ± from the m = 0 subspace also increase with γ s2 . For γ s2 = 0, the nuclear spin dissipation becomes very slow. In this case, the contribution from the processes involving the electron spin flip (not included in our analytical treatment) is no longer negligible. This leads to the discrepancy between the analytical results (black solid lines) and the numerical results (black dotted lines) in Fig. 4(d).
Finally we set the scale factor η = 1 and compare our theoretical results with the experimental measurements 2 . In this case, the strong HFI makes the NV center a highly non-Markovian bath, so our analytical theory only provides a qualitative description for the nuclear spin dissipation. Since B y , the magnetic field term γ N B y e y and the HFI contribution a g , a e ∝ B y e y dominates the average nuclear spin precession frequencyω, so the nuclear spin quantization axis e Z ∝ω is almost along the y axis. Since a g and a e (b g and b e ) are nearly parallel (perpendicular) to e Z , the m = 0 (m = ±1) subspace mainly contribute to the nuclear spin pure dephasing (relaxation), so that Γ ϕ ≈ Γ (0) ϕ increase quadratically with the magnetic field, while Γ ± ≈ Γ (1) ± is nearly independent of the magnetic field. In other words, we expect that the nuclear spin 1/T 2 to increase appreciably with B y and the nuclear spin 1/T 1 to be nearly independent of B y , as confirmed in Fig. 5(a). According to the Bloch equation Eq. (7), since the experimentally used initial states are eigenstates ofÎ z and I x , their decay time is largely determined by T 2 . Indeed, for B y ≫ 1 mT, Fig. 5(a) shows reasonable agreement between the numerically calculated 1/T 2 and the experimentally measured decay time of different initial states. Note that the twofold degenerate 1/T 2 correspond to identical decay of Î x and Î z [see Fig. 5(d)]. By contrast, for B y → 0, the average nuclear spin precession frequencyω is dominated by a small term Ŝ g,z b g + Ŝ e,z b e along the z axis (neglected in our analytical treatment). In this case, the fluctuation of b g and b e of the m = ±1 subspace mainly contribute to nuclear spin pure dephasing, while the fluctuation of a g and a e of the m = 0 subspace mainly contribute to nuclear spin relaxation. Correspondingly, in Fig. 5(a), the nuclear spin relaxation 1/T 1 ∝ B 2 y vanishes at B y = 0, while the nuclear spin 1/T 2 is two-fold degenerate, corresponding to near identical decay of Î x and Î y [see Fig. 5(a)]. Due to the switch of the nuclear spin quantization axis at intermediate magnetic field B y ∼ 1 mT, the association of the solid line with 1/T 1 and the dashed lines with 1/T 2 in Fig. 5(a) near the crossover region is meaningless.

IV. CONCLUSION
We have presented a numerical and analytical study for the nuclear spin dephasing and relaxation induced by an optically illuminated NV center at room temperature. When the NV center undergoes a single cyclic transitions, our analytical results provide a physically transparent interpretation that substantiates the previous results 16 and demonstrate the possibility to control the nuclear spin dissipation by tuning the magnetic field 16 . For general optical illumination of the NV center incorporating finite inter-system crossing, our numerical results agree with the experimental measurements 17 . Our analytical results suggests that the random hopping between the m = 0 (or m = ±1) triplet states and the corresponding remained subspace of the NV center could significantly contribute to nuclear spin dissipation. This means that increasing the spin polarization degree of NV center would effectively suppress the optical induced dissipation process. This contribution referred here is not suppressed under saturated optical pumping and provides a possible solution to the puzzling observation of nuclear spin dephasing in zero magnetic field 2 .

Appendix A: Lindblad master equation of nuclear spin
Here we derive a closed equation of motion for the nuclear spin from the coupled equation of motion Eq. (1), where L e is the Liouville superoperator of the two-level or seven-level NV model. First, we calculate the steady state density matrixP of the NV center from L eP = 0. For the two-level model, the steady state populations on |e and |g are P e = R/(2R + γ 1 ) and P g = 1 − P e , where R = 2π(Ω R /2) 2 δ ((γ 1 +γ ϕ )/2) (∆) is the optical transition rate from |g to |e and δ (γ) (x) = (γ/π)/(x 2 + γ 2 ) is the broadened δ-function. For the seven-level model at room temperature, due to the large orbital dephasing rate γ ϕ ∼ 10 7 MHz, the spin-conserving optical transition rates from the ground orbital |g to the excited orbital |e are all equal to R ≈ Ω 2 R /γ ϕ for different spin states. The steady-state population on |0 g is .
The populations on other NV levels are and P S ≈ (2γ s1 /γ s )P ±1 e . When the leakage from the m = 0 subspace to the m = ±1 subspaces are neglected by setting γ s2 = 0, we have P ±1 e = P ±1 g = P S = 0 and P 0 e = 1 − P 0 g = R/(2R + γ 1 ), which recovers the two-level fluctuator model. Second, we decompose the HFI into the mean-field part F e ·Î and the fluctuation part (F − F e ) ·Î ≡F ·Î, where F e ≡ TrFP is the average Knight field from the NV center, e.g., F e = P g ω g + P e ω e for the two-level model and F e = P 0 g a g + P 0 e a e for the seven-level model. Under this decomposition, Eq. (1) becomeṡ whereω ≡ γ N B + F e is the total magnetic field that defines the nuclear spin quantization axis. Consequently, the nuclear spin dephasing and relaxation should be defined in the cartesian frame (e X , e Y , e Z ), where e Z ≡ω/|ω|. Third, we decomposeF ·Î into the sum of the longitudional partF ZÎZ and the transverse part (F +Î− +F −Î+ )/2, where O ± ≡ O X ± iO Y . Then treatingF ·Î by the adiabatic approximation 25 up to the second order gives Eq. (2) for the nuclear spin density matrixp(t) ≡ Tr eρ (t), where is the nuclear spin pure dephasing rate due to the fluctuation ofF Z at zero frequency, and is the nuclear spin-flip rates due to the fluctuation ofF ∓ at the nuclear spin precession frequency |ω|, and is the steady-state correlation at frequency ω between the fluctuationã ≡â − TrâP and the fluctuationb ≡b − TrbP.