# Quantum embedding study of strain and charge induced Stark effects on the
NV- center in diamond

Gabriel I. López-Morales  Department of Physics, City College of the City
University of New York, New York, New York 10031, USA  Joanna M. Zajac
Instrumentation Division, Brookhaven National Laboratory, Bldg. 535 P.O. Box
5000, Upton, New York 11973-5000, USA  Johannes Flick  Department of Physics,
City College of the City University of New York, New York, New York 10031, USA
Center for Computational Quantum Physics, Flatiron Institute, 162 5th Avenue,
New York, New York 10010, USA  The Graduate Center of the City University of
New York, New York, New York 10016, USA  Carlos A. Meriles  Department of
Physics, City College of the City University of New York, New York, New York
10031, USA  The Graduate Center of the City University of New York, New York,
New York 10016, USA  Cyrus E. Dreyer  Center for Computational Quantum
Physics, Flatiron Institute, 162 5th Avenue, New York, New York 10010, USA
Department of Physics and Astronomy, Stony Brook University, Stony Brook, New
York, 11794-3800, USA

###### Abstract

The NV- color center in diamond has been demonstrated as a powerful nanosensor
for quantum metrology, due to the sensitivity of its optical and spin
properties to external electric, magnetic, and strain fields. In view of these
applications, we use quantum embedding to derive a many-body description of
strain and charge induced Stark effects on the NV- center. We quantify how
strain longitudinal to the axis of NV- shifts the excited states in energy,
while strain with a component transverse to the NV axis splits the
degeneracies of the ${}^{3}E$ and ${}^{1}E$ states. The largest effects are
for the optically relevant ${}^{3}E$ manifold, which splits into $E_{x}$ and
$E_{y}$ with transverse strain. From these responses we extract strain
susceptibilities for the $E_{x/y}$ states within the quasi-linear regime.
Additionally, we study the many-body dipole matrix elements of the NV- and
find a permanent dipole 1.64 D at zero strain, which is $\sim 30\%$ smaller
than that obtained from recent density functional theory calculations. We also
determine the transition dipole between the $E_{x}$ and $E_{y}$ and how it
evolves with strain.

##  I Introduction

The negatively-charged nitrogen vacancy (NV-) center in diamond is the most
well-developed platform for realizing color-center based quantum and nanoscale
technologies. As its name suggest, the NV- is formed through a combination of
a carbon vacancy and an adjacent nitrogen substitution, with an excess
electronic charge that typically comes from a nearby donor defect in the
diamond lattice. Extensive experimental and theoretical work on its mechanisms
of formation, charge-state stability and dynamics, as well as its interesting
spin properties, has led to a well-founded understanding about the overall
electronic structure of this color center in diamond [1, 2, 3, 4, 5, 6, 7, 8].

In particular, nanoscale sensing constitutes an area of applications where the
optoelectronic properties of the NV- have been exploited with great success
[9, 10, 11]. High-precision temperature sensing, magnetometry-based bio-
imaging, and noise spectroscopy at the nanoscale, are all recent applications
based on optical and microwave manipulation of the excited states of the NV-
center through its spin-selective transitions [12, 13, 11, 14, 15]. The NV- is
also useful in sensing local electric fields at the nanoscale with high
precision [16, 17, 18, 19, 20, 21, 22]. The lack of inversion symmetry in this
color center facilitates the emergence of static electric dipoles for some of
its electronic states [23, 4, 24], which interact with the external (e.g.,
local or applied) electric fields such that it results in spectral shifts of
optical transitions. This phenomenon is broadly known as the electric Stark
effect. The presence of local or macroscopically-applied strain fields has
somewhat equivalent effects on the NV-; i.e., they introduce strain-induced
Stark shifts or splitting of the optical lines [23, 4, 25, 24, 26, 27].
Altogether, these Stark effects in the NV- could be further harnessed to
implement novel NV-based spectroscopy and sensing protocols [25, 26, 28, 29],
including for high-energy detectors [30, 31, 32, 33, 34]

For use of NV- as a nanosensor, it is crucial to have a quantitative
understanding of the strain and electric field induced Stark shifts. The first
experimental study to estimate the electric-field susceptibility (permanent
dipole) of the NV-, provided 1.2–1.5 D [16]; more recently, this number was
estimated to have longitudinal and transverse components of 2.82 and 3.63 D,
respectively [21]. Studies employing the modern theory of polarization through
density functional theory (DFT) have found a permanent dipole of 4.33 D [24],
compared to 2.23–2.68 D obtained through direct computation of the electric-
field Stark shift in the NV- through DFT [35]. An application of group theory,
on the other hand, finds 0.79 D [23]. Other works employing few-level
Hamiltonians and molecular-orbital approaches to treat the properties of the
NV- under strain provided useful insights on the strain responses, but did not
provide explicit numbers for the dipole couplings [1, 2]. Although some of
these theoretical results are relatively close to experimental observations
and provide valuable information about Stark effects on the NV-, further work
is required to fully quantify its dipole couplings and understand the physics
underlying the evolution of the optical excited states of this color center
under minute external (e.g., strain, electric) fields.

In particular, this defect center features a non-negligible degree of
electronic correlations within its set of $sp^{3}$-hybridized localized
electrons, two of which are unpaired. It has become increasingly evident that
an adequate treatment of such correlations is required to quantitatively
describe the excited-state properties of the NV- [36, 37, 38, 39, 40]. In
these regards, a full ab-initio description of the many-body physics of this
defect center, most of which is central for applications such as nanoscale
sensing, is highly desirable and constitutes an area of ongoing research
efforts [36, 37, 38, 39, 40, 41, 42]. Recent work employing methods that
combine quantum chemistry and density functional theory (DFT) via quantum
embedding have demonstrated quantitative descriptions of the NV- excited
states, at convenient computational costs [36, 37, 39, 42]. Thus, developing a
quantitative description of strain- and charge induced Stark effects on the
optical transitions of the NV- through such methods has important implications
for NV-based nanoscale sensors.

In this paper, we a first-principles approach based on quantum embedding to
derive such a description. We study the NV- under strain to derive a many-body
description of strain susceptibilities in its optical excited states. We also
access the many-body dipole moments of the NV- center and strain induced
effects in the optical transitions of the NV center by considering many-body
transition dipoles between excited states.

The paper is organized as follows. In Sec. II we discuss the computational
methodology employed in this work and the computational details. Then, we
present the obtained results for strain susceptibilities and dipole couplings
in Secs. III.1 and III.2. This is followed by a discussion of the results in
Sec. IV, and the conclusions in Sec. V.

##  II Computational methodology

###  II.1 DFT for structure and strain

