# Nitrogen overgrowth as a catalytic mechanism during diamond chemical vapour
deposition

Lachlan M. Oberg1 Marietta Batzer2 Alastair Stacey3 and Marcus W. Doherty1,4

(1Laser Physics Center, Research School of Physics, Australian National
University, Australian Capital Territory 2601, Australia
2Department of Physics, University of Basel, Klingelbergstrasse 82,
Switzerland
3School of Science, RMIT University, Melbourne, Victoria 3001, Australia
4Quantum Brilliance Pty Ltd, 116 Daley Road, Acton, Australian Capital
Territory 2601, Australia
)

###### Abstract

Nitrogen is frequently included in chemical vapour deposition feed gases to
accelerate diamond growth. While there is no consensus for an atomistic
mechanism of this effect, existing studies have largely focused on the role of
sub-surface nitrogen and nitrogen-based adsorbates. In this work, we
demonstrate the catalytic effect of surface-embedded nitrogen in nucleating
new layers of (100) diamond. To do so we develop a model of nitrogen
overgrowth using density functional theory. Nucleation of new layers occurs
through C insertion into a C–C surface dimer. However, we find that C
insertion into a C–N dimer has substantially reduced energy requirements. In
particular, the rate of the key dimer ring-opening and closing mechanism is
increased 400-fold in the presence of nitrogen. Full incorporation of the
substitutional nitrogen defect is then facilitated through charge transfer of
an electron from the nitrogen lone pair to charge acceptors on the surface.
This work provides a compelling mechanism for the role of surface-embedded
nitrogen in enhancing (100) diamond growth through the nucleation of new
layers. Furthermore, it demonstrates a pathway for substitutional nitrogen
formation during chemical vapour deposition which can be extended to study the
creation of technologically relevant nitrogen-based defects.

## Introduction

It is well established that trace amounts of nitrogen (N) accelerate the rate
of diamond growth during chemical vapour deposition (CVD)[1, 2, 3, 4]. While
the magnitude of this enhancement varies with plasma composition, pressure,
temperature[5], and activation technique[6], a ten-fold increase in the growth
rate of (100) diamond can be achieved at high plasma intensities[5]. Despite
the widespread adoption of N-enhanced diamond growth for both academic and
commercial purposes, there does not exist a consistent and universally
accepted mechanism for the effect. Moreover, a process detailing overgrowth of
N on the (100) surface to form a bulk substitutional defect remains unknown.
This remains a critical pre-cursor for understanding and enhancing the
formation of technologically relevant N-based defects in diamond during
growth. In particular, increasing the low yield and alignment of nitrogen-
vacancy (NV) centers which has applications in quantum computing[7],
communications[8], and metrology[9, 10].

While the catalytic role of N has generated considerable theoretical interest,
existing literature has failed to identify an atomic mechanism which can
demonstrably produce order-of-magnitude enhancements to diamond growth rates.
The majority of studies have focused on the role of sub-surface N and N
adsorbates in reducing reaction and transition state energies for hydrogen-
terminated (100) diamond growth[11, 12, 13, 14, 15, 16]. These ab initio works
have primarily studied the impact of N on carbon (C) insertion into a C–C
dimer of the reconstructed (100) surface[17, 18, 19, 20, 21]. This process has
received substantial attention because it is fundamental for all (100) diamond
growth by CVD. Most studies follow or slightly develop the widely-accepted
model for C insertion presented by Garrison et al.[22]. This plasma-surface
reaction is itself built on extensive earlier work[23, 24, 25, 26] and
consists of six key steps depicted in Figure 1.

Figure 1: Conventional mechanism for C insertion into a C–C dimer of the
H-terminated (100) diamond surface[22, 17, 19], a pivotal stage of CVD diamond
growth. For simplicity we depict only the C–C dimer and neighbouring ions.
Full depictions of the surface are presented in the supplementary material. In
step $1\rightarrow 2$, a terminating H atom is abstracted by an atomic H
radical within the plasma. This creates a radical site allowing for the
adsorption of a methyl radical in step $2\rightarrow 3$. Further H abstraction
of the methyl adsorbate ($3\rightarrow 4$) creates a radical CH2 adsorbate.
This adsorbate incorporates into the diamond surface via a ring
opening/closing mechanism in which the dimer expands ($4\rightarrow 5$) to
allow formation of an energetically favourable six-membered ring
($5\rightarrow 6$). Finally, the remaining surface radical is stabilized the
adsorption of atomic H in step $6\rightarrow 7$.

Early work by Frauenheim et al. demonstrated that sub-surface N reduced the
binding energy of terminating H atoms and therefore enhanced the rate of
initial H abstraction[11]. This was attributed to charge stabilization of the
resulting surface radical (structure 2 in Figure 1) due to electron transfer
from the substitutitonal N lone pair[27]. However, this makes methyl
adsorption (step $2\rightarrow 3$) energetically unfavourable. An alternative
growth mechanism was therefore proposed involving charge transfer to the
metastable anti-bonding orbital of C–C surface dimers and subsequent
adsorption of CH2 radicals. This hypothesis is challenged by CVD plasma
experiment and modeling which demonstrated that methyl radicals (CH3) are the
majority C-based radical species proximate to the diamond surface during
growth and exceed CH2 concentrations by three orders of magnitude[28, 29].

Alternatively, more recent works have rationalised Frauenheim et al.’s
original findings by assuming that H abstraction, not methyl adsorption, is
the rate-limiting step for diamond growth. These studies identify that near-
surface N defects enhance H abstraction from methyl adsorbates ($3\rightarrow
4$) through charge transfer[13]. Additionally, near-surface N defects promote
the creation of adsorption sites ($1\rightarrow 2$) through weakening of
proximal C–C surface dimer bonds[16]. Further ab initio work has demonstrated
that sub-surface N increases activation barriers for CH2 insertion into the
C–C dimer and CH2 surface migration[15]. Other studies have found that co-
adsorbed N atoms have no bearing on H abstraction[16], reduces the adsorption
rate of CH3 on step edges, but increases the rate of adsorption (and
hypothesized to enhance migration) of CH2 on step-edges[14]. It has also been
proposed that N accelerates (111) growth by nucleating new diamond layers[30].
In summary, extensive ab initio modeling has not yet provided a clear
consensus on the catalytic role of N in CVD.

