Loughborough University

Institutional Repository

Ab initio study of point

defects in magnesium oxide
This item was submitted to Loughborough University's Institutional Repository
by the/an author.

Citation: GILBERT, C.A. ... et al., 2007. Ab initio study of point defects
in magnesium oxide. Physical Review B (Condensed Matter and Materials
Physics), 76 (18), 10pp.

Additional Information:

• This article was published in the journal, Physical Review B

(Condensed Matter and Materials Physics) [ c The Ameri-
can Physical Society]. The denitive version is available at:
http://dx.doi.org/10.1103/PhysRevB.76.184103

Metadata Record: https://dspace.lboro.ac.uk/2134/11635

Version: Published

Publisher: c The American Physical Society

Please cite the published version.

 This item was submitted to Loughborough’s Institutional Repository
(https://dspace.lboro.ac.uk/) by the author and is made available under the
following Creative Commons Licence conditions.

For the full text of this licence, please go to:

http://creativecommons.org/licenses/by-nc-nd/2.5/
 PHYSICAL REVIEW B 76, 184103 共2007兲

Ab initio study of point defects in magnesium oxide

C. A. Gilbert, S. D. Kenny,* R. Smith, and E. Sanville

Department of Mathematical Sciences, Loughborough University, Loughborough LE11 3TU, United Kingdom
共Received 4 May 2007; revised manuscript received 25 July 2007; published 6 November 2007兲
Energetics of a variety of point defects in MgO have been considered from an ab initio perspective using
density functional theory. The considered defects are isolated Schottky and Frenkel defects and interstitial
pairs, along with a number of Schottky defects and di-interstitials. Comparisons were made between the
density functional theory results and results obtained from empirical potential simulations and these generally
showed good agreement. Both methodologies predicted the first nearest neighbor Schottky defects to be the
most energetically favorable of the considered Schottky defects and that the first, second, and fifth nearest
neighbor di-interstitials were of similar energy and were favored over the other di-interstitial configurations.
Relaxed structures of the defects were analyzed, which showed that empirical potential simulations were
accurately predicting the displacements of atoms surrounding di-interstitials, but were overestimating O atom
displacement for Schottky defects. Transition barriers were computed for the defects using the nudged elastic
band method. Vacancies and Schottky defects were found to have relatively high energy barriers, the majority
of which were over 2 eV, in agreement with conclusions reached using empirical potentials. The lowest
barriers for di-interstitial transitions were found to be for migration into a first nearest neighbor configuration.
Charges were calculated using a Bader analysis and this found negligible charge transfer during the defect
transitions and only small changes in the charges on atoms surrounding defects, indicating why fixed charge
models work as well as they do.

DOI: 10.1103/PhysRevB.76.184103 PACS number共s兲: 61.72.⫺y, 31.15.Ew, 61.82.Ms

I. INTRODUCTION ergetics of a Schottky defect in MgO have already been stud-

ied using a Hartree-Fock cluster embedded within a crystal
In recent years, computer simulations that predict the be- by Grimes et al.,7 and the energy was found to be 8.2 eV.
havior of materials are becoming more and more sophisti- Isolated Schottky defects were studied using density func-
cated and therefore useful. Modeling the effects of radiation tional theory 共DFT兲 by De Vita et al.8,9 Their result using
damage in oxide materials is one such area of study that has quite a small supercell containing 32 atoms shows a predic-
been receiving a lot of attention, due to the fact that there is tion of a formation energy lying in the region between 6 and
an increasing interest in storing radioactive waste safely for 7 eV. A more recent DFT calculation was performed by Alfe
long periods of time. It is therefore important to be able to and Gillan10 in which they found the formation energy of an
predict accurately the behavior of the affected materials over isolated Schottky defect to be 6.99 eV.
long time scales. Materials that have come under scrutiny in In this paper, these isolated defects have also been stud-
the past are ceramics based on the MgO-Al2O3 system, the ied, but using a much larger supercell allowing the calcula-
simplest of which is magnesium oxide 共MgO兲. tion of accurate formation energies. We have also considered
MgO has a relatively simple crystal structure, known as not only isolated defects but also interacting defects in the
the NaCl 共or rocksalt兲 structure, and this is one of the reasons form of Schottky defects and di-interstitials. Furthermore, we
that make this oxide such an appealing material to consider. have studied a number of energy barriers involving the mi-
Its simplicity and well understood properties have resulted in gration and transformation of these defects.
MgO being utilized as a test material for both experimental
analysis and computer simulation. There are well established II. METHODOLOGY
empirical potentials existing for MgO and these have given
rise to detailed molecular dynamics studies1–6 of various The main body of work carried out in this paper has been
properties of the material. performed using the ab initio code PLATO 共package for linear
Simulations that utilize empirical potentials are extremely combination of atomic-type orbitals兲.11 This code utilizes
useful as they can be used to describe the defects in the density functional theory12,13 using localized, numerical or-
material that would result from exposure to radiation. It was bitals. Pseudopotentials for magnesium and oxygen were
shown previously2,3 that the main types of defect present taken from the work of Hartwigsen et al.14 and it is important
after a collision cascade in MgO were Frenkel defects, to note that a semicore pseudopotential was used for Mg.
Schottky defects, and di-interstitials. It is due to this obser- The local density approximation15 was used to describe ex-
vation that the study in this paper considers these specific change correlation within the system as this has been shown
defects. to give good agreement with experiment for MgO 共Refs. 8
The empirical potential model itself, although successful, and 9兲 and with higher order methods.10 Our basis sets con-
is still not necessarily creating the most accurate picture of tained 23 basis functions for Mg and 14 basis functions for
the defect energetics. It is therefore of interest to study these O. The required numerical integrations were carried out on
defects from a more accurate ab initio perspective. The en- an atom centered grid consisting of 40 points in the radial

1098-0121/2007/76共18兲/184103共10兲 184103-1 ©2007 The American Physical Society

GILBERT et al. PHYSICAL REVIEW B 76, 184103 共2007兲

direction and 35 points in the angular direction.

Tests were carried out upon the convergence of the equi- (a)
librium lattice constant and bulk modulus of MgO with a (c)
(e)
variety of Monkhorst-Pack16 k-point mesh sizes. A 6 ⫻ 6
⫻ 6 mesh was settled upon 共for the eight atom cubic unit
cell兲 as this provided sufficiently converged values for both
the bulk parameters and the total energy of the eight atom (b)
cell. For larger cells, equivalent meshes were used. The re- (d)
sulting lattice constant was predicted to be 4.17 Å and the
bulk modulus was 1.65 Mbar. De Vita et al.8,9 computed an
identical prediction for the lattice constant with a slightly
lower bulk modulus of 1.54 Mbar using a plane wave basis
set, indicating that our localized orbital basis set is perform- FIG. 1. The initial configurations of atoms for the 共a兲 first, 共b兲
ing well for bulk systems. These values compare well with second, 共c兲 third, 共d兲 fourth, and 共e兲 fifth nearest neighbor di-
the experimental parameters, which are 4.21 Å 共Ref. 17兲 for interstitials. The black circles represent the Mg atoms, while the
the lattice constant and 1.55 Mbar 共Ref. 18兲 for the bulk white circles represent the O atoms.
modulus of MgO.
For the isolated defects, two additional electrons were Periodic boundary conditions were utilized in the DFT
added to the system in the case of a Mg vacancy, and two calculations and so a supercell approach was taken. A rela-
were removed for an O vacancy. Equivalent electron addi- tively small number of atoms were considered and these
tions and subtractions were made for the isolated interstitials. were effectively repeated throughout space due to the condi-
This is to take into account the ionic nature of the material. tions imposed by the boundaries. When doing this, since any
Testing was performed by varying this integer number of defects placed in the supercell are also repeated, it was nec-
electrons between 0 and 3, and it was shown that the two essary to create a supercell of suitable size and shape such
electron insertion or removal produces the minimum forma- that the interaction between the defect and its periodic re-
tion energies. peats was minimized. We have performed a series of calcu-
We have also used an empirical potential approach to cal- lations for each defect in which we increased this separation
culate the defect energetics. A model without shells, as this is between the defect and the repeats, thus creating a set of
typically used in dynamical simulations, was implemented formation energies for each defect.
and the potential used is given in Eq. 共1兲. This is a standard In order to optimize the supercells so that the minimum
Buckingham potential. The parameters A, C, and ␳ are from number of atoms is considered for a specific separation dis-
Lewis and Catlow.19 tance, thus increasing computational efficiency, a methodol-
ogy was implemented that allowed a variety of shapes of

冉 冊
␾共r兲 = A exp −
r
␳
C
− 6 + V共r兲.
r
共1兲
supercell to be considered. The computational tool CRYSTAL
共Ref. 20兲 also implements the method. This method utilizes a
transformation matrix which repeats the two atom MgO
The most common defects found in empirical potential primitive cell to create a larger supercell of the required size
simulations of collision cascades in MgO are di-interstitials and shape. This is carried out by first considering a transfor-
and Schottky defects. Thus, we have calculated the formation mation matrix T, which contains integers tij for i , j
energies of these as a function of their separation distances. = 1 , . . . , 3, and this is used to construct three supercell vec-
It was also of interest to calculate the formation energies of tors 共s1 , s2 , s3兲 from the three primitive cell vectors
isolated Schottky and Frenkel defects and interstitial pairs. 共p1 , p2 , p3兲 as shown in Eq. 共2兲,
The isolated defects allow us to assess when the defect pairs s1 = t11p1 + t12p2 + t13p3 ,
are sufficiently separated such that there is no longer any
interaction between the two defects comprising the pair.
For the defect pairs, we have considered three different (a) (c)
Schottky defects and five different di-interstitial configura-
tions, all defined by the separation between the two defects
in the pair. Schottky defects are formed with the removal of
one Mg atom and one O atom from the system, and in a (b)
similar way the di-interstitials consist of adding an additional
Mg and an additional O to the system. The initial configura-
tions of these defects are shown in Figs. 1 and 2. It is impor-
tant to note that the distance between the vacancies in a
second nearest neighbor Schottky defect configuration is ini- FIG. 2. The initial configurations of the atoms for the 共a兲 first,
tially the same as the distance between the interstitials in a 共b兲 second, and 共c兲 third nearest neighbor Schottky defects. The
third nearest neighbor di-interstitial configuration, and the circles represent atoms, while the squares represent vacancies.
same comparison is true for a third nearest neighbor Black squares are Mg vacancies and white squares are O vacancies.
Schottky defect and a fifth nearest neighbor di-interstitial. This notation is used in all subsequent figures.

184103-2
AB INITIO STUDY OF POINT DEFECTS IN… PHYSICAL REVIEW B 76, 184103 共2007兲

TABLE I. The sizes of the supercells used in the simulations of the defects using DFT. The minimum
separation distance, in Å, between the periodic repeats of the defects is given in square brackets alongside the
number of atoms.

Number of atoms 关Separation distance 共Å兲兴

Defect共s兲 Cell 1 Cell 2 Cell 3 Cell 4

1NN di-interstitial/ 32 关5.90兴 48 关7.52兴 90 关8.85兴 128 关10.43兴

1NN Schottky defect
2NN di-interstitial 32 关5.90兴 58 关7.80兴 84 关8.85兴 132关10.22兴
3NN di-interstitial/ 32 关5.90兴 56 关7.52兴 90 关8.85兴 144 关10.43兴
2NN Schottky defect
4NN di-interstitial 32 关5.90兴 54 关6.60兴 84 关8.85兴 128 关9.33兴
5NN di-interstitial/ 40 关5.90兴 64 关7.52兴 104 关8.85兴 160 关10.43兴
3NN Schottky defect
Isolated interstitial/vacancy 8 关4.17兴 64 关8.34兴 128 关11.65兴 180 关12.86兴

s2 = t21p1 + t22p2 + t23p3 , two local minima. Forces were computed at each step, and
atoms were positioned in new locations, until the maximum
s3 = t31p1 + t32p2 + t33p3 . 共2兲 force was found to be below 0.07 eV/ Å. The maximum en-
ergy of this minimum energy path 共the saddle point兲 was
Given a variety of supercells containing the same number then found by fitting a cubic spline through the obtained data
of atoms, the most suitable one is that which maximizes the points, clamped in such a way to force the gradient to be zero
distance between a defect placed within the cell and with its at the local minima that represent the relaxed structures. An
periodic repeat generated through application of periodic example of this is shown in Fig. 3.
boundary conditions. The results of this optimization are Charges within the defective systems were computed us-
shown in Table I. From this, it can be seen that we have ing a Bader charge analysis with Bader volumes allocated
maintained approximately the same separation distances for using the methodology found in Ref. 22.
each of the different defect types. Since some of the defects
have a larger separation, it is unavoidable that larger super-
cells, containing more atoms, are needed to maintain the III. RESULTS
same separation distance between the defects and their re- A. Defect formation energies and geometrical structure
peats.
The formation energy of a defect, E form, which is the en- DFT calculations of the formation energies of isolated
ergy required to form the defect in the perfect crystal struc- Schottky and Frenkel defects and interstitial pairs for the cell
ture, is calculated using Eq. 共3兲 for the case of a pair of sizes presented in Table I give results as shown in Table II.
isolated defects. En共i , j兲 is the total cell energy for a supercell From the largest cell size for these defects, the formation
containing n lattice sites with i defects of one species and j energy of a pair of isolated vacancies appears to converge to
defects 共of the same type兲 of the other. In the parentheses, the 6.0 eV, whereas for the pair of isolated interstitials, the for-
⫾ becomes positive for the case of interstitial defects and mation energy is 16.6 eV. These results differ by less than
becomes negative for the case of vacancies. En共0 , 0兲 is the
perfect lattice energy.

冉 冊
Energy above minimum (eV)

n±1 1.0
E form = En共1,0兲 + En共0,1兲 − 2 En共0,0兲. 共3兲
n 0.75
A simplified version of this 关Eq. 共4兲兴 is used for calculating
0.5
the formation energy of di-interstitials and Schottky defects,
where now both interstitials or vacancies are placed within 0.25
the same supercell and so can interact.

E form = En共1,1兲 − 冉 冊
n±2
n
En共0,0兲. 共4兲
0.0
0 2 4 6
Intermediate State Number
8

Transition barriers have been calculated by using the FIG. 3. 共Color online兲 This graph shows how a clamped cubic
nudged elastic band method with a modified tangent spline is fitted to a series of calculated energies for the transition of
approach.21 We utilized a spring constant of k = 0.4 eV/ Å2 to an isolated interstitial in a 54 atom system. The black dots represent
compute the states of the system at 10–12 points 共depending the intermediate states between the initial and final system
on the system size兲 along the minimum energy path between configurations.

184103-3
GILBERT et al. PHYSICAL REVIEW B 76, 184103 共2007兲

TABLE II. The formation energies from the DFT calculations 7

for the pairs of isolated defects as a function of the system size. 1NN
2NN
Formation energy 共eV兲 3NN
6
Frenkel defect
Number of Schottky Isolated
atoms defect interstitial pair Magnesium Oxygen
5
8 7.23 20.04 13.55 13.72
64 5.65 16.35 10.07 11.93
128 5.88 16.74 10.33 12.31 4
180 5.97 16.58 10.35 12.17 5 7 9 11
(a) Repeated Image Separation Distance (Å)

3 eV from those formation energies predicted by an empiri- 19

cal potential simulation of a cell containing 64 000 atoms, 1NN
which utilizes the potential stated in Eq. 共1兲. The classical 2NN
3NN
potential predicts formation energies of 8.8 eV for the iso- 17 4NN
lated vacancies and 18.4 eV for isolated interstitials, both of 5NN
which are overestimates of the energies calculated using 15
DFT. The isolated Mg Frenkel defect has a formation energy
of approximately 10.3 eV compared with an empirical poten-
tial estimate of 14.1 eV, and for the O Frenkel defect, DFT 13
computes a formation energy of 12.2 eV and the empirical
potential simulation predicts 13.6 eV. The energies of the 11
isolated Frenkel defects have also been calculated using the
5 7 9 11
empirical potentials with shells. This gives a formation en- (b) Repeated Image Separation Distance (Å)
ergy of 12.6 eV for the Mg Frenkel defect and 12.3 eV for
the O Frenkel defect.
FIG. 4. Graphs that show the convergence of the formation en-
In addition to calculating the formation energies for these ergies as larger supercells are utilized. The energies in 共a兲 are for
isolated defects, the energetics of di-interstitials and Schottky Schottky defects and in 共b兲 are for di-interstitials.
defects were studied, and in this case interaction is allowed
between the two defects comprising the pair. Figure 4 shows the Schottky defects by 1 – 3 eV, but underestimate the for-
how these formation energies converge as the supercells are mation energies of the first, third, fourth, and fifth nearest
increased in size to negate any interactions between the de- neighbor di-interstitials by less than 1 eV.
fect pairs and the periodic images of the defect pairs 共which Since the empirical potential calculations allow for many
are introduced by periodic boundary conditions兲. The defect more atoms to be considered in the simulations, we contin-
energies are converged to within ±0.2 eV for sufficiently ued the formation energy calculations for a variety of di-
large supercells. interstitials having larger separation distances. This allowed
Table III shows the formation energy for each of these us to gain an insight into the point at which the isolated
defects taken from the computation for the largest cell for defect energy is reached. Figure 5 shows the formation en-
that particular defect. These are compared to the empirical ergy of the di-interstitials as a function of the separation
potential simulation predictions for the formation energies distance.
for each of the different defects. This graph shows that for a separation distance of larger
From these comparisons, it can be seen that the empirical than 7 Å, the formation energy lies close to 18 eV. This is in
potential simulations overestimate the formation energies for agreement with the converged value of 18.4 eV. For the larg-

TABLE III. The defect formation energies for the nearest neighbor Schottky defects and di-interstitials for
both the DFT and the empirical potential calculations as a function of the defect separation.

DFT formation energy 共eV兲 Emp. pot. formation energy 共eV兲

Defect
separation Schottky defect Di-interstitial Schottky defect Di-interstitial

1NN 4.6 12.4 5.7 11.6

2NN 5.9 12.1 7.7 12.4
3NN 5.7 15.4 7.2 15.3
4NN 16.6 16.5
5NN 12.1 11.5

184103-4
AB INITIO STUDY OF POINT DEFECTS IN… PHYSICAL REVIEW B 76, 184103 共2007兲

20 a) b) C
A
18
B D
16