As a basis for the quantum embedding, we employ DFT calculations to obtain the
relaxed atomic structure and the corresponding electronic ground state of a
diamond supercell containing a single NV- center, with axis along the $\langle
1,1,1\rangle$. Calculations are performed using the 216-atom
(3$\times$3$\times$3) diamond supercell, which is found to yield converged
many-body states within our methodology [39]. To study the evolution of the
electronic states of the NV- under strain, we introduce distortions to the
supercell via direct elongation (contraction) of the lattice vectors for
tensile (compressive) strain. For simplicity, we focus on two cases of strain:
(i) along the $\langle 1,1,1\rangle$ axis ($\varepsilon_{\langle
1,1,1\rangle}$), (ii) along the $\langle 1,-1,1\rangle$ axis
($\varepsilon_{\langle 1,-1,1\rangle}$). These two directions of strain fields
are also experimentally relevant because the NV- is known to respond
differently under longitudinal and transverse strains [23, 2, 3, 24]. The
amount of strain applied on each component of the lattice vectors is taken as
$|\varepsilon|$, and defined with respect to the length of the pristine
lattice vectors.

###  II.2 Quantum embedding

From the ground-state DFT calculation of NV-, we downfold onto a minimal
active space containing only the defect states of NV- via Wannierization
(using the Wannier90 code) [43]. The down-folded problem is mapped onto the
effective Hamiltonian [36, 37, 39]

$\begin{split}\hat{H}_{\text{eff}}&=-\sum_{\braket{ij},\sigma}(t_{ij}c^{\dagger}_{i\sigma}c_{j\sigma}+h.c.)\\\
&+\frac{1}{2}\sum_{\braket{ijkl},\sigma,\sigma^{\prime}}U_{ijkl}c^{\dagger}_{i\sigma}c^{\dagger}_{j\sigma^{\prime}}c_{l\sigma^{\prime}}c_{k\sigma}\\\
&-\hat{H}_{DC},\end{split}$ (1)

where $c^{\dagger}$ ($c$) represent creation (annihilation) operators, and
$i$,$j$,$k$,$l$ label different correlated defect states within the active
space, with spin $\sigma$, $\sigma^{\prime}$. The first term in Eq. (1)
represents the single-particle part of the effective Hamiltonian and is
characterized by the hopping matrix elements $t_{ij}$ in the Wannier basis.
The hopping matrix elements are determined from the Wannierized basis
functions $\psi$ as

$t_{ij}=-\braket{\psi_{i}}{\hat{H}_{\text{KS}}}{\psi_{j}}.$ (2)

In the absence of disentanglement from bulk bands, the eigenvalues of $t_{ij}$
reproduce the KS eigenvalues exactly [44]. The second term in Eq. (1) captures
the effective electron-electron interactions between correlated states within
the active space, in this case represented via the four-body screened Coulomb
tensor $U_{ijkl}$. We calculate these effective interactions using the
constrained random-phase approximation (cRPA) [45], via

$U_{ijkl}=\braket{\psi_{i}\psi_{j}}{\hat{U}}{\psi_{k}\psi_{l}}.$ (3)

To do so, we evaluate a decomposed version of the full polarizability tensor
$\hat{\Pi}_{\text{full}}$ containing all screening effects except those among
states within the correlated subspace

$\hat{U}=[1-\hat{V}(\hat{\Pi}_{\text{full}}-\hat{\Pi}_{\text{NV}})]\hat{V}.$
(4)

Here, “NV” refers to the correlated subspace comprised of the defect states
associated with the NV- center, and $\hat{V}$ is the bare Coulomb interaction.
In this way, the correlated defect manifold described by Eq. (1) couples to
the diamond matrix through the screening effects included while evaluating the
Coulomb tensor. Screening effects and correlations within the active space of
the defect are included explicitly through the configuration interaction (CI)
by diagonalization Eq. (1) exactly (using the triqs software library [46]).
The last term in Eq. (1), so-called double-counting (DC) correction, is added
to deal with unwanted Coulomb interactions partially included (from DFT) in
the single-particle parts of the full effective Hamiltonian. Here, we adopt a
DC correction of the form [39, 37, 36]

$\hat{H}_{DC}=\sum_{\braket{ij},\sigma}c^{\dagger}_{i\sigma}c_{j\sigma}\sum_{\braket{kl}}P_{kl}(U_{iljk}-\alpha
U_{ilkj}),$ (5)

where, following Ref. [36], $\alpha$ relates to the amount of Hartree-Fock
exchange mixed to construct the DFT functional. Since we use the Heyd-
Scuseria-Ernzerhof (HSE) hybrid functional in all DFT calculations, we set
$\alpha=0.25$. Employing this DC scheme has shown to yield NV- excitation
energies in good agreement with those from experiments [36, 39]. In any case,
the quantities of interest in this work and their trends with strain are
quantitatively very similar to those obtained without DC corrections.

###  II.3 Many-body dipole couplings

To describe the dipole couplings to electric fields of the NV-, we evaluate
the dipole matrix elements across the many-body spectrum via

$\bm{\mu}_{ij}=\braket{\Psi_{i}}{{\hat{\textbf{r}}}}{{\Psi_{j}}},$ (6)

where, $\Psi$ denotes many-body wave functions of states $i$ and $j$
(respectively), and $\hat{\textbf{r}}$ represents the many-body dipole
operator (along the Cartesian directions). In our formalism, we construct the
many-body dipole operator as

${\hat{\textbf{r}}}=\sum_{nm}\textbf{r}^{\text{wann}}_{nm}({\bm{0}})c^{\dagger}_{n}c_{m},$
(7)

using the matrix elements of the position operator between Wannier functions
$m$ and $n$, given by

$\textbf{r}^{\text{wan}}_{nm}(\textbf{R})=\braket{n\textbf{0}}{\textbf{r}}{m\textbf{R}},$
(8)

as implemented in Wannier90 [43]. This approach allows for permanent and
transition dipole matrix elements between explicitly _many-body_ ground and
excited states of NV-. We note that some properties of the dipole matrix
elements might depend on the specific gauge chosen for constructing out
Wannier functions. To quantify this, we evaluated the dipole matrix (in the
absence of strain) using a Wannierization procedure slightly different from
that employed for the main results (Sec. II.4), i.e., increasing the
disentanglement window slightly, but other properties of out many-body model
(i.e., energies, spin state, symmetries) are not significantly altered. We
make use of these results to comment on possible uncertainties regarding the
absolute orientation for some of the dipoles in the NV- (see Sec. III.2).

###  II.4 Wannierization

