# Interaction between dislocation and vacancies in magnesium oxide: Insights from atomistic simulations and elasticity theory

Marie Landeiro dos Reis, Philippe Carrez, Patrick Cordier

► To cite this version:

Marie Landeiro dos Reis, Philippe Carrez, Patrick Cordier. Interaction between dislocation and vacancies in magnesium oxide: Insights from atomistic simulations and elasticity theory. Physical Review Materials, 2021, Physical Review Materials, 5 (6), pp.063602. 10.1103/physrevmaterials.5.063602. hal-03255168

HAL Id: hal-03255168

https://hal.univ-lille.fr/hal-03255168

Submitted on 9 Jun 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Interaction between dislocation and vacancies in magnesium oxide: Insights from atomistic simulations and elasticity theory

Marie Landeiro Dos Reis and Philippe Carrez

Univ. Lille, CNRS, INRAE, Centrale Lille, UMR 8207 UMET Unité Matériaux et Transformations, F-59000 Lille, France

Patrick Cordier

Univ. Lille, CNRS, INRAE, Centrale Lille, UMR 8207 UMET Unité Matériaux et Transformations, F-59000 Lille, France
and Institut Universitaire de France, 75005 Paris, France

<figure><img src="image_2.png" /></figure>

(Received 27 November 2020; revised 22 February 2021; accepted 30 March 2021; published 4 June 2021)

We determine the interaction between an 1/2⟨110⟩{⟨110⟩} edge dislocation and charged vacancies in MgO, using both molecular static simulations and the elasticity theory developed in the framework of the elastic dipole approach. In this study, the confrontation of these two methods highlights the specific role of the dislocation core structure on the interaction. We thus show that in MgO, the edge dislocation core, within a region across the glide planes that expands over several Burgers vector, is strongly attractive for vacancies, especially those of oxygen. However, the resulting pinning force on the dislocation remains weak and should not contribute to a significant hardening.

DOI: 10.1103/PhysRevMaterials.5.063602

## I. INTRODUCTION

Plastic properties of crystalline materials depend not only on the nature of the defects present in the crystal but also and more substantially on their mutual interactions [1]. In metals, the interactions between point defects and dislocations are well established and known to be responsible of different mechanical behaviors. For instance, at low or moderate temperature, there is a strong influence of point defects on strain hardening or softening [2–6] and on strain ageing [7]. At higher temperature, point-defect absorption or emission along the dislocation lines allow the dislocation climb mechanism [8–11] and can thus impact creep [12].

In ionic materials, the presence of different type of ions leads to a more complex variety of interactions between dislocations and charged point defects [13,14], a well-known example being the photolytic reaction in AgBr [15–17]. As in metals, dislocations are supposed to be strong point-defect sinks and induce short circuit pathways through the pipe diffusion process [18–21] which may lead for instance to the formation of jogs by vacancy absorption. However, the specificity of ionic crystals is that dislocation cores can be charged or can carry charges [13,14,22]. The consequence of this electrically charged core is that an external electric field can influence the dislocation motions and the plastic deformation [23] and can increase the electrical conductivity [24].

The understanding of these peculiar effects in ionic materials requires a precise description of the dislocation core, and a fine analysis of the interaction between a point defect and the dislocation. As an example, whereas edge dislocations are uncharged in typical rock salt crystals, the jogs or kinks resulting from either dislocation-dislocation interaction or from the vacancy adsorption carry charges [1,14,25].

In this work, we study the interaction between vacancies and a 1/2⟨110⟩{⟨110⟩} edge dislocation in magnesium oxide. Indeed, in compounds such as MgO, the creep mechanism and more generally the mechanical behavior is often interpreted through the interactions between vacancies, impurities and dislocations (see for instance Refs. [26–28] and references therein). In particular, the experimentally observed hardening or strain ageing is assumed to be controlled either by impurities [29,30] or by intrinsic vacancies [31]. As an explanation, one usually refers to either a strong pinning of dislocations by point defects [29] or to a reduction of the glide velocity for dislocations surrounded by an atmosphere of charged impurity ions or vacancies [26]. Here we focus on the interaction energies with a charged vacancy ($V''_{Mg}$ and $V''_{O}$, respectively the cationic and anionic vacancies) relying on both molecular static simulations and the elasticity theory [32–35] using the elastic dipole approach [36–44]. As shown in the following, we extract the vacancy elastic dipole tensor, defined as the negative derivative of the defect energy against the strain on the crystal lattice, from atomistic simulations. Doing so, we show that the charged vacancy can be fully described within the framework of the elastic dipole tensor which naturally captures the local deformation resulting from the electrostatic repulsion between the firsts neighbor ions of the vacancy. Our results highlight the necessity to account for the exact dislocation core structure to properly infer the nature of the interaction dislocation-vacancy. Relying on both the elasticity theory and our atomistic simulations, we finally

<footnote>* philippe.carrez@univ-lille.fr</footnote>

<footnote>† patrick.cordier@univ-lille.fr</footnote>
<table><caption>TABLE I. Properties of MgO as computed here with molecular statics (MS), compared to density-functional theory calculations (DFT) [46,51] and experimental measurements. $a_0$ is the lattice parameter, $b = a_0\sqrt{2}/2$ is the Burgers vector, $C_{11}$, $C_{12}$, and $C_{44}$ are the elastic constants and $B = (C_{11} + 2C_{12})/3$ is the bulk modulus.</caption><thead><tr><th></th><th>MS (DFT [46,51])</th><th>Experiments</th></tr></thead><tbody><tr><td>$a_0$ (Å)</td><td>4.218 (4.237)</td><td>4.212 [52]</td></tr><tr><td>$b$ (Å)</td><td>2.983 (2.996)</td><td>2.978 [52]</td></tr><tr><td>$C_{11}$ (GPa)</td><td>283 (279)</td><td>297.8±0.4 at 293 K [53]<br>295.9 at 300 K [54]<br>306.7 at 4.2 K [54]</td></tr><tr><td>$C_{12}$ (GPa)</td><td>138 (93)</td><td>97±0.5 at 293 K [53]<br>95.40 at 300 K [54]<br>93.71 at 4.2 K [54]</td></tr><tr><td>$C_{44}$ (GPa)</td><td>138 (146)</td><td>156.3±0.2 at 293 K [53]<br>153.9 at 300 K [54]<br>157.6 at 4.2 K [54]</td></tr><tr><td>$B$ (GPa)</td><td>186 (155)</td><td>163.9±0.6 at 293 K [53]</td></tr></tbody></table>

discuss the pinning of the dislocation line by vacancies at low temperatures.

## II. METHODS AND MODELS

