Nitrogen-vacancy (NV) centers in diamond have attracted a great deal of attention because of their possible use in information processing and electromagnetic sensing technologies. We examined the atomistic generation mechanism for the NV defect aligned in the [111] direction of C(111) substrates. We found that N is incorporated in the C bilayers during the lateral growth arising from a sequence of kink propagation along the step edge down to [1¯1¯2]. As a result, the atomic configuration with the N-atom lone-pair pointing in the [111] direction is formed, which causes preferential alignment of NVs. Our model is consistent with recent experimental data for perfect NV alignment in C(111) substrates.

Nitrogen-vacancy (NV) centers in diamond1 have been recognized as representative examples of entangled S = 1 systems in solid-state materials. Because of the externally controllable entanglement between the electronic and nuclear spin states2–4 and excellent coherence properties in time5 as well as space,6 numerous technological applications of the NV centers have been demonstrated to date.7–12 

An NV center in a diamond crystal is composed of an N atom substituting for a carbon atom with an adjacent carbon vacancy (V) sitting in various orientations relative to N. For example, the NVs near the C(111) surface could have eight possible NV orientations, four with N in one sublattice in the bilayer [α layer, Fig. 1(a)] and the remaining four with N in the other sublattice [β layer, Fig. 1(b)].

FIG. 1.

NV center near the C(111):H surface. In panels (a) and (b), the gray balls adjacent to N illustrate four possible positions of V. (a) α and β represent the two triangular sublattices in one bilayer. The α and β layers for third bilayer are omitted for clarity. (a) N in the α sublattice of the second bilayer. (b) N in the β sublattice of the first bilayer.

FIG. 1.

NV center near the C(111):H surface. In panels (a) and (b), the gray balls adjacent to N illustrate four possible positions of V. (a) α and β represent the two triangular sublattices in one bilayer. The α and β layers for third bilayer are omitted for clarity. (a) N in the α sublattice of the second bilayer. (b) N in the β sublattice of the first bilayer.

Close modal

For practical applications, it is desirable to have the NV orientations aligned in one direction. In spin-dependent NV fluorescence experiments in C(001) samples grown by chemical vapor deposition (CVD), the magnetic field sensitivity was doubled relative to samples with random NV orientations.13 This was because most of NV centers were arranged in two of four possible orientations. The partial alignment of the NV centers was also reported for CVD grown (110) substrates.14 

Recently, various groups15–17 have reported the nearly perfect alignment of the NV centers along the [111] direction in C(111) samples grown by CVD (94% ± 2%,15 ∼97%,16 and ≥99% (Ref. 17) relative to the total number of NVs generated). The formation mechanism was explained in terms of the step-flow growth16 in a similar way to the (110) case.14 A first-principles study published prior to these experiments predicted that V and N are in the first and second layers, respectively [Fig. 1(b)].18 Although the theory also discussed how the N is located in the α layer [Fig. 1(a)], the atomistic reason for the selective NV alignment in [111] (Refs. 15–17) is still unknown.

In this study, we address this question by using first-principles energetics. We found that the N incorporated at the kink site of the [1¯1¯2] step edge is embedded in the subsurface in a sequence of kink-flows along this step and finally becomes the [111]-oriented NV center. Our theoretical scenario qualitatively explains the experimental results.15–17 

We assumed that the N atoms participating in the formation of NV centers are supplied from the vapor phase during CVD. This assumption is justified by the experimental observation that the NV centers are aligned in the CVD growth directions for all (100),13 (110),14 and (111) substrates.15–17 Hence, we studied the behavior of an external N atom attached to C(111) growing in the step-flow mode. The step-down direction was [1¯1¯2] in accordance with structure identification experiments for the steps on C(111).19,20 A kink was introduced as a terminator of the step edge because these experiments19,20 also show that the lateral step-flow growth of C(111) in CVD proceeds via kink flow. We did not consider the case of N incorporated at the β carbon layers, because the displacement of N to the β site, once incorporated in the α layer of the kink,17 should be negligible before the subsequent carbon deposition.

