Simulation of Atomic Diffusion in the Fcc NiAl System: A Kinetic Monte

Apr 28, 2015 - National Energy Technology Laboratory, U.S. Department of Energy, 626 Cochrans Mills Road, Pittsburgh, ... URS Corporation, P.O. Box 19...
1 downloads 0 Views 966KB Size
Download PDF
Article pubs.acs.org/JPCC

Simulation of Atomic Diffusion in the Fcc NiAl System: A Kinetic Monte Carlo Study Dominic R. Alfonso*,† and De Nyago Tafen†,

‡

†

National Energy Technology Laboratory, U.S. Department of Energy, 626 Cochrans Mills Road, Pittsburgh, Pennsylvania 15241, United States ‡ URS Corporation, P.O. Box 1959, Albany, Oregon 97321, United States S Supporting Information *

ABSTRACT: The atomic diffusion in fcc NiAl binary alloys was studied by kinetic Monte Carlo simulation. The environment dependent hopping barriers were computed using a pair interaction model whose parameters were fitted to relevant data derived from electronic structure calculations. Long time diffusivities were calculated and the effect of composition change on the tracer diffusion coefficients was analyzed. Our results indicate that this variation has noticeable impact on the atomic diffusivities. A clear reduction in the mobility of both Ni and Al is demonstrated with increasing Al content. Examination of the pair interaction between atoms was carried out for the purpose of understanding the predicted trends.

I. INTRODUCTION Driven by the urgent need to enhance the efficiency of fossil fuel fired power plants, significant research efforts have been devoted to the design of structural materials that are robust under extreme operating environment. Much progress has been accomplished, particularly in the improvement of the thermochemical properties of alloy materials from which vital power plant components are constructed. Despite these advances, technical innovations are still needed to come up with alloys possessing exceptional properties with respect to strength and plasticity. Ni-based alloys with minor Al alloying element are promising base materials because of their mechanical compliance which is crucial for elevated temperature applications.1−3 This high thermal stability makes them particularly suitable for used in turbine blades for power generation. They are also attractive for structural applications under harsh corrosive conditions since the alloying element Al promotes the formation of surface Al2O3 under reactive gaseous oxygen species, as in stainless steel.4,5 The protection imparted by the oxide layer is a result of the sluggish diffusion in the corundum Al2O3 lattice, suppressing its growth at an acceptable level. The presence of such passive layer at the interface between the material and the environment provides resistance to corrosion, thereby improving the service lifetime of the alloy. The formation of external Al2O3 layer depends critically on solid state atomic diffusion occurring in the alloy materials. Although it largely determines whether a stable film would form and at what time scale it would occur, little information is available. The atomic mobility which is quantitatively described by a key kinetic quantity called the tracer diffusion coefficient/ diffusivity can be measured, for example, by secondary ion mass © 2015 American Chemical Society

spectroscopy (SIMS). However, the multicomponent nature and the inevitable interplay of the different metal species in the alloys can make direct experimental measurement challenging. In the case of fcc NiAl alloy, which is of direct relevance to the present study, the unavailability of appropriate isotope for Al poses additional conundrum. Numerical simulations provides a valuable complement to experiment as they are well suited for systematically understanding diffusion trends intrinsically connected with the chemical composition of the alloy. In particular, the application of kinetic Monte Carlo (KMC) method has become more prevalent since it provides an efficient way to fundamentally simulate atomic diffusion process.6−9 By treating diffusion as stochastic transition between locally metastable states, it becomes a decisive toolbox in simulating the interplay of ensemble of thermally activated atomic hops at practical time scales significantly larger than those that can be achieved by molecular dynamics method. Using a replica of the actual crystal structure in which every lattice site has a different environment, it permits the prediction of diffusion behavior that is not captured in approximate analytical model expression.10,11 Additionally, input kinetic information obtained from electronic structure calculations has allowed for improved accuracy in its predictive quality. We previously explored the KMC approach to examine diffusion in Ni containing a minor Fe alloying element.12 Fundamental data obtained from ab initio calculations was used Received: January 23, 2015 Revised: April 27, 2015 Published: April 28, 2015 11809

DOI: 10.1021/acs.jpcc.5b00733 J. Phys. Chem. C 2015, 119, 11809−11817

Article

The Journal of Physical Chemistry C

results wherever possible. Finally, we summarize our results and draw conclusions from them in section 4.

as input to the simulation. In lieu of computationally expensive on-the-fly ab initio determination of diffusion barriers, a computationally efficient environment dependent energy barrier model was employed. The KMC-determined temperature dependent self-diffusion coefficients in pure Ni and tracer diffusion coefficient of dilute Fe in Ni were found to be in very good agreement with experiments. Furthermore, the model predicts that Fe diffusivity is mildly reduced when the concentration of Fe is increased. In the present paper, we report on the application of the same technique to simulate diffusion in the technologically important fcc NiAl alloy. The purpose of this paper is 2-fold. The first one is to provide tracer diffusion coefficients data in the technologically important NiAl alloy. The second purpose is to identify trends in the rates of diffusivities of the metals in the alloy and to what degree the chemical composition influences these quantities. Previous experimental determination of tracer diffusion coefficients in select NiAl alloys exist. For Ni containing dilute Al, secondary ion mass spectroscopy measurements yielded Al tracer diffusion coefficient data.13−15 But their utility could be more beneficial if each effort was not confined to limited temperature range. Extrapolation to other temperatures is typically done by fitting the measured data to an Arrhenius expression. For Ni3Al (alternatively referred to as Ni75Al25), tracer diffusion coefficients of Ni from radio tracer experiments were reported.16−19 However, measurements of the corresponding tracer diffusivity of Al are rare because of the lack of a suitable stable Al radioisotope. Computer simulations using empirical potential based molecular statics and molecular dynamics have been used to study diffusion in Ni3Al.20−24 These works are devoted mainly on fundamental understanding of the mechanisms contributing to diffusion. For example, Debiaggi et al. predicted that a simple vacancy−atom interchange jump is the most energetically favorable for Ni.23 This same mechanism was also found by Duan to be the most favorable for solute Al.24 For calculations of tracer diffusivities, the methodologies employed are however technically limited to rather high temperature range and short periods of real time. Simulation studies using Metropolis Monte Carlo and kinetic Monte Carlo were also reported.25,26 Using an efficient pair interaction model, the dynamics of atomic ordering in Ni3Al was examined as a function of temperature. Neither works, however, provides detailed information about the tracer diffusivities over a wide temperature range. Such calculations have never been performed in the literature to the best of our knowledge. We note that diffusion in NiAl alloys can be understood much more completely by examining postulated mechanisms for Ni and Al jumps in the systems such as those performed in previous studies.20−24 While such characterization should provide a more profound understanding and would be a subject of our future efforts, determination of tracer diffusion coefficients would be a meaningful initial step as it allows us to directly benchmark our computational approach with available experimental diffusivity data. The layout of the paper is as follows. Section 2 contains the technical details of the KMC model utilized in this work. It includes the discussion of procedure for the determination of relevant kinetic parameters incorporated into the model. In particular, the parametrization of the pair interaction model used to represent the environment dependence of the activation barrier was outlined. In section 3, simulation results for atomic diffusion in Ni lattice with varied amount of Al were presented. The predictions were validated with experimental