In this study, all atomistic calculations are performed using classical molecular statics simulations. The ionic interactions are based on rigid-ion pair potentials. The short range interactions are modeled with a Buckingham form while the long range interactions are classically described through Coulombic terms. We use the parametrization proposed in Ref. [45] for MgO, relying on partial ionic charges of $\pm 1.7e$. Partly parameterized against *ab initio* calculations, the potential formulation is known to fairly reproduce a large set of MgO intrinsic properties (Table I): lattice parameters, elastic constants (although the rigid ion model satisfy the Cauchy relation) [45], point-defect diffusion [45], dislocation core properties [46] or the structure of low angle grain boundaries [47]. Practically, we use the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) [48] for all simulations. As implemented in LAMMPS, the Coulomb part is computed by means of a particle-particle particle-mesh method [49] and the atomic relaxation are ensured by a Fast Inertial Relaxation Engine algorithm [50] until a force criterion of 1 meV/Å.

### A. Vacancy elastic dipole

Within the framework of elasticity theory and the elastic dipole approach, point defects are described by the local deformation they induce on the crystal lattice [43,44]. The elastic dipole tensor of point defects, quoted $P$ is then the negative derivative of the defect energy, $\delta E$, against the crystal strain, $\epsilon$:

$$P = -\frac{\partial(\delta E)}{\partial \epsilon} \quad (1)$$

Among the various methods reviewed in Ref. [36] to extract $P$ from atomistic calculations, we use the definition derived from the stress tensor. For an unstrained volume $V$, the elastic

<figure><img src="image_3.png" /><figcaption>FIG. 1. Atoms in a vicinity of a cationic vacancy, $V''_{Mg}$. Small spheres show the positions of anions, and the bigger spheres show the cations. The atoms are colored according to their elastic dipole tensor.</figcaption></figure>

dipole tensor can be deduced from the point-defect residual stress, quoted $\sigma_{ij}$:

$$P_{ij} = -V\sigma_{ij} \quad (2)$$

Practically, $P$ is determined from the computation of the stress state of a periodic crystal in which only one vacancy is introduced. Since vacancies are charged, simulation volumes are neutralized with a uniform neutralizing background charge density [42,55]. Because of the interaction between the defect and its neutralizing background, it is worth noticing that the elastic dipole tensor has to be corrected as described in Appendix A.

In MgO, as shown in Fig. 1 and as previously discussed by Leslie and Gillan [42], the symmetry of the rock-salt structure results in the elastic dipole tensors of both cationic and anionic vacancies to be merely a multiple of the identity matrix with $P = pI$. We give the elastic dipole component $p$ in Table II. The $V''_{Mg}$ elastic dipole is found to be stronger than $V''_O$ elastic dipole. For both vacancy type, the elastic dipole is positive, indicating that the relaxation volume around a vacancy is also positive. Indeed, the relaxation volume $\Omega$ is conveniently deduced from the elastic dipole tensor as $\Omega = \frac{\text{tr}(P)}{3B} = p/B$. As already discussed in Refs. [39–42], a positive dilatation volume results, in ionic systems, from the strong coulombic repulsion that prevails between the ionic neighbors of the vacancy.

<table><caption>TABLE II. Elastic dipole component and relaxation volume for an isolated vacancy in MgO computed with the empirical potential [45].</caption><thead><tr><th></th><th>V''<sub>Mg</sub></th><th>V''<sub>O</sub></th></tr></thead><tbody><tr><td>p (eV)</td><td>5.35</td><td>4.88</td></tr><tr><td>Ω (Å<sup>3</sup>)</td><td>4.61</td><td>4.20</td></tr></tbody></table>
<figure><img src="image_4.png" /><figcaption>FIG. 2. View of a simulation cell containing the dislocation dipole, composed of a positive dislocation situated in (x+, y+) and a negative one in (x-, y-), and a Mg vacancy, V''$_{Mg}$.</figcaption></figure>

## B. Simulation cells

The definition of the simulation setup follows the one proposed in Ref. [25]. The simulated volume corresponds to a dipole of edge dislocations of opposite Burgers vectors $b = \pm 1/2[110]$. The orientation of the MgO crystal is chosen to match the conventional reference frame (see for instance Ref. [1]), i.e., the Burgers vector direction [110] is along x, the dislocation line [001] is along z and the normal to the glide plane [110] is parallel to the y axis. In the following, we quote (x+, y+) and (x-, y-) the positions of the positive and negative dislocation respectively. To satisfy with the periodic boundary conditions, the dislocations are located one above the other (x+=x-) and distributed along y so that the distances between them and their periodic images are equal. All calculations are performed within a constant volume (417 312 atoms) of lengths, $L_x \approx L_y \approx 80b \approx 240$ Å and $L_z = 16a_0 \approx 67$ Å. In order to compute the interaction energy between a point defect and an edge dislocation, we introduce a single charged vacancy in the neighborhood of one of the dislocations as illustrated Fig. 2. The use of a neutralizing background does not affect the dislocation since its core is charge neutral.

# III. DISLOCATION-VACANCY INTERACTION ENERGY

## A. Atomistic simulations

From our atomistic simulations, we compute the dislocation-vacancy interaction energy, $E_{\text{int}}^{\text{at}}$, according to the following expression:

$$E_{\text{int}}^{\text{at}} = E - E_{\text{dipole}} - \delta E. \quad (3)$$

In Eq. (3), E and $E_{\text{dipole}}$ are respectively the energy of the simulation cell for a dislocation dipole interacting with a vacancy, and for the dislocation dipole only. $\delta E$ corresponds to the extraction energy of a given type of ion and is computed as the difference between a perfect 3D-periodic crystal and an equivalent system containing a single vacancy. One notes that the extraction energy is merely linked to the sum of the

<figure><img src="image_5.png" /><figcaption>FIG. 3. Color maps of the interaction energy between the dislocation dipole and (a) V''$_{Mg}$ and (b) V$_{O}$ computed by atomistic simulations [Eq. (3)]. The dislocation core position (x+, y+) is highlighted in the figure by a purple “⊥”. The first plane above the dislocation, $n_y=1$, is the most attractive plane.</figcaption></figure>

defect formation energy and the chemical potentials of ions [56,57]. $\delta E$ is thus computed with systems scaled on the size of the cell containing the dislocation dipole. According to our potential parametrization, we find $\delta E = 17.89$ eV for the cationic vacancy and $\delta E = 18.49$ eV for the anionic vacancy. Such values are very close to the point-defect extraction energy of an isolated defect embedded in an infinite volume [58]. According to Eq. (3), the interaction energy is then negative when a vacancy binds to an attractive atomic site. By definition, the binding energy of a vacancy with a dislocation corresponds to the opposite of the interaction energy.

