Journal of Computer Chemistry, Japan
Online ISSN : 1347-3824
Print ISSN : 1347-1767
ISSN-L : 1347-1767
研究論文
Frenkel Defects and Interstitial Atoms in Periclase (MgO) at High Temperature by Molecular Dynamics Simulation
Yasuhiro UEDAAkira MIYAKE
著者情報
ジャーナル フリー HTML

2015 年 14 巻 4 号 p. 117-123

詳細
Abstract

An analysis program was newly developed to explore the vacancies and interstitial atoms and to determine atomic migration, i.e., atomic diffusion. We applied this program to the results calculated by the molecular dynamics (MD) simulation of periclase (MgO) in which Schottky defects (vacancies) were not initially introduced. Generation, migration and extinction of Frenkel defects and interstitial atoms for only magnesium ion were first observed at high temperature in this MD system and they strongly corresponded to the change of the mean square distance (MSD) of magnesium ion in MD system. On the other hand, we could not observe Frenkel defects and interstitial atoms for oxygen ion and MSD value of oxygen ion had almost constant value. Generation, migration and extinction of Frenkel defects and interstitial atoms cannot be ignored for the diffusion process at high temperature.

1 Introduction

Diffusion process is a very important physical property not only in material science but also in earth and planetary sciences and controls many kinetic phenomena such as creep of solid phases and chemical mass transportation. Periclase (MgO) is thought to be one of the important component minerals of the Earth's lower mantle. Therefore, diffusion in periclase has been studied in the material science and earth science. For example, Van Orman et al. [1] measured self-diffusion coefficients of periclase at 15–25 GPa pressures. Moreover, the diffusion coefficient in periclase was also calculated by computer simulation. For example, Vočadlo et al. [2], calculated self-diffusion coefficients by computer simulations based upon lattice dynamics and Ito and Toriumi [3] calculated the pressure effect by molecular dynamics (MD) calculations. Point defects play a very important role in diffusion processes. The important point defect is not only Schottky defects but also Frenkel defects and interstitial atoms. In the previous studies such as Ito and Toriumi [3], however, the only Schottky defects such as vacancies of magnesium and oxygen pair were initially introduced into the MD structure, but Frenkel defects and interstitial atoms have not been mentioned from the simulation results. Furthermore, the migration of individual atoms and vacancies could not be analyzed in detail.

In the present study, we developed newly an analysis program to explore the vacancies and interstitial atoms and to determine the atomic migration. We first report that Frenkel defects and interstitial atoms generate and become extinct at high temperature even by simulation results from an MD crystal without Schottky defects in initial structure.

2 Method

2.1 Molecular dynamics simulation

In the present study, an interatomic potential model was employed, where the interatomic potential function (Φij) between two atoms, i-th and j-th atoms, consists of the Coulombic, the short range repulsion, the van der Waals attraction, and the Morse potential terms as given by   

φ i j ( r i j ) = z i z j e 2 r i j + f 0 ( b i + b j ) exp ( a i + a j r i j b i + b j ) c i c j r i j 6 + D i j { exp [ 2 β i j ( r i j r i j * ) ] 2 exp [ β i j ( r i j r i j * ) ] } ,
where rij is an interatomic distance, f0 ( = 6.9511 × 10−11N) a constant, e the electronic charge, z, a, b and c parameters for each atomic species, and Dij, βij and r*ij parameters for cation-anion pairs. The parameters used in the present study were empirically evaluated by Miyake [4]. MD simulations for silicate minerals such as enstatite (MgSiO3) and cordierite (Mg2Al4Si5O18) were sucessfully carried out using this parameter set in the previous studies (e.g., Miyake et al. [5]; Miyake [6],).

The present MD simulations were carried out by using a program MXDORTO and MXDORTOP (Kawamura [7],), by applying the Ewald method for the summations of Coulombic interactions, and by integrating the atomic motion by means of Verlet's algorithm (Δt = 2.0 femto seconds). Periodic boundary conditions were imposed with the MD basic cell. The MD basic cell was composed of 216 unit cells (aMD = bMD = cMD = 6apericlase) for periclase (containing 1728 atoms). Every MD simulation was carried out by an NPT ensemble, where N is the number of atoms, P pressure and T temperature in the system. Temperature and pressure were controlled by scaling particle velocities and simulating cell parameters, respectively. Temperature range of 300 − 3000 K was investigated, while pressure was kept around 1 atm (0 GPa) throughout the MD simulations. Various crystallographic and physical properties were obtained by the calculations for at least 2 nano seconds (ns) ( = 1,000,000 steps) after annealing for at least 20 ps.