In terms of the Wannierization scheme specifically implemented for the NV-
defect center, we rely on initial projections of $sp^{3}$-hybrid Wannier
functions centered around the NV to obtain the localized atomic basis
orbitals. This basis is utilized to evaluate the screened Coulomb interactions
and subsequently perform all cRPA-CI calculations. The single-particle
electronic states relevant for most of the optical excitations of the NV-
center are the $a_{1}(1)$, $a_{1}(2)$, $e_{x}$, and $e_{y}$. While the
$a_{1}(2)$, $e_{x}$, and $e_{y}$ lie isolated within the bandgap, the
$a_{1}(1)$ state is resonant with the valence band (VB) of diamond, requiring
disentanglement during the Wannierization procedure. Hence, we set our
disentanglement window to include the $a_{1}(1)$ state, while the $a_{1}(2)$,
$e_{x}$ and $e_{y}$ states are kept “frozen” such that Eq. (2) reproduces the
DFT eigenvalues exactly. All Wannierized states are maximally-localized by
minimizing the Wannier spreads. Overall, our active space is comprised of the
four $a_{1}(1)$, $a_{1}(2)$, $e_{x}$, and $e_{y}$ Wannierized orbitals and 6
electrons. This active space has demonstrated to be sufficient to yield
quantitative agreement between calculated and experimental excitation energies
for the NV- [36, 39]. The main task concerning our current work centers on
studying how this active space and its many-body states evolve under the
external influence of strain, as we discuss in Sec. III.

##  III Results

###  III.1 Strain susceptibilities of NV-

We start by considering the case of strain along the $\langle 1,1,1\rangle$.
Since this strain propagates along the high-symmetry axis of the NV- center
($A_{1}$ strain), the $C_{3v}$-symmetry of the defect is preserved. As a
result, the ${}^{3}E$ states remain degenerate, but shift in energy with
respect to the ${{}^{2}}A_{2}$ ground state as a function of
$\varepsilon_{\langle 1,1,1\rangle}$ [Fig. 1(a)]. The lack of inversion
symmetry in the NV center is reflected in the asymmetric evolution of the
excited states under $\pm$ $\varepsilon$ (elongation or compression). Overall,
these results indicate that in the scenario where local strain is only
longitudinal to the NV axis, a red-shift in the optical transition of the NV-
corresponds to local compressive strain, while tensile strain will manifest as
a blue shift in the optical transition. Lastly, Fig. 1(a) clearly illustrates
that the energy shifts are much larger for the ${}^{3}E$ manifold compared to
the singlet excited states, confirming that the triplet manifold is much more
sensitive to longitudinal strain.

Figure 1: Many-body spectrum of the NV- center in the presence of strain. (a)
$\varepsilon_{\langle 1,1,1\rangle}$ introduces energy shifts of the ${}^{3}E$
manifold with respect to the ${}^{3}A_{2}$ ground state. Energy shifts of the
other manifolds are comparatively much smaller. (b) $\varepsilon_{\langle
1,-1,1\rangle}$ introduces an energy split of the ${}^{3}E$ into the
${}^{3}E_{x/y}$ sub-manifolds. The ${}^{1}E$ manifold also splits under such
strain, while the ${}^{1}A_{1}$ state shows larger shifts compared to the
$\langle 1,1,1\rangle$-strained case.

Under $\varepsilon_{\langle 1,-1,1\rangle}$, the picture is more complicated.
Here, the strain has both longitudinal and transverse components to the local
NV- axis. The $C_{3v}$ symmetry is thus broken and the degeneracy of the $E$
states is lifted [Fig. 1(b)]. The evolution of the excited states under
$\varepsilon_{\langle 1,-1,1\rangle}$ still includes energy shifts caused by
the longitudinal strain components. Again, the ${}^{3}E$ states show the
largest strain induced effects. However, in this case the ${}^{1}E$ states
also display a relatively strong response. This splitting correlates well with
observations of split optically-detected magnetic resonance through the
singlet transition [47, 48]. In a sense, the singlet manifold could be an
interesting probe for strain sensing under specific conditions, e.g., in the
near-infrared spectral region, given its weak/strong response to
longitudinal/transverse strains [Figs. 1(a) and 1(b)]. However, the inter-
system crossing required for populating the singlet manifold is typically a
slow process, which, combined with a smaller radiative efficiency and longer
decay lifetime [47, 48], makes the photon collection in the singlet manifold
poorer compared to the triplet optical transitions. With this in mind, we make
the ${}^{3}E$ manifold the focus of the remainder of this paper.

Figure 2: Strain susceptibility ($\chi$) for the ${}^{3}E$ many-body states of
the NV- center. (a) ${}^{3}E_{y/x}$ under $\varepsilon_{\langle
1,1,1\rangle}$, covering a wide range of strain. (b) Left panel: Same as in
(a), but for the case of $\varepsilon_{\langle 1,-1,1\rangle}$. A set of
higher–resolution calculations for this case of strain (within shaded region)
is included in the right panel of (b), to better resolve the response within
the strain region most relevant to experiments ($\varepsilon<10^{-5}$). The
scatter and solid lines represent numerical data and second-degree polynomial
fits, respectively.

In Fig. 2, we make use of the strain evolution of the ${}^{3}E$ states to
extract the strain susceptibility (denoted here as $\chi$) of the optical
transitions in the NV center. As mentioned previously, $\varepsilon_{\langle
1,1,1\rangle}$ does not split the ${}^{3}E$ manifolds. In this case, $\chi$ is
identical for the ${}^{3}E_{y/x}$ states [Fig. 2(a)], and is found to follow
an almost-linear response, with linear coefficient $\chi_{\text{L}}=19.2$ eV.
The wide range of strain considered here allows to also capture the quadratic
component of the strain susceptibility $\chi_{\text{Q}}$, although the
realistic range of strain (typically $<10^{-5}$) renders this contribution
negligible within linear-strain regimes most relevant for experiments. In the
case of $\varepsilon_{\langle 1,-1,1\rangle}$, the resulting quadratic
susceptibility over the same range of strain is more pronounced [left panel of
Fig. 2(b)], due to splitting of the ${}^{3}E$ manifold.

The magnitude of strain fields experienced by NV- centers in typical
experiments is well below the strains considered up to this point, usually
below $\varepsilon=\pm 10^{-5}$ [18, 49, 50, 28, 51, 21, 22, 27]. Again,
within this range of strain, the responses to $\varepsilon_{\langle
1,1,1\rangle}$ for the ${}^{3}E_{y/x}$ states are effectively linear. Thus, to
quantify $\chi_{\text{L}}$ more reliably for transverse strain, we expand the
resolution by considering a narrower range of strains around zero
($\varepsilon<10^{-4}$) for $\varepsilon_{\langle 1,-1,1\rangle}$ [highlighted
region in the left panel of Fig. 2(b)]. The obtained results for these
calculations is presented in the right panel of Fig. 2(b). The fit to the
numerical data obtained in this case is also improved by exclusion of the
extreme points $\varepsilon>1.0\%$, which are in any case well outside any
realistic experimental situation. Compared to the left panel in Fig. 2(b), the
set of points around zero-strain resemble a clearer linear response and the
$\chi_{\text{L}}$ extracted from this fit should be closer to that expected
for calculations in the linear regime of strain.