As a supplement to ab initio calculations, extensive Monte-Carlo modeling has
been performed to understand the complexity of mesoscale diamond growth. These
works have highlighted the importance of surface migration[31] and nucleation
of new layers[32] to reproduce experimental growth rates and morphologies. For
example, the surface of CVD-grown (100) diamond often exhibits a relatively
smooth and terraced structure[33]. However, note that this is not ubiquitous
and depends on reactor conditions[34, 5]. Such morphologies are indicative of
step-flow modes, in which migration of hydrocarbon adsorbates and preferential
adsorption at step-edges is believed to contribute largely to layer growth.
Recent Monte-Carlo modeling has emphasized the importance of critical nuclei
for propagating layer growth. These are immobile surface features, such as a
lone C–C dimer, which act as nucleation points[35, 36]. Furthermore, the
inclusion of super-nucleating species in Monte-Carlo models has been found to
catalyze diamond growth. These are adsorbates or surface defects (hypothesised
to be N-based) which quickly form critical nuclei following formation. When
growth is limited by nucleation of new layers (i.e., growth is dominated by
step-flow), super-nucleating species have demonstrated ten-fold enhanced
growth rates[35]. The search for a N-based super-nucleating species is
therefore well founded and forms the primary aim of this work.

We identify a potential super-nucleating species by producing the first
atomistic model for sub-surface N formation during CVD. Our calculations
demonstrate that surface-embedded N catalyzes the nucleation of new diamond
layers during its overgrowth and subsequent encapsulation into bulk diamond.
The atomic structure of surface-embedded N is depicted in Figure 2. It
consists of a substitutional defect which maintains the structure of the
dimerised 2x1-(100) surface. The bonding of the C–N dimer resembles that of a
typical C–C surface dimer with a C-H covalent bond substituted by a N electron
lone pair. Surface-embedded N has previously been identified as the most
energetically stable form of substitutional N in the (100) diamond surface[9].
Furthermore, it is reported to be evident in NEXAFS scans of diamond (100)
surfaces following N plasma treatment[37]. Consequently, it is likely that
surface-embedded N is a common surface defect during CVD growth of the (100)
surface. It is therefore chosen as the initial point for our overgrowth
mechanism.

Further motivation for this work includes the refinement of existing ab initio
techniques used in previous CVD-growth literature. Firstly, many works employ
cluster models[15, 17, 21, 20, 19] to represent the diamond surface. While
cluster models are capable of reproducing results consistent with periodic
slab calculations[38, 39], this is not always the case[40], and careful
optimisation of cluster dimensions is always required. As has been noted in
previous works[19], the relaxed steric constraints associated with some
cluster calculations result in an underestimation of structural stabilities
during the ring opening/closing mechanism[20]. Secondly, some diamond
studies[19, 16] employ the drag method for determining the transition state
(or do not explicitly mention the method used[17, 21, 20]) which can produce
inaccurate activation barriers for complex potential energy surfaces[41, 42].
Instead, it is preferable to use chain-of-states methods when possible, such
as the nudged-elastic-band (NEB) technique, which provide greater consistency
in determining the transition state[43]. Consequently, in this work we employ
state-of-the-art density functional theory, including a fully quantum slab
model of the H-terminated diamond (100) surface, the use of hybrid
functionals, and the NEB method for determining transition states.

Figure 2: Surface-embedded N on the (100) surface forms a dimer unit with a
neighbouring C atom.

Our model for overgrowth begins with surface-embedded N as presented in Figure
2 and is separated into three distinct stages. Firstly, nucleation of a new
diamond layer through C insertion into the C–N dimer is presented in Section
2. We find that surface-embedded N reduces energy requirements for nearly
every stage of the reaction depicted in Figure 1. In particular, the rate of
the ring opening/closing mechanism ($4\rightarrow 5\rightarrow 6$) is
increased by a factor of approximately 400. This provides a compelling
atomistic mechanism for the N catalysis of diamond CVD growth through enhanced
layer nucleation. In Section 3 we present the second growth stage where we
consider layer growth emanating from the nucleation point. We demonstrate the
formation of a new C–C dimer and C bridging to form a new dimer row. The key
result is that surface-embedded N does not impede typical diamond layer
growth. The final growth stage is considered in Section 4, where we
demonstrate that full encapsulation of the N defect into the surface (i.e.,
the formation of sub-surface N) is mediated by charge transfer.

##  1 Method

Density functional theory was performed using the VASP plane-wave code[44, 45,
46, 47] using PBE[48] and verified using B3LYP[49] functionals. PBE belongs to
the generalised-gradient approximation (GGA) class of exchange-correlation
functionals. While GGA functionals are computationally efficient, they are
known to over-delocalize electrons and over-stabilize stretched bonds of
transition states. This can result in underestimations of energy barriers,
especially for gas-surface reactions[50]. Hybrid functionals, such as B3LYP,
correct the over-delocalisation of GGA functionals by including some
contribution of exact exchange energy thereby improving estimations of
reaction barriers[51]. However, this increase in accuracy requires significant
computational costs and so we perform the majority of calculations using the
PBE functional. To validate the accuracy of these PBE calculations we
calculate the key ring opening step ($4\rightarrow 5$ of Figure 1) using B3LYP
functionals as discussed in the following section.

The $2\times 1$ reconstructed H-C(100) surface has been represented using a
periodic slab model consisting of 240 C ions with dimensions $5\times 4\times
6$ (in units of bulk-diamond primitive unit cells). As presented in the
supplementary material, these dimensions have been stringently optimised
against the surface work function, reaction and transition state energies for
the ring opening/closing mechanism, and interaction energies between N defects
in adjacent supercells. Slab images in the <100> direction were separated by a
distance of 10 Å and long-range dipole interactions have been corrected for.
$\Gamma$-point sampling and real-space projection operators have been used for
all calculations. The cut-off energy for the plane-wave basis is 600 eV. The
electronic tolerance for successive iterations of the self-consistent field
method is taken as $0.1$ meV, while all geometry optimisations are performed
to a tolerance of $0.05$ eV/Å per ion. Ions composing the top two layers of
primitive unit cells (those primarily involved in surface reactions) were
allowed to relax in all directions, while the remaining ions were permitted to
relax only in the <100> direction. The effects of spin polarization on both
stable and transition states has been fully accounted for. The climbing NEB
technique has been used to calculate transition states for all reaction
barriers[52]. This extension to NEB determines the exact transition state by
forcing the highest energy image to the saddle point.

##  2 Nucleation of a new diamond layer

The first stage of growth consists of new layer nucleation through C insertion
into the C–N dimer. The results of our ab initio calculations are presented in
Table 1, which documents the reaction energies ($\Delta E$) and transition
state energies (TS) for each of the steps required for C insertion as per
Figure 1. This is performed using PBE functionals for both C–C and C–N dimers
for the purposes of comparison. We also include several B3LYP calculations
performed using both our slab geometry as well as a small adamantane cluster
(C10H16) as discussed below. Additionally, we compare our results to previous
values for C–C dimers by Cheesman et al[19]. This is a large-scale hybrid
study in which the reactive C–C dimer and neighbouring atoms were treated
using B3LYP functionals, while the remaining $\sim 1500$ atoms of the cluster
were treated using molecular mechanics (MM). We have chosen this particular
work as it reproduces or demonstrably improves on previous studies for C
insertion that employ both fully quantum[18, 17] and hybrid cluster models[20,
21].