2. COMPUTATIONAL APPROACH The KMC method adopted was described in our previous studies.12 For clarity, its application to fcc NiAl is summarized here. It is based on the residence time algorithm,27−29 where at each step, the simulation cell is scanned in order to determine all possible hopping processes and their corresponding jump rates ki. The hop j was chosen according to the following selection rule: Σji=1− 1ki ≤ σ1K < Σji=1 ki, where σ1 is a random number between 0 and 1, while K represents the sum of rates of all possible hopping processes. The time at each step was incremented by Δt = −ln(σ2)/K where σ2 is another random number. From the updated configuration, a new process table was recreated and the whole procedure starts again. The atomic configuration evolves by Ni or Al atoms exchanging their positions with a vacancy on the first neighbor site in the lattice. The associated hop rates are ⎛ E ⎞ k T q+ ⎛ E ⎞ ki = νexp exp⎜ − i ⎟ = B vib exp⎜ − i ⎟ 0 h qvib ⎝ kBT ⎠ ⎝ kBT ⎠

(1)

where kB, h, T, and Ei are the Boltzmann constant, Planck constant, temperature, and activation energy of the hop. The terms q+ and qo in the pre-exponential factor vexp represent the vibrational partition function of the transition and initial states given by +(o) qvib

=

∏ i

(

) 1 − exp(− ) hν

exp −0.5 k Ti B

hνi kBT

(2)

with νi corresponding to vibrational mode. The hopping barrier Ei depends on the local atomic configuration and in our framework, a computationally efficient Bronsted−Evans− Polanyi-type based energy barrier model30−32 was adopted to evaluate this quantity: Ei = Eoi + (Efinal − Einit)/2. In this expression, Eio is the intrinsic activation energy which corresponds to the vacancy mediated diffusion barrier of hopping atom in the pure Ni host. Here, we used the DFT calculated values of Eoi for Ni and Al hop. This energy barrier model was used previously to study the kinetic properties of Febased metallic alloys.33−35 It was employed to assess the concentration dependence of the diffusion coefficients of fcc metal alloys.36 Our previous KMC simulation of diffusivities in fcc NiFe alloys was carried out based on this model.12 Recently, this approach has been demonstrated to be a reasonable approximation for NiCr alloy.37 It should be noted that the Bronsted−Evans−Polanyi-type based energy barrier model here is a simplified scheme. In this approach, no contribution to the intrinsic barrier from the local atomic environment is included and it is accounted for only in the (Efinal − Einit)/2 term. A possible improvement would be to introduce this effect into the model. While our adopted model neglects any potential effect of the local environment on the baseline intrinsic barrier Eoi , we assume that the majority of the effect of the local environment on the total activation energy will be captured by this treatment. It may, however, has shortcoming particularly for order-preserving jump of an Al atom through a so-called window formed by 2−4 surrounding Al atoms as shown in previous work.26 But this scenario would entail a priori the formation of a hypothetical configuration where the 11810

DOI: 10.1021/acs.jpcc.5b00733 J. Phys. Chem. C 2015, 119, 11809−11817

Article

The Journal of Physical Chemistry C diffusing Al atom has at least 2−4 first neighbor Al atoms (i.e., presence of unrealistic local clustering of Al in the alloy). Thus, while assuming that a fixed value of the intrinsic barrier as a function of composition may not be an appropriate approximation for Al atom diffusing through an Al-rich local environment, there is evidence that it is reasonable for the alloys we examined here since the formation of such configuration in our simulation would occur with a large energy penalty. The terms Einit and Ef inal are the system energies before and after the jump. They depend on the local atomic configuration (i) and defined as Einit( final) = ΣiΣj>kε(i) jk , where εjk is the ith neighbor pair interaction energy for atoms i and j. Pair interaction up to second neighbor was considered. Our model (i) (i) (i) contains 10 pair interaction energies: ε(i) NiNi, εAlAl, εNiVac, εAlVac (i) (i) (i) and εNiAl; i = 1,2. For εNiNi and εNiVac, we used calculated values that were reported previously.12 The pair interaction between Al atoms ε(i) AlAl was fitted using experimental cohesive energy of Al based on the following relation 2 cohesive EAl = −∑ i

Another kind of model is the use of pair interactions for the saddle configuration.39 While this would yield a refined description, the challenge with such approach is the total number of DFT jump barriers to be calculated and included in the data set to parametrize the pair interactions. The size is determined by the number of neighbor sites included in the definition of the local atomic environment of the saddle point. In the fcc NiAl alloy, this number is 24 (=16) if only the first neighbor sites of the saddle point is considered. Inclusion of further neighbor sites would give rise to an overwhelming number of calculations. For FeCu alloy, for example, this led Le Bouar et al.39 to resort to interatomic potential to build the data set. The use of interatomic potentials, however, restricts the range of alloy materials which can be modeled with such refined scheme. The quantity of relevant interest here is the tracer diffusion coefficient D* which is a useful quantity to characterize mobility. Exploiting the random-walk theory, it can be extracted numerically from KMC simulations through ∞

ni (i) εAlAl 2

D* = lim t →∞

(3)

where terms n1 and n2 are the number of first- and secondneighbor which for fcc crystal are n1 =12 and n2 =6, while Ecohesive is the experimental cohesive energy of the Al metal. The A pair interaction energies between Al and a vacancy, ε(i) AlVac, and between Ni and Al, ε(i) NiAl, were estimated using the following relation b(i) (i) (i) (i) (i) EXY = εNiY + εNiX − εXY − εNiNi ,

i = 1, 2

(6)