An example of the kinked step edges considered is shown in Fig. 2. The kink flow along the step edge arises from the movement of the kink site by the successive attachment of carbon atoms to the kink. The step-flow growth of C(111) is in units of bilayers,19–22 which are double layers composed of α and β layers [see Fig. 1(a)]. Thus, we restricted ourselves to considering this kink-flow growth mode that proceeds in units of two carbon atoms, one in the α layer and the other in the β layer. According to the energetics, the β carbon atom should form a kink terminal, shown as xk in Fig. 2(a).

FIG. 2.

Atomistic kink flow processes on C(111) considered in this study. Blue and yellow circles represent carbon atoms in the α and β sublattice layers, respectively. Two carbon bilayers are shown from the topmost surface; the first and second bilayers are represented in dark and light colors, respectively. The H atoms are omitted for clarity. See Ref. 30 for how the surface is hydrogenated in the calculation. The dashed green line segments are visual guides to show the [1¯10] direction of the step edge along which the kink flow occurs. The step-down direction is [1¯1¯2] for both panels. (a) Process I, in which the kink terminated at xk grows to xk+1 by absorbing two C atoms. (b) Process II, in which a C atom occupies the β site at xk+2 and leaves a vacancy at xk+1.

FIG. 2.

Atomistic kink flow processes on C(111) considered in this study. Blue and yellow circles represent carbon atoms in the α and β sublattice layers, respectively. Two carbon bilayers are shown from the topmost surface; the first and second bilayers are represented in dark and light colors, respectively. The H atoms are omitted for clarity. See Ref. 30 for how the surface is hydrogenated in the calculation. The dashed green line segments are visual guides to show the [1¯10] direction of the step edge along which the kink flow occurs. The step-down direction is [1¯1¯2] for both panels. (a) Process I, in which the kink terminated at xk grows to xk+1 by absorbing two C atoms. (b) Process II, in which a C atom occupies the β site at xk+2 and leaves a vacancy at xk+1.

Close modal

We considered two elementary processes of kink flow (Fig. 2). In process I [Fig. 2(a)], two carbon atoms (shown as circles bounded by dotted curves) are attached to the kink at xk. Then the xk+1 site forms a new kink site in order to minimize the kink energy. In process II [Fig. 2(b)], a carbon atom occupies the β site at xk+2 while the other C atom is the same as that in Fig. 2(a). The latter process creates a vacancy at the β site at xk+1 in the surface bilayer, which is critical to the formation of the NV center aligned in [111].

The C(111) surface structure was modeled using a five-bilayer thick slab in a repeated supercell with 5 × 4 periodicity in the surface directions. A vacuum layer about 1.1 nm thick was inserted to decouple the spurious interactions between the slab and its self-images. All the C atoms on the terrace were hydrogenated. For the C atoms at the step edge, we investigated both the hydrogen-rich and hydrogen-poor cases. The hydrogen-rich step edge was constructed by attaching two H atoms to each C atom in the step edge. Lower-terrace C atoms nearest to the step-edge were terminated with one H atom each. In the hydrogen-poor step edge, only one H atom was attached to each C step-edge atom and the corresponding lower-terrace C atoms were not hydrogenated except the one closest to the kink position, to which one H atom was attached.

To calculate the total energies of the N impurity in the surface, we performed first-principles electronic structure calculations23–25 throughout this work based on density functional theory (DFT)26,27 within the generalized gradient approximation to the exchange-correlation energy.28 Further details of a calculation method are described in Ref. 29. Because the NV centers are found in not only the negative (NV) but also neutral (NV0) charge states,1 we performed structure optimization for the NV defects in the respective charge states. For that purpose, we constructed two sets of supercell geometries: NV in C(111) in the one and NV0 in C(111) in the other. Then, we optimized structures of the NV and NV0 centers separately to calculate their total energies.