|  | PBE |  | B3LYP |  | B3LYP |  | B3LYP
---|---|---|---|---|---|---|---|---
|  | (slab) |  | (slab) |  | (adamantane) |  | (Cheesman et al.[19])
Step |  | C–C | C–N |  | C–C | C–N |  | C–C |  | C–C
$1\rightarrow 2$ | $\Delta E$ | -0.026 | -0.135 |  |  |  |  |  |  | -0.007
| TS | 0.185 | 0.136 |  |  |  |  |  |  | 0.274
$2\rightarrow 3$ | $\Delta E$ | -4.556 | -3.671 |  |  |  |  |  |  | -3.868
| TS | 0 | 0 |  |  |  |  |  |  | 0
$3\rightarrow 4$ | $\Delta E$ | -0.228 | -0.142 |  |  |  |  |  |  | -0.312
| TS | 0.355 | 0.358 |  |  |  |  |  |  | 0.306
$4\rightarrow 5$ | $\Delta E$ |  | -0.143 |  |  | -0.235 |  | 0.019 |  | 0.290
| TS |  | 0.057 |  |  | 0.094 |  | 0.370 |  | 0.450
$5\rightarrow 6$ | $\Delta E$ |  | -0.476 |  |  | -0.329 |  | -0.540 |  | -0.813
| TS |  | 0.283 |  |  |  |  | 0.571 |  | 0.550
$4\rightarrow 6$ | $\Delta E$ | -0.350 |  |  | -0.299 |  |  |  |  | -0.513‡
| TS | 0.959 |  |  | 0.965† |  |  |  |  | $0.900^{\ddagger}$
$6\rightarrow 7$ | $\Delta E$ | -4.141 | -4.166 |  |  |  |  |  |  | -4.573
| TS | 0 | 0 |  |  |  |  |  |  | 0
Table 1: Reaction ($\Delta E$) and transition state (TS) energies for C
insertion into C–C or C–N dimers on the H-C(100) surface (c.f. Figures 1 and
3). All energies are in eV. †Transition state geometry optimised using 400 eV
cut-off energy. ‡Effective energy for the direct transition $4\rightarrow 6$
if state 5 was not stable. Figure 3: The mechanism for C insertion into a C–N
dimer on the C(100) surface. Reaction energies (in eV) are presented between
each step in the reaction, and where applicable we also include barrier
energies in square brackets. Figure 4: Transition state for the ring
opening/closing reaction ($4\rightarrow 6$) during C insertion into a C–C
dimer.

###  2.1 C insertion into C–C

Our results for C insertion into the C–C dimer are largely consistent with
existing ab initio works. We identify energy barriers and transition states
for initial H abstraction ($1\rightarrow 2$) similar to Cheesman et al. and
several other studies[17, 18]. However as expected the PBE functional
underestimates the transition state energy by approximately 100 meV. As with
other works, we identify no barrier to CH3 adsorption ($2\rightarrow 3$) or H
adsorption ($6\rightarrow 7$)[18, 19, 21, 20]. While our transition state
energy for the H abstraction in step $3\rightarrow 4$ overestimates that of
Cheesman et al. by 49 meV, it underestimates that obtained in the B3LYP
cluster studies of Kang and Musgrave by approximately 100 meV[17].

The greatest deviation between our study and existing literature is for the
ring opening/closing mechanism, $4\rightarrow 5\rightarrow 6$. We now discuss
this discrepancy in-depth due to its relevance for the validity of our results
and to a general understanding of conventional diamond growth. The majority of
studies find that state 5, the open-ring configuration, is stable with the CH2
adsorbate realizing an sp2 bonding configuration (c.f. Figure 1)[19, 18, 21,
17]. However, in agreement with the GGA cluster calculations of Oleinik et
al.[20], we find that state 5 is not stable and that C insertion occurs
through the direct transition $4\rightarrow 6$. The transition state for
$4\rightarrow 6$ does not possess the sp2 bonding configuration of state 5.
Instead, as depicted in Figure 4, the CH2 radical inserts into the ring
structure simultaneously with breakage of the dimer bond.

We have confirmed that the instability of state 5 is not a deficiency of the
PBE functional. B3LYP functionals in combination with a range of geometry
optimisation algorithms (RMM-DIIS[53], conjugate gradient, and velocity
damping) also failed to identify a stable sp2-like intermediate between states
4 and 6. As presented in Table 1, we also perform climbing NEB calculations
with B3LYP for the direct transition $4\rightarrow 6$ and obtain a transition
state consistent with the PBE functional. Note that due to computational
constraints the transition state geometry was optimised using a lower plane
wave cut-off energy of 400 eV. However, as demonstrated in Figure 1 of the
supplementary material, this does not significantly reduce the accuracy of our
calculations.

Instead, the stability of state 5 appears to be conditional on the use of
either cluster or hybrid quantum/MM models. We have calculated the energetics
for ring opening/closing using B3LYP functionals on an adamantane molecule
(C10H16, resembling the models in Figure 1 if all dangling bonds were hydrogen
terminated), the simplest possible cluster model for a surface dimer. Note
that the corresponding C10 carbon cluster is also the fully quantum region
used in several hybrid studies[19, 20, 21]. Here we identify that state 5 is
stable and that the reaction $4\rightarrow 5\rightarrow 6$ possesses
transition state and reaction energies in close agreement with Cheesman et al.
as shown in Table 1. Moreover, we find that the geometry of state 5 in our
cluster calculations is considerably relaxed with respect to our slab
calculations. In particular, the separation distance between C atoms of the
broken dimer of state 5 is large at 2.95 Å. As demonstrated in the
supplementary, more than 1.4 eV is required to achieve a separation distance
of this scale within the slab geometry. It is clear that the periodicity
enforced by slab models limits the relaxation of the reactive dimer and
surrounding ions. Clusters are therefore inadequate for modeling some aspects
of diamond surface chemistry which are better suited to periodic slabs.

Without the erroneous identification of state 5 as stable, cluster models
produce results in agreement with our slab model. To see this, consider the
effective energies for the direct transition $4\rightarrow 6$ deduced from the
work of Cheesman et al. As shown in Table 1, if state 5 is ignored we obtain a
transition-state energy of 0.900 eV and a reaction energy of -0.513 eV,
similar to respective values of 0.959 eV and -0.350 eV for our work using PBE
functionals.

###  2.2 C insertion into C–N