The structure of periclase determined by the X-ray diffraction method (Wyckoff [8],) was adopted as the initial state. Ito and Toriumi [3] prepared pair vacancies of magnesium and oxygen ions (Schottky defects) in MD system of periclase. In the present study, we simulated the perfect MD crystal with no vacancies (no Schottky defects) as the initial state to reveal the appearance of Frenkel defects and interstitial atoms.

2.2 Analysis program to explore the vacancy and interstitial atom and to determine the atomic migration

In MD simulation, we can track the motion of each atom because time variation of all atomic positions in structure can be obtained. At high temperature, an atom oscillates hardly and moves too far from the initial position. Consequently, an atom will migrate from the original position to around the neighbored atom and the generation of vacancy (Frenkel defect) and interstitial atom might be expected. In the previous studies, however, Frenkel defect and interstitial atom only from the time variation of atomic position could not be identified and they had not been mentioned. In the present study, we developed newly an analysis program by introducing the original concept as "room," to explore the vacancy and interstitial atom and to determine the atomic migration. We explain this original program in detail below;

2.2.1 Original concept, “room”

"Room" is defined as the enclosed space whose dimensions with the half distance of the first neighbor atom from the atomic position in perfect structure. Here, the first neighbor atom is defined as same element. The “room” is individually decided on each element and it is clearly distinguished between “room” of magnesium ion and “room” of oxygen ion. Figure 1 shows the 2 dimensional schematic illustration in the case of pericalse.

Figure 1.

 Schematic illustration of (a) MgO crystal, (b) “room” of O ion and (c) “room” of magnesium ion, defined in this study. magnesium ions are yellow and O ions are green.

2.2.2 Determination method of atomic “room” address

For determination of the “room” address where each atom stayed, “room” defined as above was first applied to each snapshot of MD simulation results at each step. We decided basically the address of a “room” to which an atom stays as the atomic “room” address. The schematic illustration was shown in Figure 2.

Figure 2.

 Schematic illustration of (a) one snapshot in MD simulation results and (b) the “room” of O ion. Mg ions are yellow and O ions are green. The “room” with Frenkel defect has no O ion in the "room." The “room” with two O ions has the interstitial atoms.

Here, if the structure changed significantly and the phase transformed to a new phase, we needed to redefine newly “room”. Therefore, we calculated the structure factors of MD simulated structure at each simulated condition using the method proposed by Miyake et al. [9] and identified the phase. In the present study, the phase transition could not be confirmed except for melting over 3100 K.

In addition to these procedures, we corrected the decision method of the atomic “room” address in the following methods, because we distinguish between the temporary migration of an atom to the other “room” due to the atomic thermal vibration at high temperature and the migration for the diffusion process. If the “room” to which the atom x stays at i-1 step (“room-A”) was different from the “room” at i step (“room-B”), we applied the correction as following. At arbitrary j step after i step (i < j),if the atom x did not go back to the “room-A” and the angle ABXj is larger than 90 degrees ( ABXj > 90°), the atom x stays to “room-B” after i step. Here, the center positions for “room-A” and “room-B” were defined as positions A and B, respectively, and the atomic position of atom x was defined as position Xj at j step. ABXj > 90° means that atom x moved beyond the half of “room-B” at j step. At j step after i step, on the other hand, if the atom x went back to the “room-A” without ABXj > 90°, the atom x temporarily migrated to “room-B” and kept to stayed in “room-A” after i step. This correction means that the “room” where atom x stays at i step depends on the atomic position of atom x not only at i step but also after i step. The schematic illustration of this corrected method was shown in Figure 3.

Figure 3.

 Schematic illustration of the correction method of determination of atomic "room."

2.2.3 Criteria for judgment of vacancy (Frenkel defect), interstitial atom and atomic migration

After we decided the “room” to which each atom stayed at each step according to determination method as above at first, we judged vacancy (Frenkel defect), interstitial atom, and atomic migration according to the following criteria.

Frenkel defect: if “room” where no atom stayed to exited at certain step, the “room” contained the vacancy (Frenkel defect) at the step (Figure 2b).

Interstitial atom: if “room” where two or more atoms stayed to exited at certain step, the “room” had one/some interstitial atoms at the step (Figure 2b).

Atomic migration: if a “room” where a certain atom stayed to at certain step was different from a “room” at previous step, we judged that the atom migrated at the step