The entry point for an externally supplied N atom into diamond is the α site at the kink, as shown in Fig. 3(a). The preference for the α site over the β site was theoretically identified in a previous study.17 

FIG. 3.

Growth of [111]-aligned NV center in C(111). Red circles represent the N atoms. The other conventions for the atoms are the same as in Fig. 2(a) N at the α site of kink ∼0.5 eV lower in energy than the N at β, irrespective of the hydrogen content at the step edge.17 Panels (b) and (c) are possible structures in stage (1). The N at the kink is about to be incorporated in the surface terrace via process I [(b)] and II [(c)]. Panels (d) and (e) are possible structures in stage (2). The N in the surface terrace is about to be overgrown via process II [(d)] and I [(e)].

FIG. 3.

Growth of [111]-aligned NV center in C(111). Red circles represent the N atoms. The other conventions for the atoms are the same as in Fig. 2(a) N at the α site of kink ∼0.5 eV lower in energy than the N at β, irrespective of the hydrogen content at the step edge.17 Panels (b) and (c) are possible structures in stage (1). The N at the kink is about to be incorporated in the surface terrace via process I [(b)] and II [(c)]. Panels (d) and (e) are possible structures in stage (2). The N in the surface terrace is about to be overgrown via process II [(d)] and I [(e)].

Close modal

The structures studied here are shown in Figs. 3(b)–3(e), 4(a), and 4(b). We classify these into three stages depending on the location of N with each stage having two branches. In stage (1), N is incorporated in the topmost surface terrace [Figs. 3(b) and 3(c)], and then in stage (2) the C bilayer is overgrown [Figs. 3(d) and 3(e)]. In stage (3), the NV is embedded in the second bilayer [Figs. 4(a) and 4(b)]. In stage (4) [Fig. 4(c)], the NV is separated from the step edge. It should be noted that we are considering the case where the position of N is not changed but the kink and step edge move toward [1¯10] and [1¯1¯2], respectively, through the addition of the carbon atoms traveling from the vapor phase.

FIG. 4.

(a) and (b) A [111]-aligned NV center in stage (3) near completion. Red circles represent the N atoms. The other conventions for the atoms are the same in Fig. 2(a) The NV is surrounded by the carbon atoms in the surface bilayer through process I [Fig. 2(a)]. (b) The structure that is degenerate with (a) in terms of energy when the step edge is hydrogen-poor. (c) The NV completely embedded in the C(111) surface, which is final stage (4). This configuration can be obtained from both (a) and (b) by adding two C atoms.

FIG. 4.

(a) and (b) A [111]-aligned NV center in stage (3) near completion. Red circles represent the N atoms. The other conventions for the atoms are the same in Fig. 2(a) The NV is surrounded by the carbon atoms in the surface bilayer through process I [Fig. 2(a)]. (b) The structure that is degenerate with (a) in terms of energy when the step edge is hydrogen-poor. (c) The NV completely embedded in the C(111) surface, which is final stage (4). This configuration can be obtained from both (a) and (b) by adding two C atoms.

Close modal

Comparing the total energies of the two structures in each stage (Table I) indicates that the [111]-aligned NV center can grow via the structural evolution shown in Figs. 3(a), 3(b), 3(d), and 4(a). In stage (1), the N atoms [Figs. 3(b) and 3(c)] are about to be embedded in the surface terrace as a result of adding two carbon atoms to the kink at xk via processes I [Fig. 3(b)] and II [Fig. 3(c)]. The energy of the structure in Fig. 3(b) is much lower than that of Fig. 3(c), irrespective of the hydrogenation content at the step edge (Table I). Thus, the lateral step-flow growth arising from the kink flow can embed the N atom at the step edge into the surface just like the step flow of pure carbon bilayers.

TABLE I.