The presence of surface-embedded N greatly enhances C insertion. Firstly, our
results in Table 1 show that the transition barrier for the initial H
abstraction (step $1\rightarrow 2$) is reduced by 49 meV while stability is
increased by 109 meV. This is attributable to the greater electronegativity of
N in comparison to C, where we note that the length of the C–N dimer bond is
reduced by 5.5% with respect to C–C. A similar reduction in the transition
state energy is not observed for the abstraction $3\rightarrow 4$ as this
reaction step does not directly involve bonding with surface-embedded N.

The greatest enhancement to C insertion due to surface-embedded N occurs
during the ring opening/closing mechanism. In contrast to the C–C dimer, state
5 for the C–N dimer forms a stable but shallow minimum in which the CH2
adsorbate realizes an sp2 bonding configuration. Moreover, as per our PBE
calculations the transition $4\rightarrow 5$ requires an activation energy of
approximately 5 meV; negligible at typical diamond growth temperatures. Given
the significance of transition $4\rightarrow 5$ for our model of N catalysis,
we validate the accuracy of these PBE calculations using B3LYP functionals.
Full geometry optimisation and the climbing NEB method produce $\Delta
E=-0.235$ eV and a transition state energy of $0.094$ eV. These are in good
agreement with the PBE values of $-0.143$ eV and $0.057$ eV respectively. This
indicates that the PBE functional possesses similar accuracy to B3LYP for
diamond surface reactions (potentially excluding gas-surface interactions as
discussed above), significantly reducing computational costs.

The combination of this low transition state energy and the stability of state
5 drastically enhances the rate of the ring opening/closing mechanism. This
can be quantified using kinetic rate equations to compare the reaction times
between C–C and C–N dimers. Defining the population of states 4, 5, and 6 as
$\phi_{4}$, $\phi_{5}$, and $\phi_{6}$ respectively we have that

$\displaystyle\phi_{4}^{\prime}(t)$ $\displaystyle=-\Gamma_{4\rightarrow
5}\phi_{4}+\Gamma_{5\rightarrow 4}\phi_{5}$ (1)
$\displaystyle\phi_{5}^{\prime}(t)$ $\displaystyle=-(\Gamma_{5\rightarrow
4}+\Gamma_{5\rightarrow 6})\phi_{5}+\Gamma_{4\rightarrow
5}\phi_{4}+\Gamma_{6\rightarrow 5}\phi_{6}$
$\displaystyle\phi_{6}^{\prime}(t)$ $\displaystyle=-\Gamma_{6\rightarrow
5}\phi_{6}+\Gamma_{5\rightarrow 6}\phi_{5}.$

Here $\Gamma_{i\rightarrow j}=A_{i\rightarrow j}e^{-\text{TS}_{i\rightarrow
j}/k_{B}T}$, where $A_{i\rightarrow j}$ is the attempt frequency for the
transition between states $i$ and $j$, $\text{TS}_{i\rightarrow j}$ is the
transition state energy, and $T=1200$ K is the approximate temperature of the
diamond surface during CVD[33]. For C insertion into the C–C dimer, one only
needs to consider the direct transition $4\rightarrow 6$ and the rate
equations are simply

$\displaystyle\phi_{4}^{\prime}(t)$ $\displaystyle=-\Gamma_{4\rightarrow
6}\phi_{4}+\Gamma_{6\rightarrow 4}\phi_{6}$
$\displaystyle\phi_{6}^{\prime}(t)$ $\displaystyle=-\Gamma_{6\rightarrow
4}\phi_{6}+\Gamma_{4\rightarrow 6}\phi_{4}.$

A quantitative determination of the pre-exponential factors requires
calculating the vibrational partition functions of each stable and transition
state[54]. To compare the pre-exponential rate of ring opening/closing between
C–C and C–N dimers, we may instead approximate the ratio of thermal occupation
between their respective normal modes. We find that this does not differ by
more than a factor of three for vibrational energies in the expected range of
$100-200$ meV. Therefore for simplicity we assume that $A_{i\rightarrow j}=A$
is constant for all processes involved in ring opening and closing. This
allows us to solve the above system of equations given the initial conditions
$\phi_{4}(0)=1$ and $\phi_{5}(0)=\phi_{6}(0)=0$. Processes interrupting ring
opening/closing, such as H adsorption and abstraction, are not included. These
are assumed to impact C insertion into C–C and C–N at equal rates and
therefore do not influence our relative comparison.

In Figure 5 we plot $\phi_{6}$ as a function of time for ring opening/closing
with C–C and C–N dimers. These results indicate that ring closure occurs
approximately 400 times faster in the presence of surface-embedded N. Not only
does this drastically enhance the rate of C insertion, it also reduces the
residence time of state 4. This would likely result in less etching of the C
adsorbate, believed to occur predominately via $\beta$-scission[55, 31, 33].

Given the large transition state energies required for ring opening/closing
during conventional diamond growth, step $4\rightarrow 6$ appears to be a
limiting step for new layer nucleation. Consequently, we propose that surface-
embedded N can catalyze new layer growth by substantially reducing energy
requirements for C insertion into the surface. This enhancement is expected to
be most pronounced during (100) step-flow modes when growth rates are limited
by new layer nucleation. Under such circumstances, we propose that surface-
embedded N acts as a super-nucleation species for new-layer C growth and
subsequently the formation of critical nuclei on the surface. This provides an
atomic mechanism for enhanced growth rates observed in some Monte-Carlo
studies[35]. Note, however, that the catalytic mechanism identified here
cannot be universal for all N-enhanced diamond growth, particularly that
observed for diamond surfaces other than (100).

It is difficult to quantify the enhancement this catalytic effect has on
macroscopic growth rates. This would require a mesoscopic study (such as a
Monte-Carlo simulation) which includes the effects of surface-embedded N
amongst the complex diversity of surface chemistry that occurs during CVD.
However, amongst existing ab initio studies of (100) diamond growth, this work
demonstrates the greatest enhancement to any key reaction process. For
example, previously suggested catalytic effects of N include enhanced
H-abstraction through weakening of surface bonds proximate to sub-surface
N[16]. This increases the rate of H abstraction local to the sub-surface N by
2.4. This enhancement is minor compared to the 400-fold increase in the rate
of ring opening and closing in the presence of surface-embedded N. By this
metric, our work demonstrates the most compelling atomistic mechanism for N
catalysis to date.

Figure 5: Occupation of state six as a function of time for the ring
opening/closing mechanism on a C–C dimer (orange) and C–N dimer (blue).

##  3 Layer growth