In eq 6, t is the simulation time and the vector r(⃗ t + Δt) − r(⃗ t) represents the change in position after a time interval Δt. The angular brackets indicate an expectation value which is obtained by averaging over long-time intervals. In this work, 5 × 108 KMC steps were performed using 5000 sampling passes to collect data on the trajectory of the atoms The total energies required in eq 5 were determined using first-principles spin-polarized DFT calculations as implemented in the Vienna ab initio simulation package (VASP) code.40 The local density approximation (LDA) was used to calculate the exchange-correlation energy41 and the electron−ion interaction was described by the projector-augmented wave (PAW) method.42 The uniform density based assumption of DFT breaks down in describing materials with vacancy since it introduces a steep variation of electronic density near the vacant site. Even the use of slowly varying density based generalized gradient approximation (GGA) does not address this issue. LDA, however, tends to provide a better description of internal surface due to the well-known cancellation effect.43−45 The wave functions were expanded in plane-wave basis sets with an energy cutoff of 400 eV. Our calculations yielded a lattice constant for Ni (3.42 Å) which is within ∼3% of the measured value while the vacancy formation energy for pure Ni was predicted to be 1.63 eV.12 DFT calculations by Delczeg et al.,45 Krcmar et al.,46 and Janotti et al.47 using the same LDA functional yielded values of 1.67, 1.71, and 1.70 eV which compare well with our prediction. The cohesive energies (i) (i) used in the parametrization of εNiNi and εAlAl are the experimental ones since DFT overpredicts these properties. For Ni, the LDA and GGA values are −6.27 and −5.10 eV/ atom while the experimental value is 4.44 eV/atom.48 For Al, the LDA and GGA values are −4.14 and −3.68 eV/atom versus experimental value of −3.39 eV/atom.48 The Ni host containing entities X and Y was represented by a 3 × 3 × 4 periodic fcc supercell. Integrations over the Brillouin zone were carried out using a Monkhorst−Pack 2 × 2 × 2 kpoint grid.49 The coordinates of all the atoms were allowed to relax to determine the equilibrium geometries. The intrinsic activation energy Eoi for Al hop was determined using a host Ni lattice represented by a 3 × 3 × 3 fcc supercell of 108 sites with

(4)

where Eb(i) XY is the binding energy between two entities X and Y in the Ni host. The binding energy between i-th neighbor X and Y in the Ni host in eq 4 is obtained from DFT from the following relation, b(i) EXY = ENi[(N − 2)XY ∞] − ENi[(N − 2)XY (i)]

1 ⟨∑ [ r (⃗ t + Δt ) − r (⃗ t )]2 ⟩ 6t t = 0

(5)

where ENi[(N − 2)XY∞] and ENi[(N − 2)XY(i)] are the total energies of Ni host containing two infinitely separated X and Y and i-th neighbor X and Y, respectively. The binding energy in eq 4 and eq 5 is defined as positive if the interaction between two entities X and Y is attractive. For example, the binding energy between Al and vacancy in the Ni host is attractive between first nearest neighbor (Eb(1) Alvac ∼ +0.07 eV) and repulsive between second nearest neighbor (Eb(2) Alvac ∼ −0.06 eV). The cluster expansion technique with N-body interaction terms can also be used to estimate the Einit and Efinal terms as for instance was done by Lee et al.31 In this case, the database used for numerical fitting consists of DFT predicted formation energies of a set structures with randomly chosen configurations. This more sophisticated approach was also employed in the examination of thermodynamics of interfacial defects in NiAl alloys.38 Though more accurate, cluster expansion remains quite heavy numerically. For the system under investigation in particular, it would require the construction and training of ternary Ni−Al-vacancy cluster expansion containing complicated many body contributions. Its implementation and the deployment within the diffusion model would be difficult and much more time-consuming. It should also be noted that our long-term goal is to model Ni alloys containing more than one kind of solute atoms for which realistic cluster expansion model is very tedious to build. 11811

DOI: 10.1021/acs.jpcc.5b00733 J. Phys. Chem. C 2015, 119, 11809−11817

Article

The Journal of Physical Chemistry C a 2 × 2 × 2 k-point grid. A Ni atom in the lattice was replaced by Al and the diffusion path which consists of exchange of Al atom with first neighbor vacancy along the ⟨110⟩ direction was examined. The climbing image nudged elastic band (CI-NEB) method for finding saddle points and minimum energy paths between known two states was employed.50 Five movable images connecting the initial and final positions were used and the activation energy was obtained from the energy difference between the end point and the image of highest energy. The coordinates of the Ni host and diffusing atom were allowed to fully relax. The vibrational frequencies, νi, in eq 2 were calculated by fixing the host atoms in their relaxed positions and computing the normal-mode frequencies of the diffusing atom (all the 3N degrees of freedom were considered) using the harmonic oscillator approximation.

That is, the larger Al atom has relatively more enhanced compressible directional bonds with the Ni host giving rise to the stabilization of the transition state and concomitant lowering of the diffusion barrier compared to Ni. 3.2. Dilute Al Diffusion in Ni. The KMC approach described in section 2 was now used to calculate tracer diffusivities in fcc NiAl alloy. The DFT computed kinetic parameters are given in Table 1. The tracer diffusivity of dilute Table 1. Calculated Kinetic Parameters, Vibrational Frequencies Entering the Prefactor of the Hop Rates and Pair Interaction Data for Ni and Al Diffusion Processes in the Ni Host process

Eoi (eV)

vo (eV)

v+ (eV)

Ni + Vac → Vac+ Ni 1.22 0.029, 0.029, 0.024 0.036, 0.026 Al + Vac → Vac + Al 0.75 0.043, 0.043, 0.032 0.055, 0.041 pair interaction first neighbor (eV) second neighbor (eV) from

3. RESULTS AND DISCUSSION 3.1. DFT Calculations. The intrinsic activation energy Eoi for Al hop was calculated with DFT and the result is displayed in Figure 1. The diffusing atom moves roughly linearly from the

ε(i) NiNi ε(i) AlAl ε(i) NiAl ε(i) NiVac ε(i) AlVac

−0.59 −0.45 −0.66 −0.18 −0.31

−0.30 −0.23 −0.26 −0.09 0.01

ref 12 eq 3 eq 4 ref 12 eq 4

