⇦Chinese Physics Letters, 2017, Vol. 34, No. 7, Article code 076102⇨ Anisotropic Migration of Defects under Strain Effect in BCC Iron * Ning Gao(高宁)1**, Fei Gao(高飞)2, Zhi-Guang Wang(王志光)1 Affiliations 1Laboratory of Advanced Nuclear Material, Institute of Modern Physics, Lanzhou 730000 2Department of Nuclear Engineering and Radiological Sciences, University of Michigan, Michigan 48109, USA Received 20 February 2017 ^*Supported by the National Natural Science Foundation of China under Grant Nos 11375242, 11675230 and 91426301.
^**Corresponding author. Email: ning.gao@impcas.ac.cn Citation Text: Gao N, Gao F and Wang Z G 2017 Chin. Phys. Lett. 34 076102 Abstract The basic properties of defects (self-interstitial and vacancy) in BCC iron under uniaxial tensile strain are investigated with atomic simulation methods. The formation and migration energies of them show different dependences on the directions of uniaxial tensile strain in two different computation boxes. In box-1, the uniaxial tensile strain along the $\langle 100\rangle$ direction influences the formation and migration energies of the $\langle 110 \rangle$ dumbbell but slightly affects the migration energy of a single vacancy. In box-2, the uniaxial tensile strain along the $\langle 111\rangle$ direction influences the formation and migration energies of both vacancy and interstitials. Especially, a $\langle 110 \rangle$ dumbbell has a lower migration energy when its migration direction is the same or close to the strain direction, while along these directions, a vacancy has a higher migration energy. All these results indicate that the uniaxial tensile strain can result in the anisotropic formation and migration energies of simple defects in materials. DOI:10.1088/0256-307X/34/7/076102 PACS:61.80.-x, 61.82.Bg, 61.72.J-, 61.72.Bb © 2017 Chinese Physics Society Article Text Atomic mass transport in materials through the diffusion of simple defects (e.g., self-interstitial atom and vacancy) has been studied for a long time.[1] Without external stress or strain field, the diffusion coefficient $D$ and migration energy $E_{\rm m}$ are usually considered as constants in materials. According to theoretical calculations, the migration energies of a single vacancy and interstitial in BCC iron (Fe) are close to 0.62 and 0.34 eV,[2,3] respectively. For a single vacancy, one step jump to one of eight nearest neighbors is its basic diffusion mechanism. However, for a single interstitial, the transition-rotation (TR) mechanism of a $\langle 110 \rangle$ dumbbell[4] is the main feature for its migration. The parameter $D$ can also be described as a tensor[5] $D_{ij}(\epsilon)= D(0)\delta_{ij} + d_{ijkl}\epsilon_{kl}$, where $d_{ijkl}$ is the elastodiffusion tensor, which can be theoretically calculated or experimentally measured. Therefore, the anisotropic diffusion phenomena induced by the external stress or strain field are generally related with parameter $D$ because of its anisotropic property. However, since $D$ is a function of $E_{\rm m}$, hence except $D$, understanding whether $E_{\rm m}$ can also be described as an orientation dependent variable under stress field, may provide a new physical understanding of anisotropic diffusion. When materials are performed under the extreme condition of radiation in nuclear power, the irradiation induced supersaturated point defects and clusters are inevitably affected by the coupling effect of stress and temperature, resulting in radiation damages and mechanical degradation of materials.[6] The irradiation creep, which occurs under the tensile stress lower than the yield stress, has been found to seriously affect the lifetime of structural materials used in nuclear power.[6-9] The irradiation induced stress corrosion cracking has also been recognized as a critical factor to affect the lifetime of reactors, especially for those that have reached threshold damage levels.[10-12] Since these two damage processes can be described by mass transport at atomic scale under the tensile stress effect,[6] thus to understand the anisotropic diffusion of supersaturated point defects under tensile strain would be useful to explore these irradiation induced phenomena. Although previous studies by DFT have indicated the above influences, only a few cases are considered due to the limitation of DFT.[13] In this Letter, we simulate the anisotropic diffusion of point defects under the tensile strain effect in BCC Fe and determine the orientation dependence of $E_{\rm m}$, to provide possible answers to the above question. The knowledge obtained from BCC Fe in this work is expected to provide useful understanding to these Fe-based structural materials. The empirical Ackland-2004 Fe potential[2] is used in the present work for Fe–Fe interaction, which is able to well reproduce the properties of point defects and the defect properties under the stress effect.[14] Two computational boxes are created with $X$, $Y$ and $Z$ being taken as different coordinate systems: the box-1 with three directions along [100], [010] and [001], and box-2 along [111], [11$\bar{2}$] and [$\bar{1}$10], respectively. Different uniaxial tensile strains are applied on these two boxes. In box-1, because of its symmetrical property, only the tensile strain along the $X$ direction is applied, while in box-2, the tensile strains are applied along $X$, $Y$ and $Z$, respectively. The applied strains along the given direction are 0.5%, 1.0%, 1.5%, 2.0%, 2.5% and 3.0%. In this work, box-1 and box-2 consist of 2000 and 2880 atoms, respectively, to contain a single vacancy or interstitial. After a defect is introduced to the pre-set strain box, the system is firstly relaxed with the molecular static (MS) method. The formation energy of the defect is calculated as a function of the applied strain, $$ E_{\rm f}(D, \epsilon)=E_{\rm t}(D,\epsilon) - NE_{\rm Fe}^{\rm p}(\epsilon),~~ \tag {1} $$ where $E_{\rm t}(D,\epsilon)$ is the total energy of the system containing the defect under strain $\epsilon$, $N$ is the number of atoms in the system, and $E_{\rm Fe}^{\rm p}(\epsilon)$ is the atomic energy per atom in a perfect lattice under strain $\epsilon$. After MS relaxation, the nudged elastic band (NEB) method[15] is used to calculate $E_{\rm m}$ of a single defect under different strains. The molecular dynamic (MD) simulations are continuously performed in the temperature range 0–400 K, in which the velocity scaling method is used in the NVT ensemble. The time-step is 1 fs and the total simulation time for each case is up to 100 ps. A periodic boundary condition is applied along three directions of the computational box. With MD simulations, the kinetic diffusion process of a point defect under strain is explored.