In this section we demonstrate that layer growth propagating from the
nucleation point is not hindered by the presence of surface-embedded N. Due to
the complexity of the CVD environment there are many possible processes that
can contribute to layer growth. Hence, we consider two representative
mechanisms; formation of a new surface dimer and C addition between dimer
rows. We note that these two mechanisms are not necessarily the atomic
processes that produce step-flow growth (which have not yet been conclusively
identified in existing literature). Instead, the two representative mechanisms
involve key elements ubiquitous to most layer growth processes and therefore
assumed to be highly relevant for step flow. These elements include H
abstraction, CH3 adsorption, ring opening/closing, and migration. We calculate
the reaction and transition state energies of these two mechanisms in the
presence of surface-embedded N and compare these results to existing
literature of analogous processes without N[19].

###  3.1 Dimer formation

Figure 6: Mechanism for C insertion into a surface dimer adjacent to an
incorporated C. This results in the formation of a new dimer unit on the
surface and the beginning of a new layer. Reaction energies (in eV) are
presented between each step in the reaction, and where applicable we also
include barrier energies in square brackets.

Figure 6 presents the reaction mechanism and energy requirements for the
insertion of C adjacent to the C already incorporated into the C–N dimer. The
result is the formation of a new C–C dimer on the surface. This mechanism
follows that identified by Cheesman et al. for formation of the dimer unit in
the absence of N. Our calculations have similar energy requirements with the
exception of a secondary dimer/opening closing mechanism as discussed
below[19].

In step $1\rightarrow 2$, H abstraction forms a radical site on the previously
incorporated C. This C adopts a bonding geometry similar to a free methyl
radical, indicating that the N atom does not donate an electron and maintains
its lone pair. This allows for the stable adsorption of CH3 in step
$2\rightarrow 3$ with no reaction barrier. Steps $3\rightarrow 4$ and
$4\rightarrow 5$ consist of H abstractions on the adsorbed methyl and a H on
the adjacent C–C dimer. N maintains its lone pair throughout these processes.
This permits the CH2 adsorbate to bond at the adjacent surface site. Further H
abstraction on the bridging CH2 produces a reactive radical. In step
$7\rightarrow 8$, the C–C dimer bond is severed to accommodate C insertion
into the dimer. This is also a ring opening/closing mechanism and has an
activation barrier commensurate with that observed for C insertion on a lone
C–C dimer (0.860 eV vs 0.959 eV, c.f. Table 1).

In the work of Cheesman et al., step $7\rightarrow 8$ was considered in the
absence of N and found to require a substantially lower activation barrier of
0.324 eV[19]. As this suggests that N hinders the formation of the dimer unit,
we have also calculated step $7\rightarrow 8$ in the absence of N using PBE
functionals and our slab geometry. However, we instead find a transition state
energy of 0.898 eV which indicates that N has little impact on ring
opening/closing in step $7\rightarrow 8$. As discussed previously, the
deviation between our results and that presented in Cheesman et al. is likely
attributable to over-relaxation of the transition state due to their cluster
model.

The new C–C dimer is highly stable and readily adsorbs a further H radical to
satiate the final dangling bond seen in state 8 of Figure 6 ($\Delta E=-4.179$
eV). Note in particular that N does not partake in any immediate chemistry
throughout dimer formation and maintains its lone pair. All reactions possess
similar energetic requirements for abstraction and adsorption as that for a
N-free surface (c.f. Table 1), and consequently N does not appear to benefit
or hinder dimer formation.

###  3.2 C addition between dimer rows

Following formation of a new surface dimer, layer growth may propagate between
dimer rows (often denoted the “trough”). The addition of C across the trough
is the beginning of a new surface dimer adjacent to the one previously formed,
and therefore represents the first step in establishing a new dimer row. Once
again, we identify that surface-embedded N does not hinder or benefit this
process as it does not directly participate in bonding.

The reaction mechanism for C addition across the trough including reaction and
barrier energies is presented in Figure 7. In step $1\rightarrow 2$, H is
abstracted from a ring-closed C on the surface layer. This permits the
adsorption of CH3 onto the radical site in step $2\rightarrow 3$. Two further
H abstractions may then occur on the CH3 adsorbate (step $3\rightarrow 4$) and
on H from a C on the adjacent dimer row (step $4\rightarrow 5$). The latter H
abstraction results in the immediate bridging of the trough by the CH2
adsorbate. This bridging reaction also represents the first half of CH2
migration between dimer rows. Alternatively, through a mechanism similar to
that presented in Figure 6, a new dimer can be formed across the trough.

The addition of C across the trough appears to be a facile mechanism for
propagating dimer rows, requiring only two successive H abstractions following
CH3 adsorption. It is therefore a possible candidate for an atomic process
which drives step-flow growth. Regardless, although the energies we determine
in Figure 7 are similar to those identified by Cheesman et al. for the
analogous process without N[19], they argue that the trough bridging mechanism
is not substantially faster than C insertion into a C–C dimer.

Through the three processes considered above, C insertion into a C–C dimer,
formation of a new dimer, and C addition across dimer troughs, layer growth
may propagate across the surface. When surface-embedded N is not directly
involved in the formation of new bonds and growth, it does not appear to
hinder or benefit growth processes. Instead, the N maintains its stable
electron lone pair rather than donating charge to radical acceptor sites on
the surface. The exception to this typical growth behaviour is C addition in
the trough adjacent to N itself, which we treat in the next section.

Our results indicate that surface-embedded N catalyzes new-layer nucleation
during step-flow modes but does not necessarily increase the rate of step-flow
growth itself. Consequently, it may be expected that high N concentrations
would result in a growth mode dominated by nucleation. Indeed, this effect is
observed in experimental results in which increasing the N concentration leads
to morphology changes indicating a transition from step-flow to nucleation
dominated growth[5].

Figure 7: Mechanism for C addition between two dimer rows. Also known as the
“trough bridging” mechanism, this begins the formation of a new dimer row
along the C(100) surface. Reaction energies (in eV) are presented between each
step in the reaction, and where applicable we also include barrier energies in
square brackets.

##  4 Encapsulation of N

Surface growth may therefore proceed unimpeded until the encapsulation of
surface-embedded N by the new layer. A four-fold N bonding configuration
requires C addition across the trough between the N and the C in the adjacent
dimer row. However, a trough bridging process analogous to that depicted in
Figure 7 is not viable. This scenario is depicted in Figure 8. In state 1, a
radical CH2 adsorbate attempts to form a covalent bond with surface-embedded
N. However, our calculations indicate that such a bond is not energetically
stable. Similarly, it is also not possible for a CH3 radical to adsorb
directly onto surface-embedded N. Any covalent bonding is prohibited by the
stability of the N lone pair.