Ideally, one would want to consider even weaker strain, such that the
calculated strain evolution of the ${}^{3}E_{y/x}$ states is explicitly
matched to that obtained experimentally. However, such strain fields require
very small atomic displacements, as well as the resolution of correspondingly
small energy splittings of the single-particle/many-body states which are
close to the numerical precision of our methods.

###  III.2 Dipole couplings of NV- under strain

Table 1: Summarized dipole couplings for the NV- under strain. The format of the dipole couplings presented below is “$\Delta p_{\text{x}}$, $\Delta p_{\text{y}}$, $\Delta p_{\text{z}}$, $|\Delta\textbf{p}|$” (or $\Delta\mu_{\text{x}}$, $\Delta\mu_{\text{y}}$, $\Delta\mu_{\text{z}}$, $|\Delta\bm{\mu}|$ for the transition dipole moment in the last column). The “$x$, $y$, $z$” represent directions along the supercell lattice vectors. In all cases, the applied strain is along the $\langle 1,-1,1\rangle$ direction (negative values correspond to compressive strain). $\varepsilon$ (%) |  $\braket{{}^{3}E_{y}}{{\hat{\textbf{r}}}}{{}^{3}E_{y}}$-$\braket{{}^{3}A_{2}}{\hat{\textbf{r}}}{{}^{3}A_{2}}$ |  $\braket{{}^{3}E_{x}}{{\hat{\textbf{r}}}}{{}^{3}E_{x}}$-$\braket{{}^{3}A_{2}}{{\hat{\textbf{r}}}}{{}^{3}A_{2}}$ |  $\braket{{}^{3}E}{{\hat{\textbf{r}}}}{{}^{3}E}$-$\braket{{}^{3}A_{2}}{{\hat{\textbf{r}}}}{{}^{3}A_{2}}$ |  $\braket{{}^{3}E_{y}}{{\hat{\textbf{r}}}}{{}^{3}E_{x}}$
---|---|---|---|---
-0.008 |  0.054, 2.7, 0.075, 2.7 |  1.8, -0.82, 1.8, 2.7 |  0.94, 0.94, 0.94, 1.6 |  1.5, 0.016, 1.5, 2.2
-0.004 |  1.5, -0.73, 2.1, 2.7 |  0.42, 2.6, -0.21, 2.7 |  0.94, 0.94, 0.94, 1.6 |  1.7, 0.55, 1.3, 2.2
0.00 |  1.7, 1.7, -0.55, 2.4 |  0.19, 0.21, 2.4, 2.4 |  0.94, 0.94, 0.94, 1.6 |  1.6, 1.6, 0.94, 2.4
+0.004 |  1.6, -0.80, 2.0, 2.7 |  0.27, 2.7, -0.12, 2.7 |  0.94, 0.94, 0.94, 1.6 |  1.6, 0.27, 1.4, 2.2
+0.008 |  1.8, -0.82, 1.8, 2.7 |  0.072, 2.7, 0.051, 2.7 |  0.94, 0.94, 0.94, 1.6 |  1.5, 0.023, 1.5, 2.2

Within our methodology, we can obtain a many-body description of the permanent
dipole and transition dipole moments (TDMs) of the NV-, as well as possible
strain induced effects on these dipole couplings. The diagonal elements of the
full dipole matrix described in Sec. II.3 [Eq. (6)] are used to extract the
permanent dipoles, $\textbf{p}_{i}\equiv{\bm{\mu}}_{ii}$. To quantify the true
permanent dipole character of the excited states of the NV-, we take the many-
body ground-state position as reference [24], i.e.,
$\Delta\textbf{p}_{y/x}\equiv\braket{{}^{3}E_{y/x}}{\hat{\textbf{r}}}{{}^{3}E_{y/x}}-\braket{{}^{3}A_{2}}{\hat{\textbf{r}}}{{}^{3}A_{2}}$
for transitions to individual ${}^{3}E$ states and
$\Delta\textbf{p}_{\text{avg}}\equiv\braket{{}^{3}E}{\hat{\textbf{r}}}{{}^{3}E}-\braket{{}^{3}A_{2}}{\hat{\textbf{r}}}{{}^{3}A_{2}}$
for the average. A summary of these results is provided in Table 1. First, the
individual permanent dipoles for the ${}^{3}E_{y/x}$ states are found to be
tilted with respect to the local axis of the NV- (here the $\langle
1,1,1\rangle$ axis, which we will refer to as NV${}^{-}_{z}$); the $x$ and $y$
components of $\Delta\textbf{p}_{y}$ and $\Delta\textbf{p}_{x}$ are anti-
parallel (i.e., point in opposite directions) with respect to NV${}^{-}_{z}$,
while the $z$ components of $\Delta\textbf{p}_{y}$ and $\Delta\textbf{p}_{x}$
are parallel and aligned along the NV${}^{-}_{z}$ axis. The net
$|\Delta\textbf{p}_{\text{avg}}|$ (averaged between ${}^{3}E_{x}$ and
${}^{3}E_{y}$) is thus fully aligned with the local NV${}^{-}_{z}$ (see
“$\braket{{}^{3}E}{{\hat{\textbf{r}}}}{{}^{3}E}$-$\braket{{}^{3}A_{2}}{{\hat{\textbf{r}}}}{{}^{3}A_{2}}$”
column in Table 1). In the absence of strain, the individual permanent dipoles
have magnitude of $|\Delta\textbf{p}_{y/x}|=2.4$ D, while
$|\Delta\textbf{p}_{\text{avg}}|=1.6$ D. As a result, if one considers the
unresolved ${}^{3}A_{2}\leftrightarrow$ ${}^{3}E$ optical transition, i.e.,
via off-resonant excitation and/or at room temperature, the net dipole is
effectively weaker. Furthermore, since the off-axis projections of
$|\Delta\textbf{p}|$ are lost for the unresolved excited manifold, the
possible sense of directionality from electric fields experienced by the NV-
should be substantially suppressed under such conditions. Overall, this would
translate into an enhanced electric-field sensitivity when addressing the
individual ${}^{3}A_{2}\leftrightarrow^{3}E_{y/x}$ transitions (resonantly),
which correlates very well with experimental observations [16, 17, 19, 20, 21,
22]. From our results, the “resonant enhancement” of $|\Delta\textbf{p}|$ is
$\sim 1.5$.