The total energies of the structures in Figs. 3(b), 3(d), and 4(a) relative to those in Figs. 3(c), 3(e), and 4(b), respectively. We list the values calculated for the slabs of five carbon bilayers, convergent within 0.1 eV per N in both charge states.

H-richH-poorH-richH-poor
N0NN0NNV0NVNV0NV
Stage(eV)(eV)Stage(eV)(eV)
(1) −2.6 −2.3 −1.0 −1.0 (3) −0.5 −0.3 +0.1 +0.6 
(2) −1.6 −1.6 −1.4 −1.4  
H-richH-poorH-richH-poor
N0NN0NNV0NVNV0NV
Stage(eV)(eV)Stage(eV)(eV)
(1) −2.6 −2.3 −1.0 −1.0 (3) −0.5 −0.3 +0.1 +0.6 
(2) −1.6 −1.6 −1.4 −1.4  

Stage (2) [Figs. 3(d) and 3(e)] explains how the vacancy is created on top of the N atom to produce the NV center aligned in the [111] direction. Here, process II occurs and contributes to the formation of the [111]-aligned NV center. The geometries in stage (2) illustrate the location in which the N atoms, now incorporated in the terrace, are about to be overgrown by additional carbon bilayers. The structure of Fig. 3(d) corresponds to process II in the attachment of two carbon atoms to the kink at xk. The total energy of this structure in both charge states is ∼1.5 eV lower on average than its counterpart in Fig. 3(e), generated by process I, insensitive to the hydrogenation level of the step edge (Table I). The results for stage (2) arise from the stability of the electron lone pair of the N atom in the surface terrace that is threefold coordinated.

The NV center is then grown on the C(111) surface in stage (3) [Figs. 4(a) and 4(b)]. The total energy difference between the structures in Figs. 4(a) and 4(b) shows that the completion of the NV center in process I of the kink flow may be favored when the step edge is hydrogen-rich. For the hydrogen-poor step edge, the two structures are close in energy for NV0 and the structure in Fig. 4(b) is favored for NV (Table I). Thus, a high hydrogen chemical potential is necessary to proceed to stage (4) [Fig. 4(c)], in which the NV center is created.

The atomistic mechanism that we have described above delivers several messages. (i) V is closer to the vacuum than N, as shown in Fig. 1(a). The lone-pair orbital of an N atom placed in the α layers, on which a V is created, points into the vacuum. Experimental data show that a majority of the NV centers have this structure.15 Combining the C(111) bilayer structure and the energy gain arising from the reduction of the surface dangling bonds, the most achievable stacking order of carbon layers from the surface to the bulk during CVD is αβαβ⋯, where the N atoms are located only in the α layers. This is in accordance with a suggestion by a preceding first-principles study.18 (ii) Hydrogen-rich environments are required for preferential growth of the [111]-aligned NV centers. The experiments in which the nearly perfect alignment of the NV centers occurs in the samples grown with the 1.5% CH4/H2 ratio or less15–17 imply that as many of the surface carbon atoms as possible should be terminated with H atoms. Our calculation results are consistent overall with this requirement. The exception is that, in stage (3), the NV defect at the H-poor step edge may occur, in which the V is not fully isolated in the terrace [Fig. 4(b)], depending on the chemical potential of H. The grand potential of the H-rich (H-poor) step edge is lower than that of the H-poor (H-rich) one, if we define the hydrogen chemical potential as the free energy of H atoms (H2 molecules) at a temperature of ∼1000 °C and a pressure of ∼150 Torr, similar to the experimental conditions.15–17 To refine the modeling resolution of these competing energetics further, it may be necessary to take account of the vibrational entropies of the substrate at finite temperatures and/or the deviation of the H-atom content from equilibrium that may occur in reactors by the excitation of high-density H plasmas. Although both these refinements are beyond the scope of the present study, they are interesting future areas of research. (iii) Throughout this study, our results suggest that the [111]-aligned NV centers are grown in the crystal during CVD, which is consistent with the argument based on experimental results.14 The growth scenario is further supported by a recent advanced DFT study that shows that V diffusion may induce the formation of divacancies (V2) instead of NVs, and also the immobilization of V because of the production of the V defects.31 This study also shows that the ratio [NV]/[NV0] is a decreasing function of [V2]/[Ns], where Ns is a substitutional N. The Fermi energy decreases as [V2]/[Ns] increases, too. The two relationships both originate from the fact that the V2 is more electronegative than NV. This suggests that carbon dangling bonds may behave like electron absorbers. Turning to our system, the H-poor step edge may be regarded as a reservoir of carbon dangling bonds, since the C-C distance in the edge is 13% larger than the bulk C-C counterpart so that the self-termination among the C dangling bonds is relatively weak. So we expect that the [NV0]/[NV] ratio could be smaller on the H-rich step edge than the H-poor step edge if the charge transfer properties among these defects are similar in the surface and bulk regions. On the other hand, another recent first-principles study32 proposes that nearly 90% of NV defects can be aligned in the [111] direction if the sample is under a bi-axial, 2% compressive strain normal to [111] and thermally annealed at 970 °C. This supports the post-growth process for the pre-existing NV defects. It should be noted that the question of which of the alignment mechanisms is more feasible is unresolved at present, which is another intriguing issue.