The interaction energies are summarized in Fig. 3. As expected in ionic crystals, because of the positive relaxation volume of vacancies, the interaction is globally repulsive in the compressive region above the dislocation glide plane, whereas it is attractive in the tensile region below the glide plane. Yet, we observe a strong attractive area close to the dislocation glide plane [see the marked blue horizontal atomic rows on Fig. 3]. It corresponds to vacancy sites on the two first atomic planes, i.e., planes above and below the glide plane quoted respectively $n_y=1$ and $n_y=-1$ [see arrows on Fig. 3]. The interaction energies reach values below -1 eV and are maximum at the close vicinity of the core with -1.01 eV for V''$_{Mg}$ and -1.35 eV for V$_{O}$ (Table III). For these two planes and whatever the nature of the vacancy, the attractive sites spread beyond 22 atomic rows (≈11b). Thus, the attractive region is larger than the dislocation core spreading (≈4b), estimated from the atomic disregistry.

Previous investigations of dislocation-point-defects interaction in MgO already reported the occurrence of highly attractive sites close to the core [18,33,34,59,60] emphasizing that the strong attraction does not result from the rigid ion potential used here (Table III). However, mostly dealing with
<table><caption>TABLE III. Extremum values of the interaction energy for the cationic and anionic vacancy in eV. The maximum attractive interaction is found for the configuration $M_1$ and the maximum repulsive interaction is found for the configuration $M_2$ (Fig. 7).</caption><thead><tr><th rowspan="2">Configurations:</th><th colspan="2">V''<sub>Mg</sub></th><th colspan="2">V''<sub>O</sub></th></tr><tr><th>M<sub>1</sub></th><th>M<sub>2</sub></th><th>M<sub>1</sub></th><th>M<sub>2</sub></th></tr></thead><tbody><tr><td>This works</td><td>-1.01</td><td>0.62</td><td>-1.35</td><td>0.51</td></tr><tr><td>Woo et al. [59]</td><td></td><td>1.1<sup>a</sup><br>0.77<sup>b</sup></td><td>-1.09<sup>a</sup><br>-1.40<sup>b</sup></td><td>1.16<sup>a</sup></td></tr><tr><td>Puls [33]</td><td></td><td>0.72<sup>a</sup><br>0.76<sup>b</sup><br>0.84<sup>c</sup></td><td>-1.51<sup>a</sup><br>-1.42<sup>b</sup><br>-1.61<sup>c</sup></td><td></td></tr><tr><td>Zhang et al. [18]</td><td>-1.7</td><td></td><td>-1.5</td><td></td></tr></tbody></table>

<sup>a</sup>Shell model.

<sup>b</sup>Point ion model.

<sup>c</sup>Breathing-shell model.

only few sites around the core, previous attempts [33,34,59] to characterize the attraction of point defects by the dislocation core did not capture the extension of this highly attractive region.

## B. Elasticity theory

Straight edge dislocation in MgO being uncharged, the interaction between the vacancy and the dislocation can be essentially described in the framework of the elasticity theory using the elastic dipole approach. The interaction energy of a point defect with a dislocation can be deduced from Eq. (4) giving the interaction energy, $E_{\text{int}}^{\text{th}}$, as a simple product between the elastic dipole, $P_{ij}$, and the intensity at the point-defect location of an external strain field, $\epsilon^{\text{ext}}$ [42–44]:

$$E_{\text{int}}^{\text{th}} = -P_{ij}\epsilon_{ij}^{\text{ext}} \quad (4)$$

Assuming that the external strain field $\epsilon^{\text{ext}}$ is the one created by the dislocation at the vacancy site and introducing the isotropic characteristic of $V''_O$ and $V''_{Mg}$, the interaction energy can be further simplified as the following:

$$E_{\text{int}}^{\text{th}} = -p^{\text{cell}}\Delta. \quad (5)$$

Since the vacancy elastic dipole component $p$ is sensitive to size effects, we used here the value $p^{\text{cell}}$ corresponding to the non cubic simulation cell defined in Sec. III. The vacancy elastic dipole components is then $p^{\text{cell}} = p_{xx} = p_{yy} = 5.55$ eV for $V''_{Mg}$ and 5.07 eV for $V''_O$. $\Delta = \text{tr}(\epsilon)$ refers to the dilatation field of the edge dislocation. In the following, we use the anisotropic elasticity theory to derive the dilatation field of the dislocation. Indeed, we verified in Appendix B that the dilatation field developed in the framework of anisotropy elasticity [see Eq. (B7) in Appendix B] reproduces the dilatation field extracted from our simulations. The comparison between $E_{\text{int}}^{\text{th}}$ and our atomic results $E_{\text{int}}^{\text{at}}$ is shown Fig. 4. Apart from a few atomic sites on planes above and below the glide plane, we note that the elastic description of the interaction energy matches the atomistic calculations within a few meV (i.e., 5 meV which is indeed the accuracy of the molecular static simulations given in the previous section). The largest energy differences between $E_{\text{int}}^{\text{th}}$ [Eq. (5)] and $E_{\text{int}}^{\text{at}}$ [Eq. (3)] are around

<figure><img src="image_6.png" /><figcaption>FIG. 4. Difference between atomistic simulations [Eq. (3)] and elasticity theory [Eq. (5)] for the interaction energy in the far field of the dislocation dipole for (a) $V''_{Mg}$ and (b) $V''_O$.</figcaption></figure>

0.07 eV and seem to be localized on the sides of the region of dislocation core spreading.

## IV. DISCUSSION

### A. Vacancies in the far field of the dislocation

Far from the dislocation core, the elasticity theory is expected to reproduce the interaction energies computed atomistically. According to Fig. 4 $E_{\text{int}}^{\text{th}}$ is qualitatively in reasonable agreement with $E_{\text{int}}^{\text{at}}$ although the global shape of the interaction predicted by the theory seems to be more spread than the one obtained by atomistic simulations.

To quantitatively assess this comparison, we have chosen to compute the coefficient of determination, quoted $R^2$, between $E_{\text{int}}^{\text{at}}$ and $E_{\text{int}}^{\text{th}}$ within a 2b wide disk as a function of the distance to the dislocation core (Fig. 5). $R^2$, between a set of data y (here $E_{\text{int}}^{\text{at}}$) and the theoretical predictions $\hat{y}$ (in our case $E_{\text{int}}^{\text{th}}$) writes:

$$R^2 = 1 - \frac{\sum_i (y_i - \hat{y}_i)^2}{\sum_i (y_i - \bar{y})^2} \quad (6)$$

with $\bar{y} = (\sum_{i=1}^n y_i)/n$ the average value of the data. $R^2$ in that convention may be negative, and tends to 1 when $y = \hat{y}$, i.e., when the predictions are in perfect agreement with the data.