14

12 c) d)
F H
10 E G
0 3 6 9 12 15
Distance between Interstitials (Å)

FIG. 5. Graph showing the formation energy of di-interstitials,

calculated using empirical potentials, as a function of the separation
between the two interstitial atoms. FIG. 7. The direction of displacement of the atoms surrounding
the isolated defects: 共a兲 Mg vacancy, 共b兲 O vacancy, 共c兲 Mg inter-
est considered separation of interstitials 共13– 14 Å兲, the em- stitial, and 共d兲 O interstitial. The magnitudes of the displacements of
pirical potential predicts a formation energy of 17.9 eV. Also selected atoms are shown in Table IV.
of interest is the fact that in addition to the first, second, and
fifth nearest neighbor di-interstitials having a relatively low atoms are repelled. The opposite is true for the isolated in-
formation energy of 11.5– 12.5 eV, another favorable con- terstitials. Directions of displacement are the same for both
figuration is predicted when the di-interstitials are in an ab initio and classical potential calculations. Isolated inter-
eighth nearest neighbor configuration 共separated by a dis- stitials cause atoms of the same species to move a greater
tance of 6.32 Å兲, in a 具221典 direction. This has a formation distance than atoms of the other species. The DFT results
energy of 11.5 eV. There exists another configuration with show that the Mg displacement for the Mg interstitial is
this same separation distance, but with the interstitials lying 0.45 Å compared to 0.023 Å for the O atoms and the O
in a 具100典 direction, and this is predicted to have a much movement for the O interstitial is 0.48 Å compared to
higher formation energy of 17.3 eV. 0.017 Å for the Mg atoms. For the displacements of the at-
We studied the convergence of formation energies of oms surrounding the isolated vacancies, the empirical poten-
Schottky defects in a similar manner, and the results of this tial simulation overestimates the magnitudes of the displace-
analysis are presented in Fig. 6. Only a small deviation from ments by about 50%. The agreement between DFT and the
the converged value of 8.8 eV occurs for separation dis- empirical potentials is much stronger in the case of isolated
tances greater than 6 Å. interstitials, with less than 0.05 Å difference between the
In order to study these defects further, the relaxation of magnitudes of displacement predicted by the two computa-
the structures and the defects was also analyzed. For the tional methods. An earlier comparison between empirical po-
isolated defects, the directions of displacement of the atoms tentials and ab initio calculations 共GAMESS兲 was also given in
during the relaxation are shown in Fig. 7. Ref. 25. Although the relaxation distances calculated for at-
Table IV shows the magnitude of these displacements. It oms around the vacancies were similar, using the empirical
can be seen in this table, and in Fig. 7, that the surrounding
Mg atoms are attracted toward a Mg vacancy, whereas the O TABLE IV. The magnitude of the displacements of the atoms
during relaxation by the DFT and the empirical potential methods.
10 The labels refer to the letters in Fig. 7.

