The bulk conversion depth of the NV-center in diamond: computing a charged defect in a neutral slab

The negatively charged nitrogen vacancy (NV−) center in diamond has properties that make it a promising candidate for applications such as a qubit in room temperature quantum computing, single-molecule photoluminescence and NMR sensor, and as a single photon source for quantum cryptography. For many of its uses it is desirable to have the NV-center close to the diamond surface. In this work, we use density functional theory simulations to investigate how the distance of the NV− center to a surface, and its orientation, affect its properties, including the zero-phonon-line. We study the three technologically important surfaces terminated with fluorine, oxygen/hydroxyl and nitrogen. Since the NV-center is charged it requires special measures to simulate within a slab-model. We use the recently proposed charging with a substitutional donor in the diamond lattice resulting in a neutral super-cell, which provides very satisfactory results. We have found that the NV-centers properties converge to bulk values already at 5 Å depth.


Introduction
Engineering of defects close to surfaces is an ongoing field of research, and is vital for several different applications. The negatively charged nitrogen-vacancy center in diamond has remarkable properties and is especially well suited to many applications [1]. One of those applications is nanoscale sensors of many physical quantities, such as magnetic [2][3][4] and electric fields [5,6], magnetic moments of electron and nuclear spins [7][8][9][10], strain [11], temperature [12][13][14][15][16] and pressure [17], to mention a few. Another application is as qubits in quantum information processing at room temperature [18][19][20][21][22]. It could also be used as a cellular biomarker in medicine [23,24]. The properties that make these applications possible include long spin coherence time T 2 (up to 2 ms in bulk diamond), photostablity at room temperature, magnetic field dependent spin sublevels and spinflipping during the optical excitation and decay cycle.
To create NV-centers in a nitrogen containing diamond sample, one bombards with particles such as electrons, ions, protons or neutrons, and subsequently anneal the sample at temperatures around 800°C-1000°C [25]. The bombardment creates vacancies and the annealing makes them diffuse towards the nitrogens. This is the preferred technique if one wants large ensembles of NV-centers, to use in ensemble magnetometry for example [26,27]. However, using high-density nitrogen diamonds does not allow control over single NVcenters. And diamond with high nitrogen content also contains other defects which affects the properties of NVcenters, such as lifetime, linewidth and decoherence of spin. Another approach is to use an ultra-high purity diamond sample, that contains almost no nitrogen atoms or other defects (1 ppb) [28]. One can then implant this pure sample with nitrogen atoms, followed by annealing to diffuse the created vacancies toward the nitrogen atoms, thereby creating NV-centers. This will, however, also create other defects such as interstitial nitrogens (N s ), divacancies (V 2 ), N 2 V, N 2 and NVH [29]. In this way one can create a diamond with sparsely distributed NV-centers, which can be used in applications such as qubits in quantum computers. For some applications it is necessary to reduce the amount of 13 C in the growth of the diamond sample, since they provide a mechanism for nuclear spin interactions, which in turn leads to spin decoherence for the NV-centers [30].
In addition to the issues discussed above, you must also take into consideration the surface, as NV-centers in many cases need to be placed near a surface to couple for example to optical fibers or external cavities, or for nanoscale NMR sensing with target spins on the surface. Results have shown that the vertical distribution of NVs depends strongly on the nitrogen concentration in the diamond sample [31]. There are several challenges with near-surface (<5 nm) NVs: lowered NV -/NV 0 ratios (charge state instability), significantly broadened ESR linewidths, and shortened T 2 [32]. For example, NV-centers at a depth of <10 nm was found to have a T 2 of less than 100 μs [33,34], while bulk NV-centers (in 12 C enriched diamond) have a reported T 2 up to 2 ms [30,35], because of the surface induced noise [33,34,36]. One way to lower the influence of the surface is to minimize undesirable near-surface spins by surface oxidation methods, such as annealing in a dry oxygen environment [37,38] and oxygen soft plasma treatment [39]. Also, NV-centers less than 5 nm from the surface has been found to exhibit photo-blinking; the photoluminescence disappear and reappear suddenly [40]. This occurs because of charge state instability where the negatively charged NV centers lose their electrons to the surface and get converted to neutral NV centers. The reason for the conversion are electron traps on the surface, leading to upward band bending beneath the surface [41,42]. Therefore, oxygen-terminated (high electronegativity, negative electron affinity at surface) surfaces are preferred for the charge stabilization of shallow NVcenters [41,43]. Recently, Yamano et al stabilized the charge state of very shallow NV-centers (2.6±1.1 nm from surface) created by 1.2 keV, 10 9 cm −2 fluence nitrogen ion implantation into 12 C-enriched high-purity diamond films. The implantation was followed by formation of an oxygen-terminated surface [44].
We charge the NV-center by an additional donor-N, thereby achieving an uncharged supercell, as described in our previous work. This charging method has been shown to render the correct properties of the NV -(minus) charge state, and result in zero spin density on the additional donor nitrogen [45,46]. A charged system must be treated with great care since the pseudopotential expressions for the total energy are derived in a chargeneutrality condition [47]. A charged supercell is problematic as the electrostatic energy would be divergent. This is conventionally solved by imposing an artificial background jellium charge of opposite charge density that exactly compensates for the net charge in the supercell [48]. However, the jellium charge will produce another problem due to the finite-size effect and periodic boundary conditions: the charged defect in the supercell will suffer a spurious long-range Coulomb interaction between the periodic charge images and the artificial background charges. To overcome this problem, many correction schemes have been proposed in recent years [49][50][51][52]. These correction schemes have successfully worked for materials in bulk phase. They might however be insufficient for charged defects in two-dimensional systems such as thin-films and surfaces. The compensating background charge leads to a dipole, which diverges as the vacuum between the slab images increases [53,54]. Several different methods were proposed to solve this spurious electrostatic interaction for calculating charged defects in low-dimensional structures [53,[55][56][57]. By charging the NV-center with a donor-N, resulting in a neutral slab, we avoid correction schemes and the undesired arbitrary dependence on the thickness of the vacuum.
There have been several earlier theoretical works done on the NV-center in slabs, focusing on how different surface terminations affect the properties of the NV-center [58][59][60]. Some terminations investigated are H (hydrogen), F (fluorine), O (oxygen), OH (hydroxyl), N (nitrogen), N/H (half of the surface carbon atoms terminated with N, the other half with H), H/O/OH (mixture of hydrogen, ether, and hydroxyl terminators), and also clean surfaces. In this work, we investigate how the negatively charged NV-center in diamond is affected by the close proximity to surfaces of F, H/O/OH, and N terminations, using an uncharged supercell. The H/O/ OH terminated surface remove the surface states from the bandgap, and have a predicted slight positive electron affinity of 0.5 eV, while the F-and N-terminated surfaces have shallow surface states (see supplementary information is available online at stacks.iop.org/NJP/21/053037/mmedia for details of energy level positions and distribution near the band edges and surface layers), as concluded in previous works [58,61]. The F-terminated surface have a positive electron affinity of 3.0 eV [62], while the N-terminated surface is predicted to have a positive electron affinity of 3.46 eV [61]. We look at how the zero-phonon-line (ZPL), the defect energy levels, and the spin-density are affected by the NV-center's distance and orientation to these surfaces. We investigate NV-centers at very shallow depths, between 0.1 and 1.4 nm from the surface. We determine at what distance the NV-center properties are significantly perturbed by the proximity to the surface.

Method
All our calculations were carried out with the Vienna ab initio Simulation Package (VASP5.4.1), within the density functional theory framework, using projector augmented wave potentials. Using the supercell method with periodic boundary conditions, the lattice parameter for a diamond primitive cell (2 atoms) was optimized From the unit cells of the slab models we constructed (6×6) supercells (see figures 1 and 2), resulting in a lattice constant (width) of 15.2 Å. The supercells contained 1080(72), 1080(144) and 1008(72) carbon (termination) atoms for the F-, H/O/OH-and N-terminated slabs, respectively. The electronic structure for the slabs was also calculated with the hybrid functional HSE06 (Heyd-Scuseria-Ernzerhof), using the PBEoptimized coordinates. For the large supercell slab models we have used only the gamma-point. HSE06 have been shown to give more accurate transition and band gap energies than PBE. Full HSE06 geometry optimizations would be too expensive on our large slab models, and our tests on bulk models of the defect (64 atom supercell) showed that a full HSE06 optimization has a rather small effect on the transition energies (less than 0.07 eV). Visualization of the slabs and their spin densities were done using VMD [63].
We achieve a neutral supercell by using an additional substitutional N to donate an electron to the NVcenter, thereby giving it a negative charge state. This charging method has been shown to give the correct structure, energetics and spin density of the NVcharge state, and result in zero spin density on the donor-N [45]. This donor-N was placed 7.7-8.1 Å from the NV-center, and at the same depth from the top layer as the NV-center to avoid the creation of an artificial dipole interaction across the layers of the slab. In these relaxed supercells containing NV-centers we investigate the properties the ZPL excitation energy, eigenvalues in the Kohn-Sham bandgap, and spin-density. The ZPL was calculated within the Franck-Condon approximation (see supplementary information for computational details and pictorial presentations).
We note that one of those models required about 30 000 CPU hours (including HSE06 calculations), using 40 supercomputer nodes with a total of 1120 cores (28 cores/node). This resulted in a total time of over 1 million CPU hours for all the models.

Results and discussion
We have simulated the NV-center in these diamond slab models (see figures 1 and 2), and we have varied its distance and configuration to the surfaces, in order to investigate how the surfaces and their termination affect the transition energies of the NV-center. We name the configurations as N(x)V(y), where N and V is replaced by numbers representing the carbon plane parallel to the surface the nitrogen and the vacancy of the NV-center is placed in, respectively (starting at 1 for the first plane at the surface, see figure 2). The x and y is the plane perpendicular to the surface the nitrogen and the vacancy is placed in (see figure 3), respectively (denoted as i, j, k, l for the F-and N-terminated surface, and i, j for the H/O/OH-terminated surface).
We have investigated how the NV-center ZPL excitation energy, eigenvalues in the Kohn-Sham bandgap, and spindensity, are affected as we vary its distance and configuration relative to the F-, H/O/OH-, and N-terminated surfaces. PBE-nitrogen atom closest to surface From inspection of the black points in figure 4 we immediately note that the ZPL of the NV-center, when placed more than 3 Å from the F-terminated surface, is not greatly affected regardless of configuration. For NVsurface distances (distance to the first layer of carbon atoms below the termination, see figure 2) below 3 Å, the ZPL is affected by the surface, and the configuration it is placed in determine how it changes. For the 4th layer below the surface, we see that for the 4(k)5(l) configuration the ZPL energy increases by 9.2% (from 1.694 eV (average of the 15(k)16(l) and 15(i)16(i) configurations), to 1.850 eV). However, comparing with the 4(i)5(j) configuration we instead see a decrease by 15.6% (to 1.429 eV). For the 3rd layer in the 3(i)4(i) configuration, the ZPL energy increases only by 1.9% (to 1.727 eV) from the bulk converged value. Also in the 3(k)4(k) configuration, the ZPL energy decreases by a modest 7.2% (to 1.574 eV). We, thus, find that the ZPL energies are  . PBE computed ZPL transition energy for different distances between the NV-center and the surface for our F-terminated surface, for each N(x)V(y) configuration (see figure 3(a)). The numbers in the diagram represent the layer that the N of the NV-center is placed in (see figure 2). The black and red points represent the NV-center when placed with the nitrogen atom and the vacancy closest to the surface, respectively. closer to the converged value in the 3rd layer than in the 4th layer, despite the 3rd layer being closer to the surface, because of the difference in configuration. Looking at the 2nd layer in the 2(j)3(k) configuration, we see that the ZPL energy decreases by 11.5% (to 1.500 eV). This is deviating more than for the 3rd layer, but still not as much as the largest deviation in the 4th layer. From our PBE calculations we can conclude that the ZPL is affected when the NV-center is placed less than 3 Å from the surface, and that it is especially affected in the 4th layer, with a maximum perturbance of 15.6%.

PBE-vacancy closest to surface
These configurations are represented by the red points in figure 4, and we have found that the ZPL is not greatly affected for N-surface distances larger than 6 Å. For distances closer than this, we see a decrease in the ZPL for the configurations investigated. When placing the NV-center (the N of the NV-center) in the 6th layer below the surface in the 6(l)5(l) configuration, the ZPL decreases by 6.9% (from 1.694 eV (average of the 15(k)16(k) and 15(i)16(i) configurations to 1.576 eV)). Placing the NV-center in the 5th layer in the 5(l)4(k) configuration reveals a decrease of the ZPL by 16.8%. When placed in the 4th layer in the 4(k)3(k) configuration the decrease in ZPL is 7.1%. And finally, when the NV-center is placed in the 3rd layer in the 3(k)2(j) configuration, the ZPL decreases by 19.6%. From our PBE calculations we can therefore conclude that the NV-center is perturbed by the surface significantly starting from layer 6.

HSE06-nitrogen atom closest to surface
Comparing the black points in figure 4 (PBE) and figure 5 (HSE06) we see that they have fairly similar behavior. For NV-surface distances larger than 4 Å, the ZPL is not greatly affected. For closer distances than 4 Å, we do get a considerable effect. A difference that stands out between the PBE and HSE06 computations is the ZPL of the NV-center placed in the 5(l)6(l) configuration in the 5th layer. For the PBE computations, the ZPL is almost unaffected, while for the HSE06 computations the ZPL is affected considerably, with an increase of 13.2% (from 1.961 eV (average of the 15(i)16(j) and 15(k)16(k) configurations), to 2.220 eV). For the NV-center in the 4th layer, we see a similar behavior as the PBE calculations: Placing the NV-center in the 4(k)5(l) configuration increases the ZPL energy with 10.4% (to 2.164 eV), while placing the NV-center in the 4(i)5(j) configuration instead decreases the ZPL energy with 13.8% (to 1.690 eV). For the 3rd layer in the 3(i)4(i) configuration the ZPL energy increases by 0.5% (to 1.971 eV), while it decreases by 8.7% (to 1.790 eV) in the 3(k)4(k) configuration. Comparing the 3rd and 4th layer in our HSE06 computations, we see that the 3rd layer have less diverged ZPL energies than the 4th layer, despite being closer to the surface, which is the same result we arrived at in our PBE calculations. From our HSE06 calculations we can conclude that the ZPL is affected when the NV-center is placed less than 4 Å from the surface, and that it is especially affected in the 4th and 5th layers.

HSE06-vacancy closest to surface
The corresponding comparison of the red points in figure 4 (PBE) and figure 5 (HSE06) we see that they have a similar behavior. The ZPL is not greatly affected for N-surface distances larger than 6 Å, which agrees well with our PBE calculations. Placing the NV-center in the 6th layer below the surface yields a decrease of the ZPL by 6.7% (from 1.961 eV (average of the 15(i)16(j) and 15(k)16(k) configurations) to 1.829 eV). When placed in the Figure 5. HSE06 computed ZPL transition energy for different distances between the NV-center and the surface for our F-terminated surface, for each N(x)V(y) configuration (see figure 3(a)). The numbers in the diagram represent the layer that the N of the NV-center is placed in (see figure 2). The black and red points represent the NV-center when placed with the nitrogen atom and the vacancy closest to the surface, respectively. 5th layer in the 5(l)4(k) configuration, the decrease in the ZPL is 15.8%. And when placed in the 4th layer in the 4(k)3(k) configuration, the ZPL decreases by 7.2%.
We have found that the placement of the N or the vacancy of the NV-center closest to the surface does not significantly alter the dependence of the ZPL on the N-surface distance. In addition, we have found that the bulk conversion depth of the ZPL is qualitatively described well with the PBE-functional, but that quantitative results are only achieved with the HSE06-functional. Thus, for the H/O/OH-terminated and N-terminated surfaces we probe the bulk conversion using configurations with N closest to the surface, and we only report the results obtained with the HSE06-functional.
We find that the ZPL of NV-center is not greatly affected when placed more than 2 Å from the H/O/OHterminated surface, regardless of configuration. For NV-surface distances below 2 Å, the ZPL is noticeably affected by the surface. When placing the NV-center in the 3rd layer below the surface, we see that for the 3(j)4(j) configuration the ZPL decreases by 7.5% (from 1.959 eV (for the 15(j)16(j) configuration) to 1.813 eV). We have found similar results for the PBE-functional (see figure S2 in supplementary information).
Similarly to the two previous surface terminations we have found that for NV-surface distances larger than 2 Å, the ZPL is not greatly affected. While for distances below 2 Å, the ZPL is greatly affected. When placing the NV-center in the 3rd layer below the surface in the 3(i)4(j) configuration, we see that the ZPL decreases by 21.7% (from 1.963 eV (for the 14(k)15(k) configuration) to 1.536 eV). And when placing the NV-center in the 2nd layer below the surface in the 2(k)3(k) configuration the ZPL decreases by 24.2%. So, our HSE06 computations show that for NV-centers placed below 2 Å from the N-terminated surface, the ZPL changes more than 20%, thereby resulting in a strongly perturbed NV-center. This is in agreement with our PBE computations (see figure S3 in supplementary information).
Our overall finding is that the ZPL energy does not vary appreciably for NV-center placed as close as 4-5 Å from the surface, regardless of termination. For distances closer than this, the ZPL energy starts to fluctuate, and the fluctuation is dependent on the surface-termination and configuration of the NV-center. Placing the NVcenter closer than 2 Å from the surface results in a change greater than 20% for H/O/OH-and N-terminated surfaces. For the F-terminated surface, the ZPL never changes more than 20% regardless of the configuration. It is, however, important to bear in mind that this dependence is for perfect crystals free of additional defects and faults, even with regard to the surface and its termination. Thus, the results obtained only include the effect of the presence of a surface.  figure 3(b)). The numbers in the diagram represent the layer that the N of the NV-center is placed in (see figure 2).
Comparing the bulk conversion values for the different terminations (HSE06 computations), we get a value of 1.961 eV for the F-terminated surface (average value of the two configurations in the 15th layer), for the H/O/ OH-terminated surface we get 1.959 eV (for the configuration in the 15th layer), and for the N-terminated surface we get 1.963 eV (for the configuration in the 14th layer). This can be compared to the bulk value of 1.962 eV for a 512 atom neutral supercell (containing a NV-center and an electron donor nitrogen atom placed 12.23 Å from the NV-center along the symmetry axis of the supercell in the NV-N configuration) in our previous work [45]. As can be seen, the values are very close to each other, deviating 0.004 eV at most. Such small differences are of course important for sensor applications, and individual NV-center baseline ZPL setting/verification is needed to determine its unperturbed energy.
The study of the ZPL excitation energy shows the bulk conversion for an energy difference and will not reveal if there are unison changes or shifts in energy levels caused by the surface, which could be a more long-range effect. Thus, in the next section we have investigated how the eigenvalues of the spin down e 1 -and v-levels (corresponding to the vertical transition, see figure S1 in supplementary information) in the Kohn-Sham bandgap vary as you vary the NV-center depth (large supercell slab models with gamma-point evaluation).
F-terminated surface-defect energy levels v-level (HSE06) As can be seen from figure 8, the v-level deviates by less than 0.2 eV from the bulk value (average energy of the 15(k)16(k) and 15(i)16(i) configurations) for all configurations except 3(k)2(j) and 4(k)5(l). The largest deviation is 0.45 eV found for the 3(k)2(j) configuration. We also see that for NV-center placed in the same layer, the configuration influences the changes; In layer 5, placing the NV-center in the 5(l)6(l) or the 5(l)4(k) configuration, yields changes in energy of −0.03 and 0.16 eV, respectively. Similarly in layer 4, for the NV-center in the 4(k)5(l), 4(i)5(j) and 4(k)3(k) configurations, we find energy changes of −0.23, 0.07 and 0.11 eV, respectively. For layer 3, when placing the NV-center in the 3(i)4(i), 3(k)4(k) and 3(k)2(j) configurations, the change in energy is −0.10, −0.01 and 0.45 eV, respectively. We have found very similar results using the PBEfunctional (see supplementary information). e 1 -level (HSE06) In figure 9 we see that for NV-surface distances larger than 4 Å, the deviation to the bulk value for the e 1 -level (average energy of the 15(k)16(k) and 15(i)16(i) configurations) is less than 0.12 eV. For closer distances the largest deviation is found in the 3(k)2(j) configuration, where the energy increases by 0.28 eV. Also for e 1 -level, we find that the configuration affects the change in energy when the NV-center is placed in the same layer; for the 7th layer the 7(i)8(i) and 7(k)8(k) configurations result in energy changes by −0.06 or −0.02 eV, respectively. Looking at the 5th layer configurations 5(l)4(k) and 5(l)6(l) we find changes in energy by −0.24 and −0.06 eV, respectively. For the configurations 4(i)5(j), 4(k)3(k) and 4(k)5(l) in the 4th layer, the changes are −0.24, −0.01 and 0.04 eV, respectively. Looking at the 3(k)4(k), 3(i)4(i), and 3(k)2(j) configurations in the 3rd layer, we find energy changes of −0.21, −0.02 and 0.28 eV, respectively. We have found very similar results using the PBEfunctional (see supplementary information). HSE06 computed ZPL transition energy for different distances between the NV-center and the surface for our N-terminated surface, for each N(x)V(y) configuration (see figure 3(c)). The numbers in the diagram represent the layer that the N of the NV-center is placed in (see figure 2).