Al in Ni was first examined. A cubic simulation cell with linear dimension L = Na was used, where a is the lattice parameter and N = 10. Periodic boundary conditions were then applied to simulate the infinite crystal. A vacant site and about 0.1 atomic percent of Al were introduced in the 10 × 10 × 10 periodic cell which was then equilibrated before sampling was carried out. Tracer diffusion coefficients for temperatures in the range 900− 1700 K were calculated. The vacancy concentration under the simulation condition is not necessarily similar to the one under real conditions. To make a direct comparison with experiments, the predicted tracer diffusion coefficient D* was scaled according to the following expression: D = (Cvac,real/C*vac′)D*, where C*vac is the vacancy concentration in the simulation box and Cvac,real is the real equilibrium vacancy concentration.51,52 Since Al is present in dilute amount, Cvac,real was set to the experimental vacancy concentration in pure Ni which is previously determined using differential dilatometry measurements (see inset in Figure 2).53 A straight line was used to fit the experimental data in order to obtain the vacancy concentration at a specific temperature. On the basis of this, Cvac,real =5.65 × 10−8, 4.18 × 10−7, 2.15 × 10−6, 8.40 × 10−6, 2.67 × 10−5, 7.17 × 10−5, 1.69 × 10−4, 3.58 × 10−4 and 6.95 × 10−4 for T = 900, 1000, 1100, 1200, 1300, 1400, 1500, 1600, and 1700 K, respectively. The predicted Al tracer diffusion coefficients DAl at different temperatures are shown in Figure 2. They are plotted on a logarithmic scale as a function of 1000/T. The data lies on a straight line, indicating that Al diffusivity is an Arrhenius function of the temperature. For comparison, the experimental results by Gust et al. based on single-crystal measurements obtained for temperatures ranging from 914 to 1212 K were included.13 The predicted values of DAl compare very well with the measured one. The experimental data essentially lies on a straight line over the whole investigated temperature showing that Al diffusivity is an Arrhenius function of temperature as well. The trend in DAl with temperature (i.e decreasing diffusion coefficients versus inverse temperature) is reproduced. In Figure 2, DAl measured by Swalin et al.14 and Allison et al.15 based on polycrystalline

Figure 1. DFT predicted barrier of solute Al to first-neighbor vacant site in Ni host (along the ⟨1 1 0⟩ direction). Barrier for Ni selfdiffusion is included for comparison. Diffusing atom, vacant site and host atoms are represented by unfilled circle, filled square, and filled circles, respectively.

initial state to the final state recreating a vacancy. The highest point in the diffusion pathway occurs when the atom is about midway to the target vacant site. Thus, the transition state was approximated by this split vacancy configuration, i.e., when the diffusing Al sits roughly midway between the adjacent Al− vacancy pair. Examination of Figure 1 yields a diffusion barrier of about 0.75 eV. Vibrational frequencies were calculated as outlined above and the transition state was verified to have a single imaginary frequency. We previously examined the energy landscape for Ni in the bulk12 and the DFT result was included in Figure 1 for comparison. The diffusion barriers show a noticeable variation with Eoi (Al) < Eoi (Ni). In previous DFT studies of diffusion of 4d and 5d solutes in Ni, it was found that there is a correlation between the solute atomic sizes and their hopping barriers.46,47 That is, the smaller the solute atom, the larger its diffusion barrier. Our calculations concur with this finding as the metallic radius of Al is larger than Ni. A descriptor which is the atomic compressibility was identified as the basis for this trend.46,47 11812

DOI: 10.1021/acs.jpcc.5b00733 J. Phys. Chem. C 2015, 119, 11809−11817

Article

The Journal of Physical Chemistry C

on the diffusivity of Ni. The predicted diffusion coefficients of Ni for this system are within ≤4% of the pure Ni values that we reported previously.12 On the other hand, Al was predicted to have a larger diffusivity compared to Ni. Our KMC simulations yielded D0,Al = 6.35 × 10−5 m2/s which is about five times larger than D0,Ni = 1.29 × 10−5 m2/s. The predicted QAl is slightly lower than QNi (2.73 vs 2.77 eV, see Table 1). This resulted into tracer diffusion coefficients of Al that is about 6−9 times larger compared to Ni. If QNi is expressed as a sum of the vacancy formation in pure Ni and Ni diffusion barrier while QAl is defined as the sum of the vacancy formation in pure Ni, Al diffusion barrier and Al-vacancy binding energy, then there is about 0.4 eV difference between QNi and QAl. It should be noted, however, that the value of QAl based on this description yields a value of QAl = 2.44 eV (using DFT-LDA predicted vacancy formation in pure Ni of 1.63 eV) which is not in quantitative agreement with the reported experimental values of 2.7 ± 0.14, 2.78 ± 0.04, and 2.58 eV (see Table 2). It appears that in this case, such simple characterization of QAl does not lend itself to a plausible explanation of the marked disagreement with experiments. This leads us to a conclusion that aside from the above-mentioned energetic contribution, one may have to consider in detail other underlying factors entering QAl. This would include the contribution of temperature dependent correlation effect to the activation energy that is taken into account in the KMC simulation. Previous studies have shown that solutes with noticeable size and mass misfit with the host atoms have prominent correlation effects.46,47 On the basis of this result, Al which has large size and less mass compared to its host Ni atom (1.43 Å radius and 26.98 amu versus 1.24 Å radius and 58.69 amu) could exhibit such highly correlated motion as well. Since Al was predicted to have lower hopping barrier, there is an increase probability for Al-vacancy exchange to take place in the Ni host. The correlation effects, however, influence whether this hopping event would result in the solute migrating away or returning back to its original lattice site.46 In the KMC simulation, this effect is inherently taken into account and this could be the reason why the predicted QAl with KMC

Figure 2. Tracer diffusion coefficients of dilute Al in fcc Ni computed using kinetic Monte Carlo plotted versus 1000/T. Experimental data based on single crystal (Gust (1981)13) and polycrystalline samples (Swalin (1956)14 and Allison (1959)15) are included for comparison. Experimental vacancy concentration in Ni53 plotted as a function of 1000/T is shown in the inset.

samples are included as well for reference. Each set of data also shows linearity of the Arrhenius relation for the temperatures considered. Quantitative estimates of the diffusion prefactor (D0) and the effective activation energy (QAl) can be obtained by fitting the set of data to the Arrhenius relation, DAl = Do exp(−QAl/kBT). The results are summarized in Table 2. Our KMC results are in very good agreement with that obtained from the single-crystal experimental data of Gust et al. The predicted QAl and D0 values of 2.73 eV and 6.35 × 10−5 m2/s compare quite well with the reported experimental values of 2.7 ± 0.14 eV and 1.0 ± 3.4 × 10−4 m2/s. The tracer diffusion coefficients of Ni in the presence of dilute Al are not reported in the previous experiments. In this work, the presence of dilute Al was found to have a small effect

Table 2. KMC Predicted and Available Experimentally Measured Diffusion Prefactor D0,i and Effective Activation Energy Qi for species i on NiAl alloy D0,Al (m2/s) KMC Gust (1981) Swalin (1956) Allison (1959)

