# Continuous description of grain boundaries using crystal defectfields: the example of a 3 1 0/[0 0 1] tilt boundary in MgO

Xiao-Yu Sun, Patrick Cordier, Vincent Taupin, Claude Fressengeas, Bijaya B Karki

► To cite this version:

Xiao-Yu Sun, Patrick Cordier, Vincent Taupin, Claude Fressengeas, Bijaya B Karki. Continuous description of grain boundaries using crystal defectfields: the example of a 3 1 0/[0 0 1] tilt boundary in MgO. European Journal of Mineralogy, 2017, 29 (2), pp.155-165. 10.1127/ejm/2017/0029-2609. hal-02938037

HAL Id: hal-02938037

https://hal.univ-lorraine.fr/hal-02938037

Submitted on 14 Sep 2020

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
# Continuous description of grain boundaries using crystal defect fields: the example of a {310}/[001] tilt boundary in MgO

XIAO-YU SUN<fnref n="1" />*, PATRICK CORDIER<fnref n="1" />, VINCENT TAUPIN<fnref n="2" />, CLAUDE FRESSENGEAS<fnref n="2" /> and BIJAYA B. KARKI<fnref n="3" />,<fnref n="4" />

Corresponding author: Dr. Xiao-Yu Sun

Université Lille 1

Unité Matériaux et Transformations - UMR CNRS 8207 - Bat C6

59655 Villeneuve d'Ascq Cedex, France

Tel: +33-320-434341

E-mail: xiaoyu.sun@univ-lille1.fr

## Abstract

Continuum dislocation/generalized-disclination modelling of a MgO {310}/[001] tilt boundary is proposed by using an atomistic-to-continuum crossover method. Using the relaxed and un-relaxed atomic positions obtained from first-principles simulations, calculation of the displacement, distortion (strain and rotation), curvature and second order distortion fields is performed in the boundary area. The dislocation and disclination density fields are used to smoothly describe the discontinuities in the displacement and rotation fields across the boundary. In addition, a generalized-disclination density tensor field is introduced to reflect shear and stretch discontinuities across the boundary. This continuous model provides new mechanical insights into the structure of MgO grain boundary, as well as a basis for coarse-grained representations of the polycrystalline structure.

**Keywords:** grain boundary; MgO; dislocation; generalized-disclination; continuum mechanics; atomistic simulations

<footnote n="1">¹ Unité Matériaux et Transformations, UMR 8207 CNRS/Université Lille1, 59655 Villeneuve d'Ascq, France</footnote>

<footnote n="2">² Laboratoire d'Etude des Microstructures et de Mécanique des Matériaux (LEM3), Université de Lorraine/CNRS, Ile du Saulcy, 57045 Metz Cedex, France</footnote>

<footnote n="3">³ Department of Computer Science, Louisiana State University, Baton Rouge, Louisiana 70803, U.S.A</footnote>

<footnote n="4">⁴ Department of Geology and Geophysics, Louisiana State University, Baton Rouge, Louisiana 70803, U.S.A.</footnote>
# 1. Introduction

The Earth's lower mantle is mainly composed of silicates with perovskite structure and of ferropericlase (Mg, Fe)O. The structure of crystal defects such as grain boundaries in these compounds is of primary importance for their rheology. The structure of grain boundaries (GBs) in MgO has been studied mainly at the atomic scale, on the basis of first-principles methods and molecular dynamics (MD) simulations. For instance, a number of MgO tilt boundaries were constructed at atomistic level to investigate the effects of pressure and temperature on their structure and stability (Harris et al., 1996), the related vacancy migration (Harris et al., 1997), the pressure-induced structural transition (Harris et al., 1999), and the impurity-induced structural transformation (Saito et al., 2014; Yan et al., 1998). The structure and stability of oxygen vacancies on the terraces, sub-surface and low-coordinated surface sites of MgO were studied by performing *ab initio* Hartree-Fock cluster model calculations (Pacchioni and Pescarmona, 1998). The electronic properties of oxygen vacancies, proton impurities, and their stable combinations at a MgO boundary were calculated by employing the embedded cluster method (McKenna and Shluger, 2009). The atomic and electronic structures of three {n10}/[001] MgO tilt GBs (n=2, 3, and 4) were explored as a function of pressure from first principles within the density functional theory (Verma and Karki, 2010).

Atomistic studies are useful to provide information on fundamental material properties at GBs such as defect migration, diffusion and stability of vacancies. However, the size of the systems investigated and the time scale over which simulations can be conducted are limited for computational reasons. Further, they do not give direct insight into the role played by GBs in the dissipative mechanical behaviour of the material during plastic yielding (mediated by shear, sliding or migration), because describing the latter essentially needs the stress, strain and strain-rate fields, and their relationships in the body, particularly in the GB area. Continuous modeling of the structure of GBs using smooth distortion (strain and rotation) fields at interatomic scale is attractive because it provides a basis for a mechanical description of the latter. In such a framework, it is understood that work density is defined and stress transmitted by matter at this scale (Nielsen and Martin, 1985). The attractiveness of such continuous modeling of the GB structure further derives from its ability to serve as a basis for coarse-grained representations of polycrystalline media (Fressengeas et al., 2014) to be used in more macroscopic engineering studies.

Within the recent elasto-plastic framework proposed by Fressengeas et al. (2011), it was found that dislocation/disclination density fields were appropriate objects for a consistent continuous description of the atomic structure of GBs (Cordier et al., 2014; Fressengeas et al., 2014). Based on this work, an atomistic-to-continuum crossover method has been proposed to construct the dislocation/disclination density fields and the associated elastic distortion (strain and rotation), curvature and second-distortion fields in the GB area, starting from the relaxed and unrelaxed atomic configurations of GBs (Sun et al., 2016). In the present paper, we apply this method to build a continuous model of a symmetric grain boundary in MgO, the atomistic structure of which is calculated from first-principles simulations. The strong variations of the displacement, strain and rotation tensor fields across GBs are discussed in details, and crystal defect density fields smoothly describing the structure of the GB are introduced accordingly.

