Theoretical study of hyperfine interactions and optically detected magnetic resonance spectra by simulation of the C291[NV]-H172 diamond cluster hosting nitrogen-vacancy center

Single nitrogen-vacancy (NV) centers in diamond coupled to neighboring nuclear spins are promising candidates for room-temperature applications in quantum information processing, quantum sensing and metrology. Here we report on a systematic density functional theory simulation of hyperfine coupling of the electronic spin of the NV center to individual 13C nuclear spins arbitrarily disposed in the H-terminated C291[NV]-H172 cluster hosting the NV center. For the ‘families’ of equivalent positions of the 13C atom in diamond lattices around the NV center we calculated hyperfine characteristics. For the first time the data are given for a system where the 13C atom is located on the NV center symmetry axis. Electron paramagnetic resonance transitions in the coupled electron–nuclear spin system 14NV-13C are analyzed as a function of the external magnetic field. Previously reported experimental data from Dréau et al (2012 Phys. Rev. B 85 134107) are described using simulated hyperfine coupling parameters.


Introduction
The ability to create, control and measure coherence in multi-spin systems in solids is crucial for scalable applications of quantum information processing, quantum sensing and metrology. Coupled electron-nuclear spin systems where electrons act as fast processing qubits and, additionally, form an interface with photons, while nuclei can store quantum information for a long time owing to their exceptional isolation from the environment, are especially useful for these purposes.
Along with experimental characterization of hyperfine interactions (HFIs) in different NV-13 C spin systems, a complementary approach to getting the desired information is provided by quantum chemistry modeling, which has been shown [26,[35][36][37][63][64][65][66][67] to be effective in calculating the characteristics of hyperfine interactions of NV centers with surrounding nuclear spins. In earlier works [63,64] this was done using rather small supercells [63] and Hterminated carbon clusters [64] hosting the NV centers, while subsequent works, in which larger 512-atom supercells [65,66] and C 84 [NV] -H 78 clusters [35][36][37] were studied, were focused on the simulation of HFI characteristics for NV- 13 C spin systems wherein the 13 C atom was located quite close to the NV center vacancy (in the first or third coordination sphere of the center) because these spin systems were being actively investigated experimentally at that time.
Here, for the above reasons, we will mainly pay attention to the systematic computational chemistry simulation of HFI characteristics for NV- 13 C spin systems involving more distant 13 C nuclear spins. For this purpose, we will use the much larger NV-hosting H-terminated carbon cluster C 291 [NV] -H 172 and calculate the A KL C ( ) matrices describing HFIs between the e-spin of a single NV center and 13 C n-spins taking all possible positions in the cluster. Further, the calculated HFI matrices are used in the ground-NV-state spin Hamiltonian of the NV center to simulate the signatures of these HFIs in optically detected magnetic resonance (ODMR) spectra of an arbitrary 14 NV-13 C spin system. Numerical diagonalization of respective spin Hamiltonians are accompanied by simplified analytical consideration within the secular approximation. We also show how calculated HFI characteristics correlate with spatial positions of 13 C nuclei in the cluster. The HFI characteristics are presented for the 'families' [26,27] of equivalent positions of the 13 C atom in the diamond lattice around the NV center. The effect of external magnetic field on the rates of 'allowed' and 'forbidden' EPR transitions in an arbitrary 14 NV-13 C three-spin system is studied and shown to be different in the manifolds depending on the sign of the HFI matrix element A ZZ C ( ) . General consideration is illustrated by detailed analysis of the experimental data presented in the work [27].

Model and methods
To model the NV center in bulk-or nano-sized diamond we chose the diamond cluster C 291 [NV] -H 172 (figure 1(a)) which was constructed from a piece of ideal diamond lattice by removing one C atom in its central part, substituting the neighboring C atom with the N atom and saturating the surface C atoms' dangling bonds with H atoms. With the available computer resources the size of this cluster was optimal to ensure simulation of HFI characteristics for rather distant 13 C nuclear spins, in particular, for those disposed on the NV axis.
The geometric structure of the studied C 291 [NV] -H 172 cluster was optimized by the total energy minimization using the density functional theory (DFT) method with the B3LYP hybrid functional [68,69] and the MINI basis set [70]. In the cluster, atomic relaxation of the diamond lattice around the formed NV center results in increased distances between the three C atoms that are NNs of the vacancy (from 2.527 Å to 2.649 Å) and between these atoms and the N atom (from 2.527 Å to 2.764 Å), while the distances from the N atom to the three nearest C atoms are reduced from 1.547 Å to 1.512 Å. For the other C atoms in the relaxed cluster the C-C bond lengths range from 1.55 to 1.59 Å. These bond lengths reflect the overestimated bulk bond length obtained at this level of theory [71].
For the relaxed cluster the spin density distribution was calculated using the same method and the 3-21G basis set [72]. Calculations have been performed for singly negatively charged clusters in the triplet ground state (S = 1) using the unrestricted Kohn-Sham (UKS) procedure for open shell structures. We have used the PC GAMESS (US) software package [73] for geometry optimization of the cluster and the ORCA software package [74] to calculate the spin characteristics, in particular for full HFI matrices A k C ( ) for all possible positions k = 1 ÷ 291 of the 13 C atom in the cluster. These calculations have been done in the principle axes system (PAS) of the NV center where the z-axis coincides with the C 3V symmetry axis of the center while the x-and y-axes are chosen arbitrarily. Various 13 C lattice sites exhibit different and generally anisotropic hyperfine interactions with the NV e-spin, leading to different spin (and optically detected) properties of various NV-13 C spin systems. The calculated HFI A k C ( ) matrices can be converted into their respective diagonal ones = by unitary transformations U k from the NV PAS to the 13 C PAS with elements of the U k matrix being the direction cosines between various axes of both PASs.
In this article we are going to present not only the results of our calculations of the HFI matrices A k C ( ) for all possible positions of the 13 C atom in the studied cluster, but also show how DFT calculations can be employed for the quantitative description of experimental data. Here we will focus on ground-state ODMR spectra of different 14 NV-13 C spin systems. To model the NV center we introduce the standard spin Hamiltonian H k of an arbitrary 14 NV-13 C system transition in the presence of the magnetic field B||OZ = 50 gauss, (c) energy levels of the typical 14 NV-13 C system in low external magnetic field B|| OZ numerated from 1 to 18 in accordance with the increase in their energy and transitions between them (right part) shown along with simplified fine-structure energy levels of the NV center without and with external magnetic field B.
comprising the e-spin S = 1 of the ground-state NV center coupled to the internal n-spin I N ( ) = 1 of the 14 N atom of the center ( 14 N is the most common nitrogen isotope having 99.63% natural abundance) and to the additional n-spin I k C ( ) = 1/2 of a 13 C atom disposed in the kth lattice position in the cluster. These spin Hamiltonians are of the form Here the first term is the common zero-field spin Hamiltonian of the ground-state NV center  2 of equal area and using the halfwidth Γ as a fitting parameter.
3. Arbitrary 14 NV-13 C spin system: analytical consideration in the secular approximation

ODMR spectra, energy levels and eigenstates
Typical HFI structure of an ODMR spectrum of an arbitrary 14 NV-13 C system in low magnetic field is shown in figure 1(b) for the case of rather distant 13 C n-spin for which the HFI with the NV center e-spin is weaker than that with the 14 N n-spin. Both manifolds of the spectrum exhibit three characteristic pairs of lines corresponding to possible EPR transitions in the system with their frequencies determined by the HFI with both 14 N and 13 C n-spins, and also by the interaction of these spins with the external magnetic field. Note that the respective ODMR spectrum of the 14 NV center, having no nearby 13 C n-spin, consisted of just three lines with HFI splitting of ∼2. 16 MHz between them [75][76][77]. It is the HFI with additional distant 13 C n-spin splitting each of these three lines into pairs of lines, with each pair corresponding to definite projection m I N ( ) = −1, 0, 1 of the 14 N nspin, as shown in figure 1(c).
In many practical cases (excluding those at B ∼ 1027 gauss where avoided-crossing of sublevels with m S = 0 and m S = −1 takes place) good approximation for energy levels, eigenstates and transition rates for the 14 NV-13 C spin system is provided by the secular approximation where only the terms with S Z are kept in the Hamiltonian (1). In particular, using the approximation one can find that, at fixed m I N ( ) projection (for definiteness, we choose here the case m I N ( ) = +1 [19] and use below the state numeration adopted in figure 1(c)), the energy levels of an arbitrary threespin system 14 NV-13 C in a magnetic field B, aligned along the z-axis, are while those for the m S = ±1 substates ( Δ − = − E E 10 9 and Δ − = + E E 14 13 ) are determined by the parameters Δ ± that describe the combined effect of the HFI and the external magnetic field on the 13 C n-spin being conditioned on the electronic spin projection = ± m 1 S [60]. The same is valid for the sublevels corresponding to other m I N ( ) projections. Note the different dependence of the splittings Δ ± on B which can be used to determine the sign of the HFI matrix element A ZZ , as will be discussed in more detail later. At zero B field both splittings Δ ± are The eigenfunctions Ψ α , corresponding to the eigenvalues (4), can be written in the secular approximation as ZX where the angles θ ± and φ are chosen on account of the signs of A A , one needs to take the angles θ ± within the first quarter π [0, /2], while at the angles must be taken from the second quarter π π [ /2, ].

HFI-induced dynamics of a single 13 C nuclear spin in a 14 NV-13 C spin system
One can see from (5) that in the general case of an arbitrary spin system NV- 13 with coefficients depending on the HFI matrix elements and on the magnetic field. From this it follows that the 13 C n-spin does not change its projection when the NV e-spin is in the m S = 0 state while it will oscillate between the basis states due to the HFI with the NV espin having m S = −1 or m S = +1 projections, respectively. One can show within the secular approximation that in the case of the m S = −1 e-spin state, for example, the 13 C n-spin prepared in the = ⇓ ↑ ↗ V state at the moment t = 0 can be found in the state . Analogously, in the m S = +1 manifold the probability of the 13 C n-spin flip from the initial state . Therefore, the parameters Δ − and Δ + determine the rate of n-spin flips for 13 C induced by its HFI with the NV e-spin in m s = −1 and m s = +1 states.
From the above it follows that the HFI-induced 13 C n-spin flipping will be absent if the respective HFI matrix is diagonal. This case is realized for the NV-13 C spin system with the 13 C atom located at the NV's symmetry axis. The quantization axis of such a 13 C n-spin is always parallel to the NV axis and therefore it can be completely polarized using the technique demonstrated in [77,80,81] for the 14 N nuclear spin belonging to the NV center. Looking ahead, we would like to point out here that in the simulated C 291 [NV] -H 172 cluster there are three such 'on-NV-axis' positions for the 13 C n-spin, with two of them disposed on the edge of the cluster while the third one (which is the fifth neighbor of the vacancy) is rather well inside the cluster so that its HFI characteristics are barely influenced by passivating H-atoms. In fact, it is this last position which was chosen for the simulation at B = 50 gauss of the ODMR spectrum, shown in figure 1(b). The calculated A C ( ) matrix for this specific 13  , one can get for the on-NV-axis 13 giving the result ≈187.387 kHz, which is very close to the value Δ (0) = 187.380 kHz obtained by direct computer diagonalization of the spin Hamiltonian (1) with the above diagonal A C ( ) matrix. Calculated numerically spin wave functions in this special case practically coincided with the basis states m m m , ,  figure 1(c). As was explained above, the only fitting parameter for figure 1(b) was the halfwidth of Lorentzians which was taken to be Γ = 20 kHz. One can check that in this special case the splitting of paired lines in ODMR spectra both in the m S = −1 and m S = +1 states is equal to Δ (0) = 187.38 kHz. Thus, the 14 NV-13 C spin system with the 13 C atom in the on-NV-axis position can be identified experimentally by monitoring the above characteristic splitting in its ODMR spectrum.

Allowed and forbidden transitions
Along with HFI-induced 13 C n-spin flipping, mixed eigenstates (5) The frequency difference of the pairs of allowed lines at low magnetic field is Δ (0) . In particular, the above case is realized for the considered special case of the on-NV-axis 13 C atom and respective allowed transitions are shown in figure 1(c). At A ZZ < 0 the allowed and forbidden transitions are mutually interchanged.
The increase of a magnetic field B||OZ acting on a given NV-13 C spin system results in a modification of both probabilities and frequencies of EPR transitions in the spin system, which is different for the  for the analogous initially allowed transitions 1-14 and 3-13 in the = ↔ = + m m 0 1 S S manifold will constantly decrease from W a (0) up to zero with the B growth thus resulting in the transformation of the respective transition from the allowed to the forbidden. Moreover, from (4) it follows that the difference of frequencies Obviously, at A ZZ < 0 the situation with the modification of transition probabilities with B growth will be reversed for the forbidden and allowed transitions. Note that in all cases the threshold magnetic field B at which forbidden transitions become allowed and vice versa is γ = B A / ZZ n C ( ) . The above approximate analytical consideration of an arbitrary spin system NV-13 C indicates that, in experimentally recording ODMR spectra for various individual NV centers, one can not just discover the presence of the 13 C atom in the vicinity of the NV center but also, measuring the HFI splittings Δω of paired ODMR lines for the spin system and their modification with an applied external magnetic field, extract HFI parameters Δ (0) , A ZZ and A nd , and, moreover, determine the sign of the most essential element A ZZ of the HFI matrix for the specific kth position of the 13 C n-spin in a diamond lattice. Such experimental work was systematically performed recently in the studies [26,27], in which the ODMR spectra of hundreds of NV centers were studied for their HFI structures. Many of them proved to be coupled with single nearby or distant 13 C n-spins. Such spin systems exhibited discrete possible values of HFI splittings . Evidently, if we compare such experimental data with those calculated for all possible 14 NV-13 C systems by diagonalization of their spin Hamiltonians (1) we will be able to distinguish the position of the specific 13 C nucleus among others. practically equal values of calculated HFI parameters (see also [26]). These sets of nearequivalent lattice sites can be termed as 'families' [27]. In table 1, we present data for 26 such families labeled by letters of the English alphabet, indicated in the first column. The second column contains the number N C (=3 or 6) of equivalent members for each family. The third, fourth and fifth columns of the table show the values ofĀ ZZ , A nd and zero-field splitting Δ k (0) being the averages of respective quantities A ZZ , A nd and Δ k (0) calculated for all members of each family. These data are characteristic for each family. Moreover, the next characteristic property of the above families is the absolute value of the cosine between the z-axis in the NV PAS and the z-axis of the 13 C PAS given by the element U ( ) k Zz of the respective unitary matrix U k that transforms the DFT simulated HFI A k . The averages of these directional cosines over family members are shown in the sixth column of table 1. One can see from these data that only the on-NV-axis 13 C position exhibits exactly U ( ) k Zz = Zz cos( ) k = 1 indicating that this 13 C n-spin has a quantization axis coinciding with that of the e-spin of the NV center. It should, however, be noted that our cluster simulation predicts that there are also few families (e.g. C, H, T) for which the principle z k axes of 13 C nuclear spins are aligned very closely to the NV symmetry z-axis.
Analysis of the spatial locations of 13 C positions belonging to specific families showed that all of them have near equal Z coordinates and are also situated in the plane near perpendicular to the z-axis, so that they are near equidistant from it. These calculated Z coordinates and distances from the z-axis averaged over family members are given in the seventh and eighth columns of table 1, respectively. Note that the coordinate origin for the NV PAS was set by computer after relaxation of the cluster and was disposed approximately at the position of the vacancy of the NV center. The N atom of the NV center in this PAS has the coordinates X N = 0.001 Å, Y N = 0.002 Å and Z N = 1.731 Å. In turn, the coordinates of the nearest on-NV-axis 13 C atom are as follows: X C = −0.006 Å, Y C = −0.006 Å and Z C = −4.734 Å. Finally, the ninth column in table 1 shows the averaged distances from the N atom to the positions of family members which are also characteristic of each family. Note that columns 6-9 of table 1 give the data for specification of 'shells' and 'cones' used recently in [84] to discuss the contribution to FID dynamics of various 13 C n-spins disposed differently with respect to the NV e-spin. Additionally, to demonstrate situations when a few families exhibit very close values of Δ k (0) Figure 2. Calculated values of the HFI splittings Δ k (0) for families A-N (a) versus those measured in [27] (b). 2(b) is reproduced with permission from [27]. Copyright 2012 by the American Physical Society.
but have different positions with respect to the NV center, we single out in table 1 three pairs of families K1-K2, O1-O2 and Z1-Z2 that have such properties. Also, we do not show the largest HFI splittings of ∼130 MHz for the three sites that are NNs of the vacancy, as they are well known from the literature.

Comparison with experiment
In figure 2 we compare the calculated values of HFI splittings Δ k (0) shown in table 1 with those experimentally measured in [27]. Both figures clearly demonstrate discrete values of possible HFI splittings Δ k (0) , corresponding to different families. Typically, the calculated values of HFI splittings were slightly lower than those experimentally measured. As it was pointed out in section 2 (see also [71]) the basic reason for that is the level of theory used here to simulate the HFI characteristics for the studied cluster. It should also be noted that the accord of the theoretical predictions with the experimental data could be partially improved if we took into account that the experimental data were obtained at a low magnetic field of B ∼ 20 gauss which gave a small additional contribution ∼0.04 MHz to the measured splitting.
From figure 2 one can see also that experimentally obtained values for 14 NV-13 C spin systems, attributed to a specific family in some cases (e.g. those of the families D, G and I), deviate from theoretical predictions, probably due to a short experimental data set in [27]. Nevertheless, qualitative near-coincidence of the two figures demonstrates that the HFI parameters simulated by DFT for the C 291 [NV] -H 172 cluster in conjunction with subsequent spin Hamiltonian calculations provide a reasonable fit with the experimental HFI splittings, allowing one also to identify possible positions of the 13 C atom in the diamond lattice as belonging to a definite family.
Evidently, achieving this last goal would depend on the experimental frequency resolution. Therefore, in figure 3 [27] while figures 3(b) and (c) demonstrate the distributions for the case of much higher resolution, h = 10 kHz. One can see from figure 3 that in the low-resolution case a measured value of Δ k (0) only allows the assignment of the respective 14 NV-13 C system to a definite family for rather large values of HFI splittings Δ k (0) (Δ k (0) > 1 MHz). In turn, high-resolution experiments distinguishing values of Δ k (0) with an accuracy of 10 kHz provide much more detailed information regarding the location of the 13 C atom in the diamond lattice with respect to the NV center. For example, one can see from figure 3(c) that in the case of the specific 13 C atom located on the NV axis and characterized by the value Δ k (0) = 187.4 kHz which falls within the frequency interval [180][181][182][183][184][185][186][187][188][189][190], there are three other positions providing HFI splittings Δ k (0) that also fall into the above interval. Evidently, in these cases additional experiments are needed to differentiate the 13 C atom within a definite frequency interval and determine the family it belongs to. Figure 3(c), which differs from figure 3(b) by the x-axis scale, shows that in the case of the distant 13 C atoms weakly coupled to the e-spin of the NV center there are typically rather many (10-20) different positions for the 13 C atom, demonstrating HFI splittings Δ k (0) from the same frequency interval, and making it difficult to identify their positions in the diamond lattice solely by measuring the value of HFI splitting Δ k (0) at zero magnetic field. Experiments in diamond samples with a lower 13 C density than natural samples where * T 2 can be high enough to resolve the hyperfine structure with high accuracy will be ideal for spectroscopy of 14 NV-13 C systems. It should be noted, however, that in such samples the probability of finding an NV with 13 C spin in a specific position is reduced and time-consuming systematic work is required to find the desired system 14 NV-13 C among others.
It would be instructive to demonstrate the ability to identify the position of a 13 C n-spin in a diamond lattice having measured the six-line ODMR spectrum of the 14 NV-13 C spin system. In particular, available experimental data can be found in the work [27] where the magnetic- manifold of the studied 14 NV-13 C spin system, in comparison with the experimental data of [27]. Curves 5 and 6 correspond to the allowed transitions 14-3 and 13-1 in the = ↔ = + m m 0 1 S S manifold, which were not studied in [27]. Left and right parts of (d) exhibit, respectively, the experimental and simulated ODMR spectra for the studied 14 NV-13 C spin system undergoing magnetic field B||OZ of three different strengths, with the lines indicating the theoretical parts of the transitions. See more details in text. Experimental data in (d) are courtesy of Vincent Jacques, ENS de Cachan.
field-dependent ODMR spectra of some specific 14 NV-13 C systems were carefully studied in the = ↔ = − m m 0 1 S S manifold. The low-field (20 gauss) HFI splitting for the studied spin system was 1.12 MHz. According to table 1, one can attribute the 14 NV-13 C system to the families K1 or K2, each having three equivalent positions. Additional data for the attribution gave the studied B-field dependence of the ODMR spectra, presented in figures 3 and 4 of article [27], which made it possible to determine experimentally the values of A ZZ C ( ) = 1.02 MHz ( )2 1/2 = 0.51 MHz. According to table 1, the finding A nd C ( ) = 0.51 MHz is in agreement with the K1 family but not with the K2 one. Therefore, we conclude that the center 14 NV-13 C that was carefully investigated experimentally in [27] belongs to the K1 family.
In more detail, for one representative member of the K1 family our simulation gave the following HFI matrix  (1) one can straightforwardly find exact eigenvalues (energy levels) α E and eigenfunctions Ψ α of the analyzed 14 NV-13 C spin system in the presence of an arbitrary magnetic field. Comparing them with those obtained using simple approximate analytical expressions (4), (5) one can make sure that these expressions work very well for the studied spin system. Further, we have simulated ODMR spectra and studied their modification with an applied magnetic field. The results of our simulations for the above specific 14 NV-13 C spin system are summarized in figure 4 in comparison with the available experimental data from [27].
Spatially, three equivalent 13 C n-spin positions belonging to the K1 family are fourth neighbors of the vacancy of the NV center. Their positions in the diamond lattice with respect to the NV center are illustrated by figure 4(a), which relates to the studied relaxed cluster C 291 [NV] -H 172 and shows side (upper figure 4(a)) and top (lower figure 4(a)) views of essential carbon atoms around the NV center. Here, numeration of the atoms is as follows: 1 is the N atom of the NV center; 2, 3 and 4 are NNs of the N atom; 5, 6 and 7 are NNs of the vacancy of the NV center; 8, 14, 17, 18, 21, 23, 24 and 25 are the second neighbors of the vacancy; 9, 11, 13, 16, 19 and 22 are the third neighbors of the vacancy; and finally, 10, 12 and 20 are the three positions belonging to the K1 family. One can see from the lower part of figure 4(a) that the last three positions are near-equidistant from the C 3V symmetry axis of the NV center (r XY = 2.985 Ǻ according to table 1) which passes through the N atom and is perpendicular to the plane of the lower part of figure 4(a). The distance from the N atom to the plane wherein lie the K1 family members is 3.849 Ǻ.
As follows from the above general analytical consideration of eigenvalues and eigenstates (see equations (4), (5)) for the chosen 14 NV-13 C spin system belonging to the K1 family and exhibiting negative element A ZZ C ( ) = −0888 MHz, at low magnetic field the transitions 10-1, 9-3, 14-3 and 13-1 will be allowed while the transitions 9-1, 10-3, 14-1 and 13-3 will be forbidden. In figure 4(b), which shows part of the simulated energy levels (those with the m I N ( ) = +1 nitrogen atom nuclear spin projection) of the studied 14 NV-13 C spin, these transitions are shown with thick solid vertical and thin dotted inclined lines. Note the difference between these transitions and those shown in   gauss and then further diverge so that finally the frequency distance between the respective two lines will approach A ZZ = 0.888 MHz. In turn, in the upper = ↔ = + m m 0 1 S S manifold one will observe at increased magnetic field B||OZ only two lines of the HFI structure associated with the allowed transitions 14-3 and 13-1 whose frequency difference is changed only slightly from Δ (0) at low B field to A ZZ at high B.
Finally, the left and right parts of figure 4(d) show a comparison of experimental ODMR spectra obtained in [27] for the = ↔ = − m m 0 1 S S manifold with fixed m I N ( ) = +1 nuclear spin projection at three different magnetic fields B||OZ = 363, 598 and 785 gauss with respective spectra simulated for the chosen representative of the K1 family. In theoretical spectra the only fitting parameter was the width of Lorentzians Γ which was taken to be equal to 0.12 MHz. One can see that the simulated spectra are in good coincidence with the experimental ones. Concluding this section, it should also be noted that analogous analysis done for the system 14 NV-13 C consisting of the 13 C atom in the on-NV-axis position reveals no forbidden transitions at any magnetic fields B||OZ, as should be clear from figure 1 which shows that, for this system, all combining substates are pure states with definite spin projections m I C ( ) .

Summary
In summary, we have systematically studied HFI and ODMR spectra of three-spin systems 14 NV-13 C including the NV electronic spin S = 1, the intrinsic nuclear spin I N ( ) = 1 of the 14 N nuclei belonging to the center and the nuclear spin I C ( ) = 1/2 of the isotopic 13 C atom disposed in any possible position in the relaxed H-terminated C 291 [NV] − H 172 cluster, hosting the NV center. Using DFT we have calculated full HFI matrices A KL C ( ) (in the NV PAS for coordinates) for the 14 NV-13 C spin systems which were substituted into a standard ground-electronic-state spin Hamiltonian to calculate energy levels and eigenstates of the systems and further simulate ODMR spectra using ODMR linewidths as the only fitting parameters. In particular, HFIinduced splittings of the m S = ±1 states at zero magnetic field were calculated for all positions of the 13 C nuclear spin in the cluster, along with the values of the most essential HFI parameters . We tabulated these data for 30 families of equivalent positions of the 13 C atom in a diamond lattice around the NV center, each consisting of three or six members due to the C 3V symmetry of the NV center. It was shown that the members of each family are disposed symmetrically on cones, each of which has its vertex at the N atom of the NV center (the axis coinciding with the NV symmetry axis) and is equidistant from the N atom of the NV center. For the first time HFI data are obtained for the specific system 14 NV-13 C in which the 13 C atom is disposed on the NV center symmetry axis in the nearest-to-the-vacancy position. At zero field the HFI splitting of the m S = ±1 states for this specific spin system was found to be 187.38 kHz. Simple analysis of an arbitrary 14 NV-13 C spin system in the axial magnetic field B||OZ done within the secular approximation showed that the effect of an external magnetic field on the rates of EPR transitions in the system is different in the data from the previous work [27] which is described well using simulated HFI parameters. Our results may contribute to the general understanding of NV-13 C interactions in diamond, to the realization of multipartite entangled states in multi-nuclear spin systems, to the implementation of single-shot read-out measurements and to the understanding of details of the NV electron spin decoherence under a fluctuating 13 C nuclear spin bath, especially in the presence of an external magnetic field.