We find that a stable bond between N and the bridging CH2 adsorbate is only
possible following abstraction of a nearby surface H atom. For example, H
abstraction on the newly formed dimer unit in step $1\rightarrow 2$ of Figure
8. This creates an acceptor site for a migrating charge from the N lone pair.
In step $2\rightarrow 3$, the CH2 radical is able to bridge the trough and
form a stable covalent bond with surface-embedded N. Encapsulation of surface-
embedded N defect therefore requires breaking its lone pair and charge
migration to a surface acceptor site. While state 3 is stable, reaction
$2\rightarrow 3$ is endothermic and indicates that the CH2 adsorbate
preferentially remains radical. This reflects the stability of the N lone
pair.

Figure 8: Encapsulation of surface-embedded N through C addition between two
dimer rows. In state 1, the bridging CH2 radical is unable to bond with
surface-embedded N. The N lone pair prevents formation of a C–N covalent bond.
However, abstraction of a nearby surface H (step $1\rightarrow 2$) creates an
acceptor site permitting charge migration of a N electron. This breaks the N
lone pair and allows for covalent bonding with the bridging CH2 in step
$2\rightarrow 3$. Reaction energies (in eV) are presented between each step in
the reaction, and where applicable we also include barrier energies in square
brackets.

Exothermic encapsulation of surface-embedded N can still occur through a
different reaction pathway. Consider state 1 depicted in Figure 9. In this
scenario, growth of a new layer is complete with the exception of
encapsulation of surface-embedded N. Methyl adsorption and subsequent H
abstractions result in a CH2 adsorbate attempting to bridge the final trough
between a C and the surface-embedded N. In step $1\rightarrow 2$, the bridging
C adsorbs to the sub-surface C (analogous to step $5\rightarrow 6$ in Figure
6). In step $2\rightarrow 3$, H abstraction on the bridging C creates a
radical site ready to bond with the surface-embedded N. During this mechanism
the remaining H on the bridging C rotates away from the N lone pair, greatly
reducing steric interactions. Consequently, step $2\rightarrow 3$ has a low
barrier energy of 0.085 eV and results in a 1.025 eV increase in stability.

Once again, encapsulation of surface-embedded N requires breaking the lone
pair to accommodate covalent bonding with the C radical. This is only possible
if a lone-pair electron can be donated to a suitable acceptor site. In step
$3\rightarrow 4$ this is realised through H abstraction of a surface H. In
this circumstance N donates an electron to the acceptor site and bonds with
the bridging C immediately following abstraction. This fully encapsulates the
surface-embedded N and forms a positively-charged defect and negatively-
charged acceptor site on the surface. Alternatively, given that during CVD
growth anywhere between $1-10$% of the surface consists of acceptor sites
rather than H termination (depending on growth conditions[33]), step
$3\rightarrow 4$ may occur immediately without further H abstraction. In
contrast to step $2\rightarrow 3$ in Figure 8, N encapsulation is exothermic
by 765 meV. This can be attributed to the strain of the C–C dimer bond of
state 3 in Figure 9. Not only does N encapsulation satisfy all radicals, it
results in the formation of an unstrained surface dimer.

The facile pathway from surface-embedded to substitutional N identified in
Figure 9 provides a potential explanation for the low yield of NV centers
during CVD growth. Experimentally, it has been demonstrated that at most 0.5%
of incorporated N forms NV centers[56]. Experimental results suggest that NV
centers form as a unit on (110) surfaces, with N incorporating into the
lattice first before the the vacancy is produced in the overgrown layer[57].
Ab initio calculations support this hypothesis and predict a similar mechanism
for NV formation on the (100) surface[9].

This would therefore suggest that the overgrowth mechanism in Figure 9 is a
critical point which determines either formation of substitutional N or an NV
center. If the CH2 adsorbate were to be prevented from bridging the final
trough (or CH3 initially adsorbing to the surface), it is conceivable that a
vacancy would form during overgrowth of the next layer (the second layer).
While we do not identify any such vacancy formation mechanism in this work,
one possibility could be if growth of the second layer occurs simultaneously
with the first layer. The step-edge of the secondary layer could then overgrow
the local surface depicted in Figure 9 before the reaction can proceed,
subsequently forming an NV center. However, given the low energy requirements
for N encapsulation in Figure 9 this is unlikely, thus giving rise to the low
NV to substitutional N ratios observed experimentally.

Figure 9: Encapsulation of surface-embedded N. Here we consider complete
growth of a new surface layer with the exception of a final bond between the
surface-embedded N and a bridging C. Bonding only occurs if the N lone pair is
broken through donation of an electron to an acceptor site (step $3\rightarrow
4$).

##  5 Conclusion

Three decades of ab initio research has produced no consensus on the atomistic
mechanism for N catalysis of CVD diamond growth. In this work we identify a
new catalytic effect relevant to step-flow dominated growth modes on (100)
surfaces. Specifically, surface-embedded N drastically enhances the rate of
new layer nucleation by reducing the energy barrier for C insertion. The
presence of N increases the rate of the ring opening/closing reaction by a
factor of 400, the greatest enhancement to any key (100) diamond growth
process to date. Experimental support for our nucleation model may be obtained
through correlation of surface-embedded N and nucleated layers following CVD.
One way to achieve this could be use of an NV-based quantum microscope to
locate N defects in the surface through their hyperfine structure[58].
Conventional surface imaging techniques such as STM could then be used to
identify surrounding layer growth which has propagated from the defect[59].

This work has also established the first atomic model describing overgrowth of
N during (100) diamond CVD. We determine that surface-embedded N does not
impede typical layer growth processes. N maintains its strong electron lone
pair and does not participate in bonding during formation of the dimer unit or
C bridging of non-adjacent troughs. The exception is encapsulation of surface-
embedded N to form a sub-surface defect through C bridging on the adjacent
trough. This reaction requires donation of a lone pair electron to a surface
acceptor site. Moreover, the process may be endothermic or exothermic
depending on the presence of surrounding new-layer growth. The resulting sub-
surface defect is four-fold coordinated and positively charged.

## Acknowledgments

We acknowledge funding from the Australian Research Council (DE170100169).
This research was undertaken with the assistance of resources and services
from the National Computational Infrastructure (NCI), which is supported by
the Australian Government.