Inside the dislocation core, there is strictly no correlation ($R^2 \simeq 0$). However, above 2b, we note a monotonic trend with a rapid increase of $R^2$ which ultimately reaches almost 1 at a distance of 10b whatever the vacancy type. Therefore, we define the near field of the dislocation, as being the zone for which $R^2 < 0.8$, that is at $\simeq 5.5b$ from the dislocation core (Fig. 5). Oppositely, in the far field of the dislocation (above 5.5b, corresponding to the planes $|n_y| = 12$) the agreement between the atomistic simulations and the elasticity theory
<figure><img src="image_7.png" /><figcaption>FIG. 5. Coefficient of determination, $R^2$, between the interaction energy, $E_{\text{int}}^{\text{at}}$, computed from atomistic simulations [Eq. (3)] and $E_{\text{int}}^{\text{th}}$, resulting from elasticity theory [Eq. (5)], as a function of the distance to the dislocation core.</figcaption></figure>

becomes remarkable. As an example, Fig. 6 shows that the elasticity theory is able to capture the order of magnitude of the interaction energy and its evolution as a function of the vacancy sites along [110]. Hence above 5.5b, an elastic treatment is able to reproduce both the intensity and the shape of the interaction energy with a reasonable accuracy.

## B. Vacancies in the near field of the dislocation

For the vacancy sites in the dislocation core region, i.e., below 4b, the elasticity theory is obviously not able to predict the interaction energies (see Fig. 3 and Fig. 5). In particular, above the glide plane, the first atomic layer ($n_y = 1$) expected to be repulsive corresponds atomistically to the most attractive one (see configuration $M_1$, blue squares compared to blue line, in Fig. 7).

A close inspection of the relaxed configuration of $V_O^{\ddot{}}$ bound to the core, shows that neighbors cations benefit from the local free volume accessible in the dislocation core to strongly repel [see green arrows on Fig. 8(a)]. Indeed, due to the local atomic arrangement and the core spreading of the edge dislocation investigated here, we observe a region [see region colored in red in Fig. 8(c)] for which the atomic volume is locally in excess with respect to the bulk configuration [Fig. 8(c) showing the atomic volume using the voronoi formulation]. For vacancies binding in this tensile region, the interaction energies is negative. This behavior follows the classical mechanism of relaxation around a point defect associated with a positive elastic dipole tensor in ionic crystals [42]. For instance, layers adjacent to the glide plane, $n_y = 1$ and $n_y = -1$, are the zones containing the most free space [Figs. 8(a) and 8(c)] and thus are the most attractive sites for vacancies. Surprisingly, plane $n_y = 1$ is even more attractive than plane $n_y = -1$ (Fig. 3). An explanation might result in the second neighbors relaxation capabilities. Indeed, in case

<figure><img src="image_8.png" /><figcaption>FIG. 6. A comparison between atomistic simulations [symbols, Eq. (3)] and elasticity theory [lines, Eq. (5)] at $|n_y| = 12$, the limit between the near and the far field of the dislocation, i.e., at 5.5b from the core (Fig. 5), for the interaction energy between the dislocation and (a) $V''_{\text{Mg}}$ and (b) $V_O$.</figcaption></figure>

of a vacancy in configuration $M_1$, the neighbor ions belonging to $n_y = 2$ seem to undergo slight relaxation. In such a case the compressed atoms will gain some space due to the local rearrangement (Fig. 9), so the system energy will decrease. Oppositely, for a vacancy belonging in plane $n_y = -1$, the neighbor ions belonging to the plane $n_y = -2$ are already in a tensile zone, so a gain in space will not be equally beneficial than introducing a vacancy in configuration $M_1$.

The interaction energy is thus strongly related to the free volume of atoms induced by the core structure of the dislocation which results from a complex relaxation scheme of atoms. In the near field, since the elasticity theory treats the dislocation as a singularity, such a behavior inherent of the peculiar core structure of the dislocation and its spreading can not be inferred by the elasticity theory.
<figure><img src="image_9.png" /><figcaption>FIG. 7. A comparison between atomistic simulations [symbols, Eq. (3)] and elasticity theory [continuous lines for $n_y = 1$, dashed lines for $n_y = 2$, Eq. (5)] of the interaction energy between the dislocation and (a) $V''_{Mg}$ and (b) $V''_O$. The interaction energy extrema are spot by letters: $M_1$ and $M_2$.</figcaption></figure>

This close correlation between the core structure and the interaction energy is also visible for repulsive areas. For instance, in Fig. 8(b), we show the atomic configuration $M_2$ of $V_O$. For such configuration, the vacancy is located at plane $n_y = 2$ corresponding to the layer of highest compression state [Fig. 8(c)] and therefore corresponds to the most repulsive interaction energy for a vacancy site (see yellow dots in Fig. 7). It is almost the only strong unstable vacancy configuration in close vicinity of the dislocation core (Fig. 3).

### C. Implications

Experimentally, it has been reported that vacancies, and especially ionic vacancy, can efficiently pin dislocations [31] in MgO. Such an assumption can be revisited according to

<figure><img src="image_10.png" /><figcaption>FIG. 8. [(a) and (b)] Snapshot of the atomic configurations near $V_O$ and the dislocation core ($x_+, y_+$) determined by atomistic simulations. (a) For configuration $M_1$; (b) For configuration $M_2$ (Fig. 7). The red spheres represent the oxygen ions and the green ones the magnesium ions. The arrows highlight the ion motions (yellow for oxygen ions and green for magnesium ions). (c) Magnification of the atomic configuration near the perfect dislocation core ($x_+, y_+$). The atoms are colored according to the difference between the atomic volume near the dislocation and the perfect atomic volume. The atomic volume are determined using the voronoi formulation.</figcaption></figure>

our atomistic calculations by computing the pinning force of an edge $1/2 \langle 110 \rangle \{\bar{1}10\}$ resulting from the interaction with $V''_{Mg}$ or $V''_O$.

By definition, the pinning force is the maximum of the attractive interaction force, $f_p = \max(|f_{int} - 0|)$, between the dislocation and the vacancy. $f_{int}$ is derived from the interaction
<figure><img src="image_11.png" /><figcaption>FIG. 9. Snapshot of the atomic configurations in the vicinity of the dislocation core. The O atoms are represented by small red spheres and the Mg atoms by larger green spheres. The bonds between the O atoms are shown. (a) Perfect dislocation core. [(b) and (c)] The most stable configuration, for which $V_O$ is situated in the plan $n_y = 1$. For a better readability the O atoms in vicinity of the vacancy are colored in blue.</figcaption></figure>

energy [61,62]:

$$f_{\text{int}} = -\frac{\partial E_{\text{int}}}{\partial x} \qquad (7)$$

We thus extract the interaction force from our simulations by numerically deriving the interaction energy, $E_{\text{int}}^{\text{at}}$ [Eq. (3)], in such a way $f_{\text{int}}^{\text{at}} = -\partial E_{\text{int}}^{\text{at}}/\partial x$ (see symbols in Fig. 10).

Moreover, in the dislocation far field, the elasticity theory can accurately predict the interaction energy (Fig. 5), so the interaction forces can be given by the following elastic formulation:

