Exploring the Spin Dynamics of a Room-Temperature Diamond Maser using an Extended rate Equation Model

Masers - the microwave analogue of lasers - are coherent microwave sources that can act as oscillators or quantum-limited amplifiers. Masers have historically required high vacuum and cryogenic temperatures to operate, but recently masers based on diamond have been demonstrated to operate at room temperature and pressure, opening a route to new applications as ultra-low noise microwave amplifiers. For these new applications to become feasible at a mass scale, it is important to optimise diamond masers by minimising their size and maximising their gain, as well as the maximum input power of signals that can be amplified. Here, we develop and numerically solve an extended rate equation model to present a detailed phenomenology of masing dynamics and determine the optimal properties required for the cavity, resonator and gain medium in order to develop portable maser devices. We conclude by suggesting how the material parameters of the diamond gain media and dielectric resonators used in diamond masers can be optimised and how rate equation models could be further developed to incorporate the effects of temperature and nitrogen concentration on spin lifetimes.


I. INTRODUCTION
Due to their long spin relaxation times, negatively-charged nitrogen-vacancy defect centres (NV − ) in diamonds 1,2 have been widely used as room-temperature quantum materials for quantum metrology [3][4][5] , communications and quantum information processing [6][7][8][9] .Recently, their use has been extended to the microwave analogue of lasers, known as masers.Although masers are unmatched for low-noise microwave amplification, their use has historically been limited to niche applications in deep space communications and radio astronomy due to the need for high vacuum and cryogenic temperatures required to allow them to operate.The demonstration of roomtemperature masers based on organic [10][11][12][13] and inorganic [14][15][16][17][18] gain media thus opens up a path to the widespread use of masers as ultra-low noise microwave amplifiers, which in addition to existing applications in radio astronomy and communications technology are of increasing interest for the low-power microwave signal processing essential to many quantum technologies 19 .In this latter regard, masers may offer an advantage over Josephson Parametric Amplifiers currently used for such applications, which have a limited dynamic range and are easily saturated 20 .
The NV − defect is an S = 1 spin centre that can be described using a simple Hamiltonian consisting of a ground state triplet and an excited-stated triplet.NV − defects can be photoexcited to the excited state triplet using laser pulses with wavelengths at or below 637 nm, from which they may relax either through a spin-conserving radiative decay or through a spin-selective nonradiative decay whose net effect is to preferentially transfer electrons from the excited m s = ±1 states to the m s = 0 ground state via a metastable singlet state, as seen in Fig. 1.These spinselective transitions mean that NV − centres can be optically controlled and read out.The defect exhibits a long T 1 spin-lattice relaxation time (up to 6 ms 21 ) and T 2 spin-spin relaxation time (1.8   ms in isotopically purified 12 C diamond 22 , 0.7 ms for diamonds that have a natural abundance of 13 C isotopes 23 ) at room temperature.It is these properties that allow an ensemble of NV − defects to be used as a gain medium for a room-temperature maser, as well as for lasers 24 .The maser threshold pump rate is proportional to (T 1 (C − 1)) −1 , where C = 4g 2 N/(κ c κ s ) is the cooperativity, g is the coupling strength between each spin and the cavity field, N is the number of spins, κ c is cavity loss rate, κ s is the dephasing rate of the spins and T 1 is the spin-lattice relaxation time 15 The masing frequency is determined by the energy difference between the m s = 0 and m s = −1 level, thus, the masing frequency can be tuned via an external magnetic field through the Zeeman effect 3,25 .
Exploring the Spin Dynamics of a Room-Temperature Diamond Maser using an Extended rate Equation Model Despite impressive experimental demonstrations and characterisation of room-temperature diamond masers 15,17,18 , the widespread adoption of diamond masers will remain challenging until technical difficulties such as making diamond gain medium with sufficient amount of NV − centres, high Q-factor cavity, strong and homogeneous external magnetic field and alignment with NV − axis are addressed 26 .Here we use extended rate equation simulations to investigate how maser performance depends on the properties of the diamond gain medium and the dielectric resonator.Specifically, we extend the standard first-order rate equations sometimes used to model the photophysics of NV − centres 16 to a system of coupled second-order equations that incorporate the number of photons in the cavity.The extended second-order rate equations incorporate the number of photons in the cavity and allow the time-evolution of the populations of NV − centres in different energy states and the photon number to be simulated explicitly.The rate equations are employed to model the spin dynamics of the NV − centres under optical pumping at a wavelength of 532 nm.The results show masing dynamics with various parameters such as the loaded Q-factor of the cavity and the mode volume of the resonator.As well as aiding the design of room-temperature diamond masers, the extended rate equations presented will be useful for describing microwave mode-cooling devices 11 and quantum heat engines based on NV − centres in diamond 27 .