6.35 × 10−5 (1.0 ± 3.4) × 10−4 (1.9 ± 1.4) × 10−4 1.1 × 10−4

KMC

6.25 × 10−5

KMC

1.30 × 10−4

KMC

3.68 × 10−4

KMC

1.51 × 10−2

KMC Frank (1995) Shi (1995) Hancock (1971) Bronfin (1975)

9.72 × 10−2 − − − −

D0,Ni (m2/s)

QAl (eV)

Ni + Dilute Al 2.73 2.7 ± 0.14 2.78 ± 0.04 2.58 Ni95Al5 2.81 Ni90Al10 2.96 Ni85Al15 3.15 Ni80Al20 3.69 Ni75Al25 4.07 − − (1.9 − − 11813

QNi (eV)

T (K)

1.29 × 10−5 − − −

2.77 − − −

900−1700 914−1212 1357−1538 1072−1200

1.72 × 10−5

2.82

1000−1700

2.16 × 10−5

2.87

1000−1700

2.63 × 10−5

2.92

1000−1700

3.07× 10−5

2.96

1000−1700

3.60 × 10−5 1.14 × 10−4 ± 0.5) × 10−4 3.79 × 10−4 1.0 × 10−4

3.00 3.03 3.06 ± 0.03 3.14 2.97

1000−1700 1104−1259 1169−1472 1157−1526 1166−1625

DOI: 10.1021/acs.jpcc.5b00733 J. Phys. Chem. C 2015, 119, 11809−11817

Article

The Journal of Physical Chemistry C differs with that based on simply considering energetic contributions. 3.2. Non-Dilute Al Diffusion in Ni. The calculations were extended to examine the case of nondilute Al diffusion in Ni. Four compositions ranging from 5 up to 25% Al (i.e., Ni containing 0.05, 0.1, 0.15, 0.2, and 0.25 atomic fraction of Al which are alternatively referred to as Ni95Al5, Ni90Al10, Ni85Al15, Ni80Al20, and Ni75Al25) were selected. Al atoms were randomly distributed in the lattice sites in the 10 × 10 × 10 periodic supercell with each configuration equilibrated using KMC as a preparatory step. For Ni75Al25, the respective sublattice sites of the 10 × 10 × 10 cubic cell were filled with Ni and Al so that a starting configuration with an ideal L12 structure was constructed. A single vacancy chosen at random was then introduced to the cell. For Ni75Al25, the enthalpy, ΔHVac,real, and entropy, ΔSVac,real, vacancy formation were measured from previous positron-annihilation spectroscopy experiments.54 To obtain the real vacancy concentration, the following expression was used: CVac,real = exp(−ΔGVac,real/kBT).55 The terms ΔGVac,real and kB correspond to the effective free energy of vacancy formation and the Boltzmann constant, respectively. The free energy term ΔGVac,real defined as ΔGVac,real = ΔHVac,real − TΔSVac,real. From this, the real vacancy concentration as a function of temperature can be ascertained. With no experimental information available for the other compositions, their real vacancy concentration at specific temperature was obtained from linear interpolation between the pure Ni and Ni75Al25 values. The KMC results for Ni diffusivity in Ni75Al25 were first examined. Radiotracer experimental studies diffusion in this alloy was performed by Frank et al.16 They carried out series of secondary ion mass spectroscopy measurements in which they measured the temperature dependence of the tracer diffusion coefficients of Ni in the single crystal sample. The measurements were done at T = 1104−1259 K and the resulting data can be fitted into a single Arrhenius plot, thereby serving as an additional benchmark for the KMC model here. Higher temperature data are also available from other groups but the measurements are based on polycrystalline samples.17−19 Only the tracer diffusivity of Ni was reported in view of unavailability of suitable isotope for Al. Figure 3 shows the comparison. For reference, the experimental results based on the polycrystalline samples are also included. The KMC predictions compare very well with the single crystal data. The temperature dependence of DNi follows to a very good approximation, an Arrhenius law. The single crystal data gives effective activation energy, QNi of 3.03 eV (see Table 2) which is in very good agreement with the predicted value of 3.00 eV. Additionally, the predicted D0,Ni of 3.60 × 10−5 m2/s agrees quite good with the measured value of 1.14 × 10−4 m2/s. This provides support that the KMC approach adopted here can yield acceptable diffusion coefficients at increased Al concentration. For reference, the degree of long-range order (LRO) and short ranged order (SRO) for Ni75Al25 were also determined. For LRO, a Bragg-Williams type parameter, η, was employed varying between 0 (complete disorder) and 1 (perfectly ordered):25

η=1−

Figure 3. Tracer diffusion coefficients of Ni in Ni75Al25 computed using Kinetic Monte Carlo plotted versus 1000/T. Experimental data based on single crystal (Frank (1995),16) and polycrystalline samples (Shi (1995),17 Bronfin (1975),18 and Hancock (1975)19) are included for comparison.

where N(Al) Ni is the number of Ni-antisites occupying the Alsublattice and N(Al) is the number of Al-sublattice in the simulation cell. For SRO, an antisite correlation pair type parameter χ defined as25 χ=

(Al) N Ni

(8)

where represents the number of first neighbor pairs of Ni and Al antisites. On the basis of this definition, χ = 0 for perfect L12-type order in Ni75Al25. In Figure S1 and S2 of the Supporting Information, the η and χ values at different temperatures considered here are plotted versus the KMC steps. For all plots, a qualitatively similar behavior was observed. Both η and χ eventually evolve to their predicted equilibrium level from starting values corresponding to perfectly ordered Ni75Al25. They always reached their equilibrium levels far before the 5 × 108 KMC sampling steps employed in this work. The temperature dependence of the equilibrium values of η and χ are plotted in Figure 4. The Ni75Al25 alloy retains significant long-range ordering for up to 1500 K. It is then followed by a discontinuous drop of ηeq to zero at Tc ∼ 1600 K indicative of order−disorder transformation expected for this alloy. The values of χeq, on the other hand, showed a systematic increase with essentially loss of SRO seen at 1600 K. Qualitatively similar behavior for the temperature dependence of these two order parameters was also seen in a previous work based on a Metropolis Monte Carlo approach.25 Our predicted Tc ∼ 1600 K, however, underestimates their value of ∼1720 K. It should be noted that their model is fitted to specifically yield an order−disorder transition point close to the extrapolated value reported for Ni75Al25.56 Our results indicate that Ni75Al25 alloy exhibits transition from fully L12 ordered to disordered state with some degree of SRO roughly estimated at about ∼1550−1600 K temperature range. A combination of dilatometry and high-temperature diffractomentry measurements on Ni−Al−Fe alloys a hypothetical transition temperature of about 1723 K for stoichiometric Ni75Al25 was estimated by extrapolation to zero N(Al)(Ni) NiAl