$$f_{\text{int}}^{\text{th}}(x, y) = -\frac{\partial E_{\text{int}}^{\text{th}}}{\partial x} = \frac{pkb}{\pi} \frac{2xy[\Lambda(x^4 + y^4) - 2x^2y^2]}{(x^4 - 2\Lambda x^2y^2 + y^4)^2} \quad (8)$$

with

$$\kappa = \frac{C_{12} - C_{11}}{C_{11}} \sqrt{\frac{C_{11}C_{44}}{(C_{11} - C_{12})(C_{11} + C_{12} + 2C_{44})}} \quad (9)$$

and

$$\Lambda = \frac{C_{11}^2 - 6C_{11}C_{44} - C_{12}(C_{12} + 2C_{44})}{(C_{11} - C_{12})(C_{11} + C_{12} + 2C_{44})} \quad (10)$$

Indeed, a comparison between $f_{\text{int}}^{\text{at}}$ and $f_{\text{int}}^{\text{th}}$ at the border of the dislocation far field ($|n_y| = 12$) [Figs. 10(a) and 10(b)]

<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /></figure>

<figure><img src="image_14.png" /><figcaption>FIG. 10. A comparison between atomistic simulations [symbols, Eq. (7)] and elasticity theory [continuous lines for $V''_{\text{Mg}}$, dashed lines for $V''_{\text{O}}$, Eq. (8)] of the interaction forces at plane (a) $n_y = 12$, (b) $n_y = -12$, i.e., the limit between the far field and the near field of the dislocation and at plane (c) $n_y = 1$, i.e., the most attractive plane in the dislocation core.</figcaption></figure>

shows that in the dislocation far field, Eq. (8) predicts interaction forces in fairly good agreement with $f_{\text{int}}^{\text{at}}$ derived from our atomistic calculations.
From Eq. (8), one note that $f_{\text{int}}(x, y) = f_{\text{int}}(-x, -y)$ and $f_{\text{int}}(x, y) = -f_{\text{int}}(x, -y)$. Therefore, for a vacancy situated at $n_y = 12$, the interaction is repulsive (Fig. 6), as well as the force [Fig. 10(a)] while at plane $n_y = -12$, the force is attractive and is strictly the opposite of the repulsive force occurring at $n_y = 12$. As already pointed out, the singularity of the core in the elasticity theory leads to an abrupt inversion of the interaction force as the core is crossed whereas such a discontinuity cancels in case of $f_{\text{int}}^{\text{at}}$ accounting for the spreading of the core.

Yet, for vacancies far from the dislocation glide plane [as illustrated Figs. 10(a) or 10(b)], the interaction force has practically vanished.

Since the interaction between a point defect and the dislocation is strongest in the vicinity of the glide plane, it is more relevant to focus on the atomic layers just above and below the glide plane. For the $n_y = 1$ plane, we already reported that the strongest attractive sites spread over 22b from the dislocation for both $V''_{\text{Mg}}$ and $V_O''$ (see blue squares for $x \in [-20\AA, 20\AA]$ in Fig. 7). It is thus not surprising to find that the interaction force reaches its highest values for vacancies located in the first atomic layers [as illustrated Fig. 10(c)]. According to our simulations, $f_{\text{int}}^{\text{at}}$ reaches 0.097 eV/Å for $V''_{\text{Mg}}$ and 0.142 eV/Å for $V_O''$. The interaction force reaches its maximum [Fig. 10(c)] at the border of the dislocation core region (i.e., ≈6Å, consistently with a half-width of the core of ≈2b) and evolves smoothly throughout the dislocation core. The spreading of the core region has thus an effect on the actual values of the pinning forces. Despite the dislocation core being strongly attractive for vacancies (Fig. 3), the spreading of the dislocation allows a gradual variation of $E_{\text{int}}^{\text{at}}$ in the core region. Consequently, the pinning forces as computed with our atomistic simulations for the $n_y = 1$ plane are relatively weak, at least weaker than any elastic prediction. It is worth noticing that $V_O''$ are stronger pinning points than $V''_{\text{Mg}}$ in agreement with the conclusions of Ref. [31]. However, with a few tens of eV/Å, even anionic vacancies should be considered as weak pinning points. Hence the interaction of $V''_{\text{Mg}}$ or $V_O''$ with edge dislocation might only slightly contribute to the hardening in MgO. Interaction between dislocations and tetragonal defects, such as impurity-vacancy dipole, are expected to govern the mechanism of hardening in ionic crystal. However, to our knowledge, the resulting pinning forces have only been evaluated through various elastic treatments [61,62]. Based on the previous investigation, the elasticity theory of point-defect-dislocation interaction in the far field region of the edge dislocation seems fairly robust but its prediction in the vicinity of the dislocation may be largely incorrect. As an example, the use of Eq. (8) for vacancies binding in the core leads to overestimate the pinning forces by a factor of 4 [as illustrated Fig. 10(c)]. A treatment purely based on the elasticity theory as done in Refs. [61,62] is thus questionable and may require an atomistic investigation in order to infer reasonable values for the pinning forces.

To complete the discussion of the interaction of point defects with dislocation in MgO, we shall add that one of the striking results of our atomistic simulations is the occurrence of a large region of attractive sites for vacancies. Such an occurrence of strong binding raft spread in the $x[110]$ direction may results in a sweep-up mechanism [63] from which the

dislocation might acquire and carry electric charges during its glide. This peculiar region in the vicinity of the dislocation core can also affect the diffusion by inducing short-circuit pathways leading to potentially athermal pipe diffusion.

## V. CONCLUSION

In this study, we quantified the interaction between $1/2\langle 110 \rangle\{\overline{1}10\}$ edge dislocations and vacancies in MgO using molecular static simulations and elasticity theory. Using the framework of the elastic dipole tensor to characterise a vacancy, we showed that in the far field of the dislocation (above ≈ 5.5b), the anisotropic elasticity theory accurately predicts the interaction between the edge dislocation and vacancies. Indeed, both the interaction energy and the interaction force are well captured by elastic treatments. Otherwise, in the near field of the dislocation that contains the most attractive and repulsive sites, we show that the interaction strongly depends on the dislocation core structure. Hence, the elasticity theory describing the dislocation core as a singularity is not able to reproduce the correct vacancy response in the near field. Within a region across the glide planes that expands over several Burgers vector, the edge dislocation core is strongly attractive for vacancies, especially the anionic ones. However, isolated vacancies remain weak pinning points which are not expected to contribute significantly to the hardening observed experimentally in MgO.

## ACKNOWLEDGMENTS

This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program under Grant No 787198 TimeMan. Philippe Carrez thanks the support of the French government through the Programme Investissement d'Avenir (I-SITE ULNE / ANR-16-IDEX-0004 ULNE) managed by the Agence Nationale de la Recherche, for the project named LASCO. Computational resources have been provided by the DSI-calcul intensif at Université de Lille.

## APPENDIX A: CELL SIZE DEPENDENCY ON THE ELASTIC DIPOLE TENSOR