II. THEORETICAL MODEL
The populations of the electronic energy levels of the NV − can be described by a 7-level model in Fig. 1, where Levels 1 and 4 correspond to the m s = 0 states, and levels 2, 5 and 3, 6 correspond to the m s = −1 and m s = +1 respectively.The spin dynamics of different states under an optical excitation on NV − spin defects are described by rate equations based on a 7-level model based on an extended Bloch equation formalism 28,29 .In this model, a laser pumps electrons from the ground-state sublevels 1, 2 and 3 respectively to the excited-state sublevels 4, 5 and 6 all at the same optical pumping strength Ω l .The laser pump strength is defined below.
where E is the electric field of the laser and D ij is the dipole moment of the induced transition between level i and level j, such as the ground triplet states and the excited triplets.When a green laser at 532 nm is used to excite the NV − spin defects from the 1, 2 and 3 spin sub-levels of the ground state, the defects "overshoot" the transition to the corresponding excited-state sublevels The Ω m indicates microwave excitation.i indicate the levels, the κ i j indicates the depopulation rate from level i to level j, γ hB is the Zeeman splitting, W st and W sp is the stimulated emission rate and the spontaneous emission rate.The rates in the dashed line are larger than the rates in the dotted lines, and the more dashes or dots in the lines the higher the rate.
before relaxing non-radiatively to the ground phonon state of the excited state indicated by levels 4, 5 and 6.Since the phonon relaxation is generally fast compared to the optical pump rate, the transitions 1 → 4, 2 → 5 and 3 → 6 are given a common rate defined by Ω l .From the excited states, the NV − centres can either decay radiatively through spontaneous emission of a photon at 637 nm (indicated by red lines with respective decay rates) or non-radiatively through the intersystem crossing, which we represent with a singlet state 7.
The populations of excited states then decay back to the corresponding ground-state sublevels with decay rates, κ 41 , κ 52 and κ 63 .Additionally, the population of the excited levels also decay to where ω i is the frequency of the ith energy level, |i , of the NV − centres, ω m is the microwave excitation strength, Ω l is the coupling strength of the optical field, which is the pump strength that is determined by the electric field amplitude of the excitation laser at the NV − and the transition dipole moment of the ground state to the excited states of the NV − centres.
The rate equations can be derived by taking the matrix elements of the Hamiltonian.The decay rates between the levels are added in diagonal terms, representing the populations, whereas the dephasing terms are added in the non-diagonal terms, representing the coherences.The photon number is treated as a second-order rate equation that depends on the number of thermal photons, and the spontaneous and stimulated transition rates.These additional terms allow the rate equations to be extendable, such as the operating frequency through the effect of the spontaneous and stimulated emission rate, as well as the dephasing and detuning terms in the coherence terms that describe the coherence between the spins and the cavity fields and the spectral properties of the model.
The full matrix equations are shown in Appendix A. Using rates reported in 27 , they were solved numerically using the RK45 method and an implicit Backward Differentiation method 30 to simulate the time evolution of a diamond maser.
In general, the dipole moment for the NV − ground and excited states can depend on the amount Exploring the Spin Dynamics of a Room-Temperature Diamond Maser using an Extended rate Equation Model of strain in particular diamond samples, but for these simulations, we use a dipole moment of 5.2 C•m, estimated from the fluorescence lifetime 31 .The electric field amplitude is determined by the laser pump power available at the NV − centres with respect to the absorption dipole direction, laser polarisation and laser spot size.