DFT Empirical potential

displacements displacements
8
Defect Label 共Å兲 共Å兲

Mg vacancy A 0.08 0.14

6 B 0.12 0.21

O vacancy C 0.07 0.11

4 D 0.16 0.21
0 3 6 9 12 Mg interstitial E 0.023 −0.023
Distance between Vacancies (Å)
F 0.45 0.49

FIG. 6. Graph showing the formation energy of Schottky de- O interstitial G 0.017 −0.024
fects, calculated using empirical potentials, as a function of the H 0.48 0.47
separation between the two vacant lattice sites.

184103-5
GILBERT et al. PHYSICAL REVIEW B 76, 184103 共2007兲

E A C E
F
A B D
C D
C F D
E
B F (a)
(a)
F H
G
A J F
D E
B I
C
B A
B
(b) (b)

G
FIG. 8. The direction of displacement, during relaxation, of at-
E
oms surrounding the 共a兲 first and 共b兲 second nearest neighbor F
B
Schottky defects. Magnitudes of these displacements are given in
Table V.
I F G
D F
potentials, to those given here, the relaxation distances cal-
H A
culated using GAMESS were about double those from our ab
initio results. C H I
For the interacting defect pairs, a similar analysis was H
carried out. For the Schottky defects, the directions of dis- (c)
placements of atoms surrounding the first and second nearest
neighbor Schottky defects are shown in Fig. 8, with the mag- A B C D
nitudes for these displacements given in Table V. E F H
It can be seen in Table V that the displacements surround- I G H
F G J
E
ing the first nearest neighbor Schottky defect follow a trend A
similar to that of the isolated vacancies. Here, however, the B C D
magnitude of the oxygen displacement is predicted to be (d)
greater by the empirical potential simulations than the mag-
nitudes calculated using DFT. This is also observed in the FIG. 9. The direction of displacement, from perfect lattice sites,
of atoms surrounding the 共a兲 first, 共b兲 second, 共c兲 third, and 共d兲
TABLE V. The magnitude of the displacements of the atoms fourth nearest neighbor di-interstitials. Magnitudes of these dis-
during relaxation of the first and second nearest neighbor Schottky placements are given in Table VI.
defects by the DFT and the empirical potential methods. The labels
refer to the letters for each on the atoms in Fig. 8. second nearest neighbor Schottky defect, with the empirical
potential simulation overestimating the displacement of the
DFT Empirical potential O atoms by 54% but only overestimating the movement of
displacements displacements surrounding Mg atoms by 12%. The third nearest neighbor
Defect Label 共Å兲 共Å兲 Schottky defect follows a similar trend, in that the O atoms
are predicted to move 35%–57% further than the DFT simu-
1NN Schottky defect A 0.12 0.19 lation shows. The displacements of the Mg atoms are again
B 0.17 0.18 in good agreement, with only a 9%–11% overestimate by the
C 0.12 0.19 empirical potentials.
D 0.16 0.18 Information about the atomic displacement of the di-
E 0.09 0.10 interstitial configurations is shown in Fig. 9 and Table VI in
the same way as for the Schottky defects. The first and fourth
F 0.09 0.10
nearest neighbor di-interstitials are aligned in a 具100典 direc-
2NN Schottky defect A 0.13 0.20 tion and the interstitial atoms remain an equidistance from
B 0.17 0.19 the atoms at each of the cell sides but are attracted to each
other. The second nearest neighbor and third nearest neigh-