The outline of the paper is as follows. In Section 2, atomic scale modelling of a {310}/[001] GB of MgO, and details of the calculation techniques bridging atomic and continuum representations are briefly recalled. In Section 3, the continuous displacement, distortion, curvature and second-order distortion fields of the GB are calculated from discrete atomic positions. In Section 4, the strong variations of the displacement and rotation fields in the GB area are discussed, and dislocation and disclination density fields are introduced. In Section 5, generalized disclinations are presented to further describe the strong variations in the strain field. Discussion and conclusions follow in Section 6.
## 2. Modelling methods

### 2.1 {310}/[001] GB in MgO from atomistic simulations

The MgO grain boundary studied here is a {310}/[001] symmetric tilt boundary of misorientation $\theta = 36.8^\circ$, which was recently modeled by using first principles methods (Verma and Karki, 2010). The initial positions of atoms were obtained by creating an interface between two MgO single crystals, as described in the following two steps. At first, a bulk crystal of MgO was rotated to align the (310) plane with the GB plane. Then, the crystal was cut in the (310) plane and one part was rotated around [001] by 180°. Because the positions of atoms are mirrored across the GB plane, a symmetric boundary was formed. Details of the method for constructing the boundary are presented in (Verma and Karki, 2010).

After the initial configuration of the GB was generated, the structure was relaxed using the projector-augmented-wave (PAW) method to describe the interactions between valence and core electrons. Periodic boundary conditions were applied in all three directions. Both the atomic positions and the cell parameters were completely optimized for the boundary at zero pressure. Computations were performed by using the Vienna Ab-initio Simulation Package (VASP) (Kresse and Furthmüller, 1996).

The relaxed configuration of MgO {310}/[001] GB is shown in Fig. 1a. Visualization is carried out by using the Visual Molecular Dynamics (VMD) tool (Humphrey *et al.*, 1996). Mg and O atoms are represented by yellow and red spheres, respectively. It is found that Mg/O atoms on one side of the boundary are directly facing O/Mg atoms on the other side, which shows a quasi-continuous structure at the atomic scale. The relaxed structure of the GB can be considered as an array of pipes forming a channel structure. Both the average coordination number and bond lengths decrease at the boundary, as compared to the bulk region. The initial and relaxed configurations are superimposed in Fig. 1b, with their atoms marked as red crosses and blue circles, respectively. As can be seen in Fig. 1a, the GB structure is composed of a periodic sequence of atomic structural units. Atoms near the structural units exhibit the largest displacements while the displacements of atoms far away from the boundary are relatively small.

### 2.2 Linking atomic scale and continuum

The atomic structures of the relaxed and unrelaxed configurations calculated in Section 2.1 are used to build the continuous dislocation/disclination model for the GB. At first, the discrete values of the transformation gradient tensor (Fig. 2) are calculated at each atomic position from finite difference approximations of $F = \partial_x/\partial X$, where the vectors **X** and **x** are the positions in the reference state and deformed state, respectively. In the present case, the reference configuration is chosen to be the relaxed structure and the current configuration is the initial un-relaxed configuration.

Once the transformation tensor is calculated at each atomic position, the elastic strain, rotation, curvature, dislocation, disclination and generalized-disclination density fields are calculated in the boundary area by using the standard relationships of continuum mechanics. Finally, all the atoms are projected onto plane $(\vec{e}_2, \vec{e}_3)$, and two-dimensional linear interpolation is used to generate spatial field distributions in between atoms. It has been checked that the choice of a particular interpolation scheme does not interfere with the results qualitatively. The fundamentals of this atomistic-to-continuum method were presented in (Sun *et al.*, 2016). The method was applied to a 18.9° copper symmetrical tilt boundary as a benchmark test and its accuracy was validated by comparison with a similar recent technique (Gullett *et al.*, 2008; Zimmerman *et al.*, 2009). Additional details are presented in the following.

## 3. Elastic fields
### 3.1. Transformation and distortion tensors

To describe the variations in the shape, size and orientation of a continuous crystalline body containing defects, a fixed Cartesian coordinate system is used, as shown in Fig. 2. For each point in the material, the positions of an element in the reference and current state are described by the Lagrange and Euler coordinates **X** and **x**, respectively. Then the changes in positions are described by the total (elastic) displacement field **u** = **x** − **X**. The components of the displacement vector field **u** in the MgO {310}/[001] boundary area are shown in Fig. 3. The black dashed line in the middle of each component field represents the interface or symmetry plane. The third component $u_3$ is null everywhere. The maximum value of displacements is 0.119 Å, which is 2.8% of MgO lattice constant (4.212Å). The displacement field $u_1$ in $x_1$ direction (Fig. 3a) is symmetric about the interface while that in $x_2$ direction, $u_2$, (Fig. 3b) shows positive/negative extrema about 4 Å away across the interface. Hence, whereas the $u_1$ field is stationary across the interface, the $u_2$ field displays a strong variation by jumping from its maximum negative to maximum positive value in less than 5 Å. To an observer equipped with “mesoscopic” rulers with divisions several times larger than the typical inter-atomic distance, such a strong variation appears as a discontinuity, whereas using a resolution length scale smaller than 1 Å entails regularization and smoothness, as in Fig. 3.

In regions of the body where continuity and differentiability are assumed, the transformation tensor of **x**(**X**) is the second-order Jacobian tensor (Love, 1920; Nadeau, 1964; Nabarro, 1967):

$$ \mathbf{F} = \mathrm{grad} \ \mathbf{x} \text{ or } F_{ij} = \frac{\partial x_i}{\partial X_j} = x_{i,j}. \qquad (1) $$

Similarly, the distortion tensor, can be defined as

$$ \mathbf{U} = \mathrm{grad} \ \mathbf{u} = \mathrm{grad} \ \mathbf{x} - \mathbf{I} = \mathbf{F} - \mathbf{I} \text{ or } U_{ij} = \frac{\partial u_i}{\partial X_j} = \frac{\partial x_i}{\partial X_j} - \delta_{ij} = x_{i,j} - \delta_{ij}, \qquad (2) $$