H/O/OH-terminated surface-defect energy levels v-level (HSE06)
For the v-level in the H/O/OH-terminated surface (see figure 10), we find that there is not a large deviation in energy for any of the configurations. The largest deviation was found for the 2(i)3(j) configuration, with a 0.2 eV decrease in energy. For the NV-center in the 3(j)4(j) and 4(j)5(i) configurations in layers 3 and 4, the deviation is only −0.05 and −0.09 eV, respectively. N-terminated surface-defect energy levels v-level (HSE06) Figure 12 reveals that the deviation in energy is less than 0.12 eV for all the configurations except the 3(k)4(l) configuration in the 3rd layer, where it decreases by 0.25 eV. We also see that when NV-centers is placed in the    figure 3(b)). The numbers in the diagram represent the layer that the N of the NV-center is placed in (see figure 2). The zero-level in energy is set to the value when the NV-center is placed farthest from the surface (15(j)16(j) configuration). Figure 11. HSE06 computed energy change for the spin-down e 1 -level in the Kohn-Sham bandgap as you vary the distance between the NV-center and the surface for our H/O/OH-terminated surface, for each N(x)V(y) configuration (see figure 3(b)). The numbers in the diagram represent the layer that the N of the NV-center is placed in (see figure 2). The zero-level in energy is set to the value when the NV-center is placed farthest from the surface (15(j)16(j) configuration). Figure 12. HSE06 computed energy change for the spin-down v-level in the Kohn-Sham bandgap as you vary the distance between the NV-center and the surface for our N-terminated surface, for each N(x)V(y) configuration (see figure 3(c)). The numbers in the diagram represent the layer that the N of the NV-center is placed in (see figure 2). The zero-level in energy is set to the value when the NV-center is placed farthest from the surface (14(k)15(k) configuration). same layer, the configuration affects the energy change, as it is 0.09 eV for the other 3rd layer configuration 3(i)4(j). e 1 -level (HSE06) We find that the deviation in energy is not large for any of the configurations except the 3(i)4(j) and 2(k)3(k) configurations, where it decreases by 0.46 and 0.67 eV, respectively (see figure 13). And for NV-centers placed in the same layer the configuration matters also for this energy level, as the 3(k)4(l) configuration only changes by 0.04 eV.
In summary, we find that the energy of the v-and e-levels does not vary appreciably (<0.1 eV) for NVsurface distances greater than 4 Å for any of the surface terminations. For closer distances, the energy-levels start to fluctuate, and the fluctuation is dependent on the surface-termination and the configuration of the NVcenter. The H/O/OH-terminated surface have the best characteristics, with very small energy deviations down to and including the 3rd layer (although the deviations is quite large for an NV-center in the 2nd layer). The Fand N-terminated surfaces have energy deviations (although not very large) for the 3rd and 4th layers. Most Figure 13. HSE06 computed energy change for the spin-down e 1 -level in the Kohn-Sham bandgap as you vary the distance between the NV-center and the surface for our N-terminated surface, for each N(x)V(y) configuration (see figure 3(c)). The numbers in the diagram represent the layer that the N of the NV-center is placed in (see figure 2). The zero-level in energy is set to the value when the NV-center is placed farthest from the surface (14(k)15(k) configuration). Figure 14. Spin densities for the 3 A 2 ground state and the 3 E excited state in the F-terminated slab when placing an NV-center at layers 3 and 15 (see figure 2) in the 3(k)4(k) and 15(k)16(k) configurations, respectively (see figure 3(a)). The spin densities for the N-and H/O/OH-terminated surfaces look very similar to the spin density for the F-terminated surface, and are therefore not shown. Cyan spheres represent carbon atoms, pink spheres represent fluoride termination atoms, while blue and red spheres represent nitrogen atoms and vacancies, respectively. All the carbon atoms in the slab (except at the surfaces) are not shown to make it easier to see the charge densities (isovalue=0.05 1/Å 3 ). Miller index for plane of paper: (1,1,0). importantly we have found that the closeness to the surface does not in itself introduce any drift in the position of the defect's energy-levels in the diamond band gap.