184103-6
AB INITIO STUDY OF POINT DEFECTS IN… PHYSICAL REVIEW B 76, 184103 共2007兲

TABLE VI. The magnitude of the displacements of the atoms method indicates that for the second nearest neighbor case,
during relaxation of the di-interstitial defects by DFT and empirical the Mg interstitial atom moves very little and it is only the O
potential methods. The labels refer to letters in Fig. 9 and agreeing atom that should form a split interstitial with atom labeled D
displacements of less than 0.3 Å are not shown. in the second nearest neighbor structure shown in Fig. 9. The
empirical potential simulations predict that the Mg interstitial
DFT Empirical potential atom will form a second split interstitial with the Mg atom
displacements displacements labeled E.
Defect Label 共Å兲 共Å兲

1NN di-interstitial C 0.42 0.43 B. Transition energy barriers

D 0.43 0.42 An important factor when considering point defects in a
E 0.58 0.54 material is the mobility of such a defect. It is of interest to
F 0.57 0.54 know if it is likely that a defect can move freely within the
material or whether it becomes trapped. Transition barrier
2NN di-interstitial A 0.56 0.19
energies are used to determine the likelihood of the system of
B 0.56 0.41 atoms moving from one given configuration to another. This
C 0.40 0.24 can therefore be used to analyze potential defect configura-
D 1.72 0.80 tion changes, or simply a movement of the defect through the
E 0.44 1.04 material. We have studied a variety of transition barriers us-
F 0.16 0.41 ing the nudged elastic band method.
G 0.10 0.24 Previous empirical potential calculations2,5 indicate that
interstitials are far more mobile than vacancies in MgO. We
I 0.13 0.58
have looked at the possibility of a simple isolated vacancy
J 1.23 0.57
migration and also considered the possibility of each of the
3NN di-interstitial A 0.88 0.89 three different Schottky defect configurations transforming
B 0.45 0.39 into each other. We found that the isolated Mg vacancy has a
migration barrier of 2.20 eV, whereas the O vacancy has a
C 1.37 1.38
barrier of 2.31 eV. These values agree with those calculated
D 0.87 0.88 using DFT by De Vita et al.8 and using the semiempirical
G 0.37 0.36 approach by Kotomin and Popov23 who utilized the method
4NN di-interstitial A 0.32 0.37
of intermediate neglect of the differential overlap. Results
from these sources were 2.39 and 2.43 eV, respectively, for
D 0.34 0.35
the isolated Mg vacancy migration energy and 2.48 and
F 0.59 0.63 2.50 eV, respectively, for the isolated O vacancy migration
G 0.58 0.62 energy.
I 0.30 0.29 Table VII shows transition energies for Schottky defect
J 0.29 0.31 configuration changes. The six possible transitions between
first, second, and third nearest neighbor Schottky defect con-
figurations all have energy barriers greater than 1.3 eV. The
bor configuration di-interstitials form split interstitials. In two lowest barriers of 1.30 and 1.35 eV both correspond to
fact, the third nearest neighbor defect forms a crowdion transformations to the first nearest neighbor Schottky defect
along 具111典 consisting of four interstitial atoms along with state 共from third and second nearest neighbor configurations,
two vacant sites. respectively兲.
With the exception of the second nearest neighbor di- In considering the motion of an isolated interstitial atom,
interstitial, there is good agreement between the empirical both Mg and O separately, we have studied a number of
potential simulation predictions and the DFT results. The dif- different mechanisms for the displacement of the isolated
ferences in predicted magnitudes of displacement vary be- interstitial atom. The three possible directions considered
tween methods by no more than 0.06 Å, with the majority of were the 具100典-, 具110典-, and 具111典-type directions. Motion is
displacements agreeing to within 0.02 Å. The ab initio possible through either a direct hop of an interstitial or an