where **I** is the second order identity tensor and $\delta_{ij}$ is Kronecker delta

$$ \delta_{ij} = \begin{cases} 1, & \text{if } i=j, \\ 0, & \text{otherwise.} \end{cases} \qquad (3) $$

However, the gradient form of the tensors **U** and **F** shown above does not hold across interfaces in the body where the displacement field encounters a discontinuity. Thus, according to the above discussion, the distortion tensor **U** cannot be computed as a gradient tensor across the interface in a “mesoscopic” perspective, that is, if the resolution length scale is larger than the typical inter-atomic distance across the interface. It is only when access to the smallest length scales is allowed in the derivative of the displacement field with respect to $x_2$ that it can be taken as a gradient (“nanoscopic perspective”). The nine components of the distortion tensor field **U** in the boundary area plotted in Fig. 4 are obtained from the method described in Section 2.2 by using interpolation functions and pushing the limit in the $x_2$ derivatives below the inter-atomic distance. According to Eq. (2), only the four components $U_{11}$, $U_{12}$, $U_{21}$ and $U_{22}$ are non-zero (Note that **U** is non-symmetric). As discussed above for the displacement field, the fast sign changes of the $U_{12}$ and $U_{21}$ components across the interface are again to be understood as discontinuities of these fields in the “mesoscopic” perspective, although they show up as smooth fields in the “nanoscopic” perspective.

### 3.2. Strain, rotation and curvature fields
Strain is a measure of the deformation of the body with respect to a reference configuration (Love, 1920; Nadeau, 1964; Nabarro, 1967). In a small deformation approximation, the strain tensor $\boldsymbol{\varepsilon}$ is the symmetric part of the distortion tensor $\mathbf{U}$ and it reads

$$ \boldsymbol{\varepsilon} = \frac{1}{2}(\mathbf{U} + \mathbf{U}^T) \quad \text{or} \quad \varepsilon_{ij} = \frac{1}{2} \left( \frac{\partial u_i}{\partial X_j} + \frac{\partial u_j}{\partial X_i} \right) = \frac{1}{2} (u_{i,j} + u_{j,i}), \qquad (4) $$

when $\mathbf{U}$ is a gradient tensor, whereas the large-strain tensor $\mathbf{E}$, known as the Green-Lagrange strain tensor, is

$$ \mathbf{E} = \frac{1}{2}(\mathbf{F}^T \Box \mathbf{F} - \mathbf{I}) \quad \text{or} \quad E_{ij} = \frac{1}{2}(F_{si}F_{sj} - \delta_{ij}). \qquad (5) $$

In this work, we use the small-strain tensor $\boldsymbol{\varepsilon}$. Of course, the strain tensor cannot be calculated along surfaces where the distortion tensor itself does not exist, because the displacements are discontinuous. The nine components of strain $\boldsymbol{\varepsilon}$ in the boundary area are plotted in Fig. 5 on the basis of the distortion field $\mathbf{U}$ shown in Fig. 4. Like for $U_{11}$ and $U_{22}$, the extremum values of the contractions/ extensions $\varepsilon_{11}$ and $\varepsilon_{22}$ are essentially located on the interface, at the ends of the structural units, while the shear strain $\varepsilon_{12}$ is found all along the boundary with alternate signs across the interface and extrema located extremely close to the interface. Again, such strong variations are to be interpreted as discontinuities in the mesoscopic perspective. Note also that the strength of the contraction/ dilatation strains $\varepsilon_{11}$ and $\varepsilon_{22}$ is larger than that of the shear strain $\varepsilon_{12}$. The trace of the strain tensor is shown in Fig. 5e. It is observed that the left ends of the structural units are squeezed while their right parts are stretched, which indicates that the structural unit tend to close at their left and open at right.

The rotation tensor $\boldsymbol{\omega}$ is defined as the skew-symmetric part of the distortion tensor $\mathbf{U}$. When the latter is a gradient, it is

$$ \boldsymbol{\omega} = \frac{1}{2}(\mathbf{U} - \mathbf{U}^T) \quad \text{or} \quad \omega_{ij} = \frac{1}{2} \left( \frac{\partial u_i}{\partial X_j} - \frac{\partial u_j}{\partial X_i} \right) = \frac{1}{2} (u_{i,j} - u_{j,i}), \qquad (6) $$

and the associated rotation vector $\boldsymbol{\bar{\omega}}$ reads

$$ \bar{\boldsymbol{\omega}} = -\frac{1}{2} \boldsymbol{\omega} : \mathbf{X} = \frac{1}{2} \text{curl } \mathbf{u} \quad \text{or} \quad \omega_i = \frac{1}{2} \sum_{j=1}^{3} \sum_{k=1}^{3} \left( e_{ijk} \frac{\partial u_k}{\partial X_j} \right) = \frac{1}{2} e_{ijk} u_{k,j}, \qquad (7) $$

where $\mathbf{X}$ is the alternating third order Levi-Civita tensor with components $e_{ijk}$.

$$ e_{ijk} = \begin{cases} +1, & \text{if } i, j, k \text{ are an even permutation of } 1, 2, 3, \\ -1, & \text{if } i, j, k \text{ are an odd permutation of } 1, 2, 3, \\ 0, & \text{if any of } i, j, k \text{ are the same.} \end{cases} \qquad (8) $$

The rotation vector field $\bar{\boldsymbol{\omega}}$ in the boundary area obtained from the distortion map in Fig. 4 is shown in Fig. 6. The components $\omega_1$ and $\omega_2$ are zero while, similar to the shear strain field, dipoles in the $\omega_3$ component field appear on both sides of the boundary, located within a 5 Å thick layer along the interface. The maximum local magnitude of $\omega_3$ is found to be the order of 5° (i.e., 0.0886 rad). Hence, the $\omega_3$ field appears as discontinuous across the interface in the mesoscopic perspective, with a 10° discontinuity, although it is seen as a smooth field in the refined description of Fig. 6. When the rotation field $\bar{\boldsymbol{\omega}}$ is continuous and differentiable, the curvature tensor $\mathbf{K}$ is defined as