Bulk conversion depth for the spin-density
In addition to the ZPL transitions and the levels in the Kohn-Sham bandgap, we look at how the spin-density changes when varying the NV-surface distance for the F-, N-and H/O/OH-terminated surfaces. For the F-terminated surface we have found that the spin-density does not change appreciably when varying the NV-surface distance (see figure 14). The general appearance of the spin-density does not change considerably between the 3(k)4(k) and 15(k)16(k) configuration for either the ground state 3 A 2 (see figures 14(a) and (c)) nor the 3 E excited state (see figures 14(b) and (d)). Although the ZPE changed with −0.171 eV for the former. The same was concluded for the N-and H/O/OH-terminated surfaces. Since we have found no hybridization with the surface it won't change the overall shape, or localization, of the states of the optically active NV-center, even when the excitation energy is changing due to the presence of the surface. At least when the surfaces are defect free and without its own spin-density. It is also worth reiterating that the donor-N has no spin density at all [45], which is also verification of the success of the charging of the defect-center resulting in a neutral supercell.

Conclusion
We have used density functional theory to study the bulk conversion depth of the NV-centers ZPL excitation energy, energy levels in the Kohn-Sham bandgap, and spin-density for three different surface terminations (F-, H/O/OH-, and N-terminated). We conclude that all these properties converge already at 5 Å from the surface, which means that NV-centers placed in the 7th layer from the surface, or deeper, has very similar properties as a NV-center in the bulk. With an individual ZPL, v and e 1 deviation of less than 0.05, 0.05 and 0.12 eV, respectively. For NV-centers closer than 4 Å to the surface, the ZPL starts to fluctuate, and the configurational direction the NV-center is placed in relative to the surface also starts to play an important role. Placing the NVcenter closer than 2 Å from the surface results in a change in the ZPL greater than 20% for H/O/OH-and N-terminated surfaces. For the F-terminated surface, it never changes more than 20% regardless of the distance to the surface. We have found that the changes in ZPL for NV-centers close to the surface could go both up and down in energy, which is configuration dependent. No trends in these fluctuations could be found.
Our results are very encouraging for technological applications that require the NV-center to be close to the diamond surface. It should, however, be pointed out that our results have isolated the convergence to depend only on the surface in models that contain no other defects or flaws.
We have used a substitutional donor (acceptor) to simulate a charged defect in a slab-model to avoid introducing the ambiguities of using an added (subtracted) electron, and the accompanying counter-charge compensation. In our models the donor nitrogen is separated from the NV-center by 7.7-8.1 Å, resulting in a neutral slab (supercell) and properties in excellent agreement with experiment for the negatively charged NVcenter. Our way of charging the defect with a substitutional donor (acceptor) opens up for bulk conversion depth electronic structure computations for charged defects in any material using neutral slabs.