TABLE VII. The transition energies of configuration changes of the Schottky defects studied.

Final configuration

Initial configuration 1NN Schottky defect 2NN Schottky defect 3NN Schottky defect

1NN Schottky defect 2.61 2.35

2NN Schottky defect 1.35 2.23
3NN Schottky defect 1.30 2.23

184103-7
GILBERT et al. PHYSICAL REVIEW B 76, 184103 共2007兲

TABLE VIII. The transition energies of the displacement of isolated magnesium and oxygen isolated
interstitial atoms.

54 atom barrier 共eV兲 128 atom barrier 共eV兲

Transition mechanism Magnesium Oxygen Magnesium Oxygen

具100典 direct hop 1.21 1.26 0.99 0.99

具110典 direct hop 2.59 2.71
具100典 atom exchange 2.78 3.00
具110典 atom exchange 1.17 1.25 0.94 0.96
具111典 atom exchange 0.83 0.59 0.71 0.44

exchange with a nearest neighbor atom of the same species, di-interstitial that moves. If a mechanism does not appear in
with the initial interstitial moving to occupy the lattice site of the table 共e.g., 2NN-3NN exchange兲, it indicates that the
the neighboring atom, while the neighboring atom moves to transition occurs by passing through some intermediate
become the interstitial in its new location. We first consid- phase first.
ered simplifying the calculations at the cost of accuracy in From these results, it is observed that from either the third
order to speed up the computations. Considering a 54 atom or fourth nearest neighbor position, it is likely that the con-
cell allowed us to see which of the transitions have the low- figuration will change such that the interstitial atoms become
est energy. Following this, we then considered larger, 128 first nearest neighbors. The barriers for these transitions are
atom cells, in order to calculate these lowest energy barriers 0.30 and 0.10 eV, respectively. Recalling that the first and
to a much higher degree of accuracy. second nearest neighbor di-interstitials were observed to
Table VIII shows the resulting energy barriers for isolated have similar formation energies, it is important to note here
interstitial transitions, calculated using 54 atom cells. These that the 128 atom cell used here only produces an estimate to
values show that the 具110典 direct hop and 具100典 atom ex- the actual formation energies shown in Table III. For ex-
change are extremely unlikely, having energy barriers be- ample, the 128 atom cell predicts identical formation ener-
tween 2.5 and 3 eV. Due to the fact that the 具100典 direct hop gies for first and second nearest neighbor di-interstitials. The
and 具110典, and 具111典 atom exchanges have much lower en- transition energy for a first to second nearest neighbor tran-
ergy barriers ranging from 0.83 to no more than 1.26 eV, we sition, 0.74 eV, is therefore identical to that of a second to
studied these transitions using the larger 128 atom cell. The first nearest neighbor transition due to this approximation.
barrier energies computed using a 128 atom supercell were
found to be between 0.1 and 0.3 eV smaller than those found
using the 54 atom cells. The results are shown in Table VIII TABLE IX. The transition energies of configuration changes of
the di-interstitials studied.
and indicate that the 具111典 atom exchange is the most likely
transition to occur. Previous empirical potential simulations24
Barrier
have estimated a 0.32 eV energy barrier for the isolated Mg
From To Mechanism 共eV兲
interstitial and 0.40 eV for the O interstitial, incorrectly in-
dicating that motion of the Mg atom is the more favorable 1NN 2NN Exchange 共Mg兲 0.74
transition. 1NN 2NN Exchange 共O兲 0.51
It is also of interest to study not only the displacement of
1NN 2NN Direct 共O兲 1.43
interstitial atoms when they are isolated but also the dis-
1NN 3NN Exchange 共Mg兲 3.24
placement when the atom constitutes one-half of a di-
interstitial. By considering displacements in this fashion, it 1NN 4NN Direct 共O兲 4.36
allowed us to understand how di-interstitials move through 1NN 5NN Exchange 共Mg兲 1.48
the lattice. We computed the transition barriers for direct
2NN 1NN Exchange 共Mg兲 0.74
hops of the oxygen atom which would convert one di-
interstitial configuration to another. Since we found the ex- 2NN 1NN Exchange 共O兲 0.51
change mechanism to be favorable, we also computed the 2NN 1NN Direct 共O兲 1.43
barriers for various exchanges of magnesium atoms. Magne- 2NN 3NN Direct 共O兲 6.05
sium atoms were considered in the exchanges because in the
3NN 1NN Exchange 共Mg兲 0.30
case of an isolated interstitial, for two of the three exchanges
共in the 具100典 and 具110典 directions兲, the Mg atom exchange 3NN 2NN Direct 共O兲 3.11
had a lower energy barrier than an O atom exchange. 3NN 5NN Exchange 共Mg兲 4.04
Table IX shows the transition barriers that were calculated
4NN 1NN Direct 共O兲 0.10
for various di-interstitial configuration changes. The table
presents the configuration change, the mechanism for the 5NN 1NN Exchange 共Mg兲 1.48
transition, and the transition energy. These energies are ap- 5NN 3NN Exchange 共Mg兲 6.98
proximately the same whether it is the Mg or O atom in the