$$ \mathbf{K} = \text{grad } \bar{\boldsymbol{\omega}} \quad \text{or} \quad \kappa_{ij} = \frac{\partial \omega_i}{\partial X_j} = \omega_{i,j}. \qquad (9) $$
The curvature components $\kappa_{31}$ and $\kappa_{32}$ obtained from Eq. (9) when the limit in the $x_2$ derivative is pushed down to nanoscale are displayed in Fig. 7. This field does not allow rendering the mesoscopic discontinuity of the rotation field because it implies that the rotation gradient exists across the interface. To capture the curvature component $\kappa_{32}$ in the mesoscopic limit, when the rotation field is seen as discontinuous across the interface, we prolong to the interface the last discrete atomistic value $\omega_3$ available in each crystal. In this manner, $\omega_3$ takes two distinct values at the interface depending on which crystal is considered. The curvature $\kappa_{32} = \partial\omega_3/\partial x_2$ at the interface is then obtained by calculating the $x_2$ derivative of $\omega_3$ in each crystal as in Eq. (9), but only in the limit to the top in the bottom crystal, and to the bottom in the top crystal. Since the last discrete atomistic value $\omega_3$ is prolonged to the interface, the changes of $\kappa_{32}$ in $x_2$ direction is zero. Thus, $\kappa_{32}$ turns out to be zero and is continuous at the interface. The corresponding curvature field is shown in Fig. 7c. Note the differences with Fig. 7b.

# 4. Describing the GB structure with Volterra's defects: dislocations and disclinations

## 4.1. Dislocation and disclination density fields

As elaborated in Sections 2 and 3, there are discontinuities across the interface in the fields of displacement $u_2$ (Fig. 3b) and rotation $\omega_3$ (Fig. 6) in the mesoscopic limit. This is illustrated for example in Fig. 8, which provides links between the rotation field and the motion of segments in the structural units. In Fig. 8a, the atoms in the current (e.g., atoms $A_c$ and $B_c$) and reference (e.g., atoms $A_r$ and $B_r$) configurations are marked as red crosses and blue circles, respectively. During the relaxation process, motions including not only translations but also rotations occur, especially in the boundary area. For instance, the segment $A_rB_r$ between atoms $A_r$ and $B_r$, rotates clockwise, from $114.4^\circ$ with respect to the interface in the reference configuration, down to $108.4^\circ$ ($A_cB_c$) in the current deformation. The segment $C_rD_r$ has the opposite anti-clockwise rotation because of the symmetry with respect to the interface. The rotations in Fig. 8a can be quantified by the rotation vector component $\omega_3$, the distribution of which is plotted in Fig. 8b. It is observed that the largest rotations are mostly localized at the boundary, particularly near the structural units. The maximum value of the rotation field is 0.0886 radian (i.e., $5^\circ$), and is located at atoms A and C. As the distance between an atom and the boundary increases to values larger than about 5 Å, the rotation component at this atom nearly vanishes, consistent with Fig. 8a. The section of the field of the rotation component $\omega_3$ along the line AC as a function of the $x_2$ coordinate is plotted in Fig. 8c. When $x_2$ increases from -6 to -4 Å, the value of $\omega_3$ is close to zero and only has slight variations. However, as $x_2$ further increases across the boundary, $\omega_3$ first shows a strong increase and then a sharp drop, which points to discontinuity of the field $\omega_3$, as shown in Figs. 8a and 8b, if the mesoscopic perspective explained above is adopted. Like the rotation $\omega_3$, the displacement $u_2$ in Fig. 8c also shows a rapid jump-drop variation, which substantiates the discontinuity of $u_2$ in the mesoscopic limit.

The discontinuities of the elastic displacement and rotation fields along the interface can be described by Volterra's model defects, i.e., dislocations and disclinations, respectively. Figure 9 illustrates these defects. In Fig. 9a, a straight cut is made along a plane containing the axis of a cylinder in an elastic medium (Fig. 9a) to create a surface defect. Then the surfaces of the cut are shifted. After removing matter to avoid overlaps or adding material to fill the gaps, the cut surfaces are re-welded and different deformation states are produced. According to Volterra (Volterra, 1907), there are only six possible independent cuts and shifts, shown in Figs.
9b-9g, if the cut surfaces are left undeformed. All other cuts can be expressed as a linear superposition of these elementary cuts. If the shift is a translation, the defect is referred to as a dislocation (Figs. 9b-9d), the strength of which, defined as the Burgers vector, is the translation vector. If the shift is a pure rotation (Figs. 9e-9g), the defect is referred to as a dislocation, and the strength of the dislocation is the rotation vector, defined as the Frank vector.

Thus, a discontinuity of the elastic displacement across a bounded patch in the interface is introduced by dislocations. It can be shown that such a discontinuity is reflected smoothly by a continuous dislocation density field, referred to as Nye's dislocation density tensor field (Nye, 1953; Kröner, 1980). Nye's tensor is defined as the curl of the inverse tensor of elastic transformation

$$ \alpha = -\mathrm{curl} \mathbf{F}^{-1} \quad \text{or} \quad \alpha_{ij} = -\sum_{l=1}^{3} \sum_{k=1}^{3} \left( e_{jkl} \frac{\partial D_{il}}{\partial X_k} \right) = -e_{jkl} D_{il,k}, \qquad (10) $$

where $D_{ij}$ are the components of the inverse matrix of $\mathbf{F}$ in a Cartesian reference frame frame ($\vec{e}_1, \vec{e}_2, \vec{e}_3$), Eq. (10) implies that $\mathbf{F}^{-1}$ contains an incompatible<fnref n="1" /> non-gradient part (the curl of this non-gradient is non-zero). Similarly, if standard disclinations are present in the body (deWit, 1970; Fressengeas et al., 2011), the elastic rotation encounters discontinuities. As a result, the curvature tensor contains an incompatible non-gradient part. Its curl yields deWit's smooth dislocation density tensor:

$$ \theta = \mathrm{curl} \mathbf{K} \quad \text{or} \quad \theta_{ij} = \sum_{l=1}^{3} \sum_{k=1}^{3} \left( e_{jkl} \frac{\partial \mathbf{K}_{il}}{\partial X_k} \right) = e_{jkl} \mathbf{K}_{il,k}. \qquad (11) $$

In the reference Cartesian frame, the continuous components $\alpha_{ij}$ (respectively $\theta_{ij}$) of the tensor $\alpha$ (respectively $\theta$) reflect a discontinuity of the elastic displacement $u_i$ (respectively $\omega_i$) per unit surface introduced by a dislocation (respectively dislocation) line having line direction along $\vec{e}_j$. Hence, the tensors $\alpha$ and $\theta$ are continuous renditions of the discontinuities of the elastic displacement and rotation introduced by dislocations and disclinations.

The dislocation and dislocation density fields obtained from the fields in Figs. (4,7) in the MgO boundary area are displayed in Fig. 10. In this figure, the dislocation densities are presented through the fields of Burgers vectors represented by black arrows, and the dislocation density is shown by color-coded contours. Such a density field reflects disclinations with lines along the axis $x_3$ and rotation discontinuities around the same axis. Wedge dislocation dipoles are found at the edges of the structural units. The positive wedge disclinations tend to open the material, while the negative ones tend to close it down, which is consistent with the observation made in Fig. 5e. Both dislocations and disclinations density fields are located within a nanometer thick layer along the boundary.

## 4.2. Frank and Burgers vectors

In order to check the validity of the present dislocation/dislocation model of the boundary, the Frank and Burgers vectors ($\Omega, \mathbf{b}$) are computed. The Frank vector $\Omega$ is defined as the angular closure defect, which can be obtained by integrating the incompatible elastic curvatures along a circuit $C$

$$ \Omega = \int_C \mathbf{K}_e \cdot \mathbf{d}\mathbf{r} = \int_S \theta \cdot \mathbf{n} \, dS \qquad (12) $$

<footnote n="1">¹ In continuum mechanics, (first order) compatibility describes the situation where continuously differentiable strain and curvature fields integrated along paths in the body yield single-valued displacement and rotation fields. Such is not the case along paths encircling dislocation (resp. dislocation) lines (as in Fig. 9), where displacement (resp. rotation) discontinuities are obtained upon integrating the strain and curvature fields. Second-order compatibility will be discussed in the next section.</footnote>
where S is the surface of unit normal **n** delimited by C, **r** is a position vector. In the present MgO case, the Frank vector resulting from the distribution of disclination density θ₃₃ over a surface S in the plane (e₁, e₂) is

$$ \Omega = \int_S \theta \cdot e_3 dS = \Omega_3 e_3 = \int_S \theta_{33} dS e_3. \qquad (13) $$

The above equations lead to scale-dependent measures of lattice incompatibility. The Frank vector Ω has a dependence on the circuit C. By choosing an appropriate surface S, one can measure the Frank vectors of isolated defects as well as those of an ensemble of defects. Consider the circuit delineated by the dashed box in Fig. 1. The length of the box in the x₁ direction should correspond to the period of the lattice, but its height h in the x₂ direction may be varied. The ⟨001⟩ Frank vector component arising from the curvature field is obtained by performing the integration shown in Eq. (13). It is plotted in Fig. 11 as function of the box height h. It first increases but turns out to saturate and to be relatively stable above 2.5 Å. When the box height is 8 Å, the Frank vector component is Ω = 0.132 radian (i.e., 7.56°). Note that the integral ∫∫_[-∞, +∞] κ₃₂ dx₂ = Δω₃ = 36.8° provides the tilt angle of the boundary.

The Burgers vectors for a close circuit bounding a surface S is defined as:

$$ \mathbf{b} = \int (\alpha' - (\theta' \times \mathbf{r})' \cdot \mathbf{n}) dS. \qquad (14) $$

In the plane (e₁, e₂), Eq. (14) is

$$ \mathbf{b} = \int_S (\alpha_{13} - \theta_{33} x_2) dS e_1 + \int_S (\alpha_{23} + \theta_{33} x_1) dS e_2 = b_1 e_1 + b_2 e_2. \qquad (15) $$

Like Frank's vector, the Burgers vector has a dependence on circuit C. From Fig. 11, we find that the Burgers vector components converge to b₁=0.001 Å and b₂= -0.88 Å above 8 Å. The resulting orientation of the Burgers vector, normal to the interface, is consistent with that involved in conventional models of low angle tilt boundaries as edge dislocation walls. However the present model of the MgO tilt boundary is distinct from simple dislocation walls, as it also involves a periodic array of wedge disclination dipoles. The next section suggests that even more complexity is actually involved in such a boundary.

## 5. Beyond Volterra's defects: generalized-disclination density fields

As we have seen above, the boundary exhibits discontinuities in the elastic displacement and rotation fields across the interface. These can be reflected by a continuous distribution of densities of dislocations and disclinations, whose lines are the terminating curves of the surface of elastic displacement u₂ and rotation ω₃ discontinuity, respectively. However, we see here that the {310}/[001] tilt grain boundary also exhibits discontinuities of the shear strain component ε₁₂ and displacement gradient component U₁₂ as shown in Fig. 12. The values of ε₁₂ and U₁₂ along the profile DE in Figs. 12a and 12b are plotted as functions of the coordinate x₂ in Fig. 12c. The fast changes of sign of the components ε₁₂ and U₁₂ across the structural units herald the mesoscale discontinuity of the entire elastic distortion. This cannot be described within the usual framework of Volterra's defects and calls for a generalization of the concept.

In the sole presence of dislocations, the elastic rotation is single valued, the elastic curvature tensor is purely compatible, and only the elastic displacement is a multi-valued function. In the presence of conventional disclinations, the elastic rotation becomes a multi-valued function. In the presence of generalized disclinations (g-disclinations), the discontinuity of the elastic distortion can be described by the incompatibility of the elastic 2-distortion (i.e., the second gradient of displacement in gradient elasticity), as recently discussed by Acharya
and Fressengeas (2011). Figures 9h to 9j sketch some examples of elementary g-disclinations, which are an extension to Volterra's defects (Berbenni et al., 2014). In their presence, the entire elastic distortion tensor (strain and rotation) is multi-valued along a grain boundary. If the discontinuity in elastic strain vanishes, the g-disclinations reduce to conventional disclinations.