In conclusion, we have obtained an atomistic mechanism for the production of NV centers selectively aligned in the [111] direction of the C(111) surface. Starting from the N-substitution of the α site of the kink, N is incorporated in the surface during the lateral bilayer growth, leaving a vacancy attached to top of the N site. Our model is consistent with current experimental data.

We used the VESTA ver.3.1.1 software33 to draw the atomic models in Fig. 1.

1.
M. W.
Doherty
,
N. B.
Manson
,
P.
Delaney
,
F.
Jelezko
,
J.
Wrachtrup
, and
L. C. L.
Hollenberg
,
Phys. Rep.
528
,
1
(
2013
).
2.
F.
Jelezko
,
T.
Gaebel
,
I.
Popa
,
M.
Domhan
,
A.
Gruber
, and
J.
Wrachtrup
,
Phys. Rev. Lett.
93
,
130501
(
2004
).
3.
M. V.
Gurudev Dutt
,
L.
Childress
,
L.
Jiang
,
E.
Togan
,
J.
Maze
,
F.
Jelezko
,
A. S.
Zibrov
,
P. R.
Hemmer
, and
M. D.
Lukin
,
Science
316
,
1312
(
2007
).
4.
P.
Neumann
,
N.
Mizuochi
,
F.
Rempp
,
P.
Hemmer
,
H.
Watanabe
,
S.
Yamasaki
,
V.
Jacques
,
T.
Gaebel
,
F.
Jelezko
, and
J.
Wrachtrup
,
Science
320
,
1326
(
2008
).
5.
G.
Balasubramanian
,
P.
Neumann
,
D.
Twitchen
,
M.
Markham
,
R.
Kolesov
,
N.
Mizuochi
,
J.
Isoya
,
J.
Achard
,
J.
Beck
,
J.
Tissler
,
V.
Jacques
,
P. R.
Hemmer
,
F.
Jelezko
, and
J.
Wrachtrup
,
Nat. Mater.
8
,
383
(
2009
).
6.
H.
Bernien
,
B.
Hensen
,
W.
Pfaff
,
G.
Koolstra
,
M. S.
Blok
,
L.
Robledo
,
T. H.
Taminiau
,
M.
Markham
,
D. J.
Twitchen
,
L.
Childress
, and
R.
Hanson
,
Nature (London)
497
,
86
(
2013
).
7.
G.
Kucsko
,
P. C.
Maurer
,
N. Y.
Yao
,
M.
Kubo
,
H. J.
Noh
,
P. K.
Lo
,
H.
Park
, and
M. D.
Lukin
,
Nature (London)
500
,
54
(
2013
).
8.
J. R.
Maze
,
P. L.
Stanwix
,
J. S.
Hodges
,
S.
Hong
,
J. M.
Taylor
,
P.
Cappellaro
,
L.
Jiang
,
M. V.
Gurudev Dutt
,
E.
Togan
,
A. S.
Zibrov
,
A.
Yacoby
,
R. L.
Walsworth
, and
M. D.
Lukin
,
Nature (London)
455
,
644
(
2008
).
9.
G.
Balasubramanian
,
I. Y.
Chan
,
R.
Kolesov
,
M.
Al-Hmoud
,
J.
Tisler
,
C.
Shin
,
C.
Kim
,
A.
Wojcik
,
P. R.
Hemmer
,
A.
Krueger
,
T.
Hanke
,
A.
Leitenstorfer
,
R.
Bratschitsch
,
F.
Jelezko
, and
J.
Wrachtrup
,
Nature (London)
455
,
648
(
2008
).
10.
N.
Mizuochi
,
T.
Makino
,
H.
Kato
,
D.
Takeuchi
,
M.
Ogura
,
H.
Okushi
,
M.
Nothaft
,
P.
Neumann
,
A.
Gali
,
F.
Jelezko
,
J.
Wrachtrup
, and
S.
Yamasaki
,
Nat. Photonics
6
,
299
(
2012
).
11.
A.
Lohrmann
,
S.
Pezzagna
,
I.
Dobrinets
,
P.
Spinicelli
,
V.
Jacques
,
J.-F.
Roch
,
J.
Meijer
, and
A. M.
Zaitsev
,
Appl. Phys. Lett.
99
,
251106
(
2011
).
12.
D.
Le Sage
,
K.
Arai
,
D. R.
Glenn
,
S. J.
DeVience
,
L. M.
Pham
,
L.
Rahn-Lee
,
M. D.
Lukin
,
A.
Yacoby
,
A.
Komeili
, and
R. L.
Walsworth
,
Nature (London)
496
,
486
(
2013
).
13.
L. M.
Pham
,
N.
Bar-Gill
,
D.
Le Sage
,
C.
Belthangady
,
A.
Stacey
,
M.
Markham
,
D. J.
Twitchen
,
M. D.
Lukin
, and
R. L.
Walsworth
,
Phys. Rev. B
86
,
121202
(
2012
).
14.
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
,
Phys. Rev. B
86
,
035201
(
2012
).
15.
J.
Michl
,
T.
Teraji
,
S.
Zaiser
,
I.
Jakobi
,
G.
Waldherr
,
F.
Dolde
,
P.
Neumann
,
M. W.
Doherty
,
N. B.
Manson
,
J.
Isoya
, and
J.
Wrachtrup
,
Appl. Phys. Lett.
104
,
102407
(
2014
).
16.
M.
Lesik
,
J.-P.
Tetienne
,
A.
Tallaire
,
J.
Achard
,
V.
Mille
,
A.
Gicquel
,
J.-F.
Roch
, and
V.
Jacques
,
Appl. Phys. Lett.
104
,
113107
(
2014
).
17.
T.
Fukui
,
Y.
Doi
,
T.
Miyazaki
,
Y.
Miyamoto
,
H.
Kato
,
T.
Matsumoto
,
T.
Makino
,
S.
Yamasaki
,
R.
Morimoto
,
N.
Tokuda
,
M.
Hatano
,
Y.
Sakagawa
,
H.
Morishita
,
T.
Tashima
,
S.
Miwa
,
Y.
Suzuki
, and
N.
Mizuochi
,
Appl. Phys. Express
7
,
055201
(
2014
).
18.
M. K.
Atumi
,
J. P.
Goss
,
P. R.
Briddon
, and
M. J.
Rayson
,
Phys. Rev. B
88
,
245301
(
2013
).
19.
T.
Tsuno
,
T.
Tomikawa
,
S.
Shikata
, and
N.
Fujimori
,
J. Appl. Phys.
75
,
1526
(
1994
).
20.
N.
Tokuda
,
M.
Ogura
,
S.
Yamsaki
, and
T.
Inokuma
,
Jpn. J. Appl. Phys., Part 1
53
,
04EH04
(
2014
).
21.
C. C.
Battaile
,
D. J.
Srolovitz
, and
J. E.
Butler
,
Diamond Relat. Mater.
6
,
1198
(
1997
).
22.
C. C.
Battaile
,
D. J.
Srolovitz
, and
J. E.
Butler
,
J. Cryst. Growth
194
,
353
(
1998
).
23.
P.
Giannozzi
,
S.
Baroni
,
N.
Bonini
,
M.
Calandra
,
R.
Car
,
C.
Cavazzoni
,
D.
Ceresoli
,
G. L.
Chiarotti
,
M.
Cococcioni
,
I.
Dabo
,
A.
Dal Corso
,
S.
Fabris
,
G.
Fratesi
,
S.
de Gironcoli
,
R.
Gebauer
,
U.
Gerstmann
,
C.
Gougoussis
,
A.
Kokalj
,
M.
Lazzeri
,
L.
Martin-Samos
,
N.
Marzari
,
F.
Mauri
,
R.
Mazzarello
,
S.
Paolini
,
A.
Pasquarello
,
L.
Paulatto
,
C.
Sbraccia
,
S.
Scandolo
,
G.
Sclauzero
,
A. P.
Seitsonen
,
A.
Smogunov
,
P.
Umari
, and
R. M.
Wentzcovitch
,
J. Phys.: Condens. Matter
21
,
395502
(
2009
), see http://www.quantum-espresso.org.
24.
D.
Vanderbilt
,
Phys. Rev. B
41
,
7892
(
1990
).
25.
K.
Laasonen
,
A.
Pasquarello
,
R.
Car
,
C.
Lee
, and
D.
Vanderbilt
,
Phys. Rev. B
47
,
10142
(
1993
).
26.
P.
Hohenberg
and
W.
Kohn
,
Phys. Rev.
136
,
B864
(
1964
).
27.
W.
Kohn
and
L. J.
Sham
,
Phys. Rev.
140
,
A1133
(
1965
).
28.
J. P.
Perdew
,
K.
Burke
, and
M.
Ernzerhof
,
Phys. Rev. Lett.
77
,
3865
(
1996
).
29.