(Al) N Ni

0.75 × N (Al)

(Al)(Ni) N NiAl

(7) 11814

DOI: 10.1021/acs.jpcc.5b00733 J. Phys. Chem. C 2015, 119, 11809−11817

Article

The Journal of Physical Chemistry C

Al content induces small increase in the effective activation energy of Ni.

Figure 4. Equilibrium values of LRO and SRO parameters in Ni75Al25 as a function of temperature.

Fe content.56 Our transition temperature estimate of ∼1550− 1600 K deviates from the extrapolated measured value by no more than an order of magnitude. This agreement is fairly reasonable and makes us believe that our simulations provide an adequate description of diffusion in Ni75Al25 considering the fact that we used a relatively simplified pair interaction based energy barrier model. The plot of the temperature dependence of DNi for the various alloys is shown in Figure 5a while the corresponding

Figure 6. Composition dependent values of effective activation energy barrier Qi (i = Ni, Al).

This demonstrates the influence of pair interactions and the observed reduction in Ni mobility can be explained by a (i) comparison of the calculated values of ε(i) NiNi and εNiAl in Table 1. The first-neighbor attraction between Ni and Al (ε(1) NiAl = −0.66 eV) is larger by 0.07 eV compared to that between mutual Ni (1) (εNiNi = −0.59 eV). Though, the corresponding second (2) neighbor ε(2) NiAl= −0.26 eV is smaller than εNiNi = −0.30 eV, the difference is not enough to offset the first-neighbor one. Thus, the overall interaction between Ni and Al atoms is slightly more favorable than between Ni atoms. This has important consequence for the mobility of Ni especially with increasing Al concentration. For example, vacancy mediated exchange in which the hopping Ni moves from an Al-rich to a Ni-rich environment occurs with an energy penalty. The reduction in the number of positions through which Ni can diffuse then affects their ability to hop over long distances, thereby decreasing the diffusivity. Figure 5b shows the temperature dependence of the DAl for the various alloys with the corresponding D0,Al and QAl listed in Table 2. Here again, DAl deviates from dilute values, decreasing with increasing Al content. This time, however, the reduction is found to be more distinct compared to DNi. The values for Ni95Al5 deviate from the dilute case by factors of 0.42−0.61. For Ni90Al10, Ni85Al15, Ni80Al20 and Ni75Al25, the reduction factors are 0.14−0.43, 0.04−0.30, 3.4 × 10−3 −0.28 and 4.4 × 10−4 −0.22 respectively. This effect is likewise reflected in the composition dependence of predicted QAl plotted in Figure 6 which shows larger enhancement in QAl as Al concentration increases. As mentioned in Sec. 3.2, the tracer diffusivity of Al is larger than that of Ni for the dilute case. This does not hold, however, with increasing concentration of Al. In fact, a reversal of the trend is seen for Ni80Al20 and Ni75Al25. For Ni75Al25, in particular, the overall diffusivity of the Al component is relatively smaller with the predicted QAl larger than QNi by about ∼1 eV. Experimental measurement of diffusion properties for Al in this system is not available due to the lack of suitable isotope for use in the SIMS studies. Our work

Figure 5. KMC predicted diffusion coefficients of (a) Ni and (b) the corresponding solute Al in fcc Ni at various concentration plotted versus 1000/T. Results for the dilute Al case (Ni = 99%; Al = 0.1%) are included for comparison.

values of diffusion prefactor and effective activation energy are tabulated in Table 2. A mild decrease in the Ni diffusivity with increasing Al content was found. Compared with the dilute Al case, the values for Ni95Al5 are reduced by a factor of 0.70−0.91 for the temperature range considered. For Ni90Al10, Ni85Al15, Ni80Al20 and Ni75Al25, the reduction factors are 0.48−0.80, 0.34−0.71, 0.26−0.65 and 0.19−0.58 respectively. This effect is likewise reflected in the composition dependence of predicted QNi which is shown in Figure 6. It is found that increase in the 11815

DOI: 10.1021/acs.jpcc.5b00733 J. Phys. Chem. C 2015, 119, 11809−11817

Article

The Journal of Physical Chemistry C represents the first attempt to extract such information from KMC simulation. The trend observed here can be explained by (i) comparing the calculated values of ε(i) AlAl and εNiAl in Table 1. It shows that the former is clearly more repulsive than the latter (2) (1) (2) (ε(1) AlAl = −0.45 eV, εAlAl = −0.23 eV vs εNiAl = −0.66 eV, εNiAl = −0.26 eV). These values indicate that diffusing Al is much more likely to encounter sites that are unfavorable energetically and this effect becomes more pronounced as the content of Al increases. This means that Al hop becomes more prohibitive from an energetic standpoint, and the prospect of further long distance diffusion is somewhat less likely.

any legal liability or responsibility for the accuracy, completeness, or usefulness of any information, apparatus, product, or process disclosed, or represents that its use would not infringe privately owned rights. Reference herein to any specific commercial product, process, or service by trade name, trademark, manufacturer, or otherwise does not necessarily constitute or imply its endorsement, recommendation, or favoring by the United States Government or any agency thereof. The views and opinions of author(s) expressed herein do not necessarily state or reflect those of the United States Government or any agency thereof.

■

4. CONCLUSIONS A first-principles density functional theory based kinetic Monte Carlo framework was implemented to study atomic diffusivity in fcc NiAl alloys. Ni and Al atoms diffusion is considered to be vacancy mediated which occur through first neighbor hop. A fundamental input to the simulation is the diffusion barriers which depend on the local configuration in the vicinity of the diffusing atoms. An efficient pair interaction based energy barrier model was used to determine these barriers. Experimental data and first-principles database were used to fit the parameters. Using the KMC framework, good agreement was demonstrated for those cases where comparison with experimental data was possible. For the case where Al is present in dilute amount, the predicted Al tracer diffusion coefficients results compare well with that obtained from the single-crystal experiments. About the same level of agreement was found for Ni tracer diffusivity in alloy containing 25% Al. Our simulation uncovered the following trend. (i) Increasing the Al content retards Ni diffusivity. This was attributed to the Ni−Al interactions which are overall stronger compared to their Ni− Ni counterpart. This means that the diffusion of Ni occurs with energy penalty. (ii) The retardation in Al diffusivity was found to be more significant with composition variation. Here the more pronounced reduction was attributed to the much more energetically unfavorable Al−Al interaction.