The elastic second-distortion tensor **G** is defined as the third-order gradient tensor:

$$ \mathbf{G} = \mathbf{grad} \ \mathbf{U} \quad \text{or} \quad G_{ijk} = \frac{\partial U_{ij}}{\partial X_k} = U_{ij,k} \qquad (16) $$

when the first-order distortion **U** is continuously differentiable. However, the tensor **G** contains an additional incompatible part if discontinuities of the elastic distortion are present. Hence a third-order g-disclination density tensor is defined as:

$$ \pi = \mathbf{curl} \ \mathbf{G} \quad \text{or} \quad \pi_{ijk} = \sum_{l=1}^{3} \sum_{m=1}^{3} \left( e_{klm} \frac{\partial G_{ijm}}{\partial X_l} \right) = e_{klm} G_{ijm,l} \qquad (17) $$

$\pi$ is an areal continuous rendition of the discontinuity of the elastic distortion in a body. In the reference Cartesian frame ($\vec{e}_1, \vec{e}_2, \vec{e}_3$), the components $\pi_{ijk}$ of the tensor $\pi$ express a discontinuity of the distortion $U_{ij}$ per unit surface introduced by a g-disclination line having line direction along $\vec{e}_k$. The associated point-wise measure of incompatibility is the jump of distortion $\Pi$:

$$ \Pi = \int \pi . n dS \qquad (18) $$

Similarly to Eqs. (12) and (14), Eq. (18) is a scale-dependent measure of the incompatibility introduced by g-disclinations. The components of the third-order tensor $\pi$ in the MgO boundary are shown in Fig. 13. From Eq. (17), only two components $\pi_{123}$ and $\pi_{213}$ of $\pi$, shown in Figs. 13a and 13b respectively, are non-zero. They reflect the discontinuities of the in-plane distortions $U_{12}$ and $U_{21}$, respectively. By using Eq. (18) and the same surface $S$, the strength of the g-disclination $\pi_{213}$ dipole is found to be approximately 0.214 radian. The discontinuity of the distortion $U_{12}$ is found to be 0.131 radian, only half of the discontinuity of the distortion $U_{21}$. These discontinuities cover both the rotation discontinuity reflected by the standard wedge disclination distribution shown in Fig. 10 and the shear strain discontinuity.

The g-disclination density tensor $\pi$ complements the standard disclination density $\theta$ by reflecting a discontinuity of strain in addition to the rotation discontinuity, such that both the symmetric and skew-symmetric parts of the distortion are discontinuous. The following relation links both tensors:

$$ \pi = \xi - X.\theta \qquad (19) $$

where the third-order tensor $\xi$ is associated with discontinuities of the elastic strain tensor. According to Eq. (19), the components $\pi_{213}$, $\pi_{123}$, $\xi_{213} = \xi_{123}$ and the standard wedge disclination density $\theta_{33}$ are linked through the relations:

$$ \pi_{213} = \xi_{213} + \theta_{33} \qquad (20) $$

$$ \pi_{123} = \xi_{123} - \theta_{33} \qquad (21) $$

Hence, in addition to the tilt discontinuity of the rotation, the boundary induces a discontinuity of the shear strain $\epsilon_{12}$. To describe the strain discontinuity, the third-order tensor $\xi$ is defined as

$$ \xi = \mathbf{curl} (\mathbf{grad} \ \boldsymbol{\epsilon}) \quad \text{or} \quad \xi_{ijk} = \sum_{l=1}^{3} \sum_{m=1}^{3} \left( e_{klm} \frac{\partial g_{ijm}}{\partial X_l} \right) = e_{klm} g_{ijm,l} \qquad (22) $$
where $g_{ijk} = \partial\epsilon_{ij} / \partial X_k = \epsilon_{ij,k}$, in the manner of tensor $\pi$ in Eqs. (16) and (17). Figure 14 shows the non-zero components of $\xi : \xi_{213} = \xi_{123}$ characterizing the discontinuity of the elastic strain $\epsilon_{12} = \epsilon_{21}$. We note that the discontinuity of the shear strain $\epsilon_{12}$ is found to be locally as high as 0.096, which is far from being negligible and should contribute significantly to the elastic energy of the boundary.

**6. Discussion and concluding remarks**

In this work, we describe a {310}/{[001]} tilt grain boundary in MgO by using an atomistic-to-continuum crossover method. The unit cell of the rock-salt MgO structure can be described as a simple face-centered-cubic (FCC) lattice of Mg ions with O ions occupying all the octahedral sites, or vice versa. The two interpenetrating FCC lattices of Mg and O are displaced from each other by half the basic cell diagonal.

From the atomic positions before and after relaxation, we calculate the continuous displacement and rotation elastic fields. We find that some components of these fields exhibit abrupt jumps across the boundary, from positive to negative values. Such jumps must be interpreted as discontinuities if the spatial resolution is limited to equal-or-above-inter-atomic distances. As such, they can nevertheless be described within the smooth framework of an elasto-plastic theory of dislocation and disclination fields. From the displacement field, the distortion, strain and rotation tensor fields are calculated. The incompatible, non-gradient, part of the distortion yields Nye's dislocation density tensor, which describes the dislocation content of the boundary. If we had used Frank's relation to calculate the Burgers vector, by assuming that the tilt angle is obtained from a linear array of edge dislocations, the relationship between the Burgers vector component $b_2$ and the misorientation $\theta$ would have been $\theta = b_2 / d$, if $d$ is the spacing between two adjacent dislocations. Using the present values, $d = 6.53$ Å and $\theta = 0.642$ radian (i.e., 36.8°), the result would have been $b_2 = d \cdot \theta = 4.2$ Å, a value at odds with the present result. This discrepancy highlights the fact that, at large misorientations, tilt boundaries cannot be consistently described by pure edge dislocation arrays only. Indeed, discontinuities in the rotation field suggest the presence of disclinations. In a similar way, the curvature field is derived from the rotation field, and its incompatible part yields a disclination field made of dipoles at the apex of the structural units. This dipolar structure of the disclination distribution screens the diverging part of the long-range elastic fields of individual disclinations, and provides a low energy configuration for the boundary. Indeed, the continuous fields represented here show that the grain boundary affects only a ca. 1 nm thick layer. Beyond that distance, the crystal is mostly unaffected by the defect.

