ARTICLE pubs.acs.org/JPCC
Atomistic Modeling of Ultrathin Surface Oxide Growth on a Ternary Alloy: Oxidation of AlNiFe Byoungseon Jeon,*,† Subramanian K.R.S. Sankaranarayanan,‡ and Shriram Ramanathan† † ‡
School of Engineering and Applied Sciences, Harvard University, Cambridge, Massachusetts 02138, United States Center for Nanoscale Materials, Argonne National Laboratory, Argonne, Illinois 60439, United States ABSTRACT: By employing variable-charge molecular dynamics, surface oxide film growth on aluminumnickeliron alloys has been studied at 300 and 600 K. The dynamics of oxidation and oxide growth is strongly dependent on the composition of the initial alloy and the ambient temperature. Higher content of Ni and Fe in Al alloys is found to reduce the oxide growth kinetics; 15% Ni þ 15% Fe Al alloy yielded 3040% less growth at 400 ps oxygen exposure compared to pure Al. We observe dopant segregation, which disrupts the interaction between O atoms and Al atoms in the alloy, leading to a nonlinear oxide growth profile in the case of ternary AlNiFe alloy. Compared to oxidation at 300 K, 30% more oxide layer was yielded at 600 K, due to the elevated temperature. The simulated oxide kinetics indicates that the growth rate of anion surpasses the cation rate with higher sensitivity to the stoichiometry of the base metal substrate. Charge state analysis provides insights into the evolution of cation and anion species as the oxide layer grows. In particular, due to higher correlation, Fe shows a high rate of oxidation when the content is high, whereas the rate of Ni oxidation is consistently low. Density profile analysis suggests the segregation of dopant atoms below the growing ultrathin oxide layer, showing the presence of a layer-by-layer mode of oxide layer even with disordered structure. Coordination number (Z, the number of oxygen atoms around an aluminum atom) of aluminum oxide has been used to identify how the initial oxidation transitions into equilibrated states. Z = 3 is dominant in the early stages of oxidation and at the interface between oxide and bulk substrate, but it transitions quickly to Z = 4 (∼45%) and 5 (∼35%) as the oxide equilibrates and approaches its self-limiting thickness. Even though growth kinetics is dependent on the base metal stoichiometry, the composition of the oxide microstructure is not significantly affected, primarily segregating dopant elements, i.e., Ni and Fe outside of the oxide layer.
’ INTRODUCTION Ultrathin oxide formation on metal alloys is an important topic in both science and technology.1,2 The science is complicated by thermodynamics and kinetics of oxidation of different constituents and electronic properties of the oxides formed, while the physical and chemical nature of the oxide films determines several functional properties, such as resistance to corrosion or breakdown in aqueous media. In order to understand this early stage growth behavior and the effect of ambient conditions, recently there have been theoretical and computational works with aluminum or aluminum alloy system. Streitz and Mintmire3 developed a potential model for metal oxide, employing the embedded atom method4 and charge equilibration approach.5 They developed variable charge molecular dynamics simulations to model reactive systems such as aluminum surface oxidation, which yielded structural data of fcc- and Ralumina. Campbell et al.6 studied the oxide dynamics of aluminum nanoclusters, confirming the oxide thickness and atomic valence charges. This work is one of the large-scale parallel variable charge molecular dynamics, with the system size involving about a million particle interactions. A similar nanoparticle system has been studied by Alavi et al.7 as well. Zhou et al.8,9 investigated the r 2011 American Chemical Society
charged states and growth of metal oxide with planar structure, developing metallic potential with a charge transfer feature. Their approach is also employed in this work. Hasnaoui et al.10 reported that growth kinetics is independent of crystal orientation of aluminum lattices, with the growth occurring by a layerby-layer mode. They developed a computational model to localize a planar metal substrate for controlled oxide growth. Upon simulating the oxidation of the [100] orientation of crystalline Al, the AlO bond length was found as 1.8 Å and AlO4 was reported as the most abundant structure. Oxidation in polycrystalline Al has been studied by Perron et al.,11 providing the correlation between the crystallinity of the oxide film and temperature. Sankaranarayanan et al.12,13 investigated the early growth of oxide film, in particular, with AlNi binary alloy and different ambient conditions, such as temperature, pressure, and electric field. It turned out that the oxide in AlNi alloy is entirely comprised of alumina, while Ni oxide is formed only at high Received: November 8, 2010 Revised: February 16, 2011 Published: March 17, 2011 6571
dx.doi.org/10.1021/jp1106845 | J. Phys. Chem. C 2011, 115, 6571–6580
The Journal of Physical Chemistry C
ARTICLE
pressures. Georgieva et al.14 analyzed the oxide thin film growth of MgxAlyOz, but they employed fixed charge molecular dynamics, using formal and partial charges. On the other hand, density functional theory (DFT) has been employed to investigate the energetics of oxide film on binary Al metal alloys. Initial stages of oxidation on NiAl(110) surface has been analyzed by Lozovoi et al.,15 showing that the expulsion of Ni from oxide layer is favored for the surface stability. Kresse et al.16 investigated the favorable structure of oxide on NiAl(110) and reported that Al10O13 is a more stable structure than Al2O3. Islam et al.17 analyzed the energetics of TiAl(111) oxide interface and calculated vacancy formation energy. However, those employed systems are very small ( i φijðrij Þ þ ∑i Fi ðFi Þ
ð1Þ
Embedded atom potential describes the interaction of particles as a sum of pair interaction (φ) and embedding energy (F). Pair
interaction is described using pair distance (rij) as rij rij A exp R 1 B exp β 1 re re φðrij Þ ¼ 20 20 rij rij 1þ 1þ k λ re re ð2Þ Embedding energy is a function of electron density (F) and the local density arisen by particle j at i is rij fe exp β 1 re fj ðrij Þ ¼ ð3Þ 20 rij 1þ λ re Then the electron density at particle i is summed as Fi ¼
∑i j6∑¼i fjðrijÞ
ð4Þ
Finally, the embedding energy is calculated as !m Fi FðFi Þ ¼ Fm 1 Fe m
∑
ð5Þ
The embedding function has been fitted and splined over a wide 9 range of electron density. A, B, R, β, κ, λ, re, and Fe are fitting parameters, and they are enlisted in Tables 1 and 2. For the implementation of variable charge and charge transfer, we use the method by Streitz and Mintmire.3 In addition to the basic formula, we use penalty energy terms which has been suggested by Zhou et al.,8 in order to constrain the variations of charges. N
1 2 a 2 χi qi þ Ji qi kc pffiffiffi qi þ 2 π
∑i ∑G ΨðGÞ SðGÞ Sð GÞ þ ∑ ∑ ðkc qi Zj f½jjfi ½fi jfj g þ kc qj Zi f½ijfj ½fi jfj gÞ i j>i
Ees ¼
þ
∑i ∑
(
1 kc qi qj erfcðarij Þ þ rij j>i
)!
1 ½fi jfj rij
∑i ωU½ ðqi qi, min Þðqi qi, min Þ2 þ ∑ ωU½qi qi, max ðqi qi, max Þ2 i
þ
ð6Þ
Electrostatic energy is described by the interaction of charged particles, same as constant charge calculations, but there are additional components, such as electronegativity (χ), atomic hardness (J), effective core charge (Z), and Coulomb integration ([j|f],[f|f]). Coupled with periodic boundary conditions, the calculation of long-range interaction requires a specific numerical scheme like Ewald summation. In this study, we use smooth particle mesh Ewald (sPME)18 to expedite the calculation of reciprocal space components [∑ψ(G) S(G) S(G)]. U is a unitstep function and constrains the range of charge variations. ω is the energy penalty component, and 20 eV/e2 is employed.19 To find the equilibrated charges, iterative conjugate gradient20,21 is employed with the Brent method.22 Charge update is the most 6572
dx.doi.org/10.1021/jp1106845 |J. Phys. Chem. C 2011, 115, 6571–6580
The Journal of Physical Chemistry C
ARTICLE
Table 1. EAM Potential Parameters of Simulated Elements element
re (Å)
fe
Fe
Fs
R
β
A (eV)
Al
2.863 92
1.203 78
17.517 47
19.900 41
6.613 17
3.527 02
0.314 87
Ni
2.488 75
2.211 49
30.370 03
30.371 37
8.383 45
4.471 17
0.429 05
Fe
2.481 99
2.314 53
24.595 73
24.595 73
9.818 27
5.236 41
0.392 81
element
B (eV)
κ
λ
Fn0 (eV)
Fn1 (eV)
Fn2 (eV)
Fn3 (eV)
Al
0.365 55
0.379 85
0.759 69
2.807 60
0.301 44
1.258 56
1.247 60
Ni
0.633 53
0.443 60
0.820 66
2.693 51
0.076 44
0.241 44
2.375 63
Fe
0.646 24
0.170 31
0.340 61
2.534 99
0.059 60
0.193 06
2.282 32
element
F0 (eV)
F1 (eV)
F2 (eV)
F3 (eV)
F3þ (eV)
η (eV)
Fe (eV)
Al
2.83
0.0
0.622 25
2.488 24
2.488 24
0.785 91
2.824 53
Ni
2.70
0.0
0.265 39
0.152 86
4.585 68
1.013 18
2.708 39
Fe
2.54
0.0
0.200 27
0.148 77
6.694 65
1.182 90
2.551 87
Table 2. Paired Interaction Parameters of Embedded Atom Potential pair
re (Å)
R
β
A (eV)
B (eV)
κ
λ
OO OAl
3.648 57 2.985 20
5.440 72 8.497 41
3.597 46 4.521 14
0.349 00 0.097 38
0.574 38 0.381 21
0.080 07 0.189 67
0.393 10 0.952 34
ONi
2.957 32
7.965 28
4.424 11
0.135 21
0.253 32
0.470 77
0.655 24
OFe
3.079 92
7.523 09
4.133 30
0.171 08
0.398 69
0.223 35
0.343 80
AlNi
2.715 79
8.004 43
4.759 70
0.442 54
0.683 49
0.632 79
0.817 77
Table 3. CTIP Parameters of Simulated Elements element
χ (eV/e)
J (eV/e2) ξ (Å1)
Z (e)
qmin (e) qmax (e)
2.144 0.000 00
2
0
Al
1.47 9 14
9.072 22
0.968 1.075 14
0
3
Ni
1.708 04
9.109 54
1.087 1.444 50
0
2
Fe
1.905 87
8.998 18
1.024 1.286 12
0
3
O
2.000 00 14.995 23
time-consuming routine and is done every 10 time steps, for efficient computing. Most of potential parameters are from the work of Zhou et al., and a detailed description of potential components can be found in the references.9,19 We briefly list the relevant parameters in Table 3. Implementing those features, a variable-charge molecular dynamics code has been developed. Velocity Verlet time integration is employed with a stochastic or Berendsen thermostat, and an ffte2328 library is used as an external library for Fourier transform.
Figure 1. Configuration of the simulation unit cell. Vacuum layers are located at the top and bottom of the metal substrate. New oxygen atoms are inserted from the top of the unit cell.
’ SIMULATION OF METAL OXIDES Basic orientation of oxide surface is [100] comprised of fcc aluminum lattice. Because of planar structure, reduced periodicity is required in the molecular dynamics simulations of metal oxides. Hasnaoui et al.10 introduced a vacuum layer to prohibit the effect of periodic image along the nonperiodic direction, and this method has been proved to be effective in previously reported studies.12,13,29 This work employs the same approach and the basic configuration of the unit cell is shown in Figure 1. The planar structure is repeated along the x/y direction, but vacuum layers are located at the top and the bottom of substrate. For the long-range interactions, we use standard 3D Ewald
summation with smooth particle mesh Ewald for reciprocal summation, but the employed vacuum layer will allow us to simulate it as 2D periodic systems. All of the calculations employed 9.0 Å cutoff radius for both electrostatic and EAM potentials. The simulated alloys are summarized in Table 4: mono refers to pure aluminum and will be a reference, compared to other binary or ternary alloy systems; binA and binB are binary alloys made of Al and Ni atoms, using 95/5% and 85/15% stoichiometry in atomic percentage; ternA and ternB are ternary alloys containing Fe as well, and compositions are 90/5/5% and 70/ 15/15%, respectively. The size of the simulation unit cell is 20 6573
dx.doi.org/10.1021/jp1106845 |J. Phys. Chem. C 2011, 115, 6571–6580
The Journal of Physical Chemistry C
ARTICLE
Table 4. Stoichiometry and Configurations of Metal Oxide Simulations element
mono
binA
binB
ternA
ternB
Fractions Al (atom %)
100
95
85
90
70
Ni (atom %) Fe (atom %)
0 0
5 0
15 0
5 5
15 15
Number of Atoms Al (ea)
1000
950
850
900
700
Ni (ea)
0
50
150
50
150
Fe (ea)
0
0
0
50
150
Figure 3. The growth of oxide layer (left) and the number of oxygen atoms (right) for 600 K simulation.
Figure 2. The growth of oxide layer (left) and the number of oxygen atoms (right) for 300 K simulation.
20 80 Å3, and the smooth particle mesh Ewald employed 24 24 96 grids for the charge interpolation of reciprocal summation. Between the new oxygen insertion, oxidation is simulated as NVT, with 300 or 600 K temperature. Procedure of the Oxidation Simulations. An initial configuration is produced from the perfect lattice of fcc structure. Simulation of pure aluminum (mono) substrate may begin from this configuration, but binary or ternary alloy needs an additional treatment in order to reach minimum energy configuration. For a given alloy composition, it is important to predict the distribution of constituent atoms. In order to achieve better statistics of rare events and equilibrated initial conditions, Monte Carlo (MC) simulations30 have been employed with site exchange strategy. Every MC step tries a randomized displacement of a randomly selected atom or site exchange of two atoms of different elements31 with EAM potential only. Before the actual oxidation simulation is initiated, preheating is done with a stochastic thermostat at 300 or 600 K. Without CTIP (charge transfer ionic potential), the unit cell is relaxed 10 ps using EAM only. The employed time step is 1 fs, and the other calculations use the same value. After this initial equilibration, another relaxation is done for 10 ps, with the CTIP and Berendsen thermostat. This 20 ps relaxation will produce enough thermal perturbation to the atomic positions before the actual oxidation simulation is performed. For oxidation, a new oxygen dimer is inserted from the top of the unit cell. To prevent any overheating or overshoot of the pressure inside of the unit cell, a new oxygen dimer is added when the previously inserted oxygen atoms are consumed by oxidation or meet the topmost Al atom. Those methods have been
employed in other oxide simulations.10,12,13,29 The initial planar position of new oxygen atoms is randomly determined, while the velocity is scaled to the given temperature, going normal to the top surface of the metal substrate. This procedure results in a controlled configuration, and it will give us a guide on the early stages of oxide growth and oxidation kinetics of alloy metal substrates of different stoichiometry. Analysis of the Results. Oxidation simulation has been done for 400 ps with the above-mentioned configurations, at temperatures of 300 and 600 K, to investigate the effect of ambient conditions. In this study, we will analyze the early oxide growth rate and the resulting oxide structure for different base metal stoichiometries. Therefore, layer growth and its rate will be analyzed for various chemical compositions and at ambient temperature. In order to confirm how the stoichiometry affects the oxide growth, several inspections will be made. By analyzing charge evolution and density profile, the ionic behavior and segregation on nonoxidized element will be observed. With pairdistribution and coordination number study, the microstructure will be investigated as well. Oxide Layer Growth. The thickness of the oxide layer was determined as the distance between the z-positions of the lowest oxygen and the topmost aluminum or a charged oxygen atom at the top surface. The oxide growth curves are shown on the left of Figures 2 and 3, and the number of oxygen atoms is indicated on the right. For all the alloy compositions investigated, the oxide growth curves show a quadratic increase up to 80 ps, with a slight plateau after that. It turns out that up to 80 ps, the oxide layer grows nonuniformly, yielding oxide islands, as evident from visual inspection. The plateau at 80110 ps shows that there is equilibration of the nonuniform oxide layer, filling the surface area uniformly. Beyond 110 ps, the oxide layer grows uniformly. The initial build-up of the oxide layer (nonuniform growth) shows a slight difference between the various simulated configurations, yielding higher growth for mono, at both 300 and 600 K. The uniform growth phase shows a distinct difference between the various compositions, being linear in the case of mono but quadratic for ternB, whereas binA, ternA, and binB are quadratic at 300 K but almost linear at 600 K. By 400 ps of oxidation simulation, the oxide attains a limiting thickness, and the results of oxide film growth are summarized in Table 5. For quantitative comparison, we calculated growth rates using fitting over a log scale of uniform growth range. Higher temperature facilitates more oxidation and yielded ∼30% 6574
dx.doi.org/10.1021/jp1106845 |J. Phys. Chem. C 2011, 115, 6571–6580
The Journal of Physical Chemistry C
ARTICLE
Table 5. Results of Thin Layer Growth and Growth Rate thickness (Å)
growth rate (Å/fs)
configuration
300 K
600 K
300 K
600 K
mono
35.6
48.4
19.8
27.1
binA
31.2
42.7
15.0
22.4
ternA
31.3
43.8
16.2
23.3
binB
28.3
39.2
12.9
19.7
ternB
23.3
30.4
8.2
12.8
Figure 5. Side view of the oxide layer of the ternB configuration. Red spheres are oxygen atoms and gray ones are Al atoms. Green is for Ni and orange is for Fe. Blue dotted lines show the range of Ni/Fe atoms at the snapshot. As coalescence of the oxide islands is completed, the range is split into the top and bottom of the oxide layer.
Figure 4. The evolution curve of anion (red marks) minimum and cation (black marks) maximum positions (left, 300 K; right, 600 K).
thicker oxide layer and 3556% higher growth rate for all the configurations. In terms of stoichiometry, mono shows the highest growth rate and the largest amount of oxide. The growth curves are almost linear at both temperatures. In contrast, binA shows the smallest difference, less thick oxide layer, and smaller growth rate compared to mono. Around 300 ps at 300 K, oxide growth shows a quadratic curve, showing resistance against oxidation. But oxidation at 600 K shows almost linear growth up to 400 ps. The growth kinetics of ternA is very similar, showing quadratic behavior around 300 ps at 300 K and monotonic increase at 600 K, even though it has 5% Fe atoms rather than Ni and Al. This will be discussed below; however, this implies that the effect of stoichiometry is not linear. Similarly, binB, which has 15% Ni atoms, shows early quadratic growth, by 220 ps at 300 K and 300 ps even at 600 K. When the content of Ni is intermediate, the effect on the oxide growth is very distinct, producing 28 and 39 Å thickness of the oxide layer at 300 and 600 K, respectively, being fourth order in the simulated configurations. With the highest content of Ni and Fe among the simulation sets, ternB shows significant resistance against oxidation, yielding the least oxide growth and growth rate. The quadratic growth is initiated at ∼200 ps for both 300 and 600 K, and the oxide layer grows up to two-thirds of that of mono by 400 ps. To summarize, mono configuration shows the highest growth rate and oxide film thickness, whereas ternB shows the smallest oxide thickness, showing that higher content of Ni and Fe lowers the oxide growth rate and yields less amount of oxide. Comparing the results at two temperature conditions of 300 and 600 K, higher temperature is found to yield higher growth rate and larger oxide layer for each configuration. For further analysis of the oxide layer growth, the inward growth of anion and the outward growth of cation have been
investigated, as shown in Figure 4. At 300 and 600 K, cation (charged metal atom) growth rate results in a slight effect of stoichiometry, whereas anion (charged oxygen atom) growth curves show distinct differences. In particular, the growth rate of anion is much higher than that of cation for all configurations, matching with the MSD (mean square displacement) results of the previous study.12 Compared to the oxide layer growth shown in Figures 2 and 3, the results show that equilibration has more effect on cationic growth, yielding relative long plateau around 100 ps, whereas anion curves keep growing without equilibrating period, just showing a monotonically changing growth rate. As for in the layer growth curves, mono shows the highest growth of anions, while ternB shows the least. As the temperature increases, the growth of cations is more affected, showing more than 40% increase, while anions show around 20% increase. Visual Inspection of Segregation. Our simulations indicate that ternB shows the highest resistance to oxidation among the simulated configurations. Hence, particle motion in and around the growing oxide layer has been observed, as shown in Figure 5 (visualized by Jmol32). As oxidation begins, an island of AlO oxide is formed, enforcing the migration of neighboring Ni and Fe atoms (75 ps). By 100 ps, the oxide island layer covers the entire unit cell surface, and nonuniform oxide layer is equilibrating. This yields a uniform layer of oxide but also forces Ni and Fe atoms to move out of the oxide layer. At 116 ps, some of Ni and Fe atoms are segregated on the top of the oxide layer, and the others are pushed below. By 150 ps, a distinct group of Ni and Fe atoms is found just below the AlO oxide layer. As described above, the Ni and Fe atoms will be segregated at the top and bottom of oxide layer, and this impedes the oxidation. However, as shown in the oxide growth curve in Figures 2 and 3, the resistance at the nonuniform growth (∼80 ps) is slight, compared to the mono configuration. On the other hand, the effect is very distinct after 200 ps, during the uniform oxide growth stage. Beyond 200 ps, it turns out that the segregation at the bottom of the oxide layer disrupts the interaction between O atoms and Al atoms from metal substrates and acts as a physical barrier. For further oxidation, O atoms need to meet more Al atoms, but the migration of Al atoms from metal substrate is resisted by this segregation. As the fraction of Ni and Fe increases, the segregation will be larger and the resistance to oxidation will increase as well. In other words, binB has fewer 6575
dx.doi.org/10.1021/jp1106845 |J. Phys. Chem. C 2011, 115, 6571–6580
The Journal of Physical Chemistry C
ARTICLE
Figure 6. Radial distribution function [g(r)] of AlAl (top left), AlO (top right), NiO (bottom left), and FeO (bottom right). The pair distribution is averaged over 5000 time steps at 300 ps. Ambient temperature is 300 K.
Figure 8. Charge distribution at 300 ps of oxide simulation at 300 K. A circle shows a single Al and a box is for O atom. A plus sign is for Ni and a triangle is for Fe. Two out of three Al and O atom plots are skipped for better visualization, but Ni and Fe are shown entirely.
Figure 7. Radial distribution function [g(r)] of NiNi and FeFe in ternA/ternB (left) and NiNi pair distribution in binary/ternary alloys (right). Results are averaged over 5000 time steps at 300 ps. Ambient temperature is 300 K.
Figure 9. Charge distribution at 300 ps of oxide simulation with 600 K. The marking and plotting method is similar to that of Figure 8.
non-Al atoms (Ni) and is expected to have lower resistance, compared to ternB. The higher growth rate of binA and ternA, compared to ternB, can be attributed to the occurrence of similar phenomena, which explain the little resistance for smaller fraction of Ni and Fe. This result indicates the role of alloy segregation in disrupting the oxide growth kinetics. Pair Distribution. The above results show the influence of segregation on the alloy oxidation kinetics. However, it is not entirely clear why binA and ternA show similar oxidation kinetics, even with the different content of Fe. We investigate the differences in the microstructure using the radial distribution of paired particles, as shown in Figure 6. The bond length of AlO is determined as 1.7 Å consistently, and it is slightly shorter than the previous study,12 even though distribution shapes are still same. In terms of stoichiometry, the higher content of Ni and Fe is found to reduce the AlAl pair interaction but does not affect AlO. As shown above, most of Ni and Fe atoms are segregated out of the oxide layer, and the structure of the oxide layer may not be affected by the base metal stoichiometry. NiO interaction shows quite a distinct effect on the simulated configuration. In particular, binA shows higher interaction of NiO than ternA, even though both sets have the same Ni content. In other words, adding Fe may weaken the NiO interaction. Between binB and ternB, which have same 15% Ni
Figure 10. Evolution of the average charge of each element at 300 K. Markers are the same as in Figure 8.
content, the same result is confirmed, showing less interaction of NiO at ternB. The similar consequences are found in 600 K simulations as well. The FeO pair distribution shows very different results by the simulated models. For example, ternA, which has 5% Ni and 5% Fe, shows almost negligible results of FeO interactions, while ternB shows relatively high interaction, showing that the effect of stoichiometry is nonlinear. For further analysis, the pair distribution of NiNi and FeFe is detailed in Figure 7. Ternary systems of this work have same contents of Ni 6576
dx.doi.org/10.1021/jp1106845 |J. Phys. Chem. C 2011, 115, 6571–6580
The Journal of Physical Chemistry C
ARTICLE
Figure 11. Spatial and temporal evolution of density plot (left, 300 K; right, 600 K; unit of density, amu/Å3).
and Fe, respectively, and the pair distribution functions of NiNi and FeFe will reveal how much their correlation is different, as shown in the left plots of the figure. The ternA configuration shows much higher correlation of FeFe than NiNi, even though there are small amounts of Fe and Ni atoms. As their content increases, at ternB, the difference of the correlation decreases, but still FeFe shows much higher correlation than the NiNi pair. One interesting point is that adding Fe atoms may disrupt NiNi correlation, in addition to NiO, as shown in Figure 6. The right side of Figure 7 shows the NiNi correlation, with and without Fe atoms. Both binA and ternA have same amount of Ni, while ternA has Fe. As Fe atoms are added, the magnitude of NiNi interaction decreases. This result is repeated between binB and ternB, which have the same amount of Ni atoms. From detailed analysis of the pair distribution function, we could confirm that the stoichiometry effect of Fe is nonlinear, because it has very high correlation compared to the NiNi pair, yielding Fe-enriched islands for the high-Fe-content alloy. For further analysis of microstructure, the charge states of each configuration will be investigated. Charge States. A snapshot of typical charged particle states is shown in Figures 8 and 9. From the mono configuration at 300 K, oxide layer covers from 10 to 25 of z-position. In the middle of the oxide layer, oxygen and aluminum atoms are fully oxidized, but they are partially charged at the interface or top surface (free surface), at both 300 and 600 K. Results at 600 K show similar charge distribution with a larger amount of oxide layer. From both 300 and 600 K results, Ni yields fewer charged particles, consistently. However, the charged state of Fe is very abrupt, as
indicated by the pair distribution results. At 300 K, ternA has 0% oxidized Fe but ternB has 9%. On the other hand, simulations at 600 K yielded 4% oxidized Fe with ternA and 13% with ternB. As another approach to investigate charged states, we calculated the average charge of each element as shown in Figure 10. The Al atom charge curve shows a monotonic increase, while the O atom charge shows saturation by 120 ps for all the configurations, when the uniform oxide growth begins. At 600 K, we could see similar results, while saturation begins around 100 ps, slightly earlier than for 300 K. As mentioned above, Ni shows reduced charged states, from 150 ps, where uniform layer growth begins. However, Fe shows a strong dependence on the initial alloy content. The ternA configuration shows noncharged states of Fe atoms, whereas they are relatively highly charged for ternB. Although one might expect an effect of base metal stoichiometry on the evolution of the charged states of the oxide layer, our simulations indicate that it has a minimal effect. Density Profile. From the charge state curves and visual inspection, we could observe how the segregation of Ni and Fe leads to increased oxidation. To analyze the evolution of dopant segregation, density profiles based on spatial and temporal trajectory have been calculated. Ni and Fe have much higher mass than Al and O, and the density plot will give us a spatial and temporal segregation profile. Figure 11 provides a 1D view of the density profile, averaged along z-coordinates, for all the configurations at 300 and 600 K. Base metal substrate is [100] of fcc structure, and we can see the clear orientation at the beginning of oxidation. As oxide layer grows, it smears the density distribution, mixing metal alloy with oxygen atoms. However, in the middle of oxide layer, we could find that high- and low-density states are 6577
dx.doi.org/10.1021/jp1106845 |J. Phys. Chem. C 2011, 115, 6571–6580
The Journal of Physical Chemistry C
ARTICLE
Figure 13. Time evolution of coordination number fractions at 300 K simulations. Figure 12. Evolution and distribution of coordination numbers of mono configuration at 300 K. The sum of the coordination Z = 36 is 1.0 at any point on the oxide layer.
layered, recursively. They are the layer of Al and O atoms, and it shows that the oxide layer is composed in a layer-by-layer mode, even though they are disordered. Yielding 2D planar structures, this result is consistent with the results of Hasnaoui et al.10 and is found in all of the simulated configurations. For chemical compositions other than pure aluminum, they show less amount of oxide layer and lower rate of oxidation, as discussed above. However, the density state inside the oxide layer is very consistent, regardless of the employed stoichiometry, implying the segregation of heavy dopants. We could find the high concentration of heavy elements just below the oxide layer, as observed from other analysis. This heavy layer is formed when the coalescence of oxide islands is completed (∼120 ps) and gets thicker as the oxide layer grows. Also as the content of Ni and Fe increases, the layer of heavy elements becomes more distinct and thicker, showing clear segregation. This heavy layer is not adjacent to the oxide layer but separated by oxygen atoms, showing that the heavy elements are neighboring to oxygen atoms. Partial oxidation of Fe atoms in ternB might be attributed by this close interaction with oxygen atoms. Coordination Number. Using the bond-length data from pair distribution results, the coordination number (Z) of Al has been investigated. Counting the number of nearby oxygen atoms, we analyze the structural states of Al atoms in spatial and temporal trajectories. We are primarily concerned with Z = 36, as the basic structure of aluminum oxide (AlOZ), and the results are normalized with their sum. Figure 12 shows the evolution of the coordination number of the mono configuration at 300 K. Z = 3 is very dominant at early stages of oxidation, when ions exhibit partially charged states but transitions into Z = 4 or 5 quickly as oxidation equilibrates. Z = 6 is found uniformly in the fully oxidized region, but the fraction is very low. As the temperature increases, we could observe the distinct increase of Z = 3 around the oxide interface or top surface, even though it disappears at equilibration and when oxidation is nearly complete. For the other stoichiometry sets, the characteristics of coordination number evolution are consistent. By averaging the spatial distribution, Figure 13 provides the time evolution of coordination number of tested chemical compositions at 300 K. As discussed above, Z = 4 has the highest fraction (∼45%) while Z = 5 is very close (∼35%), showing AlO4 is the dominant structure in the oxide layer in all configurations. This means that the AlO
structure is not affected by other elements like Ni and Fe. From the above results, like charge/density profile, and visual inspection, Ni and Fe atoms are segregated out of oxide layer toward the oxidemetal interface, and the interior structure of the oxide thin film may be homogeneous, independent from the initial metal alloy stoichiometry. Surface segregation may affect only the surface microstructure, but still the amount of dopant utilized is very small, and therefore, the corresponding effect is minimal.
’ DISCUSSION From our simulations on the alloys, early stages of ultrathin oxide layer growth were found to have several phases, such as quadratic growth by oxide island formation, equilibration plateau, and linear growth by uniform build-up. Yang et al.,33 who studied Cu oxidation studies by ultrahigh vacuum (UHV) transmission electron microscopy (TEM), reported that the coalescence of oxide islands reduces the growth rate by diffusion through the layer. The reduction of growth rate matches with our modeling results, showing equilibration after quadratic growth consistently for all the tested configurations. Coalescence of oxide islands has been confirmed through visual inspection as well. As an effect of base alloy metal stoichiometry, dopant elements like Ni and Fe disrupt the coalescence, but the effect seems to be negligible, because those elements are not enriched yet between the islands. However, they are segregated below the oxide layer and produce significant effect on oxide growth when their content is very high, disrupting the Al atom supply from the base metal substrate. As oxide layer grows, dopant elements move away and segregate at the top and bottom of the layer. Ni shows slightly charged states consistently, but Fe has highly charged states at only ternB, which has higher content (15%). As discussed by Kizilkaya et al.,34 alumina (ΔHAl2O3 = 1,676 kJ/mol) is more favored than Fe-based oxides (ΔHFe2O3 = 826 kJ/mol, ΔHFeO = 278 kJ/mol, ΔHFe3O4 = 1118.4 kJ/mol) in terms of formation enthalpy. Therefore, the result of ternA is typical, but the sudden difference for ternB is unexpected. As addressed in the pair distribution functions, Fe atom pairs are more highly correlated than Ni pairs, and the nonlinear increase of charged states might be resolved by the correlation. Graupner et al.35 tested FeAl surface oxidation using LEED and XPS, and derived the oxidation process: when the nearby Al atoms are consumed and depleted for oxidation, Fe-enriched islands will be exposed to the incoming oxygen atoms without competition from 6578
dx.doi.org/10.1021/jp1106845 |J. Phys. Chem. C 2011, 115, 6571–6580
The Journal of Physical Chemistry C neighboring Al atoms, because the Al depletion rate is much faster than the bulk diffusion rate around room temperature. Our ambient temperature condition is on the same order, and the charge difference of Fe atoms for ternA and ternB can be described in the same way. In particular, from the pair distribution curve, we can infer that FeFe has very high correlation compared to the NiNi pair, showing higher probability of building Fe-enriched islands. As found from visual inspection, Fe and Ni build a physical barrier to resist the growth of an oxide layer. Then the supply of Al is limited and the larger number of Fe atoms will be swamped by O atoms, reducing the competition with neighboring Al atoms. The O-swamped Fe atoms are found in the density profile of Figure 11, showing oxygen atoms between the oxide layer and heavy elements. Ni has lower correlation than Fe and it may not build its own islands, lowering its probability to be oxidized. Besides, Ni oxide is the least favored (ΔHNiO = 489.5 kJ/mol36) among the employed elements, and the less oxidation of Ni is confirmed by the work of Sakiyama et al.,37 testing FeNiAl alloys with X-ray diffraction at 600800 °C and reporting that Fe2O3 and Al2O3 are dominant scales.
’ CONCLUSIONS Employing variable-charge molecular dynamics, we have investigated the oxidation and surface oxide growth on ternary AlNiFe alloys. Among the studied configurations, pure Al shows the highest and fastest oxide layer growth at 300 and 600 K. As Ni and Fe content increases, greater segregation is found, which results in a corresponding increase in the barrier to oxidation, leading to slower kinetics. At a higher temperature of 600 K, oxidation is accelerated and 30% more oxide was found for all of the simulated configurations. Investigation of the atomic scale oxide layer growth suggests that the inward migration rate of O atom surpasses the outward migration of Al atom for all the configurations and temperature conditions. When the content of Fe is very high (such as ternB), Fe atoms are oxidized more than Ni, as indicated by the high correlation in the pair distribution function, even though the nominal concentration of Fe and Ni are similar in the alloy. The oxide layer is disordered, but it shows a layer-by-layer growth mode, yielding a 2D planar structure, as confirmed by the density profile. This is also consistent with previous studies on oxidation of pure Al.10 On the basis of the bond length from pair distribution, the coordination number has been investigated. The local atomic coordination shows a transition from Z = 3 during initial stages of oxide growth to Z = 4 and 5 at longer exposure times for all the simulated alloy compositions, demonstrating a negligible effect of base metal stoichiometry. The following summarizes the analysis and available results: • The kinetics of oxide growth is dependent on the metal alloy stoichiometry and ambient temperature, while pure aluminum shows the highest growth rate among simulated configurations. • The strong correlation among Fe atoms produces a nonlinear effect on oxide growth, yielding Fe-enriched islands at high Fe content in the aluminum alloy. • Even though oxide growth is affected by the metal substrate stoichiometry, the structure of the oxide layer is not significantly influenced, due to the segregation of the dopant elements.
ARTICLE
’ AUTHOR INFORMATION Corresponding Author
*E-mail:
[email protected].
’ ACKNOWLEDGMENT This work has been supported by the Office of Naval Research through contract No. N00014-10-1-0346. The computational facilities have been provided by CNS (Center for Nanoscale Systems)NNIN (National Nanotechnology Infrastructure Network) at Harvard University. ’ REFERENCES (1) Freund, H.-J.; Kuhlenbeck, H.; Staemmler, V. Rep. Prog. Phys. 1996, 59, 283. (2) Tsuchiya, M.; Sankaranarayanan, S. K.; Ramanathan, S. Prog. Mater. Sci. 2009, 54, 981–1057. (3) Streitz, F. H.; Mintmire, J. W. Phys. Rev. B 1994, 50, 11996–12003. (4) Daw, M. S.; Baskes, M. I. Phys. Rev. Lett. 1983, 50, 1285–1288. (5) Rappe, A. K.; Goddard, W. A. J. Phys. Chem. 1991, 95, 3358–3363. (6) Campbell, T.; Kalia, R. K.; Nakano, A.; Vashishta, P.; Ogata, S.; Rodgers, S. Phys. Rev. Lett. 1999, 82, 4866–4869. (7) Alavi, S.; Mintmire, J. W.; Thompson, D. L. J. Phys. Chem. B 2005, 109, 209–214. (8) Zhou, X. W.; Wadley, H. N. G.; Filhol, J.-S.; Neurock, M. N. Phys. Rev. B 2004, 69, 035402. (9) Zhou, X. W.; Wadley, H. N. G. Phys. Rev. B 2005, 71, 054418. (10) Hasnaoui, A.; Politano, O.; Salazar, J.; Aral, G.; Kalia, R.; Nakano, A.; Vashishta, P. Surf. Sci. 2005, 579, 47–57. (11) Perron, A.; Garruchet, S.; Politano, O.; Aral, G.; Vignal, V. J. Phys. Chem. Solids 2010, 71, 119–124. (12) Sankaranarayanan, S. K. R. S.; Ramanathan, S. Phys. Rev. B 2008, 78, 085420. (13) Chang, C.-L.; Sankaranarayanan, S. K. R. S.; Ruzmetov, D.; Engelhard, M. H.; Kaxiras, E.; Ramanathan, S. Phys. Rev. B 2010, 81, 085406. (14) Georgieva, V.; Todorov, I. T.; Bogaerts, A. Chem. Phys. Lett. 2010, 485, 315–319. (15) Lozovoi, A. Y.; Alavi, A.; Finnis, M. W. Phys. Rev. Lett. 2000, 85, 610–613. (16) Kresse, G.; Schmid, M.; Napetschnig, E.; Shishkin, M.; Kohler, L.; Varga, P. Science 2005, 308, 1440–1442. (17) Islam, M. M.; Diawara, B.; Maurice, V.; Marcus, P. J. Phys. Chem. C 2009, 113, 9978–9981. (18) Essmann, U.; Perera, L.; Berkowitz, M.; Darden, T.; Lee, H.; Pedersen, L. G. J. Chem. Phys. 1995, 103, 8577–8593. (19) Zhou, X. W.; Wadley, H. N. G. J. Phys.: Condens. Matter 2005, 17, 3619. (20) Lazic, I.; Thijsse, B. Mater. Res. Soc. Symp. Proc. 2007, 978E, 0978GG08-05. (21) Lazic, I. Atomic scale simulation of oxide and metal film growth. Ph.D. thesis, Delft University of Technology, 2009. (22) Press, W. H.; Flannery, B. P.; Teukolsky, S. A.; Vetterling, W. T. Numerical Recipes in Fortran: the Art of Scientific Computing, 2nd ed.; Cambridge University Press: New York, 1992. (23) Takahashi, D. Lect. Notes Comput. Sci. 2001, 2110/2001, 551–554. (24) Takahashi, D. Lect. Notes Comput. Sci. 2006, 2367/2002, 759. (25) Takahashi, D. Comput. Phys. Commun. 2003, 152, 144–150. (26) Takahashi, D. Parallel Comput. 2003, 29, 679–690. (27) Takahashi, D. Lect. Notes Comput. Sci. 2006, 3911/ 2006, 970–977. (28) Takahashi, D. Lect. Notes Comput. Sci. 2007, 4699/2007 1178–1187. (See also http://www.ffte.jp.) (29) Sankaranarayanan, S. K. R. S.; Kaxiras, E.; Ramanathan, S. Phys. Rev. Lett. 2009, 102, 095504. 6579
dx.doi.org/10.1021/jp1106845 |J. Phys. Chem. C 2011, 115, 6571–6580
The Journal of Physical Chemistry C
ARTICLE
(30) Metropolis, N.; Rosenbluth, A. W.; Rosenbluth, M. N.; Teller, A. H.; Teller, E. J. Chem. Phys. 1953, 21, 1087–1092. (31) Wang, G.; Hove, M. V.; Ross, P.; Baskes, M. Prog. Surf. Sci. 2005, 79, 28–45. (32) Jmol: An open-source Java viewer for chemical structures in 3D; http://www.jmol.org/. (33) Yang, J. C.; Evan, D.; Tropia, L. Appl. Phys. Lett. 2002, 81 241–243. (34) Kizilkaya, O.; Hite, D. A.; Zehner, D. M.; Sprunger, P. T. Surf. Sci. 2003, 529, 223–230. (35) Graupner, H.; Hammer, L.; Heinz, K.; Zehner, D. M. Surf. Sci. 1997, 380, 335–351. (36) Finnis, M.; Lozovoi, A.; Alavi, A. Annu. Rev. Mater. Res. 2005, 35, 167–207. (37) Sakiyama, M.; Tomaszewicz, P.; Wallwork, G. R. Oxid. Met. 1979, 13, 311–330.
6580
dx.doi.org/10.1021/jp1106845 |J. Phys. Chem. C 2011, 115, 6571–6580