■

(1) Pollock, T. M.; Tin, S. Nickel-Based Superalloys for Advanced Turbine Engines: Chemistry, Microstructure and Properties. J. Propul. Power 2006, 22, 361−374. (2) Sims, S. T.; Stollof, N. S.; Hagel, W. C. SuperAlloys II. Wiley: New York, 1986. (3) Durand-Charre, M. The Microstructure of Superalloys; Gordon and Breach: Amsterdam, 1997. (4) Brady, M. P.; Gleeson, B.; Wright, I. G. Alloy Design Strategies for Promoting Oxide Scale Formation. J. Metals 2000, 52, 16. (5) Heauer, A. H.; Hovis, D. V.; Smialek, J. L.; Gleeson, B. Alumina Scale Formation: A New Perspective. J. Am. Ceram. Soc. 2011, 94, S146. (6) Fichthorn, K. A.; Scheffler, M. Island Nucleation in Thin-Film Epitaxy: A First-Principles Investigations. Phys. Rev. Lett. 2000, 84, 5371−5374. (7) Krishnamurthy, R.; Yoon, Y.-G.; Srolovitz, D. J.; Car, R. Oxygen Diffusion in Ytria-Stabilized Zirconia: A New Simulation Model. J. Am. Ceram. Soc. 2004, 87, 1821−1830. (8) Chen, R.; Dunham, S. T. Kinetic Lattice Monte Carlo Simulations of Interdiffusion in Strained Silicon Germanium Alloys. J. Vac. Sci. Technol. B 2010, 28, C1G18−C1G23. (9) Alfonso, D. Kinetic Monte Carlo Simulation of CO Adsorption on Sulfur-Covered Pd(100). J. Phys. Chem. A 2014, 118, 7306−7313. (10) LeClaire, A. D.; Lidiard, A. B. Correlation Effects in Diffusion in Crystals. Philos. Mag. 1956, 1, 518−527. (11) LeClaire, A. D. Self Diffusion in Dilute Alloys. J. Nucl. Mater. 1978, 69−7, 70−96. (12) Alfonso, D.; Tafen, D. N. Simulation of Diffusion in FCC NiFe Binary Alloys Using Kinetic Monte Carlo Method. J. Phys. Chem. C 2014, 118, 22221−22228. (13) Gust, W.; Hintz, M. B.; Loddwg, A.; Odelius, H.; Predel, B. Impurity Diffusion of Al in Ni Single Crystals Studied by Secondary Ion Mass Spectroscopy (SIMS). Phys. Status Solidi A 1981, 64, 187− 194. (14) Swalin, R. A.; Martin, A. Solute Diffusion in Nickel-Base Substitutional Solid Solutions. J. Metals 1956, 8, 567−572. (15) Allison, H. W.; Samuelson, H. Diffusion of Aluminum, Magnesium, Silicon and Zirconium in Nickel. J. Appl. Phys. 1959, 30, 1419−1424. (16) Frank, S.; Sodervall, U.; Herzig, C. Self-Diffusion if Ni in Single and Polycrystals of Ni3Al. Phys. Status Solidi A 1995, 191, 45−55. (17) Shi, Y.; Frohberg, G.; Wever, H. Diffusion of 63Ni and 114mIn in the Gamma-Phase Ni3Al. Phys. Status Solidi A 1995, 152, 361−375. (18) Bronfin, M. B.; Bulatov, G. S.; Drugova, L. A. Study of Nickel Self Diffusion in Ni3Al and Pure Nickel. Fiz. Met. Metalloved. 1975, 40, 363−366. (19) Hancock, G. F. Diffusion of Nickel in Alloys Based on the Intermetallic Compound Ni3Al. Phys. Status Solidi A 1971, 7, 535−540. (20) Belova, I. V.; Murch, G. E. A Theory of Diffusion in the L12 Structure. I. Tracer Correlation Effects. J. Phys. Chem. Solids 1997, 58, 301−309. (21) Murch, G. E.; Belova, I. V. Diffusion in Nickel-Based Intermetallic Compounds Taking the L12 Structure. J. Mater. Proc. Technol. 2001, 118, 82−87.

ASSOCIATED CONTENT

S Supporting Information *

Time evolution of the LRO and SRO parameters. The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcc.5b00733.

■

REFERENCES

AUTHOR INFORMATION

Corresponding Author

*(D.R.A.) E-mail: [email protected]. Telephone: +1 (412) 386-4113. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work was funded by the Cross-Cutting Technologies Program at the National Energy Technology Laboratory, managed by Susan Maley (Technology Manager) and Charlie Miller (Technology Monitor). The research was executed through NETL’s Office of Research and Development’s Innovative Process Technologies Field Work Proposal. This report was prepared as an account of work sponsored by an agency of the United States Government. Neither the United States Government nor any agency thereof, nor any of their employees, makes any warranty, express or implied, or assumes 11816

DOI: 10.1021/acs.jpcc.5b00733 J. Phys. Chem. C 2015, 119, 11809−11817

Article