Fig. 1. Formation energy of a single defect under different strains. [(a), (b)] $E_{\rm f}^{\rm v}$ of a single vacancy and $E_{\rm f}^{\rm i}$ of a $\langle110 \rangle$ dumbbell in boxes 1 and 2, respectively. In the figure legend, the integer means the box number. The second parameter is the direction along which the strain is applied. The last value in the legend of (b) indicates the angle between the axial direction of the $\langle110 \rangle$ dumbbell and the strain applied direction. The angle is in units of degree.

Under the uniaxial strain, the formation energies ($E_{\rm f}$) are calculated and shown in Fig. 1. Under the uniaxial strain, $E_{\rm f}^{\rm v}$ increases while $E_{\rm f}^{\rm i}$ decreases with increasing strain in both boxes. One possible reason for the above dependence of $E_{\rm f}$ on $\epsilon$ is that the application of tensile strain increases the volume to form a vacancy and decreases the elastic interaction of two atoms in one dumbbell.[16] Thus it needs more energy to form a larger vacancy and less energy to form one dumbbell. In box-2, $E_{\rm f}^{\rm v}$ has the largest value when the strain is along the $X$ direction while it has the smallest value when the strain is along the $Y$ direction. The value of $E_{\rm f}^{\rm v}$ under the strain along the $Z$ direction is between the values along $X$ and $Y$. Thus $E_{\rm f}^{\rm v}$ is orientation-dependent under the strain effect. For the interstitial, it is different from the vacancy that the orientation-dependent strain field created by a $\langle110 \rangle$ dumbbell could interact with the applied uniaxial strain. Because of its unsymmetrical configuration, the strain field of the $\langle110 \rangle$ dumbbell along its axial direction is higher than along the directions perpendicular to its axis. The value of $E_{\rm f}^{\rm i}$ of the $\langle110 \rangle$ dumbbell with different angles $\theta$, between the axial direction of the dumbbell and the applied strain direction, in two boxes under the strain effect is shown in Fig. 1(b). From these results, it can be found that the difference of $E_{\rm f}^{\rm i}$ does show its anisotropic behavior under the strain effect, which is related with the strain induced change of the elastic field and the interaction angle under given strains along different directions.

Fig. 2. Migration energy of a single vacancy in box-1 under (a) uniaxial strain of 2.0% and (b) uniaxial strain from 0.0 to 3.0%.