184103-8
AB INITIO STUDY OF POINT DEFECTS IN… PHYSICAL REVIEW B 76, 184103 共2007兲

The only exchange mechanism in a 具111典-type direction ergies calculated using DFT. The relaxed structures for these
was that of a second to first nearest neighbor exchange. Due defects show that the empirical potential simulations overes-
to the differences in migration energies of the isolated O and timate the displacement of the surrounding oxygen atoms
Mg interstitials, we also considered an O atom exchange during relaxation. The most energetically favorable Schottky
mechanism for this transition and we found this to have a defect structure is the first nearest neighbor configuration,
transition barrier of 0.51 eV 共compared to the 0.74 eV for and this is apparent in both methodologies. The computation
the Mg atom exchange兲. This was the only transition where of transition barriers using DFT showed that Schottky de-
there was a significant difference between the O and Mg fects were unlikely to change configuration; thus, it is likely
movement. that first nearest neighbor vacancies would be formed and
Since the fifth nearest neighbor di-interstitial has a simi- would remain as such. The difference in the formation en-
larly low formation energy to the first and second nearest ergy means that it is highly likely that empirical potential
neighbor configurations, potential transitions between these calculations are underestimating the number of defects
states were investigated. The transition barrier between the formed.
fifth and first nearest neighbor was found to be 1.48 eV, Analysis of the di-interstitials showed that the formation
which is relatively high, and we did not observe a single energies generated from the empirical potential simulations
stage process between second and fifth nearest neighbors. differ from DFT results by less than 1 eV, with lower forma-
The barriers to move from a third nearest neighbor to the tion energies generally predicted by the empirical potential
much lower in energy second or fifth nearest neighbor state calculations. Unlike for the Schottky defects, the relaxed
are surprisingly high 共3.11 and 4.04 eV, respectively兲, with structures for these di-interstitial defects were in good agree-
the only remaining single stage process being that of a tran- ment between simulation methodologies, both in terms of the
sition to first nearest neighbor with a far lower barrier energy directions and magnitudes of the displacements. This would
of 0.30 eV. account for the closer agreement between empirical potential
The results in Table IX show the single stage processes. and DFT formation energies of the di-interstitials. The only
Remaining transitions were all found to be multiple stage case where this is not true is the second nearest neighbor
processes. Transitions from the fourth nearest neighbor con- di-interstitial, where the formation energy was predicted to
figurations are observed to pass through the first nearest be 0.3 eV higher than that calculated using DFT, and the
neighbor configuration 共with an energy barrier of 0.10 eV兲 structure varied from that of the DFT result. The empirical
before they potentially form either a second or third nearest potential simulation predicted that two split interstitials
neighbor di-interstitial, with energy barriers of 0.74 or
would form, whereas the DFT results show that only the
3.24 eV, respectively. Our calculations indicate that the fifth
oxygen atom forms a split interstitial. Both methodologies
nearest neighbor defect, both with a direct and an exchange
mechanism, had to pass through multiple configurations to predicted the first, second, and fifth nearest neighbor di-
reach the fourth nearest neighbor state, which is unlikely interstitials to be the most energetically favorable of the stud-
given the large energy barrier required to reach the fourth ied di-interstitial defects. Empirical potential simulations
nearest neighbor configuration. also predict that there exists another low energy state in an
An exchange mechanism was also observed in the empiri- eighth nearest neighbor configuration 共this state has not been
cal potential simulations for which a first nearest neighbor studied using ab initio methods兲.
di-interstitial would move through the bulk material. Here, The study of the transition barriers for the di-interstitials
both the magnesium and oxygen atoms would exchange with indicates that should a third or fourth nearest neighbor di-
atoms at the corner of their local cells and these displaced interstitial form, the very low energy barriers dictate that this
atoms would become the new interstitials. Empirical poten- would quickly form a first nearest neighbor di-interstitial in
tial simulations predict a barrier energy of 0.75 eV 共Ref. 24兲 favor of the other defects studied. There is then a small
for this displacement, whereas the ab initio simulations com- chance that the configuration would switch between that of
puted this energy to be 0.85 eV, which is in good agreement. the first and second nearest neighbor di-interstitials. The bar-
A Bader analysis was performed on the 128 atom cell, and rier to move from a fifth nearest neighbor configuration ap-
this showed that in a perfect MgO lattice, the charge on a Mg pears to be quite high, and so it is expected that if a di-
atom is 1.83 and on an O atom is −1.83. The charges on interstitial of this separation were to form, it would be
atoms in the systems containing di-interstitials varied by no unlikely to change configurations.
more than ±0.1 with charges on the actual interstitial atoms Isolated vacancies were found to have relatively high mi-
varying by this largest amount with the majority of the re- gration barriers of 2.20 and 2.31 eV for Mg and O vacancies,
maining atoms varying by no more that ±0.04. During the respectively, from the DFT calculations, thus to be practi-
second to first nearest neighbor direct transition, the change cally immobile once formed. The DFT calculations predict
in charge on the interstitial atom involved in the motion was much lower barriers for the isolated interstitials of 0.71 and
negligible. We also calculated the charges on atoms in the 0.44 eV for an isolated Mg atom and O atom, respectively.
systems containing single isolated vacancies. The results The most favorable diffusion mechanism for the isolated in-
showed that charge was distributed throughout the atoms in terstitials was via an exchange with another atom of the same
the system rather than localizing on the vacant site. species in the 具111典 direction. Previous empirical potential
simulations24 have calculated a 0.32 eV energy barrier for
IV. CONCLUSIONS
the isolated Mg interstitial and 0.40 eV for the O interstitial.
The Schottky defect formation energies predicted by em- Thus, the ab initio calculations predict that the motion on the
pirical potential simulations were 1 – 3 eV larger than the en- O sublattice will be much faster than that on the Mg sublat-