Of course, other effects, such as phonon scattering and local strain
environment, will likely impact the effective enhancement of the
$|\Delta\textbf{p}|$ in resonant experiments. In fact, by considering strain
effects explicitly, we find that $|\Delta\textbf{p}_{y/x}|$ increases for
$\varepsilon_{\langle 1,-1,1\rangle}$, but decreases with
$\varepsilon_{\langle 1,1,1\rangle}$ (not shown). Naturally, the amount of
enhancement/suppression of the permanent dipole depends on the magnitude of
local applied strain, but from the results in Table 1, it corresponds to
$<12.5\%$ for $|\varepsilon|\leq 0.008\%$. Since $\varepsilon_{\langle
1,-1,1\rangle}$ is not fully perpendicular to the local NV${}^{-}_{z}$ axis,
and introduces a transverse strain component, these results suggest that
$\varepsilon\perp$ NV${}^{-}_{z}$ could slightly enhance
$\Delta\textbf{p}_{y/x}$ [52]. Additionally, the average permanent dipole
remains constant with strain, which means that the relative enhancement of the
$\Delta\textbf{p}_{y}$ components with strain is anti-correlated with that of
$\Delta\textbf{p}_{x}$, and if at all present, will only be quantifiable via
resonant-excitation experiments. In terms of directionality, the obtained
results show that while the individual dipoles change directions with
$\varepsilon_{\langle 1,-1,1\rangle}$, they always maintain anti-parallel
alignment of their off-axis component, such that
$\Delta\textbf{p}_{\text{avg}}$ is always along the $\langle 1,1,1\rangle$
(Table 1). It is worth noting that the $|\Delta\textbf{p}_{\text{avg}}|$
magnitudes obtained herein are in excellent agreement with experimental
estimations [16, 21] (similar to that for the NV0 [53]), also falling close to
those derived through DFT [24, 35].

Since strain splits otherwise degenerate ${}^{3}E_{y/x}$ states, one other
dipole coupling of potential interest for strain sensing with the NV- is the
TDM between these states (denoted here $\bm{\mu}_{y/x}$), represented by the
off-diagonal elements of the full many-body dipole matrix. In the last column
of Table 1, we present such TDMs and how they evolve with strain. First, we
see that in the $\varepsilon=0$ case this TDM has the same magnitude as the
permanent dipole of the individual ${}^{3}E_{y/x}$ states. This TDM is also
comparable to that of the optical transition (${}^{3}A_{2}$ $\leftrightarrow$
${}^{3}E$), which we find to be $\sim$ 2.7 D (in the absence of strain). Of
course $\bm{\mu}_{y/x}$ does not have much physical relevance here since with
no strain to break the $C_{3v}$ symmetry, ${}^{3}E_{x}$ and ${}^{3}E_{y}$ are
degenerate to there will be no optical transitions.

Under the influence of strain, $|\bm{\mu}_{y/x}|$ shows the opposite trend
compared to the permanent dipoles, i.e., $|\bm{\mu}_{y/x}|$ is suppressed by
$\varepsilon_{\langle 1,-1,1\rangle}$ and enhanced by $\varepsilon_{\langle
1,1,1\rangle}$ (results for the latter case are not shown for brevity). The
maximum change in $|\bm{\mu}_{y/x}|$ for the strain range considered herein is
$\sim 8\%$. Interestingly, the direction of $\bm{\mu}_{y/x}$ is initially
($\varepsilon=0$) almost fully aligned with the $\langle 1,1,1\rangle$ (with
$\sim 5$-deg misalignment); under strain, one component of the TDM is
significantly suppressed. Though we caution that the details of the dipole
direction may be effected by the gauge choice for our Wannier functions (see
Sec. IV for more discussion), this corresponds to the expected behavior from
the symmetry of the system. The transverse strain breaks the $C_{3}$ rotation
symmetry of NV-, but keeps one mirror plane. Under the new $C_{s}$ point
group, states with $E$ symmetry split into $A^{\prime}$ and
$A^{\prime\prime}$. Based on the selection rules for optical transitions under
this symmetry, we expect the TDM perpendicular to the mirror plane to be
suppressed, which is approximately the $y$ direction in Table 1.

These last set of results hint at the ${}^{3}E_{y}\leftrightarrow^{3}E_{x}$
transition as an alternative candidate for strain sensing using the NV-. For
instance, provided there is enough strain to split the ${}^{3}E_{y/x}$ states
(which is true for the majority of NV centers in real experiments), the
response in $|\bm{\mu}_{y/x}|$ coupling with transverse and longitudinal
strain (opposite to that of $\Delta\textbf{p}_{y/x}$) could help identify how
much of each component is locally contributing to the overall strain
experienced by the NV-. Since typical ${}^{3}E_{y/x}$ splitting is on the
order of few GHz, the relatively large $\bm{\mu}_{y/x}$ should make it
straightforward to sense strain by monitoring the amplitude of this transition
under resonant excitation. Additionally, since this transition appears almost
fully linearly polarized, and its directionality reacts to the direction of
the overall strain, it could help decompose the local strain in its vector
components. Lastly, the measurements required to probe $|\bm{\mu}_{y/x}|$ are,
in principle, independent from those required to monitor the Stark shifts on
the ZPL through the permanent dipole. Thus, combining information from these
two complimentary dipole couplings should provide an enhanced probe for strain
sensing using the NV- center.

##  IV Discussion

The results described in Sec. III provide useful insights on the intricate
response of the NV- to local strain and electric fields. For instance, the
asymmetric evolution of the excited states with strain appears to be more
prominent when the strain has both transverse and parallel components to the
NV${}_{\text{z}}$ axis. This behavior, emerging from the absence of inversion
symmetry on the NV centre, should provide a means of differentiating local
compressive versus tensile strains. On the other hand, from the relative
energy shifts and splitting under $\varepsilon_{\langle 1,1,1\rangle}$ and
$\varepsilon_{\langle 1,-1,1\rangle}$ shown in Fig. 1, one could argue that
the singlet transition is more reliable at differentiating transverse versus
longitudinal strains than the triplet transition. This could be an interesting
application of the singlet transitions in the NV- for local strain sensing:
the triplet manifold reacts to both longitudinal and transverse strains (via
shifts and splitting of the optical line, respectively), while the singlet
reacts (splits) most prominently only to transverse components of strain. In
essence, correlated high-resolution measurements in which both the triplet and
singlet transitions are interrogated through resonant excitations should
provide useful information about strain propagation through diamond at the
nanoscale. Further, given the increasing interest in studying closely-
positioned NV- centers in diamond [54, 21, 22, 52], these proposed
measurements could serve as useful probes for possible NV–NV interactions
mediated via strain.