The migration energy $E_{\rm m}$ of a single defect is then calculated with the NEB method. For a single vacancy in box-1, Fig. 2(a) shows $E_{\rm m}^{\rm v}$ under the uniaxial strain of 2.0% in box-1 as an example. The value of $E_{\rm f}^{\rm v}$ along the 4 equivalent $\langle 111 \rangle$ directions are the same. The value of $E_{\rm m}^{\rm v}$ is around 0.62 eV, which is the same as the case without strain. With increasing strain, the change of $E_{\rm m}^{\rm v}$ is also small as shown in Fig. 2(b). For a single vacancy in box-2, $E_{\rm m}^{\rm v}$ under different uniaxial strains shows different trends since the 4 $\langle 111 \rangle$ directions are not equivalent under the stain condition of box-2, as shown in Fig. 3. To show clearly the above anisotropic properties, the relative migration energy $E_{\rm m,r}^{\rm v}=E_{\rm m}^{\rm v}(\epsilon)-E_{\rm m}^{\rm v}(\epsilon=0)$ along different directions is shown in Fig. 3. When the uniaxial strain is along the $X$ direction, there are two possible independent directions for a vacancy to jump to. When the migration direction of the vacancy is along the strain or other directions, $E_{\rm m}^{\rm v}$ increases and decreases with the strain, respectively. To understand the above results, the angles $\phi$ between the 4 equivalent $\langle 111 \rangle$ directions and applied strain directions are calculated, which shows that there are two different values of $\phi$, i.e., 0.00$^{\circ}$ and 70.53$^{\circ}$, the same as the dependence of $E_{\rm m}^{\rm v}$ on the jump direction under the uniaxial strain along the $X$ direction. Thus the angle $\phi$ may be used as one of the possible factors to account for the total number of $E_{\rm m}^f$ under the uniaxial strain in box-2. This result is also approved by $E_{\rm m}^{\rm v}$ under the uniaxial strains along $Y$ and $Z$ directions in box-2. When the uniaxial strain is along the $Y$ direction, there are three directions relating with different $E_{\rm m}^{\rm v}$ as shown in Fig. 3, corresponding to three independent values of $\phi$ 90.00$^{\circ}$, 61.87$^{\circ}$ and 19.47$^{\circ}$. When the uniaxial strain is along the $Z$ direction, the similar dependence can also be explored. Two values of $E_{\rm m}^{\rm v}$ are related with two independent values of $\phi$. The underlying physics for the above dependence is related with the effect of anisotropic distribution of the elastic strain fields on the migration path, as described by the continuum elasticity theory[17,18] that if the induced effect becomes stronger, the energy barrier increases, and vice versa.

Fig. 3. Relative migration energy $E_{\rm m,r}^{\rm v}=E_{\rm m}^{\rm v}(\epsilon)-E_{\rm m}^{\rm v}(\epsilon=0)$ of a single vacancy in box-2 under the uniaxial strains conditions. In the figure legend, the first parameter is the direction along which the strain is applied. The second parameter is the site (as shown in Fig. S1 in supplementary materials) to which the single vacancy jumps, indicating the possible migration direction.

Except the single vacancy, the possible anisotropic migration property of a $\langle 110 \rangle$ dumbbell is also calculated. In box-1, since there are two states of the $\langle 110 \rangle$ dumbbell, the migrations from the lower energy state to the lower energy state (path 1) and from the lower energy state to the higher energy state (path 2) are examined. Path 2 is the normal TR mechanism, while path 1 is the direct transition (DT) mechanism. The NEB results of $E_{\rm m}^{\rm i}$ along these two paths as a function of the uniaxial tensile strain are shown in Fig. 4. The results without applied strain are also included, that is, $E_{\rm m}^{\rm i}= 0.34$ eV with $\epsilon= 0$, for comparison. When the uniaxial strain is applied, it is clear that the migration behavior of the $\langle 110 \rangle$ dumbbell has been changed along both paths. The first feature is that $E_{\rm m}^{\rm i}$ increases under the strain, from around 0.42 eV to 0.45 eV along path 1 and around 0.47 to 0.54 eV along path 2, implying that the migration mechanism has been changed from TR to DT. The second feature is the appearance of two peaks on the curves along both paths, especially along path 2. Between two peaks on each curve, a local minimum is also observed. Along path 1, the energy difference between this minimum and the peak is around 9–57 meV with increasing strain, while along path 2, such a difference is around 0.34–0.38 eV. Therefore, the DT is preferred again for migration of the $\langle 110 \rangle$ dumbbell under the uniaxial tensile strain. From Fig. 4, it can also be found that $E_{\rm m}^{\rm i}$ under the uniaxial strain along the $X$ direction is higher than $E_{\rm m}^{\rm i}$ without strain, which is the same as the DFT results,[13] that is, the volumetric expansion increases the energy barrier when the uniaixal tensile strain is applied along the $\langle 100\rangle$ direction.