3 Result and discussion

3.1 Thermal expansion at 0 GPa and compression at 2000 K

Figure 4a shows that the thermal expansion derived from the present MD calculation at 0 GPa was compared with the experiments by an X-ray diffraction method (Dubrovinsky and Saxena [10],). Figure 4b shows that the present compression under high temperature (2000 K) was compared with those determined by previous MD calculation (Ito and Toriumi [3],) and by an X-ray diffraction method (Dewael et al. [11]). Both thermal expansion and compression of the molar volumes are in good agreement with those of experiments. Therefore we consider that reproduction of the present simulation is good and the self-diffusion process simulated by the present simulation is reliable enough to provide the base data compared with the previous study.

Figure 4.

 (a) Thermal expansions of molar volumes at 0 GPa. (b) Compressions of molar volumes at 2000 K.

3.2 Generation, migration and extinction of Frenkel defect and interstitial atom of magnesium ion

Figure 5 shows a sequence of snapshots of a part of the simulated MgO structure at 0 GPa and 3000 K. magnesium ions oscillated around the atomic position in lattice (Figure 5a). Then an magnesium ion abruptly migrated to a neighbor position at certain time. As the result, Frenkel defect generated in the former “room” and the interstitial atom formed in the migrated “room” (Figure 5b). After the interstitial atom repeated to migrate to the other neighbor “room” several times (Figure 5b − 5e), the interstitial atom migrated to “room” with the Frenkel defect and finally both became extinct (Figure 5e).

Figure 5.

 A sequence of snapshots from (a) to (e) of simulated MgO structure at 0 GPa and 3000 K. Mg ions are yellow and O ions are green. The interval between snapshots is each 80 MD steps. The arrows indicate the motion of the Mg ion and the dotted circles approximately show the Frenkel defect.

These results indicated that both the Frenkel defect and interstitial atom generated in the MD structure and that the diffusion process occurred despite the fact that the initial MD cell has no point defects (no vacancy). This phenomenon in MD simulation was the first reported and this proves the usefulness of the analysis program of atomic migration developed in the present study.

3.3 Comparison between mean square distance (MSD) and appearance frequency of Frenkel defects

MSD value of ion calculated by MD simulations has the constant value in the case where the ion merely oscillates around the original atomic position in the lattice and the diffusion process does not occur. On the other hand, MSD value increases linearly with time in the case where the diffusion process occurs. MSD in melt increases precipitously and has the distinctly different gradient. If Frenkel defects and interstitial atoms generated, they migrated respectively and the diffusion occurred, that is the MSD values should increase linearly. When the defects and the interstitial atoms met each other, both became extinct and the diffusion process stopped which meant that the MSD became to have constant value. Figure 6 shows the time dependences of MSDs of magnesium and oxygen ions at 0 GPa and 3000 K and the durations which the Frenkel defects of magnesium ion were identified in the present simulation. MSD of magnesium ion showed that a repeat of the monotonous increasing and the keeping the constant value and had the stepwise change. Frenkel defects and interstitial ions observed by the present program also showed a repeat of the existences and the extinctions (Figure 6). The monotonous increasing durations of MSD of magnesium ion were in excellent agreement with the duration that Frenkel defects and interstitial atoms generated and migrated. Moreover, the durations with almost constant MSD values also were in agreement with the duration that defects and interstitial atoms could not be observed. These results indicated that the stepwise change of MSD value corresponded to the generation, the migration and the extinction of Frenkel defects and interstitial atoms. In the present study, on the other hand, Frenkel defects and interstitial atoms of oxygen ions did not generate at all simulated conditions and therefore the duration was not shown in Figure 6. This got into alignment with that the MSD of oxygen ion having constant value.

Figure 6.

 Time dependences of MSD of Mg and O ions at 0 GPa and 3000 K. Mg ions are yellow and O ions are green. The indicated areas in green show the MD step ranges of the appearance of Frenkel defect for Mg ion.

3.4 Behvior of Frenkel defect and MSD variations with temperature