In this study, the determination of the elastic dipole tensor $P$ have been carried out using standard cubic system orientation. Once one ion has been removed from the system, we used the method described in Ref. [42] to account for the introduction of an isolated charged vacancy.

As explained in Sec. II, the elastic dipole tensor is computed through the evaluation of the stress tensor induced by the vacancy in the volume $V$. As for the vacancy formation energy, using periodic boundary conditions, the stress tensor resulting from the insertion of a vacancy is therefore system size dependent. Figure 11 shows the evolution of $P$ as a function of the cell size, $L$ (see circles in Fig. 11). In case of non neutral periodic system, the sensitivity of $P$ to the system size results from two major contribution, (1) the elastic interaction between the defect and its periodic images [37], (2) a coulombic interaction between the vacancy and the uniform background charge introduced to ensure a proper convergence
<figure><img src="image_15.png" /><figcaption>FIG. 11. Elastic dipole component plotted versus the inverse of the cell size, L, for (a) the cationic vacancy and (b) the an-ionic vacancy. The dots correspond to the uncorrected elastic dipole component, directly extracted from our simulations, the small stars correspond to the elastic dipole component corrected from the Coulomb interaction and the large stars correspond to the elastic dipole component corrected from the Coulomb and from the elastic interaction.</figcaption></figure>

of the electrostatic summation [42]. Indeed, both contribute to a correction term to the energy given by:

$$ \delta E^{\text{elas}} \propto \frac{1}{L^3}, \qquad (\text{A1}) $$

$$ \delta E^{\text{coul}} = -\frac{\alpha}{2} \frac{q^2}{4\pi\epsilon_{\text{mgo}}L}, \qquad (\text{A2}) $$

with here a Madelung constant $\alpha = 2.8373$ and $\epsilon_{\text{mgo}} = 8.5\epsilon_0$ the permittivity of MgO according to the potential parametrization.

The elastic dipole tensor scaling with the strain derivative of the energy, a similar correction has to be applied in order to compute the corrected elastic dipole component $p(\mathbf{P} = p\mathbf{I})$. Here we used thus, $\delta p^{\text{coul}} \propto \frac{1}{L}$ and $\delta p^{\text{elas}} \propto \frac{1}{L^3}$.

In Fig. 11, we distinguished the elastic dipole tensor corrected only from the Coulomb effects (small purple stars in Fig. 11) or the one obtained when all the interactions were taken into account (see big yellow stars in Fig. 11). From these results we can see that the Coulomb interaction is once again predominant over the elastic interaction. The corrected value, $p$, for the cationic and anionic vacancy are quoted in Table II.

## APPENDIX B: DILATATION FIELD, A COMPARISON BETWEEN ATOMISTIC SIMULATIONS AND ELASTICITY THEORY

In order to determine the atomic dilatation field, $\Delta^{\text{at}}$, from our minimized dislocation dipole system, we first compute a per-atom-stress tensor, using a virial development for the force acting on atom, divided by a atomic volume, expressed as a voronoi volume. The reverse Hook law gives then:

$$ \Delta^{\text{at}} = (s_{11} + 2s_{12})(\sigma_{xx} + \sigma_{yy} + \sigma_{zz}) \qquad (\text{B1}) $$

with the $s_{11}$ and $s_{12}$ the compliances in the standard orientation [64], linked to the stiffnesses (Table I) by:

$$ s_{11} = \frac{C_{11} + C_{12}}{(C_{11} - C_{12})(C_{11} + 2C_{12})} \qquad (\text{B2}) $$

and

$$ s_{12} = - \frac{C_{12}}{(C_{11} - C_{12})(C_{11} + 2C_{12})} \qquad (\text{B3}) $$

The corresponding atomic dilatation field, Eq. (B1), is plotted in Fig. 12(a).

In the frame of the elasticity theory, the dilatation of an edge dislocation, $\Delta_i$, writes [64]:

$$ \Delta_i(x, y) = \frac{\kappa b}{\pi} \frac{y(y^2 - \Delta x^2)}{x^4 - 2\Delta x^2 y^2 + y^4} \qquad (\text{B4}) $$

with

$$ \kappa = \frac{C_{12} - C_{11}}{C_{11}} \sqrt{\frac{C_{11}C_{44}}{(C_{11} - C_{12})(C_{11} + C_{12} + 2C_{44})}} \qquad (\text{B5}) $$

and

$$ \Lambda = \frac{C_{11}^2 - 6C_{11}C_{44} - C_{12}(C_{12} + 2C_{44})}{(C_{11} - C_{12})(C_{11} + C_{12} + 2C_{44})} \qquad (\text{B6}) $$

In case of isotropic elasticity, $\Lambda = -1$ and $\kappa = (C_{12} - C_{11})/2C_{11}$, and the isotropic dilatation field, scaling with the Poisson ratio $\nu = C_{12}/(C_{11} + C_{12})$, becomes symmetric as given in numerous textbooks as for instance in Ref. [1].

For a proper comparison between the elasticity theory and the atomistic results, it is necessary to consider not only the dislocation dipole but also its periodic images. In the following, this is achieved by a simple summation at the position $(x', y')$:

$$ \Delta^{\text{th}} = \sum_{n=-R}^{R} \sum_{m=-R}^{R} (\Delta_{1} + \Delta_{2}) \qquad (\text{B7}) $$

with $\Delta_1 = \Delta_i(x' - (x_+ + nL_x), y' - (y_+ + mL_y))$ and $\Delta_2 = -\Delta_i(x' - (x_- + nL_x), y' - (y_- + mL_y))$. We further verified that a summation over 10 periodic images is sufficient within the accuracy of our atomistic simulations.
<figure><img src="image_16.png" /><figcaption>FIG. 12. Dilatation around the dislocation situated in $(x_+, y_+)$ (a) $\Delta^{\text{at}}$ computed from the atomistic simulations stress tensor [Eq. (B1)], (b) $\Delta^{\text{th}}$ developed in anisotropic elasticity theory, or (c) in isotropic elasticity theory [Eq. (B7)]. The dislocation core position $(x_+, y_+)$ is highlighted in the figure by a purple $\perp$.</figcaption></figure>

A typical solution in case of anisotropic elasticity [Eq. (B7)] is given as a color map in Fig. 12(b). Qualitatively, one can see on the color maps that there are only mild differences between the dilatation determined with atomistic simulations [Fig. 12(a)] and the dilatation predicted by anisotropic elasticity theory [Fig. 12(b)]. The shape of the dilatation seem slightly sharper on the theory and the atomistic results are more spread out. This discrepancy could be related

<figure><img src="image_17.png" /><figcaption>FIG. 13. Coefficient of determination, $R^2$, between dilatation from atomistic simulations and from elasticity theory, computed on a disk of $2b$ wide, as a function of the distance to the dislocation core.</figcaption></figure>