Fig. 4. Migration energy $E_{\rm m}$ of a single $\langle 110 \rangle$ dumbbell in box-1 under the uniaxial strains conditions. In the figure legend, the strain value and the path for the $\langle 110 \rangle$ dumbbell to migrate to are marked. Paths 1 and 2 are from the lower energy state to the lower energy state and to the higher energy state, respectively.

When the $\langle 110 \rangle$ dumbbell is in box-2, different mechanisms are also examined following the method used in box-1. When the uniaxial tensile strain is applied along the $X$ direction, if the angle $\theta$ between the axial direction of the initial $\langle 110 \rangle$ dumbbell and the $X$ direction is 90.00$^{\circ}$, the $\langle 110 \rangle$ dumbbell would rotate to the $\langle 110 \rangle$ dumbbell with $\theta=35.26^{\circ}$, which has been observed with MD relaxation at 400 K. Such a rotation has changed the mass center of the $\langle 110 \rangle$ dumbbell by following the TR mechanism from the higher energy state to the lower energy state. For example, the [1$\bar{1}$0] dumbbell ($\theta=90.00^{\circ}$) will undergo the TR mechanism to the [011] dumbbell ($\theta= 35.26^{\circ}$) under the uniaxial tensile strain along the [111] direction. The rotation energy barrier $E_{\rm r}^{\rm barrier}$ is also calculated with the NEB method. The results of $E_{\rm r}^{\rm barrier}$ as a function of strains are shown in Fig. 5(a). It is clear that $E_{\rm r}^{\rm barrier}$ decreases with increasing the uniaxial tensile strain, from around 0.36 eV to 0.25 eV, which means that the higher energy $\langle 110 \rangle$ dumbbell can translate to the lower energy state, especially with increasing strain. Thus in the following we mainly examine the possible migrations between two lower energy states. Under the effect of uniaxial tensile strain along the [111] direction, the TR occurs between [110] dumbbell and [011] dumbbell, and the DT occurs between the two [110] dumbbells. The results of $E_{\rm m}^{\rm i}$ with increasing strain are shown in Fig. 5(b) for two mechanisms. For both mechanisms, $E_{\rm m}^{\rm i}$ decreases with increasing the strain. Before $\epsilon=2.5\%$, $E_{\rm m}^{\rm i}$ of the DT mechanism is higher than $E_{\rm m}^{\rm i}$ of the TR mechanism, that is, the TR mechanism is preferred. When $\epsilon=2.5\%$, these two mechanisms have an equal probability. When $\epsilon=3.0\%$, the DT mechanism is preferred. It is also observed from Fig. 5(b) that under the uniaxial tensile strain along the $X$ direction, $E_{\rm m}^{\rm i}$ with the TR mechanism is lower than the value ($\sim$0.34 eV) without the strain, which means that the uniaxial tensile strain along the $X$ direction accelerates the migration of the $\langle 110 \rangle$ dumbbell along the strain direction. The conclusion can then be reached that when the system is under the uniaxial tensile strain along the Burgers vector, the interstitial and vacancy would migrate with lower and higher migration energies along the strain direction respectively. Therefore, when more interstitials and vacancies are formed under the irradiation environment, the irradiation creep strain rate along the Burgers vector will increase under the effect of the uniaxial tensile strain along the same direction. When the uniaxial tensile strain is along $Y$ and $Z$ directions in box-2, the similar orientation-dependent migration and related mass transport process have been observed along the direction close to $Y$ and $Z$ directions, respectively. The rotation and migration energies also decrease with increasing the strain, the same as the case shown in Fig. 5. The underlying physics of the above results in both boxes is considered to be related with the change of elastic interaction of atoms in one dumbbell when the tensile strain is applied. When the applied strains weaken the interaction between two atoms in one dumbbell, the lower energy barrier is obtained, and vice versa.

Fig. 5. (a) Rotation energy barrier $E_{\rm r}^{\rm barrier}$ from higher energy $\langle 110 \rangle$ dumbbell to lower energy $\langle 110 \rangle$ dumbbell under the uniaxial tensile strain along the $X$ direction in box-2. (b) Migration energy $E_{\rm m}^{\rm i}$ between two lower energy $\langle 110 \rangle$ dumbbells in box-2 under the uniaxial strains. The values of $E_{\rm m}^{\rm i}$ of the transition-rotation (TR) mechanism and direction transition (DT) mechanism are shown with black and red lines embedded by different symbols.