Another interesting aspect of the results in Sec. III is the non-negligible
susceptibility of some of the dipole couplings in the NV- to local strain. For
example, since the permanent dipoles dictate how the optical transitions of
the NV- couple to local electric fields, strain emerges as a natural means for
potentially enhancing its electric-field sensitivity. Furthermore, the
magnitude of the transition dipoles between ${}^{3}E_{y/x}$ states also appear
slightly susceptible to local strain, with a seemingly strong change in the
individual components (polarization). However, as mentioned in Sec. II.3,
these TDMs are likely to have a degree of gauge-dependence within our
methodology, resulting in possible uncertainties in their absolute
orientations. We can estimate such uncertainties by changing slightly the
choice of Wannier functions that construct our basis. For instance, increasing
the disentanglement window by 350 meV (in the absence of strain), we find that
the net permanent dipole remains almost unchanged in both direction and
magnitude (fully gauge-invariant). Its individual components
($\Delta\textbf{p}_{y/x}$) maintain the same angle with respect to the
NV${}^{-}_{z}$ axis, with an increase of only $\sim$ 6 deg in the angle
between their Cartesian $x$ and $y$ components. It should be noted that these
slight changes might be somewhat convoluted, since there could be a basis
dependence for the $\Delta\textbf{p}_{y/x}$ since the ${}^{3}E$ states are
degenerate in the absence of strain, (full C3v symmetry). Regarding the
optical TDM (${}^{3}A_{2}$ $\leftrightarrow$ ${}^{3}E$), the changes in its
Cartesian components amount to a tilt of $\sim$ 3 deg between the two Wannier
bases, while its magnitude changes by $<2\%$. Lastly, the quantity that seems
to vary the most from a slight change in the Wannier basis is the
$\bm{\mu}_{y/x}$, with changes in its Cartesian components that amount to a
tilt of $\sim$ 15 deg with respect to the $\bm{\mu}_{y/x}$ from the main
results (maximally-localized Wannier functions). Its magnitude changes by
$<7\%$. In all cases, the individual angles between each dipole vector and the
local NV${}^{-}_{z}$ axis change by $\leq$ 3 deg between the two bases. While
these differences in the individual components of the dipole vectors translate
into possible gauge-dependence of their absolute orientation, the key figures
of merit for the sensing applications proposed herein (strength of permanent
dipoles) are qualitatively unchanged. A more exhaustive study of these TDMs
under strain is required to derive their relative orientations with full
certainty. But, given that the current uncertainties are within experimental
margins, we consider this a matter for future work.

##  V Conclusions

We have employed first-principles calculations based on density-functional
theory and quantum embedding to derive a many-body description of strain and
charge induced Stark effects in the NV- defect center in diamond. The relaxed
atomic structures of a diamond supercell containing a single NV- under the
influence of strain have been obtained from DFT. Using these as starting
points for embedding calculations, we extract the strain susceptibility of the
many-body states comprising the ${}^{3}E$ excited-state manifold, as well as
permanent and transition dipoles under strain. The results derived in this
work are in excellent agreement with recent experiments, where the dipole
couplings of the NV- are estimated and/or used for sensing applications [16,
21, 22, 27, 52]. This work thus broadens our fundamental understanding on the
optoelectronic properties of the NV- center, and also highlights the utility
of including many-body effects to describe the properties of defect-related
color centres in solids for sensing applications.

###### Acknowledgements.

CED acknowledges support from the National Science Foundation under Grant No.
DMR-2237674. GILM acknowledges funding from grant NSF-2208863. JZ and CED
acknowledge support from a joint SBU-BNL seed grant. The Flatiron Institute is
a division of the Simons Foundation.