to the spreading of the edge dislocation core. Indeed, all the elastic treatment performed assumes a Volterra dislocation, i.e., the core is a singularity in the elastic media and consequently do not account for a width of the shear distribution in the vicinity of the dislocation core. In contrast, we can see that the agreement is less pronounced between isotropic elasticity [Fig. 12(c)] and atomistic simulations [Fig. 12(a)].

To be quantitative, we measured the coefficient of determination $R^2$ between Eq. (B1) and Eq. (B7) at a given distance $r$ from the dislocation core (Fig. 13). Isotropic solution being easily tractable, it is a common temptation to refer to its expression. However, the isotropic solution, with a correlation $R^2$ below 0.8 (see purple squares in Fig. 13), fails to reproduce the atomistic calculations even at large distance from the dislocation core. In contrast, except in the dislocation core region, the anisotropic dilatation field matches our simulation results with a coefficient of determination reaching almost 1 at distance of a few $b$ (see yellow dots in Fig. 13). One could possibly increase the fidelity of the anisotropic dilatation field by accounting for the spreading of the core. Nevertheless, if one defines an accuracy criterion $R^2 = 0.8$, the elasticity theory provides a sufficient prediction ($R^2 \approx 0.9$) even in an area close to the dislocation core of half-width of nearly $2b$.

[1] J. P. Hirth, J. Lothe, and T. Mura, *Theory of Dislocations* (American Society of Mechanical Engineers Digital Collection, New York, 1983).

[2] G. Hachet, A. Oudriss, A. Barnoush, R. Milet, D. Wan, A. Metsue, and X. Feaugas, The influence of hydrogen on cyclic plasticity of <001> oriented nickel single crystal. Part I: Dislocation organisations and internal stresses, *Int. J. Plast.* **126**, 102611 (2020).

[3] G. Hachet, A. Metsue, A. Oudriss, and X. Feaugas, The influence of hydrogen on cyclic plasticity of <001> oriented nickel single crystal. Part II: Stability of edge dislocation dipoles, *Int. J. Plast.* **129**, 102667 (2020).

[4] S. Patinet and L. Proville, Dislocation pinning by substitutional impurities in an atomic-scale model for Al(Mg) solid solutions, *Philos. Mag.* **91**, 1581 (2011).
[5] E. Rodary, D. Rodney, L. Proville, Y. Bréchet, and G. Martin, Dislocation glide in model Ni(Al) solid solutions by molecular dynamics, *Phys. Rev. B* **70**, 054111 (2004).

[6] R. Fleischer, Substitutional solution hardening, *Acta Metall.* **11**, 203 (1963).

[7] J. Baird, The effects of strain-ageing due to interstitial solutes on the mechanical properties of metals, *Metall. Rev.* **16**, 1 (1971).

[8] A. Abu-Odeh, M. Cottura, and M. Asta, Insights into dislocation climb efficiency in fcc metals from atomistic simulations, *Acta Mater.* **193**, 172 (2020).

[9] C. McElfresh, Y. Cui, S. Dudarev, G. Po, and J. Marian, Discrete stochastic model of point defect-dislocation interaction for simulating dislocation climb, *Int. J. Plast.* **136**, 102848 (2020).

[10] T. D. Swinburne, K. Arakawa, H. Mori, H. Yasuda, M. Ishiki, K. Mimura, M. Uchikoshi, and S. Dudarev, Fast, vacancy-free climb of prismatic dislocation loops in bcc metals, *Sci. Rep.* **6**, 30596 (2016).

[11] J. Lothe, Theory of dislocation climb in metals, *J. Appl. Phys.* **31**, 1077 (1960).

[12] J. Weertman, Theory of steady-state creep based on dislocation climb, *J. Appl. Phys.* **26**, 1213 (1955).

[13] P. Hirel, P. Carrez, E. Clouet, and P. Cordier, The electric charge and climb of edge dislocations in perovskite oxides: The case of high-pressure MgSiO3 bridgmanite, *Acta Mater.* **106**, 313 (2016).

[14] R. Whitworth, Charged dislocations in ionic crystals, *Adv. Phys.* **24**, 203 (1975).

[15] H. Layer, M. Miller, and L. Slifkin, Effects of quenching and annealing on ionic conductivity and on dislocation decorability in silver chloride, *J. Appl. Phys.* **33**, 478 (1962).

[16] J. Hedges and J. Mitchell, XLI. Some experiments on photographic sensitivity, *Lond. Edinb. Dubl. Philos. Mag.* **44**, 357 (1953).

[17] J. Hedges and J. Mitchell, The observation of polyhedral substructures in crystals of silver bromide, *Lond. Edinb. Dubl. Philos. Mag.* **44**, 223 (1953).

[18] F. Zhang, A. Walker, K. Wright, and J. Gale, Defects and dislocations in MgO: Atomic scale models of impurity segregation and fast pipe diffusion, *J. Mater. Chem.* **20**, 10445 (2010).

[19] J. Rabier and M. Puls, Atomistic model calculations of pipe-diffusion mechanisms in MgO, *Philos. Mag. A* **52**, 461 (1985).

[20] J. Narayan and J. Washburn, Self-climb of dislocation loops in magnesium oxide, *Philos. Mag.* **26**, 1179 (1972).

[21] G. Groves and A. Kelly, Climb of dislocations in magnesium oxide, *J. Appl. Phys.* **33**, 456 (1962).

[22] J. Li, Charged dislocations and plasto-electric effect in ionic crystals, *Mater. Sci. Eng., A* **287**, 265 (2000).

[23] D. Yang and H. Conrad, Influence of an electric field on the plastic deformation of polycrystalline NaCl at elevated temperatures, *Acta Mater.* **46**, 1963 (1998).

[24] J. Narayan, R. Weeks, and E. Sonder, Aggregation of defects and thermal-electric breakdown in MgO, *J. Appl. Phys.* **49**, 5977 (1978).

[25] J. Zhai, P. Hirel, and P. Carrez, Atomic-scale properties of jogs along 1/2< 110> {110} edge dislocations in MgO, *Scr. Mater.* **181**, 66 (2020).

[26] A. Clauer and B. Wilcox, High temperature tensile creep of magnesium oxide single crystals, *J. Am. Ceram. Soc.* **59**, 89 (1976).

[27] W. Rothwell and A. Neiman, Creep in vacuum of MgO single crystals and the electric field effect, *J. Appl. Phys.* **36**, 2309 (1965).

[28] C. Hulse and J. Pask, Mechanical properties of magnesia single crystals corn pression, *J. Am. Ceram. Soc.* **43**, 373 (1960).

[29] A. Foitzik, W. Skrotzki, and P. Haasen, Correlation between microstructure, dislocation dissociation and plastic anisotropy in ionic crystals, *Mater. Sci. Eng., A* **113**, 399 (1989).