184103-9
GILBERT et al. PHYSICAL REVIEW B 76, 184103 共2007兲

tice, in contrast to the empirical potential results, which will better for the Frenkel defects, with the O Frenkel defect in
have a large effect on the way that defects are annealed after good agreement but the formation energy of the Mg Frenkel
the collision cascade. defect still differs by 1.5 eV. The transition barriers for iso-
The Bader charge analysis showed that the atoms that lated interstitials, which dominate the annealing of the de-
surround defects have charges very similar to those in the fects, are good for the migration of O interstitials but poor
bulk. Further, there is very little change in the charge during for the Mg interstitials. This has important implications for
transitions; therefore, a fixed charge model appears to be an both how the annealing takes place and the time scale over
appropriate approximation to implement when utilizing em- which annealing would occur.
pirical potentials.
The overall agreement between the two methodologies is
reasonable, which indicates that the empirical potentials in ACKNOWLEDGMENTS
the main give an acceptable description of MgO. The empiri-
cal potential employed performs well when describing One of us 共C.A.G.兲 would like to thank Loughborough
interstitial-type defects but far less well when studying va- University for providing funding. The calculations were per-
cancies, where the description is poor. This also introduces formed using the facilities of the Loughborough High Per-
errors into the calculations of the energies of the Frenkel formance Computing Centre. We would like to acknowledge
defects. The use of the shell model makes the agreement useful discussions with Blas Uberuaga.