## References

  * [1] R. Locher, C. Wild, N. Herres, D. Behr, and P. Koidl.  Nitrogen stabilized ¡100¿ texture in chemical vapor deposited diamond films.  Applied Physics Letters, 65(1):34–36, 1994.
  * [2] S. Jin and T. D. Moustakas.  Effect of nitrogen on the growth of diamond films.  Applied Physics Letters, 65(4):403–405, 1994.
  * [3] R. Samlenski, C. Haug, R. Brenn, C. Wild, R. Locher, and P. Koidl.  Incorporation of nitrogen in chemical vapor deposition diamond.  Applied Physics Letters, 67(August 1995):2798, 1995.
  * [4] Michael N.R. Ashfold, Jonathan P. Goss, Ben L. Green, Paul W. May, Mark E. Newton, and Chloe V. Peaker.  Nitrogen in Diamond.  Chemical Reviews, 2020.
  * [5] J. Achard, F. Silva, O. Brinza, A. Tallaire, and A. Gicquel.  Coupled effect of nitrogen addition and surface temperature on the morphology and the kinetics of thick CVD diamond single crystals.  Diamond and Related Materials, 16(4-7 SPEC. ISS.):685–689, 2007\.
  * [6] S. Bohr, R. Haubner, and B. Lux.  Influence of nitrogen additions on hot-filament chemical vapor deposition of diamond.  Applied Physics Letters, 68(8):1075–1077, 1996.
  * [7] T. H. Taminiau, J. Cramer, T. Van Der Sar, V. V. Dobrovitski, and R. Hanson.  Universal control and error correction in multi-qubit spin registers in diamond.  Nature Nanotechnology, 9(3):171–176, 2014.
  * [8] Stephanie Wehner, David Elkouss, and Ronald Hanson.  Quantum internet: A vision for the road ahead.  Science, 362(6412), 2018.
  * [9] M K Atumi, J P Goss, P R Briddon, and M J Rayson.  Atomistic modeling of the polarization of nitrogen centers in diamond due to growth surface orientation.  Phys. Rev. B, 88(24):245301, dec 2013.
  * [10] Marcus W. Doherty, Neil B. Manson, Paul Delaney, Fedor Jelezko, Jörg Wrachtrup, and Lloyd C.L. Hollenberg.  The nitrogen-vacancy colour centre in diamond.  Physics Reports, 528(1):1–45, 2013.
  * [11] Th. Frauenheim, G. Jungnickel, P. Sitch, M. Kaukonen, F. Weich, J. Widany, and D. Porezag.  A molecular dynamics study of N-incorporation into carbon systems: Doping, diamond growth and nitride formation.  Diamond and Related Materials, 7(2-5):348–355, 1998.
  * [12] M. Kaukonen, P. Sitch, and G. Jungnickel.  Effect of N and B doping on the growth of CVD diamond surfaces.  Physical Review B - Condensed Matter and Materials Physics, 57(16):9965–9970, 1998.
  * [13] Tanguy Van Regemorter and Karin Larsson.  A theoretical study of nitrogen-induced effects on initial steps of diamond CVD growth.  Chemical Vapor Deposition, 14(7-8 SPEC. ISS.):224–231, 2008.
  * [14] T. Van Regemorter and K. Larsson.  Effect of a NH coadsorbate on the CH3 (or CH2) adsorption to a surface step on diamond (100).  Journal of Physical Chemistry C, 113(46):19891–19896, 2009.
  * [15] T. Van Regemorter and K. Larsson.  Effect of substitutional N on important chemical vapor deposition diamond growth steps.  Journal of Physical Chemistry A, 113(13):3274–3284, 2009.
  * [16] Z. Yiming, F. Larsson, and K. Larsson.  Effect of CVD diamond growth by doping with nitrogen.  Theoretical Chemistry Accounts, 133(2):1–12, 2014.
  * [17] Jeung Ku Kang and Charles B. Musgrave.  A theoretical study of the chemical vapor deposition of (100) diamond: An explanation for the slow growth of the (100) surface.  The Journal of Chemical Physics, 113(17):7582–7587, 2000.
  * [18] Andrew Cheesman.  Investigations into the fundamentals of gas-phase and gas-surface chemistry prevalent in growth of Chemical Vapour Deposited diamond films.  Phd, University of Bristol, 2006.
  * [19] Andrew Cheesman, Jeremy N. Harvey, and Michael N. R. Ashfold.  Studies of Carbon Incorporation on the Diamond {100} Surface during Chemical Vapor Deposition using Density Functional Theory.  The Journal of Physical Chemistry A, 112(45):11436–11448, 2008\.
  * [20] I. I. Oleinik, D. G. Pettifor, A. P. Sutton, and J. E. Butler.  Theoretical study of chemical reactions on CVD diamond surfaces.  Diamond and Related Materials, 9(3-6):241–245, 2000.
  * [21] Hiroyuki Tamura and Mark S. Gordon.  Ab initio study of nucleation on the diamond (1 0 0) surface during chemical vapor deposition with methyl and H radicals.  Chemical Physics Letters, 406(1-3):197–201, 2005.
  * [22] Barbara J Garrison, Erik J Dawnkaski, Deepak Srivastava, and Donald W Brenner.  Molecular Dynamics Simulations of Dimer Opening on a Diamond {001}(2x1) Surface.  Science, 255(5046):835–838, 1992.
  * [23] Stephen J. Harris.  Mechanism for diamond growth from methyl radicals.  Applied Physics Letters, 56(23):2298–2300, 1990.
  * [24] Michael Frenklach and Hai Wang.  Detailed surface and gas-phase chemical kinetics of diamond deposition.  43(2):1520–1544, 1991.
  * [25] D. Huang, M. Frenklach, and M. Maroncelli.  Energetics of acetylene-addition mechanism of diamond growth.  Journal of physical chemistry, 92(22):6379–6381, 1988.
  * [26] Michael Frenklach and Karl E Spear.  Growth mechanism of vapor-deposited diamond.  Journal of Materials Research, 3(1):133–140, 1988.
  * [27] R. Samlenski, C. Haug, R. Brenn, C. Wild, R. Locher, and P. Koidl.  Characterisation and lattice location of nitrogen and boron in homoepitaxial CVD diamond.  Diamond and Related Materials, 5(9):947–951, 1996.
  * [28] P. W. May, M. N.R. Ashfold, and Yu A. Mankelevich.  Microcrystalline, nanocrystalline, and ultrananocrystalline diamond chemical vapor deposition: Experiment and modeling of the factors controlling growth rate, nucleation, and crystal size.  Journal of Applied Physics, 101(5), 2007.
  * [29] Yuri A. Mankelevich, Michael N.R. Ashfold, and Jie Ma.  Plasma-chemical processes in microwave plasma-enhanced chemical vapor deposition reactors operating with C/H/Ar gas mixtures.  Journal of Applied Physics, 104(11), 2008.
  * [30] James E. Butler and Ivan Oleynik.  A mechanism for crystal twinning in the growth of diamond by chemical vapour deposition.  Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences, 366(1863):295–311, 2008.
  * [31] Armando Netto and Michael Frenklach.  Kinetic Monte Carlo simulations of CVD diamond growth—Interlay among growth, etching, and migration.  Diamond and Related Materials, 14(10):1630–1646, 2005.
  * [32] C. C. Battaile, D. J. Srolovitz, and J. E. Butler.  Morphologies of diamond films from atomic-scale simulations of chemical vapor deposition.  Diamond and Related Materials, 6(9):1198–1206, 1997.
  * [33] J E Butler, Y A Mankelevich, A Cheesman, Ma Jie, and M N R Ashfold.  Understanding the chemical vapor deposition of diamond: recent progress.  Journal of Physics: Condensed Matter, 21(36):364201, 2009.
  * [34] Xianglin Li, James Perkins, Ramon Collazo, Robert J. Nemanich, and Zlatko Sitar.  Investigation of the effect of the total pressure and methane concentration on the growth rate and quality of diamond thin films grown by MPCVD.  Diamond and Related Materials, 15(11-12 SPEC. ISS.):1784–1788, 2006\.
  * [35] Paul W. May, Neil L. Allan, Michael N.R. Ashfold, James C. Richley, and Yuri A. Mankelevich.  Simplified Monte Carlo simulations of chemical vapour deposition diamond growth.  Journal of Physics Condensed Matter, 21(36), 2009.
  * [36] W. J. Rodgers, P. W. May, N. L. Allan, and J. N. Harvey.  Three-dimensional kinetic Monte Carlo simulations of diamond chemical vapor deposition.  Journal of Chemical Physics, 142(21), 2015.
  * [37] Alastair Stacey, Kane M O’Donnell, Jyh-Pin Chou, Alex Schenk, Anton Tadich, Nikolai Dontschuk, Jiri Cervenka, Chris Pakes, Adam Gali, Alon Hoffman, and Steven Prawer.  Nitrogen Terminated Diamond.  Advanced Materials Interfaces, 2(10):1500079, 2015.
  * [38] J. A. Steckel, T. Phung, K. D. Jordan, and P. Nachtigall.  Concerted Use of Slab and Cluster Models in an ab Initio Study of Hydrogen Desorption from the Si(100) Surface.  Journal of Physical Chemistry B, 105(18):4031–4038, 2001.
  * [39] D. F. Tracey, B. Delley, D. R. McKenzie, and O. Warschkow.  Molecular adsorption on silicon (001): A systematic evaluation of size effects in slab and cluster models.  AIP Advances, 3(4), 2013.
  * [40] Sandeep Nigam and Chiranjib Majumder.  Comparison between cluster and slab model for Pt-group atom adsorption on gold and silver substrate.  Surface Science, 630:78–84, 2014.
  * [41] Thomas A. Halgren and William N. Lipscomb.  The synchronous-transit method for determining reaction pathways and locating molecular transition states.  Chemical Physics Letters, 49(2):225–232, 1977.
  * [42] Li Xie, Haiyan Liu, and Weitao Yang.  Adapting the nudged elastic band method for determining minimum-energy paths of chemical reactions in enzymes.  Journal of Chemical Physics, 120(17):8039–8052, 2004.
  * [43] Graeme Henkelman, Gísli Jóhannesson, and Hannes Jónsson.  Methods for Finding Saddle Points and Minimum Energy Paths, pages 269–302.  Springer Netherlands, Dordrecht, 2002.
  * [44] G. Kresse and J. Furthmüller.  Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set.  Computational Materials Science, 6(1):15–50, 1996.
  * [45] G. Kresse and J. Hafner.  Ab initio molecular-dynamics simulation of the liquid-metalamorphous- semiconductor transition in germanium.  Physical Review B, 49(20):14251–14269, 1994.
  * [46] G. Kresse and J. Furthmüller.  Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set.  Physical Review B - Condensed Matter and Materials Physics, 54(16):11169–11186, 1996.
  * [47] D. Joubert.  From ultrasoft pseudopotentials to the projector augmented-wave method.  Physical Review B - Condensed Matter and Materials Physics, 59(3):1758–1775, 1999.
  * [48] John P. Perdew, Kieron Burke, and Matthias Ernzerhof.  Generalized gradient approximation made simple.  Physical Review Letters, 1996.
  * [49] Chengteh Lee, Weitao Yang, and Robert G. Parr.  Development of the Colle-Salvetti correlation-energy formula into a functional of the electron density.  Physical Review B, 1988.
  * [50] Benjamin G Janesko.  Density Functional Theory Beyond the Generalized Gradient Approximation for Surface Chemistry, pages 25–51.  Springer International Publishing, Cham, 2015.
  * [51] Gábor I Csonka, John P Perdew, and Adrienn Ruzsinszky.  Global Hybrid Functionals: A Look at the Engine under the Hood.  Journal of Chemical Theory and Computation, 6(12):3688–3703, 2010\.
  * [52] Graeme Henkelman, Blas P Uberuaga, and Hannes Jónsson.  A climbing image nudged elastic band method for finding saddle points and minimum energy paths.  The Journal of Chemical Physics, 113(22):9901–9904, 2000.
  * [53] Péter Pulay.  Convergence acceleration of iterative sequences. the case of scf iteration.  Chemical Physics Letters, 73(2):393–398, 1980.
  * [54] Junwei Lucas Bao and Donald G Truhlar.  Variational transition state theory: theoretical framework and recent developments.  Chemical Society Reviews, 46(24):7548–7596, 2017.
  * [55] David G. Goodwin and James E. Butler.  Theory of Diamond Chemical Vapor Deposition.  In Mark Prelas, Galina Popovici, and Louis Bigelow, editors, Handbook of Industrial Diamonds and Diamond Films, pages 527–581. Taylor & Francis Inc, Bosa Roca, 1st edition, 1997.
  * [56] Mikhail A Lobaev, Alexei M Gorbachev, Sergey A Bogdanov, Anatoly L Vikharev, Dmitry B Radishev, Vladimir A Isaev, and Mikhail N Drozdov.  NV-Center Formation in Single Crystal Diamond at Different CVD Growth Conditions.  physica status solidi (a), 215(22):1800205, nov 2018.
  * [57] A M Edmonds, U F S D’Haenens-Johansson, R J Cruddace, M E Newton, K.-M. C Fu, C Santori, R G Beausoleil, D J Twitchen, and M L Markham.  Production of oriented nitrogen-vacancy color centers in synthetic diamond.  Phys. Rev. B, 86(3):35201, jul 2012.
  * [58] L Thiel, Z Wang, M A Tschudin, D Rohner, I Gutiérrez-Lezama, N Ubrig, M Gibertini, E Giannini, A F Morpurgo, and P Maletinsky.  Probing magnetism in 2D materials at the nanoscale with single-spin microscopy.  Science, 364(6444):973–976, 2019.
  * [59] Kirill Bobrov, Andrew J Mayne, and Gérald Dujardin.  Atomic-scale imaging of insulating diamond through resonant electron injection.  Nature, 413(6856):616–619, 2001.