[30] F. Sato and K. Sumino, The yield strength and dynamic behaviour of dislocations in MgO crystals at high temperatures, *J. Mater. Sci.* **15**, 1625 (1980).

[31] D. Cáceres, I. Vergara, R. González, and Y. Chen, Effect of neutron irradiation on hardening in MgO crystals, *Phys. Rev. B* **66**, 024111 (2002).

[32] E. Clouet, The vacancy-edge dislocation interaction in fcc metals: A comparison between atomic simulations and elasticity theory, *Acta Mater.* **54**, 3543 (2006).

[33] M. Puls, Vacancy-dislocation interaction energies in MgO a re-analysis, *Philos. Mag. A* **47**, 497 (1983).

[34] M. Puls, Vacancy-dislocation interaction energies in MgO, *Philos. Mag. A* **41**, 353 (1980).

[35] R. Bullough and R. Newman, The interaction of vacancies with dislocations, *Philos. Mag.* **7**, 529 (1962).

[36] E. Clouet, C. Varvenne, and T. Jourdan, Elastic modeling of point-defects and their interaction, *Comput. Mater. Sci.* **147**, 49 (2018).

[37] C. Varvenne and E. Clouet, Elastic dipoles of point defects from atomistic simulations, *Phys. Rev. B* **96**, 224103 (2017).

[38] D. Carpentier, T. Jourdan, Y. Le Bouar, and M.-C. Marinica, Effect of saddle point anisotropy of point defects on their absorption by dislocations and cavities, *Acta Mater.* **136**, 323 (2017).

[39] R. Nazarov, J. S. Majevadia, M. Patel, M. R. Wenman, D. S. Balint, J. Neugebauer, and A. P. Sutton, First-principles calculation of the elastic dipole tensor of a point defect: Application to hydrogen in α-zirconium, *Phys. Rev. B* **94**, 241112(R) (2016).

[40] A. Goyal, S. R. Phillpot, G. Subramanian, D. A. Andersson, C. R. Stanek, and B. P. Uberuaga, Impact of homogeneous strain on uranium vacancy diffusion in uranium dioxide, *Phys. Rev. B* **91**, 094103 (2015).

[41] D. A. Freedman, D. Roundy, and T. A. Arias, Elastic effects of vacancies in strontium titanate: Short-and long-range strain fields, elastic dipole tensors, and chemical strain, *Phys. Rev. B* **80**, 064108 (2009).

[42] M. Leslie and N. Gillan, The energy and elastic dipole tensor of defects in ionic crystals calculated by the supercell method, *J. Phys. C: Solid State Phys.* **18**, 973 (1985).

[43] M. Gillan, The elastic dipole tensor for point defects in ionic crystals, *J. Phys. C: Solid State Phys.* **17**, 1473 (1984).

[44] E. Kröner, *Continuum Theory of Defects, Physics of Defects*, Vol. 35, edited by J. P. Poirier, R. Balian, and M. Kléman (North-Holland, Amsterdam, 1981).

[45] G. Henkelman, B. P. Uberuaga, D. J. Harris, J. H. Harding, and N. L. Allan, MgO addimer diffusion on MgO (100): A comparison of ab initio and empirical models, *Phys. Rev. B* **72**, 115437 (2005).

[46] P. Carrez, J. Godet, and P. Cordier, Atomistic simulations of 1/2< 110> screw dislocation core in magnesium oxide, *Comput. Mater. Sci.* **103**, 250 (2015).
[47] P. Hirel, G. Moladje, P. Carrez, and P. Cordier, Systematic theoretical study of [001] symmetric tilt grain boundaries in MgO from 0 to 120 GPa, Phys. Chem. Miner. **46**, 37 (2019).

[48] S. Plimpton, Fast parallel algorithms for short-range molecular dynamics, J. Comput. Phys. **117**, 1 (1995).

[49] J. Eastwood, R. Hockney, and D. Lawrence, P3M3DPTThe three-dimensional periodic particle-particle/particle-mesh program, Comput. Phys. Commun. **19**, 215 (1980).

[50] E. Bitzek, P. Koskinen, F. Gähler, M. Moseler, and P. Gumbsch, Structural Relaxation Made Simple, Phys. Rev. Lett. **97**, 170201 (2006).

[51] P. Carrez, D. Ferré, and P. Cordier, Peierls-nabarro modeling of dislocations in MgO from ambient pressure to 100 GPa, Modell. Simul. Mater. Sci. Eng. **17**, 035010 (2009).

[52] R. Reeber, K. Goessel, and K. Wang, Thermal expansion and molar volume of MgO, periclase, from 5 to 2900 K, Eur. J. Mineral. **7**, 1039 (1995).

[53] Y. Sumino, I. Ohno, T. Goto, and M. Kumazawa, Measurement of elastic constants and internal frictions on single-crystal MgO by rectangular parallelepiped resonance, J. Phys. Earth **24**, 263 (1976).

[54] K. Marklund and S. Mahmoud, Elastic constants of magnesium oxide, Phys. Scr. **3**, 75 (1971).

[55] J. Hub, B. de Groot, H. Grubmuller, and G. Groenhof, Quantifying artifacts in Ewald simulations of inhomogeneous

systems with a net charge, J. Chem. Theory Comput. **10**, 381 (2014).

[56] S. E. Taylor and F. Bruneval, Understanding and correcting the spurious interactions in charged supercells, Phys. Rev. B **84**, 075155 (2011).

[57] S. B. Zhang and J. E. Northrup, Chemical Potential Dependence of Defect Formation Energies in GaAs: Application to Ga Self-Diffusion, Phys. Rev. Lett. **67**, 2339 (1991).

[58] S. Mahmoud, P. Carrez, M. Landeiro Dos Reis, N. Mousseau, and P. Cordier, Diffusion mechanism of bound Schottky defect in magnesium oxide, Phys. Rev. Materials **5**, 033609 (2021).

[59] C. Woo, M. Puls, and M. Norgett, Energies of formation of point defects in perfect and dislocated ionic crystals, J. Phys. Colloq. **37**, C7 (1976).

[60] R. Skelton and A. Walker, Interactions between bare and prototyped Mg vacancies and dislocation cores in MgO, Phys. Chem. Miner. **46**, 471 (2019).

[61] J. Soullard, The interaction force between the tetragonal defects and edge dislocations in the rock salt structure, Philos. Mag. A **64**, 39 (1991).

[62] D. Barnett and W. Nix, The interaction force between tetragonal defects and screw dislocations in cubic crystals, Acta Metall. **21**, 1157 (1973).

[63] J. Rabier and M. Puls, Atomistic calculations of point-defect interaction and migration energies in the core of an edge dislocation in NaCl, Philos. Mag. A **59**, 533 (1989).

[64] J. W. Steeds, *Introduction to Anisotropic Elasticity Theory of Dislocations* (Oxford University Press, New York, 1973), p. 89.