Our results demonstrate that reconstructing a field description of lattice incompatibility in a tilt boundary from atomic positions is possible. Even for a simple symmetrical tilt boundary, the structure of incompatibility is complex and requires more refined descriptions than expected on the basis of Volterra's defects theory. Indeed, we find that some components of the strain tensor also exhibit discontinuity (in Volterra's theory, the strain tensor is fully continuous around dislocations and disclinations). This result is consistent with recent experimental observations of elastic strain fields along a Σ9(122) symmetric tilt boundary with 38.90° misorientation in Si using high resolution scanning transmission electron microscopy images and geometric phase analysis (Couillard et al., 2013). In this work, the shear strain field also exhibits sign changes across the boundary. Admittedly, the agreement with the present results is only qualitative, because the material and tilt angle are different, but the distributions of positive/negative shear strain $\epsilon_{12}$ and extensions ($\epsilon_{11}, \epsilon_{22}$) across the interface are strikingly similar, as the comparison between Fig. 5 and the experimental shear strain field in Fig. 4 of (Couillard et al., 2013) shows. Therefore, the present results call for independent experimental
investigation of the MgO tilt boundary using similar methods. They also suggest that the
theoretical tools recently proposed by (Acharya and Fressengeas, 2011; Berbenni et al., 2014;
Upadhyay et al., 2013) are needed to account properly for lattice incompatibility across grain
boundaries.

The continuum dislocation/generalized-disclination modeling of GBs can help understand their physical and mechanical properties. For example, Taupin et al. (Taupin et al., 2013) extended the static and discrete disclination structural unit model into a continuous and dynamic description of GBs, and described the plasticity mechanisms through coupled transport of dislocation and disclination densities. The extended disclination model were used to analyze the shear-coupled boundary migration mechanisms in a continuous manner (Taupin et al., 2014), and provided tools for the studies of GB-mediated plasticity. For instance, disclinations were found to decorate GBs in olivine samples and provide the missing mechanism for deforming olivine-rich rocks in the mantle (Cordier et al., 2014). In our future work, the results in this paper will be used as initial inputs for the continuous simulations to study the shear-coupled grain boundary migration.

**Acknowledgments:** This work was supported by funding from the European Research Council under the ERC Grant No. 290424 – RheoMan.

**References**

Acharya, A., Fressengeas, C. (2011): Coupled phase transformations and plasticity as a field theory of deformation incompatibility. *Int. J. Fract.*, **174**, 87–94.

Berbenni, S., Taupin, V., Djaka, K.S., Fressengeas, C. (2014). A numerical spectral approach for solving elasto-static field dislocation and g-disclination mechanics. *Int. J. Solids Struct.*, **51**, 4157–4175.

Cordier, P., Demouchy, S., Beausir, B., Taupin, V., Barou, F., Fressengeas, C. (2014): Disclinations provide the missing mechanism for deforming olivine-rich rocks in the mantle. *Nature*, **507**, 51–56.

Couillard M., Radke G., Botton G.A. (2013). Strain fields around dislocation arrays in a Σ9 silicon bicrystal measured by scanning transmission electron microscopy, *Phil. Mag.*, **93**, 1250-1267.

Duffy, D.M. (2000): Grain boundaries in ionic crystals. *J. Phys. C: Solid State Phys.*, **19**, 4393–4412.

Faul, U.H., Fitz Gerald, J.D., Farla, R.J.M., Ahlefeldt, R., Jackson, I. (2011): Dislocation creep of fine-grained olivine. *J. Geophys. Res.*, **116**, B01203.

Fressengeas, C., Taupin, V., Capolungo, L. (2011): An elasto-plastic theory of dislocation and disclination fields. *Int. J. Solids Struct.*, **48**, 3499–3509.

Fressengeas, C., Taupin, V., Capolungo, L. (2014): Continuous modeling of the structure of symmetric tilt boundaries. *Int. J. Solids Struct.*, **51**, 1434–1441.

Gullett, P.M., Horstemeyer, M.F., Baskes, M.I., Fang, H. (2008): A deformation gradient tensor and strain tensors for atomistic simulations. *Model. Simul. Mater. Sci. Eng.*, **16**, 015001.

Harris, D.J., Watson, G.W., Parker, S.C. (1996): Atomistic simulation of the effect of temperature and pressure on the [001] symmetric tilt grain boundaries of MgO. *Philos. Mag. A*, **74**, 407–418.

Harris, D.J., Watson, G.W., Parker, S.C. (1997): Vacancy migration at the {410}/{001} symmetric tilt grain boundary of MgO: An atomistic simulation study. *Phys. Rev. B*, **56**, 11477–11484.

Harris, D.J., Watson, G.W., Parker, S.C. (1999): Computer simulation of pressure-induced structural transitions in MgO [001] tilt grain boundaries. *Am. Mineral.*, **84**, 138–143.

Humphrey, W., Dalke, A., Schulten, K. (1996): VMD: Visual Molecular Dynamics. *J. Mol. Graph.*, **14**, 33–38.
Jackson, I. (1998): The Earth's Mantle: Composition, Structure, and Evolution. Edinburgh, Cambridge, 592p.

Kresse, G., Furthmüller, J. (1996): Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set. *Comput. Mater. Sci.*, **6**, 15–50.

Love, A.E.H. (1920): A treatise on the mathematical theory of elasticity. Cambridge University Press.

McKenna, K.P., Shluger, A.L. (2009): First-principles calculations of defects near a grain boundary in MgO. *Phys. Rev. B*, **79**, 224116.

Morris, S.J.S., Jackson, I. (2009): Diffusionally assisted grain-boundary sliding and viscoelasticity of polycrystals. *J. Mech. Phys. Solids*, **57**, 744–761.