Figures 7 and 8 show temperature dependences of MSD and the appearance frequency of Frenkel defect for magnesium ion at 0 GPa, respectively. From 1200 to 2500 K, the MSD of magnesium ion gradually increases with temperature. Around 2600 K, however, the value increases precipitously and jumps to higher value around 2900 K (Figure 7). Although Frenkel defects cannot be detected in the present study below 2500 K, the defects appear above 2600 K. The frequency increases exponentially with temperature and jumps to higher value around 2900 K (Figure 8). On the other hand, MSD of oxygen ion gradually increases with temperature from 1200 to 3000 K and Frenkel defects cannot be observed in the present study within the simulated temperature range. Because ionic size of oxygen ion is larger than that of magnesium ion, the formation energy of Frenkel defect and interstitial atom pair for oxygen ion is supposed to be higher than that for magnesium ion. These results show that the MSD value has a strong correlation with the number of Frenkel defects and interstitial atoms. Because Frenkel defects and interstitial atoms interact strongly the diffusion process, the increase of MSD values means the increase of the diffusion coefficient and the formation, the migration and the disappearance of Frenkel defect and interstitial atom cannot be ignored for the diffusion process at high temperature.

Figure 7.

 Temperature dependences of MSD at 0 GPa.

Figure 8.

 Temperature dependences of the appearance frequency of Frenkel defect for Mg ion at 0 GPa.

3.5 Self-diffusion at atmospheric pressure

Self-diffusion coefficient for magnesium ion was estimated from MSD results during the whole simulated time, although the MSD of magnesium ion had the stepwise change (Figure 6). The estimated values, about 3.5 × 10−13 m2/s, at 2700 K are smaller than those calculated by the experimental studies,3.0 × 10−13 m2/s (Wuensch et al. [12]) and 2.6 × 10−12 m2/s (Harding and Price [13]), and the other MD simulation study, 8.5 × 10−12 m2/s (Ito and Toriumi [3]). At 3000 K, moreover, 4.1 × 10−12 m2/s obtained in this study is much smaller than 2.1 × 10−11 m2/s estimated by Ito and Toriumi [3]. Because in these experimental studies and the past MD study, the structure included the Schottky defects as initial state, the self-diffusion processes were easy to occur. In the present system, in contrast, the generation process of vacancies which needs the large activation energy are necessary, self-diffusion coefficient estimated from MSD has a small value. Although the self-diffusion coefficient for oxygen ion could be estimated at 2.5 × 10−13 m2/s from MSD, on the other hand, the diffusion process could not be detected and it could not be compared with past experiments.

From the monotonous increasing period in Figure 6, e.g., from about 19,000 to about 20,000 MD steps, the local self-diffusion coefficient for magnesium ion could be estimated roughly at about 4 × 10−10 m2/s at 3000 K. This period shows that Frenkel defects and interstitial atoms exited and migrated. This estimated value was larger than that calculated by Ito and Toriumi [3], 2.1 × 10−11 m2/s. This suggested that diffusion process readily occurred because both Frenkel defects and interstitial atoms exited in the crystal.

Acknowledgment

We also thank N. Shimobayashi and S. Ohi for helpful comments. This work was partially supported by a Grant-in-Aid for Scientific Research from the Ministry of Education, Science, Sports, and Culture in Japan.

References
  • 1   A. James, Geophys. Res. Lett., 30, 1056 (2003).
  • 2   A. Lidunka Vočadlo, Phys. Earth Planet. Inter., 88, 193 (1995).
  • 3   Yousuke Ito and Mitsuhiro Toriumi, J. Geophys. Res. Solid Earth, 112, B4, B04206 (2007).
  • 4   A. Miyake, Mineral. J., 20, 189 (1998).
  • 5   A. Miyake, N. Shimobayashi, M. Kitamura, Am. Mineral., 89, 1667 (2004).
  • 6   A. Miyake, J. Am. Ceram. Soc., 88, 121 (2005).
  • 7   Katsuyuki Kawamura, "MXDORTO" Japan Chemistry Program Exchange Newsletter, 7(4), 71, (in Japanese) (1996-)
  • 8   Ralph Walter Graystone Wyckoff, Crystal Structures. Vol. 1. (pp. 87). R.E. Krieger Pub Co. (1982)
  • 9   A. Miyake, H. Hasegawa, K. Kawamura, M. Kitamura, Acta Crystallogr. A, 54, 330 (1998).
  • 10   L. Dubrovinsky, S. Saxena, Phys. Chem. Miner., 24, 547 (1997).
  • 11   A. Dewaele, G. Fiquet, D. Andrault, D. Hausermann, J. Geophys. Res., 105, B2, 2869 (2000).
  • 12   B. J. Wuensch, W. C. Steele, T. Vasilos, J. Chem. Phys., 58, 5258 (1973).
  • 13   B. C. Harding, D. M. Price, Philos. Mag., 26, 253 (1972).
 
© 2015 Society of Computer Chemistry, Japan
feedback
Top