The Journal of Physical Chemistry C (22) Belova, I. V.; Murch, G. E. Percolation and the Anti-Structural Bridge Mechanism for Diffusion in Ordered Alloys of the L12 Type. Intermetallics 1998, 6, 403−411. (23) DeBiaggi, S. B.; DeCorte, P. M.; Monti, A. M. Diffusion by Vacancy Mechanism in Ni, Al and Ni3Al: Calculations Based on ManyBody Potentials. Phys. Status Solidi B 1996, 195, 37−54. (24) Duan, J. Atomistic Simulations of Diffusion Mechanisms in Stoichiometric Ni3Al. J. Phys.: Condens. Matter 2006, 18, 1381−1394. (25) Orasmus, P.; Kozubski, R.; Pierron-Bohnes, V.; Cadeville, M. C.; Pfeiler, W. Monte Carlo Simulation of Order-Order Kinetics in the L12 - ordered Ni3Al Binary System. Phys. Rev. B 2001, 63, 174109−1− 174109−14. (26) Leitner, M.; Vogtenhuber, D.; Pfeiler, W.; Puschl, W. Monte Carlo Simulation of Atom Kinetics in Intermetallics: Correcting the Jump Rates in Ni3Al. Intermetallics 2010, 18, 1091−1098. (27) Bortz, A. B.; Kalos, M. H.; Lebowitz, J. L. New Algorithm for Monte Carlo Simulation of Ising Spin Systems. J. Comput. Phys. 1975, 17, 10−18. (28) Voter, A. F. Introduction to the Kinetic Monte Carlo Method. In Radiation Effects in Solids; Sickafus, K. E.; Kotomin, E. A.; Uberuaga, B. P., Eds.; Springer, Berlin, 2007; Vol. 235, pp 1−23. (29) Battaile, C. C. The Kinetic Monte Carlo Method: Foundation, Implementation and Application. Comput. Methods Appl. Mech. Eng. 2008, 197, 3386−3398. (30) Simonovic, D.; Ande, K. Diffusion of Carbon in bcc Fe in the Presence of Si. Phys. Rev. B 2010, 81, 054116−1−054116−7. (31) Lee, E.; Prinz, F. B.; Cai, W. Enhancing Ionic Conductivity of Bulk Single-Crystal Ytria-Stabilized Zirconia by Tailoring Dopant Distribution. Phys. Rev. B 2011, 83, 052301−1−052301−4. (32) Boreskov, G. K. Heterogeneous Catalysis. Nova Science Publishers Inc.: New York, 2003. (33) Vincent, E.; Becquart, C. S.; Domain, C. Microstructural Evolution Under High Flux Irradiation of Dilute Fe-CuNiMnSi Alloys Studied by an Atomic Kinetic Monte Carlo Model Accounting for Both Vacancies and Interstitials. J. Nucl. Mater. 2008, 382, 154−159. (34) Vincent, E.; Becquart, C. S.; Pareige, C.; Domain, C. Precipitation of the FeCu System: A Critical Review of Atomic Kinetic Monte Carlo Simulations. J. Nucl. Mater. 2008, 373, 387−401. (35) Becquart, C. S.; Domain, C. Introducing Chemistry in Atomistic Kinetic Monte Carlo Simulations of Fe Alloys Under Irradiation. Phys. Status Solidi B 2010, 247, 9−22. (36) Swoboda, B.; Van der Ven, A.; Morgan, D. Assessing Concentration Dependence of FCC Metal Alloy Diffusion Coefficients Using Kinetic Monte Carlo. J. Phase Equilib. Diff. 2010, 31, 250−259. (37) Barnard, L.; Young, G. A.; Swoboda, B.; Choudhury, S.; Van der Ven, A.; Morgan, D.; Tucker, J. D. Atomistic Modeling of the OrderDisorder Phase Transformation in the Ni2Cr Model Alloy. Acta Mater. 2014, 81, 258−271. (38) Woodward, C.; van de Walle, A.; Asta, M.; Trinkle, D. R. Firstprinciples Study of Interfacial Boundaries in Ni-Ni3Al. Acta Mater. 2014, 75, 60−70. (39) Le Bouar, Y.; Soisson, F. Kinetic Pathways from EmbeddedAtom-Method Potentials: Influence of the Activation Barriers. Phys. Rev. B 2002, 65, 094103−1−094103−14. (40) Kresse, G.; Furthmüller, J. Effeciency of Ab-Initio Total Energy Calculations for Metals and Semiconductors Using a Plane-Wave Basis Set. Comput. Mater. Sci. 1996, 6, 15−50. (41) Payne, M. C.; Teter, M. P.; Allan, D. C.; Arias, A. T.; Joannopoulus, J. D. Iterative Minimization Techniques for Ab Initio Total-Energy Calculations: Molecular Dynamics and Conjugate Gradients. Rev. Mod. Phys. 1992, 64, 1045−1097. (42) Kresse, G.; Joubert, D. From Ultrasoft Pseudopotentials to the Projector Augmented Wave Method. Phys. Rev. B 1999, 59, 1758− 1775. (43) Carling, K.; Wahnstrom, G.; Mattsson, T. R.; Mattsson, A. E.; Sanberg, N.; Grimvall, G. Vacancies in Metals: From First-Principles Calculations to Experimental Data. Phys. Rev. Lett. 2000, 85, 3862− 3865.

(44) Mattsson, T. R.; Mattsson, A. E. Calculating the Vacancy Formation Energy in Metals; Pt, Pd and Mo. Phys. Rev. B 2002, 66, No. 214110. (45) Delczeg, L.; Delczeg-Czirjak, E. K.; Johansson, B.; Vitos, L. Assessing Common Density Functional Approximations for Ab Initio Description of Monovacancies in Metals. Phys. Rev. B 2009, 80, No. 205121. (46) Krcmar, M.; Fu, C. L.; Janotti, A.; Reed, R. C. Diffusion Rates of 3d Transition Metal Solutes in Nickel by First-principles Calculations. Acta. Mater. 2005, 53, 2369−2376. (47) Janotti, A.; Krcmar, M.; Fu, C. L.; Reed, R. C. Solute Diffusion in Metals: Larger Atoms Can Move Faster. Phys. Rev. Lett. 2004, 92, 085901−1−085901−4. (48) Kittel, C. Introduction to Solid State Physics, 6th ed.; Wiley: New York, 1986. (49) Monkhorst, H. J.; Pack, J. D. Special Points for Brillouin-Zone Integrations. Phys. Rev. B 1976, 13, 5188−5192. (50) Henkelman, G.; Uberuaga, B. P.; Jónsson, H. A Climbing Nudged Elastic Band Method for Finding Saddle Points and Minimum Energy Paths. J. Chem. Phys. 2000, 113, 9901−9904. (51) Terentyev, D. Molecular Dynamics Study of Oxygen Transport and Thermal Properties of Mixed Oxide Fuel Cells. Comput. Mater. Sci. 2007, 40, 319−326. (52) Shewmon, G. Diffusion in Solids. McGraw Hill: New York, 1963. (53) Scholz, H. P. Messungen Der Absoluten Leerstellenkonzentration in Nickel Und Geordneten Intermetallischen Nickel-legierungen Mit Einem Differentialdilatometer. University of Gottingen: Goettingen, Germany, 2001. (54) Badura-Gergen, K.; Schaefer, H.-E. Thermal Formation of Atomic Vacancies in Ni3Al. Phys. Rev. B 1997, 56, 3032−3037. (55) De Koning, M.; Miranda, C. R.; Antonelli, A. Atomic Prediction of Equilibrium Vacancy Concentration in Ni3Al. Phys. Rev. B 2002, 66, 104110−1−104110−8. (56) Cahn, R. W.; Siemers, P. A.; Geiger, J. E.; Bardhan, P. The Order-Disorder Transformation in Ni3Al and Ni3Al Fe Alloys - I. Determination of the Transition Temperatures and their Relation to Ductility. Acta Metall. 1987, 11, 2737−2751.

11817

DOI: 10.1021/acs.jpcc.5b00733 J. Phys. Chem. C 2015, 119, 11809−11817