III. MASER SIMULATION
The 7-level rate equations can be used to predict masing operation when a magnetic field above about 102.5 mT, is applied such that the m s = −1 state is below the m s = 0 state.The optical pumping polarises the m s = 0 state, thus, population inversion is achieved.The population will interchange between m s = 0 to m s = −1 via spin-lattice relaxation to maintain thermal equilibrium, at a rate κ 12 , and via spontaneous emission and stimulated emission which results in photon emission processes.The population will be excited to the m s = 0 state from m s = −1 via spinlattice relaxation and stimulated absorption.In the microwave (MW) spectrum, the spin-lattice relaxation rate, on the order of 10 2 Hz, is much higher than the spontaneous emission rate in free space, which is on the order of 10 −13 Hz 32 .The excited spins preferentially relax through the spin-relaxation mechanism by transferring energy to the lattice.The spontaneously emitted MW photons are negligible in free space.The spontaneous emission rate is related to the Einstein coefficients A and B, which phenomenologically describe the rate of the spontaneous and stimulated emission processes and can be boosted through the Purcell effect.The magnitude of this boost is proportional to the loaded Q factor and inversely proportional to the mode volume Q/V m , which Exploring the Spin Dynamics of a Room-Temperature Diamond Maser using an Extended rate Equation Model appears in the figure of merit of maser systems known as the cooperativity C, where C > 1 can be taken as a necessary but not sufficient condition for masing 15 .
The transition rate from one state i to another state j is determined 33 by the magnetic dipole matrix element of the corresponding two states, d ij = i|H • S| j , where H is the magnetic field and S is the spin vector, and the density of states of the electromagnetic field, ρ(ω If the emitter is placed in a cavity, the density of state will alter, and the emission rate will either increase or decrease.The change in rate is determined by the Purcell factor, F = 3Qc 3 2πV m ω 3 n 3 , where Q is the cavity loaded Q factor, c is the speed of light, ω is the frequency, V m is the mode volume and n is the refractive index in the cavity 32 .The spontaneous and stimulated emission rate is directly proportional to the Purcell factor.To enhance the spontaneous emission rate in a given frequency, a high Q factor and small mode volume are desired.The Purcell factor is sensitive to the frequency, F ∼ ω −3 , hence it might be beneficial to operate a maser at low frequency if the mode volume is not scaled up rapidly compared to the frequency.The stimulated emission and absorption rate are related to the spontaneous emission rate according to the Einstein coefficients. If the emitters are placed in a cavity, the density of states of the electromagnetic field is altered by reducing the number of states, thus increasing the spontaneous emission rate.The stimulated emission rate is κ st = κ sp u(ν)c 3 π 2 /hω 3 , where u(ν) is the energy density of the MW field and κ sp is the spontaneous emission rate.In the masing simulation, the initial photon number in the cavity is determined by the thermal photon number at thermal equilibrium, n th = (e hω/k B T − 1) −1 .At 9.22 GHz and room temperature, the n th is approximately 661.MW photons are generated via spontaneous emission and trigger stimulated emission and absorption.For a realistic emission, dephasing and decoherence will contribute to the rate, which replaces the density of state by the linewidth function, g(ω).Given that the magnetic field of the MW is H,the stimulated transition rate can be expressed 34 : where can be interpreted as the normalised imaginary part of the magnetic susceptibility.For a two-level system interacting with a microwave field with positive circular polarisation, the maximum value of σ 2 is 0.5, while a microwave field with negative circular polarisation gives 0 as it will oppose the spin precession 35 .Provided that the magnetic field energy , where V m is the mode volume and n is the number of photons in the mode and W = Bn, the Einstein B coefficient at resonance is given by 34 : where the filling factor η and the mode volume V m are respectively defined by The higher the filling factor, the more strongly the sample interacts with the cavity mode.
For the results presented in the following section, the parameters from Table I and Table II are used.The finite-element software CST Microwave Studio is used to estimate the mode volume and the filling factor, using equations ( 4) and ( 5).In the CST simulation, we model a sapphire resonator with an inner diameter of 5 mm and an outer diameter of 10 mm on a sapphire post of 12 mm high in a copper cavity with an inner diameter of 31 mm and a height of 30 mm.The energy density of the magnetic field of the TE 01δ mode of the cavity is shown in Fig. 2. The B coefficient is approximately 1.98×10 −6 s −1 .The spontaneous rate is given by 35 : At 9.22 GHz, the spontaneous emission rate is about 1.57 × 10 −12 s −1 , which is much slower than the stimulated transition rate and the spin relaxation rate.Knowing both the Einstein A and B coefficients, the stimulated emission rate from level i is simply W st = Bρ ii N(ω), where N(ω) is the number of photons at frequency ω, and the spontaneous emission rate from level i is W sp = Aρ ii .
The optical pumping interaction strength is 10 Hz.The initial population is set by the Boltzmann distribution, e hω m /k B T , and the initial number of photons is (e hω m /k B T − 1) −1 , where k B is the Boltzmann constant, T is the temperature and the ω m is the frequency of the cavity field, which is the frequency difference between the m = 0 and m = −1 states.
Looking at the intracavity photon number as a function of the loaded cavity Q-factor (Fig. 3), with other parameters fixed, we may identify a masing threshold as the point where the number of photons in the cavity is an order of magnitude greater than the thermal population.The maser output power in dBm can be estimated using 10 log 10 ((nhω 2 m /Q m )/1mW), where ω m is the operating maser frequency, h is the Planck constant, n is the intracavity photon number and Q m is the Q factor for the external load, which can be estimated by the loaded and unloaded Q factors.Assuming that the loaded Q, which is the Q factor on the x-axis is 30000 and the unloaded Q is 50000, which is under coupled that the coupling constant is 0.67, the maser output is without any amplifiers is -108.5 dBm, calculated from the figure.The plot in Fig. 3 starts with a low Q factor, the number of photons is close to the thermal photon number, which suggests that there is no masing.When the Q factor increases, the number of photons also increases.At a certain threshold Q factor, the number of photons increases significantly, where the slope of the plot is the greatest.At this point, the photon emission rate is higher than the cavity decay rate, thus, the number of photons in the cavity increase, and after this point, the maser can be maintained because a higher Q factor means that the cavity has a longer photon lifetime.The emission is predominantly stimulated in this case.In Figs. 4 and 5, the performance of the maser at different pump strengths is studied at a fixed Q factor of 30,000.It is clear that the maser operation becomes saturated when the pump strength is about 1 kHz, where the maser operation does not benefit from having a higher pump strength.If the pump power is increased further as in Fig. 5, the number of photons decreases after 100 kHz of excitation pump strength.At such high pump strength, the populations from the ground states are pumped to the excited state and accumulate in the singlet state |7 , as its decay rate to the ground state κ 71 is slow compared to κ 63 , κ 52 , κ 41 .Overall, the number of spins that participate in maser reduces at high pump strength and leads to a reduction in the number of photons in the cavity.
The threshold behaviour of the mode volume, which is inversely proportional to the Purcell factor and the cooperativity, is shown in Fig. 6.As expected, smaller mode volumes increase the number of masing photons in the cavity.For the parameters in Tables 1 and 2, the threshold is roughly at 0.3 cm 3 , and for larger mode volumes masing cannot occur.Fig. 7 shows the number of masing photons as a function of the pump strength and the loaded Q factor.The dark blue region shows that the masing effect does not occur due to the high losses and low pump strength.The second region, in dark green, shows the "incipient masing" region, where stimulated emission starts to build up and the number of photons is higher than the number of thermal photons.In this region, an external MW signal fed into the maser would be amplified but the maser cannot be run as an oscillator.The maser is in operation in the light green region, where stimulated emission is maintained and the device can undergo self-sustaining maser oscillations.
The number of intracavity photons is shown as a function of the gain medium filling factor (defined by the ratio of the illuminated gain medium volume to the mode volume) in Fig. 8, again displaying a clear threshold for masing.The results show that the threshold is about 0.03, beyond which the number of photons maintained in the cavity significantly increases beyond the thermal population.The number of cavity photons increases monotonically as the filling factor approaches unity.
It is assumed that the number of NV − centres is constant during masing.This is reasonable since the pump power required to achieve maser operation is low enough that the ionisation and recombination of NV centres due to laser excitation is weak 24,36 .Similarly, all the rates are assumed to be constant.In reality, the diamond will be heated by the laser, which will lead to a reduction in T 1 time 18 as well as the T 2 time 37 .Both coherence times vary with temperature as T 5 , whereas the T * 2 caused by the field inhomogeneity is only very weakly dependent on temperature 37 .The zero-field splitting is also temperature dependent 38 , but since this effect is also weak it is ignored in this model.All these temperature effects will tend to broaden the maser linewidth.A recent experiment has shown that the high power rate will significantly increase the temperature, especially the T 1 time, which makes maser action impossible 18 by affecting the ability of the gain medium to maintain a population inversion.It would be interesting to develop the rate equation approach presented here to include such effects.
A further interesting extension to our rate equation model would be to explicitly model the effect of NV concentration on the T 2 time.To first order, changes in NV − concentration do not affect the cooperativity since the spin decay rate κ s is proportional to the NV − concentration 39 so that the factor of N/κ s in the cooperativity is constant.However, since the concentration of NV − in diamond samples depends not just on the concentration of nitrogen impurities but on the efficiency with which such impurities can be converted to nitrogen-vacancy centres and on the relative populations of the two possible charge states (NV − and NV 0 ), each of which can have (1), at high pump strengths.For the parameters listed in Table I, the number of microwave photons begins to decrease if the pump strength is increased beyond around 100 kHz because the NV − centres begin to accumulate in the |7 state due to its slow decay rate.This tends to reduce the population inversion and inhibit maser action.
an independent effect on the spin decay rate, the cooperativity is not independent of the NV − for real systems.The ratio of NV:NV − is typically 15 around 14:1, but a ratio of 6.95:1 has been reported for an optimised system 40 .Finally, we expect a higher concentration to lead to increased self-absorption 41 , which is not taken into account in our model but will again tend to increase the threshold pump power by increasing the effective cavity decay rate.
An important question is how field inhomogeneity affects diamond masers and to develop scalable room-temperature diamond masers it will be important to determine the minimum field inhomogeneity required.The homogeneous linewidth of the NV − resonance is about 1 MHz 42 and in these simulations, it is assumed that the spatial inhomogeneity of the applied magnetic field over the volume of the diamond is small enough that the resulting Zeeman shifts are within the homogeneous linewidths of the NV − centres.For an intrinsic linewidth of 1 MHz, this means the field inhomogeneity should be < 0.03 mT.In real applications, the smaller the field inhomogeneity the better the maser performance.

IV. CONCLUSION
We have used extended second-order rate equations which incorporate the intracavity photon number and coherences between different NV energy states to simulate the spin dynamics of an ensemble of NV − centres in a diamond-based maser.The model predicts the behaviour of the maser as a function of material parameters of the gain medium and resonator in good agreement with recent experimental results 18 .
Cooperativity and the related Purcell factor are the key quantities which, if increased, will reduce the threshold pump power for masing.On the one hand, this requires a high Q factor and a small mode volume, but intrinsic material parameters of the diamond gain media such as the concentration of NV − and the T * 2 play a significant role, where a long T * 2 and more excited spins, which are determined by the filling factor and the concentration of the sample, are favourable.High pump strength will eventually reduce the number of maser photons, as the populations end up in the singlet state and the population inversion is reduced due to a kinetic bottleneck leading to the accumulation of NV − centres in the |7 state.Increasing the transition rate from the |7 state through Purcell enhancement or by resonantly pumping transitions from the |7 state could thus be a pathway to higher maser power outputs.
Attention to these considerations will aid the design of diamond-based masers suitable for mass production.
We further suggested improvements in the modelling of NV − ensembles for maser applica- The number of photons becomes saturated above η = 0.3.The inset shows the onset of masing at a filling factor of 0.03.
tions.If the NV − to NV 0 charge conversion and effects of temperature and NV − concentration on the spin lifetimes are considered, the trade-offs between maximising the coherence times and increasing output power and cooperativity by maximising the number of NV − centres and the pump power will become clearer.Such improvements will require further experimental input, particularly on the relationship between charge-state ratios and spin lifetimes.
The method presented in this paper, which computes the time evolution of the maser output explicitly, will be of particular use in the "incipient masing" regime where the stimulated emission is dominant over the spontaneous emission, but insufficient for self-sustaining maser oscillation and thus less amenable to solutions based on steady-state considerations.The full time-dependent response of diamond masers will also be useful for technological applications where masers are Exploring the Spin Dynamics of a Room-Temperature Diamond Maser using an Extended rate Equation Model used to amplify arbitrary time-dependent microwave signals.
More broadly, the extended rate equation model presented in this paper will be of interest for modelling other quantum technologies based on ensembles of spin defects, such as diamond magnetometers 3 and room-temperature masers based on alternative gain media such as siliconvacancy defects in SiC [43][44][45] .
Exploring the Spin Dynamics of a Room-Temperature Diamond Maser using an Extended rate Equation Model

Appendix A: Matrix equations
We present the rate equations solved numerically to generate the results of this work in matrix form.The overall equation of motion for the density matrix ρ is given by ρ where γ m is the dephasing rate of the MW-coupled spin states, γ l is the decoherence rate of the laser-coupled states and ∆ l is the detuning of the laser from the transition frequencies for 1 → 4, 2 → 5 and 3 → 6.For the results presented in this work, the laser was modelled as resonantly exciting the ground states to the excited states so that the detuning, ∆ l , is zero.∆ m is the detuning of the microwave from the maser operating frequency.
In the extended rate equations, an additional second-order term, representing the number of photons, is added to the first-order rate equations.where κ sp is the spontaneous emission rate and κ st is the stimulated emission rate.

FIG. 1 .
FIG. 1.This diagram shows the energy levels, decay rate and pumping rate of the 7-level model.The excited phonon bands are shown in a green area in levels 4, 5 and 6.The green dashed line represents the 532 nm laser that pumps the population to the phonon bands, whereas the green solid line is the transition represented by the pump strength Ω l .The Ω m indicates microwave excitation.i indicate the levels, the

FIG. 2 .
FIG. 2. The energy density of the magnetic field of the TE 01δ mode of a cylindrical copper cavity simulated using CST Microwave Studio.The red outlines show the boundaries of a sapphire resonator and post.The walls of the copper cavity are shown as the enclosure that bounds the image.

6 FIG. 3 .
FIG.3.The steady-state number of photons as a function of the loaded cavity Q factor.The photon number shows threshold behaviour with respect to the Q factor, with the onset of masing at around Q = 22000 for the parameters listed in TableI.The inset shows the onset on a logarithmic scale.

FIG. 4 .
FIG.4.The steady-state number of photons as a function of the optical pump strength, defined by Eq (1), at low pump strengths.For the parameters listed in TableI, the photon number saturates at a pump strength of around 1.5 kHz.

FIG. 5 .
FIG.5.The steady-state number of photons as a function of the optical pump strength, defined by Eq.

FIG. 6 .
FIG.6.The steady-state number of microwave photons as a function of the cavity mode volume for the parameters listed in TableI.The masing threshold occurs at a mode volume of about 0.3 cm 3 .The inset shows the onset on a logarithmic scale.

5 FIG. 8 .
FIG.8.The number of microwave photons in the cavity as a function of the filling factor, defined by Eq. 4.

M 21 M 22 
where M ij are sub-matrices of M used to clarify the equation by separating different contributions to the time evolution.M 11 is the 7 × 7 population sub-matrix, describing the spin dynamics of the 7-level model.M 12 and M 21 are the population-coupling sub-matrices, respectively 7 × 10 and 10 × 7, that describe the coupling of the laser field and the MW field with the 7-level model.Finally M 22 is the coherence sub-matrix that describes the coherence of the coupled fields.This 17 × 17 matrix contains 7 non-zero diagonal terms and 10 non-diagonal terms.The sub-matrices Exploring the Spin Dynamics of a Room-Temperature Diamond Maser using an Extended rate Equation Model are given explicitly as 71 + κ 72 + κ 73 )

da dt = κ sp ρ 11 +
κ st (ρ 11 − ρ 22 )a − κ c (n th − a) Exploring the Spin Dynamics of a Room-Temperature Diamond Maser using an Extended rate Equation Model the singlet state, level 7, with rates, κ 47 , κ 57 and κ 67 .The populations of level 5 and level 6 are much more likely to cross to the singlet state, κ 57 , κ 67 κ 47 .The population of level 7 decays back to level 1 with a rate, κ 71 , and to level 2 and level 3 with rates, κ 72 and κ 73 .Furthermore, the

TABLE I .
27ansition rates used to describe the time evolution of the ensemble of NV − centres in a cavity to produce the results presented in Figures3-8.The values are based on experimental results from Ref.27.

TABLE II .
15riable parameters used in the rate equations to describe the time evolution of the ensemble of NV − centres in a cavity.The number of available NV − , the maser frequency, ω, the dephasing time, T * 2 and the Q factor is from15.The other parameters are obtained via CST Microwave studio simulations.