*Corresponding author. FAX: ⫹44-15-09-223-969. 4899 共2000兲.

1 E. A. Kotomin, P. W. M. Jacobs, N. E. Christensen, T. Brudevoll, 12 P. Hohenberg and W. Kohn, Phys. Rev. 136, B864 共1964兲.
M. M. Kaklja, and A. L. Dopov, Defect Diffus. Forum 143, 13
W. Kohn and L. J. Sham, Phys. Rev. 140, A1133 共1965兲.
1231 共1997兲. 14 C. Hartwigsen, S. Goedecker, and J. Hutter, Phys. Rev. B 58,
2
B. P. Uberuaga, R. Smith, A. R. Cleave, F. Montalenti, G. Hen- 3641 共1998兲.
15 S. Goedecker, M. Teter, and J. Hutter, Phys. Rev. B 54, 1703
kelman, R. W. Grimes, A. F. Voter, and K. E. Sickafus, Phys.
Rev. Lett. 92, 115505 共2004兲. 共1996兲.
3 B. P. Uberuaga, R. Smith, A. R. Cleave, G. Henkelman, R. W. 16 H. J. Monkhorst and J. D. Pack, Phys. Rev. B 13, 5188 共1976兲.

Grimes, A. F. Voter, and K. E. Sickafus, Phys. Rev. B 71, 17

R. W. G. Wyckoff, Crystal Structures 共Wiley, New York and Lon-
104102 共2005兲. don, 1964兲.
4 G. Henkelman, B. P. Uberuaga, D. J. Harris, J. H. Harding, and N. 18 M. J. L. Sangster, G. Peckham, and D. H. Saunderson, J. Phys. C

L. Allan, Phys. Rev. B 72, 115437 共2005兲. 3, 1026 共1970兲.

5 B. P. Uberuaga, R. Smith, A. R. Cleave, R. W. Grimes, A. F. 19 G. V. Lewis and C. R. A. Catlow, J. Phys. C 18, 1149 共1985兲.
20
Voter, and K. E. Sickafus, Nucl. Instrum. Methods Phys. Res. B V. R. Saunders, R. Dovesi, C. Roetti, R. Orlando, C. M. Zicovich-
250, 12 共2006兲. Wilson, N. M. Harrison, K. Doll, B. Civalleri, I. Bush, and P.
6
S. Y. Li, Z. L. Liu, Y. G. Nan, and Z. R. Zhang, Chin. J. Chem. D’Arco, CRYSTAL2003 User’s Manual 共University of Torino,
Phys. 19, 315 共2006兲. Torino, 2003兲.
7 R. W. Grimes, C. R. A. Catlow, A. L. Shluger, R. Pandey, R. 21 G. Henkelman and H. Jónsson, J. Chem. Phys. 113, 9978 共2000兲.
22 E. Sanville, S. D. Kenny, R. Smith, and G. Henkelman, J. Com-
Baetzold, and A. H. Harker, Defects in materials, MRS Sympo-
sia Proceedings 209 共Materials Research Society, Pittsburgh, put. Chem. 28, 899 共2007兲.
23
1991兲, p. 257. E. A. Kotomin and A. I. Popov, Nucl. Instrum. Methods Phys.
8 A. De Vita, M. J. Gillan, J. S. Lin, M. C. Payne, I. Stich, and L.
Res. B 141, 1 共1998兲.
J. Clarke, Phys. Rev. Lett. 68, 3319 共1992兲. 24
B. P. Uberuaga, R. Smith, A. R. Cleave, G. Henkelman, R. W.
9 A. De Vita, M. J. Gillan, J. S. Lin, M. C. Payne, I. Stich, and L.
Grimes, A. F. Voter, and K. E. Sickafus, Nucl. Instrum. Methods
J. Clarke, Phys. Rev. B 46, 12964 共1992兲. Phys. Res. B 228, 260 共2005兲.
10 D. Alfe and M. J. Gillan, Phys. Rev. B 71, 220101共R兲 共2005兲. 25 R. W. Grimes, C. R. A. Catlow, and A. M. Stoneham, J. Phys.:
11
S. D. Kenny, A. P. Horsfield, and H. Fujitani, Phys. Rev. B 62, Condens. Matter 1, 7367 共1989兲.

184103-10