In summary, the effect of uniaxial tensile strain on the formation ($E_{\rm f}$) and migration ($E_{\rm m}$) energies of a single vacancy and an interstitial in BCC Fe has been calculated with atomic methods. When the uniaxial tensile strain is along the [100] direction, the weak effect on $E_{\rm f}$ and $E_{\rm m}$ of a single vacancy is confirmed, while for a single $\langle 110 \rangle$ dumbbell, both $E_{\rm f}$ and $E_{\rm m}$ show the orientation-dependent character. The $\langle 110 \rangle$ dumbbell with its axial direction perpendicular to the uniaxial tensile strain direction has the lowest formation energy. The direct transition mechanism is a main mechanism for the migration of the $\langle 110 \rangle$ dumbbells but with a higher migration energy than the case without the strains. When the uniaxial tensile strains along the $X$ [111], $Y$ [11$\bar{2}$] or $Z$ [$\bar{1}$10] directions, the complicated effects on the properties of a single vacancy and interstitial have been obtained. The formation energy is related with the angle between the axial direction of the $\langle 110 \rangle$ dumbbell and the uniaxial tensile strain direction. For migration behavior, the $\langle 110 \rangle$ dumbbell has a lower migration energy when the migration direction is along or close to the strain direction, while along these directions, the vacancy has a higher migration energy, resulting in orientation-dependent mass transport. The change of elastic interaction of atoms in one dumbbell induced by the application of tensile strain is considered as a possible reason to result in the anisotropic migration properties. These results indicate that when the materials are operated under the irradiation condition, the formation of supersaturated vacancies and interstitials may result in anisotropic irradiation creep under the different strain or stress fields. References Ueber DiffusionDevelopment of an interatomic potential for phosphorus impurities in -ironStability and Mobility of Mono- and Di-Interstitials in

?

-FeInterstitials and Vacancies in

?

IronA new mechanism of loop formation and transformation in bcc iron without dislocation reactionDislocation loops and bubbles in oxide dispersion strengthened ferritic steel after helium implantation under stressIrradiation creep models — an overviewCorrosion and stress corrosion cracking in supercritical waterStrain-field effects on the formation and migration energies of self interstitials in

? -Fe

from first principlesEffects of applied strain on radiation damage generation in body-centered cubic ironA climbing image nudged elastic band method for finding saddle points and minimum energy pathsAnisotropic diffusion in stress fieldsThermodynamics of diffusion under pressure and stress: Relation to point defect mechanisms

[1]	Fick A 1855 Ann. Phys. 170 59
[2]	Ackland G J, Mendelev M I, Srolovitz, D J, Han S and Barashev A V 2004 J. Phys.: Condens. Matter 16 S2629
[3]	Fu C C, Willaime F, and Ordejón P 2004 Phys. Rev. Lett. 92 175503
[4]	Johnson R A 1964 Phys. Rev. 134 A1329
[5]	Flynn C P 1972 Point Defects and Diffusion (Oxford: Oxford University Press)
[6]	Was G 2010 Fundamentals of Radiation Materials Science (Berlin: Springer-Verlag)
[7]	Chen J, Gao N, Jung P and Sauvage T 2013 J. Nucl. Mater. 441 216
[8]	Chen J, Jung P, Hoffelner W and Ullmaier H 2008 Acta Mater. 56 250
[9]	Matthews J R and Finnis M W 1988 J. Nucl. Mater. 159 257
[10]	Was G S, Ampornrat P, Gupta G, Teysseyre S, West E A, Allen T R, Sridharan K, Tan L, Chen Y, Ren X and Pister C 2007 J. Nucl. Mater. 371 176
[11]	Andresen P L and Was G S 2012 Compr. Nucl. Mater. 5 177
[12]	Kain V 2011 Stress Corrosion Cracking: Theory Practice (Cambridge: Woodhead Publishing Ltd) p 199
[13]	Chen Z Z, Kioussis N, Ghoniem N and Seif D 2010 Phys. Rev. B 81 094102
[14]	Beeler B, Asta M, Hosemann P and Grønbech-Jensen N 2015 J. Nucl. Mater. 459 159
[15]	Henkelman G, Uberuaga B P and Jonsson H 2000 J. Chem. Phys. 113 9901
[16]	Ullmaier H and Schilling W 1980 Physics of Modern Materials (Vienna: IAEA publication)
[17]	Dederichs P H and Schroeder K 1978 Phys. Rev. B 17 2524
[18]	Aziz M 1997 Appl. Phys. Lett. 70 2810