The ultrasoft pseudopotentials24,25 were adopted to represent the interactions between the valence electrons and ionic cores of all C, N, and H atoms. The wavefunctions were expanded in the plane-wave basis set up to the cut-off energy of 50 Ry. The augmentation charges of the ultrasoft pseudopotentials were also expanded in plane waves up to the cut-off energy of 320 Ry. The self-consistent-field (SCF) loops were iterated until the total energies of consecutive SCF steps was less than 1.5 × 10−6 eV/supercell, where the forces were calculated. The structure optimization was performed until the absolute values of all the components of all forces were smaller than 3 × 10−3 eV/Å.

30.

We assumed that the upper and lower terrace C atoms in the α layers were mono-hydrogenated. In panel (a), one H atom was attached to the C atom at xk+1. In panel (b), we assumed that two H atoms were attached to the C atom at xk+2, whereas the one below the vacancy at xk+1 had an unpaired dangling bond.

31.
P.
Deák
,
B.
Aradi
,
M.
Kaviani
,
Th.
Frauenheim
, and
A.
Gali
,
Phys. Rev. B
89
,
075203
(
2014
).
32.
T.
Karin
,
S.
Dunham
, and
K.-M.
Fu
,
Appl. Phys. Lett.
105
,
053106
(
2014
).
33.
K.
Momma
and
F.
Izumi
,
J. Appl. Crystallogr.
44
,
1272
(
2011
).