Nabarro, F.R.N. (1967): Theory of crystal dislocations, Clarendon Press, Oxford.

Nadeau, G. (1964): Introduction to elasticity, Holt, Rinchart and Winston, New-York.

Nielsen O.H., Martin R.M. (1985): Quantum-mechanical theory of stress and force. *Phys. Rev. B*, **32**, 3780–3791.

Pacchioni, G., Pescarmona, P. (1998): Structure and stability of oxygen vacancies on sub-surface, terraces, and low-coordinated surface sites of MgO: an ab initio study. *Surf. Sci.*, **412/413**, 657–671.

Saito, M., Wang, Z., Ikuhara, Y. (2014): Selective impurity segregation at a near-Σ5 grain boundary in MgO. *J. Mater. Sci.*, **49**, 3956–3961.

Sun, X.Y., Taupin, V., Fressengeas, C., Cordier, P. (2016): Continuous description of the atomic structure of grain boundaries using dislocation and generalized-disclination density fields. *Int. J. Plast.*, **77**, 75–89.

Taupin, V., Capolungo, L., Fressengeas, C. (2014): Disclination mediated plasticity in shear-coupled boundary migration. *Int. J. Plast.*, **53**, 179–192.

Taupin, V., Capolungo, L., Fressengeas, C., Das, A., Upadhyay, M. (2013): Grain boundary modeling using an elasto-plastic theory of dislocation and disclination fields. *J. Mech. Phys. Solids*, **61**, 370–384.

Upadhyay, M.V., Capolungo, L., Taupin, V., Fressengeas, C. (2013): Elastic constitutive laws for incompatible crystalline media: the contributions of dislocations, disclinations and G-disclinations. *Philos. Mag.*, **93**, 794–832.

Verma, A.K., Karki, B.B. (2010): First-principles simulations of MgO tilt grain boundary: Structure and vacancy formation at high pressure. *Am. Mineral.*, **95**, 1035–1041.

Volterra, S. (1907): Sur l'équilibre des corps élastiques multiplement connexes. *Ann. Sci. Ecol. Norm. Sup. III*, **24**, 401–517.

Yan, Y., Chisholm, M., Duscher, G., Maiti, A., Pennycook, S., Pantelides, S. (1998): Impurity-Induced Structural Transformation of aMgO Grain Boundary. *Phys. Rev. Lett.*, **81**, 3675–3678.

Zimmerman, J.A., Bammann, D.J., Gao, H. (2009): Deformation gradients for continuum mechanical analysis of atomistic simulations. *Int. J. Solids Struct.*, **46**, 238–253.
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** (Color online) Configuration of the MgO {310}/[001] symmetrical tilt boundary, as obtained from first principles method. (a) Relaxed configuration. Yellow and red balls represent Mg and O atoms, respectively. One atomic structural unit is delineated by blue lines. (b) Reference low-energy relaxed configuration (blue circle) and current high-energy un-relaxed configuration (red cross).</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>**Fig. 2.** (Color online) Illustration of the displacement and transformation gradient associated with the motion of material particles from the reference state to the current state.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>**Fig. 3.** (Color online) Close-up showing the tilt displacement vector **u** on top of the relaxed atomic structure of the simulated {310}/[001] MgO tilt boundary of misorientation 36.8°.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Fig. 4.** (Color online) Close-up showing the displacement gradient tensor **U** along the boundary.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>Fig. 5. (Color online) Close-up showing the strain tensor $\epsilon$ along the boundary.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 6. (Color online) Close-up showing the rotation vector $\vec{\omega}$ along the boundary.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>**Fig. 7.** (Color online) Close-up showing the curvature tensor **κ** along the boundary. (a) **κ**<sub>31</sub>, (b) **κ**<sub>32</sub> and (c) Improved **κ**<sub>32</sub> by prolonging the last discrete atomistic value ω<sub>3</sub> available to the interface in each crystal.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Fig. 8.** (Color online) Continuity/discontinuity of displacement and rotation fields along the boundary. (a) Reference relaxed configuration (blue circle) and current un-relaxed configuration (red cross). (b) Close-up showing the rotation component ω<sub>3</sub> along the boundary. (c) Values of rotation component ω<sub>3</sub> and displacement component u<sub>2</sub> in the profile AC as a function of the coordinate x<sub>2</sub>. Note that, with a small resolution (e.g., less than 1 Å), all 'discontinuities' of the displacement and rotation are smeared out as smoothness of the curves show. At larger resolution (e.g., 3 Å), the rotation jumps discontinuously from its smallest value above the boundary line to largest value below the boundary. Similarly, the displacement jumps discontinuously from its lowest to largest value in a 6 Å distance across the boundary.</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>**Fig. 9.** (Color online) Volterra's distortions and examples of g-disclinations. (a) Defect line and cut surface in a reference cylinder. (b)(c) Edge and (d) screw dislocations show the discontinuities of displacement. (e)(f) Twist and (g) wedge disclinations indicate the discontinuities of rotation. (h) Shear, (i) stretch and (j) strech+rotation g-disclinations.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Fig. 10.** (Color online) Disclination density field $\theta_{33}$ and Burgers vector fields. The arrows represent the local Burgers vector, whose components are the edge dislocation densities ($\alpha_{13}$ and $\alpha_{23}$) per unit surface.</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>**Fig. 11.** (Color online) Values of Frank/Burgers vectors as a function of the height of the surface S shown in Fig. 1.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Fig. 12.** (Color online) Discontinuities of the shear strain and displacement gradient. Close-up showing the (a) shear strain component $\varepsilon_{12}$ and (b) displacement gradient component $U_{12}$. (c) Values of $\varepsilon_{12}$ and $U_{12}$ in the profile DE as a function of the coordinate $x_2$.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>**Fig. 13.** (Color online) Close-up showing the third-order g-disclination density tensor $\pi$ along the boundary.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>**Fig. 14.** (Color online) Close-up showing the third-order tensor $\xi$ along the boundary.</figcaption></figure>