## References

  * Rogers  _et al._ [2009] L. J. Rogers, R. L. McMurtrie, M. J. Sellars, and N. B. Manson, Time-averaging within the excited state of the nitrogen-vacancy centre in diamond, New J. Phys. 11 (2009).
  * Doherty  _et al._ [2011] M. W. Doherty, N. B. Manson, P. Delaney, and L. C. L. Hollenberg, The negatively charged nitrogen-vacancy centre in diamond: the electronic solution, New J. Phys. 13 (2011).
  * Doherty  _et al._ [2012] M. W. Doherty, F. Dolde, H. Fedder, F. Jelezko, J. Wrachtrup, N. B. Manson, and L. C. L. Hollenberg, Theory of the ground-state spin of the $\text{NV}^{-}$ center in diamond, Phys. Rev. B 85 (2012).
  * Doherty  _et al._ [2013] M. W. Doherty, N. B. Manson, P. Delaney, F. Jelezko, J. Wrachtrup, and L. C. L. Hollenberg, The nitrogen-vacancy colour centre in diamond, Physics Reports 528, 1 (2013).
  * Deák  _et al._ [2014] P. Deák, B. Aradi, M. Kaviani, T. Frauenheim, and A. Gali, Formation of NV centers in diamond: A theoretical study based on calculated transitions and migration of nitrogen and vacancy related defects, Phys. Rev. B 89, 075203 (2014).
  * Thiering and Gali [2018] G. Thiering and A. Gali, Theory of the optical spin-polarization loop of the nitrogen-vacancy center in diamond, Phys. Rev. B 98, 085207 (2018).
  * Manson  _et al._ [2018] N. B. Manson, M. Hedges, M. S. J. Barson, R. Ahlefeldt, M. W. Doherty, H. Abe, T. Ohshima, and M. J. Sellars, $NV^{-}$-$N^{+}$ pair centre in 1b diamond, New J. Phys. 20, 113037 (2018).
  * Razinkovas  _et al._ [2021] L. Razinkovas, M. W. Doherty, N. B. Manson, C. G. V. de Walle, and A. Alkauskas, Vibrational and vibronic structure of isolated point defects: The nitrogen-vacancy center in diamond, Phys. Rev. B 104, 045303 (2021).
  * Schirhagl  _et al._ [2014] R. Schirhagl, K. Chang, M. Loretz, and C. L. Degen, Nitrogen-vacancy centers in diamond: Nanoscale sensors for physics and biology, Nano Lett. 65, 83 (2014).
  * Balasubramanian  _et al._ [2014] G. Balasubramanian, A. Lazariev, S. R. Arumugam, and D. wen Duan, Nitrogen-vacancy color center in diamond - emerging nanoscale applications in bioimaging and biosensing, Current Opinion in Chemical Biology 20, 69 (2014).
  * Zhang  _et al._ [2021] T. Zhang, G. Pramanik, K. Zhang, M. Gulka, L. Wang, J. Jing, F. Xu, Z. Li, Q. Wei, P. Cigler, and Z. Chu, Toward quantitative bio-sensing with nitrogen–vacancy center in diamond, Nat. Commun. 6, 2077–2107 (2021).
  * Neumann  _et al._ [2013] P. Neumann, I. Jakobi, F. Dolde, C. Burk, R. Reuter, G. Waldherr, J. Honert, T. Wolf, A. Brunner, J. H. Shim, D. Suter, H. Sumiya, J. Isoya, and J. Wrachtrup, High-precision nanoscale temperature sensing using single defects in diamond, Nano Lett. 13, 2738–2742 (2013).
  * Laraoui  _et al._ [2015] A. Laraoui, H. Aycock-Rizzo, Y. Gao, X. Lu, E. Riedo, and C. A. Meriles, Imaging thermal conductivity with nanoscale resolution using a scanning spin probe, Nat. Commun. 6 (2015).
  * Rovny  _et al._ [2022] J. Rovny, Z. Yuan, M. Fitzpatrick, A. I. Abdalla, L. Futamura, C. Fox, M. C. Cambria, S. Kolkowitz, and N. P. de Leon, Nanoscale covariance magnetometry with diamond quantum sensors, Science 378, 1301 (2022).
  * Monge  _et al._ [2023] R. Monge, T. Delord, N. V. Proscia, Z. Shotan, H. Jayakumar, J. Henshaw, P. R. Zangara, A. Lozovoi, D. Pagliero, P. D. Esquinazi, T. An, I. Sodemann, V. M. Menon, and C. A. Meriles, Spin dynamics of a solid-state qubit in proximity to a superconductor, Nano Lett. 23, 422–428 (2023).
  * Tamarat  _et al._ [2006] P. Tamarat, T. Gaebel, J. R. Rabeau, M. Khan, A. D. Greentree, H. Wilson, L. C. L. Hollenberg, S. Prawer, P. Hemmer, F. Jelezko, and J. Wrachtrup, Stark shift control of single optical centers in diamond, Phys. Rev. Lett. 97, 083002 (2006).
  * Dolde  _et al._ [2011] F. Dolde, H. Fedder, M. W. Doherty, T. Nöbauer, F. Rempp, G. Balasubramanian, T. Wolf, F. Reinhard, L. C. L. Hollenberg, F. Jelezko, and J. Wrachtrup, Electric-field sensing using single diamond spins, Nat. Phys. 7, 459 (2011).
  * Bassett  _et al._ [2011] L. C. Bassett, F. J. Heremans, C. G. Yale, B. B. Buckley, and D. D. Awschalom, Electrical tuning of single nitrogen-vacancy center optical transitions enhanced by photoinduced fields, Phys. Rev. Lett. 107, 266403 (2011).
  * Mittiga  _et al._ [2018] T. Mittiga, S. Hsieh, C. Zu, B. Kobrin, F. Machado, P. Bhattacharyya, N. Z. Rui, A. Jarmola, S. Choi, D. Budker, and N. Y. Yao, Imaging the local charge environment of nitrogen-vacancy centers in diamond, Phys. Rev. Lett. 121, 246402 (2018).
  * Bian  _et al._ [2021] K. Bian, W. Zheng, X. Zeng, X. Chen, R. Stöhr, A. Denisenko, S. Yang, J. Wrachtrup, and Y. Jiang, Nanoscale electric-field imaging based on a quantum sensor and its charge-state control under ambient condition, Nat. Commun. 12 (2021).
  * Ji  _et al._ [2024] W. Ji, Z. Liu, Y. Guo, Z. Hu, J. Zhou, S. Dai, Y. Chen, P. Yu, M. Wang, K. Xia, F. Shi, Y. Wang, and J. Du, Correlated sensing with a solid-state quantum multisensor system for atomic-scale structural analysis, Nat. Photon. 18, 230–235 (2024).
  * Delord  _et al._ [2024] T. Delord, R. Monge, and C. A. Meriles, Correlated spectroscopy of electric noise with color center clusters, Nano Lett. (in press) (2024).
  * Maze  _et al._ [2011] J. R. Maze, A. Gali, E. Togan, Y. Chu, A. Trifonov, E. Kaxiras, and M. D. Lukin, Properties of nitrogen-vacancy centers in diamond: the group theoretic approach, New J. Phys. 13, 025025 (2011).
  * Gali [2019] A. Gali, Ab initio theory of the nitrogen-vacancy center in diamond, Nanophotonics 8, 1907–1943 (2019).
  * Olivero  _et al._ [2013] P. Olivero, F. Bosia, B. A. Fairchild, B. C. Gibson, A. D. Greentree, P. Spizzirri, and S. Prawer, Splitting of photoluminescent emission from nitrogen–vacancy centers in diamond induced by ion-damage-induced stress, New J. Phys. 15 (2013).
  * McCullian  _et al._ [2022] B. A. McCullian, H. F. H. Cheung, H. Y. Chen, and G. D. Fuchs, Quantifying the spectral diffusion of $\text{N}-\text{V}$ centers by symmetry, Phys. Rev. Applied 18 (2022).
  * McCullian  _et al._ [2024] B. A. McCullian, V. Sharma, H. Y. Chen, J. C. Crossman, E. J. Mueller, and G. D. Fuchs, Coherent acoustic control of defect orbital states in the strong-driving limit, arXiv arXiv:2403.10989 (2024).
  * Trusheim and Englund [2016] M. E. Trusheim and D. Englund, Wide-field strain imaging with preferentially aligned nitrogen-vacancy centers in polycrystalline, New J. Phys. 18 (2016).
  * Kehayias  _et al._ [2019] P. Kehayias, M. J. Turner, R. Trubko, J. M. Schloss, C. A. Hart, M. Wesson, D. R. Glenn, and R. L. Walsworth, Imaging crystal stress in diamond using ensembles of nitrogen-vacancy centers, Phys. Rev. B 100, 174103 (2019).
  * Budnik  _et al._ [2018] R. Budnik, O. Cheshnovsky, O. Slone, and T. Volansky, Direct detection of light dark matter and solar neutrinos via color center production in crystals, Phys. Lett. B 782, 242 (2018).
  * Rajendran  _et al._ [2017] S. Rajendran, N. Zobrist, A. O. Sushkov, R. Walsworth, and M. Lukin, A method for directional detection of dark matter using spectroscopy of crystal defects, Phys. Rev. D 96 (2017).
  * Kirkpatrick  _et al._ [2023] A. R. Kirkpatrick, G. Chen, H. Witkowska, J. Brixey, B. L. Green, M. J. Booth, P. S. Salter, and J. M. Smith, Ab initio study of defect interactions between the negatively charged nitrogen vacancy centre and the carbon self-interstitial in diamond, Phil. Trans. R. Soc. A 382, 0174 (2023).
  * Ebadi  _et al._ [2022] R. Ebadi, M. C. Marshall, D. F. Phillips, J. Cremer, T. Zhou, M. Titze, P. Kehayias, M. S. Ziabari, N. Delegan, and S. R. et al., Directional detection of dark matter using solid-state quantum sensing, AVS Quantum Sci. 4 (2022).
  * Jaskula  _et al._ [2017] J.-C. Jaskula, E. Bauch, S. Arroyo-Camejo, M. D. Lukin, S. W. Hell, A. S. Trifonov, and R. L. Walsworth, Superresolution optical magnetic imaging and spectroscopy using individual electronic spins in diamond, Opt. Express 25, 11048 (2017).
  * Alaerts  _et al._ [2024] L. Alaerts, Y. Xiong, S. Griffin, and G. Hautier, A first principles study of the stark shift effect on the zero-phonon line of the NV center in diamond, arXiv arXiv:2403.07771v1 (2024).
  * Bockstedte  _et al._ [2018] M. Bockstedte, F. Schütz, T. Garratt, V. Ivády, and A. Gali, Ab initio description of highly correlated states in defects for realizing quantum bits, npj Quantum Materials 3, 31 (2018).
  * Ma  _et al._ [2021] H. Ma, N. Sheng, M. Govoni, and G. Galli, Quantum embedding theory for strongly correlated states in materials, Journal of Chemical Theory and Computation 17, 2116 (2021).
  * Bhandari  _et al._ [2021] C. Bhandari, A. L. Wysocki, S. E. Economou, P. Dev, and K. Park, Multiconfigurational study of the negatively charged nitrogen-vacancy center in diamond, Phys. Rev. B 103, 014115 (2021).
  * Muechler  _et al._ [2022] L. Muechler, D. I. Badrtdinov, A. Hampel, J. Cano, M. Rösner, and C. E. Dreyer, Quantum embedding methods for correlated excited states of point defects: Case studies and challenges, Phys. Rev. B 105, 235104 (2022).
  * Chen  _et al._ [2023] Y. Chen, T. Jiang, H. Chen, E. Han, A. Alavi, K. Yu, E. Wang, and J. Chen, Multiconfigurational nature of electron correlation within nitrogen vacancy centers in diamond, Phys. Rev. B 108, 045111 (2023).
  * Gali [2023] A. Gali, Recent advances in the ab initio theory of solid-state defect qubits, Nanophotonics 12, 1 (2023).
  * Haldar  _et al._ [2023] S. Haldar, A. Mitra, M. R. Hermes, and L. Gagliardi, Local excitations of a charged nitrogen vacancy in diamond with multireference density matrix embedding theory, J. Phys. Chem. Lett. 14, 4273 (2023).
  * Pizzi  _et al._ [2020] G. Pizzi, V. Vitale, R. Arita, S. Blügel, F. Freimuth, G. Géranton, M. Gibertini, D. Gresch, C. Johnson, and T. Koretsune, et al., Wannier90 as a community code: new features and applications, Journal of Physics: Condensed Matter 32, 165902 (2020).
  * Souza  _et al._ [2001] I. Souza, N. Marzari, and D. Vanderbilt, Maximally localized Wannier functions for entangled energy bands, Phys. Rev. B 65, 035109 (2001).
  * Aryasetiawan  _et al._ [2004] F. Aryasetiawan, M. Imada, A. Georges, G. Kotliar, S. Biermann, and A. I. Lichtenstein, Frequency-dependent local interactions and low-energy effective models from electronic structure calculations, Phys. Rev. B 70, 195104 (2004).
  * Parcollet  _et al._ [2015] O. Parcollet, M. Ferrero, T. Ayral, H. Hafermann, I. Krivenko, L. Messio, and P. Seth, TRIQS: a toolbox for research on interacting quantum systems, Comp. Phys. Comm. 196, 398 (2015).
  * Rogers  _et al._ [2008] L. J. Rogers, S. Armstrong, M. J. Sellars, and N. B. Manson, Infrared emission of the NV centre in diamond: Zeeman and uniaxial stress studies, New J. Phys. 10 (2008).
  * Acosta  _et al._ [2010] V. M. Acosta, A. Jarmola, E. Bauch, and D. Budker, Optical properties of the nitrogen-vacancy singlet levels in diamond, Phys. Rev. B 82, 201202(R) (2010).
  * Crisci  _et al._ [2011] A. Crisci, F. Baillet, M. Mermoux, G. Bogdan, M. Nesládek, and K. Haenen, Residual strain around grown-in defects in CVD diamond single crystals: A 2D and 3D Raman imaging study, Phys. Status Solidi A 208, 2038 (2011).
  * Acosta  _et al._ [2012] V. M. Acosta, C. Santori, A. Faraon, Z. Huang, K.-M. C. Fu, A. Stacey, D. A. Simpson, K. Ganesan, S. Tomljenovic-Hanic, A. D. Greentree, S. Prawer, and R. G. Beausoleil, Dynamic stabilization of the optical resonances of single nitrogen-vacancy centers in diamond, Phys. Rev. Lett. 108, 206401 (2012).
  * Happacher  _et al._ [2022] J. Happacher, D. A. Broadway, J. Bocquel, P. Reiser, A. Jimenez, M. A. Tschudin, L. Thiel, D. Rohner, M. li Grimau Puigibert, B. Shields, J. R. Maze, V. Jacques, and P. Maletinsky, Low-temperature photophysics of single nitrogen-vacancy centers in diamond, PHys. Rev. Lett. 128, 177401 (2022).
  * [52] T. Delord, R. Monge, G. I. López-Morales, J. Flick, C. E. Dreyer, and C. A. Meriles, in preparation .
  * Kurokawa  _et al._ [2024] H. Kurokawa, K. Wakamatsu, S. Nakazato, T. Makino, H. Kato, Y. Sekiguchi, and H. Kosaka, Coherent electric field control of orbital state of a neutral nitrogen-vacancy center, Nat. Commun. 15, 4039 (2024).
  * Lozovoi  _et al._ [2021] A. Lozovoi, H. Jayakumar, D. Daw, G. Vizkelethy, E. Bielejec, M. W. Doherty, J. Flick, and C. A. Meriles, Optical activation and detection of charge transport between individual colour centres in diamond, Nat. Electron 4, 717 (2021).