Size Effect on the Oxidation of Aluminum Nanoparticle: Multimillion-atom Reactive Molecular Dynamics Simulations

The size effect in the oxidation of aluminum nanoparticles (Al-NPs) has been observed experimentally; however, the mechano-chemistry and the atomistic mechanism of the oxidation dynamics remain elusive. We have performed multimillion atom reactive molecular dynamics simulations to investigate the oxidation dynamics of Al-NPs (diameters, D ¼ 26, 36, and 46 nm) with the same shell thickness (3 nm). Analysis of alumina shell structure reveals that the shell of Al-NPs does not break or shatter, but only deforms during the oxidation process. The deformation depends slightly on the size of Al-NP. This reaction from the oxidation heats the Al-NP to a temperature of T > 5000 K. Ejection of Al atoms from shell starts earlier in small Al-NPs—at t 0 ¼ 0.18, 0.28 and 0.42 ns for D ¼ 26, 36 and 46 nm, when they all have the same shell temperature of 2900 K. As the oxidation dynamics proceeds, the total system temperature (including the environmental oxygen) increases monotonically; however, the time derivative of the total temperature, (dT system /dt), reaches a maximum at t 1 ¼ 0.20, 0.32 and 0.51 ns for D ¼ 26, 36 and 46 nm. At this peak value of (dT system /dt), the shell temperature for the three Al-NPs are 3100 K, 3300 K, and 3500 K, respectively. The time lag between t 1 and t 0 is 0.02, 0.04 and 0.09 ns for D ¼ 26, 36 and 46 nm clearly indicates the size effect. V C 2013 Author(s). All article content, except where otherwise noted, is licensed under a Creative Commons Attribution 3.0 Unported License.


I. INTRODUCTION
Combustion of aluminum is a subject of great interest in science, engineering, and technology.It has been studied since the last century due to its wide use in the field of solid propellants and explosives.For instance, combustion of micron size Al particles was studied by Dreizin 1 and Jordan et al. 2 using different experimental approaches.The effect of temperature on oxidation [3][4][5][6][7] and the effect of oxide shell on the surface [8][9][10][11][12][13][14][15] have been studied extensively, especially for nanosize Al particles.][18][19] By adding Al-NPs into micron size particles, numerous experiments have observed the enhanced oxidation reactivity in terms of burning rate, 20,21 flame speed, 22,23 activation energy, 24 ignition sensitivity, 23,[25][26][27] combustion velocity, 25 and agglomeration. 3Nanosize Al particles have lower melting temperatures than that of micron size particles, because of large surface-to-volume ratio.The combustion rate of Al nano-composites has been observed to increase significantly over micron size Al composites.Several investigators have specifically investigated the size effect on the combustion of Al-NPs. 28,29Most of them concluded that the smaller the size of Al-NPs is, the higher the reactivity becomes.However, the effect of particle size on reactivity is still controversial.For example, Gan and Qiao 30 have conducted experiments on fuel droplets with nano ($80 nm) and micron ($5 mm and $25 mm) size Al particles.They found that the combustion is longer and less complete for large agglomerate of nanosuspensions due to the formation of oxide shell on the Al-NPs surface.In a review article, Yetter et al. 31 stated that for smaller Al particles, the energy loss per unit volume due to the presence of the same oxide layer thickness is significantly larger.
Besides the above-mentioned experimental studies, there have been numerical simulations of oxidation dynamics of Al-NPs.Campbell et al. 32,33 studied the oxidation of aluminum nanoclusters with a parallel molecular dynamics (MD) approach based on dynamic charge transfer among atoms in both microcanonical and canonical ensembles and for atomic and molecular oxygen environments.Alavi et al. 34 have simulated the oxidation of Al-NPs using MD with Streitz-Mintmire electrostatic plus (ESþ) potential. 35,36uri and Yang 37 performed MD simulation to study the thermo-mechanical behavior of nano-Al particles coated with crystalline and amorphous oxide layers.Perron et al. 38 have studied the oxidation of multi-grain nanocrystalline (mean grain size ¼ 5 nm) Al surfaces in the temperature range 300-750 K using variable charge molecular dynamics simulations.Structures of cand amorphous Al 2 O 3 have been investigated using MD and density function theory (DFT) by Gutierrez et al. 39,40 Vashishta et al. 41 have developed an interaction potential consisting of two-and three-body terms for alumina to simulate the amorphous and liquid phases.Combining Vashishta's potential 41 and embedded atom method (EAM) potential, 42 Wang et al. 43 have studied atomistic mechanisms of oxidation in a laser flash heated core (aluminum)-shell (alumina) nanoparticle.In another paper by Wang et al., 44 they studied effects of the crystalline and amorphous structures of alumina shells on the dynamics of oxidation of an Al-NP.
Despite the extensive experimental and simulation research efforts mentioned above, several key questions remain unanswered.Most important of them are the collective reaction behaviors among particles; and spatially and temporally resolved oxidation dynamics of individual particles at atomic resolution.It is essential to delineate these issues unambiguously and study each factor independently. Regarding the collective reaction, Shekhar et al. 45 have investigated the mechano-chemial reaction in the oxidation of three adjacent Al-NPs by initiating the oxidation in the middle nanoparticle and studying the nature and speed of the oxidation front.However, there is no systematic study on the atomistic mechanisms underlying size-dependent oxidation dynamics of individual Al-NPs.This paper focuses on the size dependence of oxidation dynamics rate and that of reaction temperature, as well as the role of oxide shell in Al-NPs oxidation.Here, we report the results of large-scale, parallel MD simulations of the size effect on single Al-NP of diameters D ¼ 26, 36, and 46 nm.The rest of the paper is organized as follows.The interatomic interaction potential and the system setup are described in Secs.II and III, respectively.Simulation results are presented in Sec.IV, and conclusions are in Sec.V.

II. INTERATOMIC INERACTION POTENTIAL
The foundation of MD simulations is interatomic potential.The modeling of Al 2 O 3 is based on a many-body potential that incorporates ionic and covalent effect through a combination of two-and three-body terms 41 V ¼ where N is the number of atoms in system.The two-body interatomic interactions include effects of steric repulsion, Coulomb interaction, charge-dipole interaction, and Van der Waals interaction, as given below: The three-body interatomic interaction reflects the covalent effects through bond-bending and bond-stretching terms In Eq. ( 2), H ij is the steric repulsion strength; Z i is the effective charge in units of the electronic charge, e; D ij represents the charge-dipole strength; and w ij is the strength of Van der Waals attractions.g ij is the exponent of the steric repulsion, and r 1s and r 4s are the screening lengths for the Coulomb and charge-dipole interactions, respectively.r ij is the distance between any two atoms in the system.In Eq. ( 3), the three-body interaction, B ijk is the strength of the three-body interaction, h ijk is the angle formed by vectors r ij and r ik , and r 0 is the cutoff distance for the three-body interaction.This potential has been validated by comparing the MD results for structural and mechanical properties of both crystalline and amorphous alumina with the experimental data. 41or aluminum in the core of the nanoparticle, we use the EAM form of the potential 42 E ¼ where E i is the energy of atom i The first term in Eq. ( 5) describes the electrostatic interaction between atom i and its neighboring atoms j.The second term describes the attractive interaction which models placing a positively charged atom in the electron density due to the free valence sea of electrons created by the host system of atoms.The embedding function F i ( q i ) is a function of superimposed charge densities q i , due to the charge density distribution of neighboring atoms where q(r ij ) is the pair-wise electronic density as a function of the distance r ij between atom i and atom j, but without angular dependency.
In the case of Voter-Chen's EAM expression, the pairwise electronic interaction is taken to be a Morse potential where D M and R M define the depth and position of the minimum, respectively.a M is a measure of the curvature at the minimum.The density function is a hydrogenic 4s orbital with a relative normalization factor added to ensure the monotonicity where b is an adjustable parameter.
In order to implement the EAM in MD with O(N) scaling, both /(r) and q(r) are truncated at cutoff r cut using where h(r) can be /(r) or q(r), and m ¼ 20.
To describe the oxidation of Al, a bond-order based interpolation scheme is proposed to smoothly interpolate the two potentials in order to guarantee the stability of the interface and the energy conservation. 42Mathematically, the total potential energy of the system can be written as the sum of the potential energy of all aluminum atoms and that of all oxygen atoms The potential energy of each Al atom i is a weighted average between Vashishta's potential for alumina and the EAM where f is the weighting factor as a function of the oxidation degree n i of atom i The oxidation degree n i of atom i is where H(r ik ) is a continuous function that counts the degree of oxidation according to the distance of the considered aluminum atom to a neighboring oxygen atom here we use R a ¼ 3.0 A ˚and D a ¼ 0.5 A ˚.The interpolated potential was validated by comparing the MD results for the structure of several Al x O y small clusters with the corresponding quantum mechanical results. 46e have also checked the structural and dynamical stability of interface between the aluminum core and the alumina shell once they are put together in the initial configuration of the aluminum nano-particles.

III. SYSTEM SETUP FOR REACTIVE MOLECULAR DYNAMICS SIMULATIONS
The initial setup for the oxidation simulation is as follows (see Table I).The Al-NP diameters are chosen to be D ¼ 26, 36 and 46 nm, with 3 nm thick amorphous alumina shell.The metallic Al core is cut from crystalline Al at 300 K.The amorphous shell is prepared by first heating the crystalline alumina, and then slowly quenching the molten alumina into amorphous alumina at 300 K, and then, from the amorphous alumina, removing the outer and inner concentric parts accordingly to retain an oxide layer of 3 nm thick.After combining the crystalline aluminum core and amorphous alumina shell structure together, the single Al-NP is relaxed and thermalized before it is put into an oxygen environment.
The following relaxation procedure is performed for the Al-NP assembled after putting the core (aluminum) and shell (amorphous alumina) together.The reactive MD was run on the Al-NP for 30 ps, in microcanonical ensemble, to form and stabilize the aluminum/alumina interface.Figure 1(a) displays the radial temperature distribution during preparation of 26 nm Al-NP system.We plot temperature within a 10 A ˚thick spherical shell as a function of radial distance R. The red curve in Fig. 1(a) shows the temperature profile at 30 ps, which is the end of interface-formation process.The shell temperature is higher than core and the environmental oxygen because of the shell reconstruction.The Al-NP is then quenched to 5 K (blue curve in Fig. 1(a)) and the whole system is gradually heated to 300 K (as shown in time sequence of temperature profiles from cyan to green to orange in Fig. 1(a)).
In order to start the oxidation dynamics in the Al-NP, velocity scaling is performed to preheat nanoparticle from 300 K to 1100 K in steps of 100 K (see Fig. 1(b)).In each step, the Al-NP is heated at a constant heating rate, followed by thermalization for 0.1 ns during each interval.During the Al-NP heating, the oxygen environment is kept at 300 K as shown in Fig. 1(b).After the preheating, systems undergo oxidation for 1 ns in microcanonical ensemble.Exactly the same procedure is applied to D ¼ 36 and 46 nm Al-NPs.

IV. RESULTS FROM THE ANALYSIS REACTIVE MOLECULAR DYNAMICS SIMULATIONS A. Oxidation dynamics of three Al-NPs
To visualize the reaction process of the three Al-NPs, Figs.2-4 show snapshots of the central slice for D ¼ 26, 36 and 46 nm Al-NPs at 0.0, 0.2, 0.4, 0.6, 0.8 and 1.0 ns, respectively.Blue color represents environmental oxygen.Green and black colors indicate the oxygen and aluminum atoms, respectively.In the alumina shell, the red color represents the core aluminum in the initial configuration.The oxidation process starts after the initial conditions are set up (core and shell preheated to 1100 K).As the oxidation begins, due to the molten core expansion at 1100 K, the size of the Al-NP grows larger than its initial size.
No shell breaking or shattering of the shell is observed in the oxidation process; only the deformation of Al-NP is observed.From the extent of the core and shell deformation, reaction is apparently more intense for D ¼ 26 nm Al-NP than D ¼ 36 and 46 nm Al-NPs.

B. Temperature profile of three Al-NP systems
To determine the size dependence on oxidation reactivity of Al-NPs, we study the heat release for all three systems by plotting temperature of whole systems as functions of time showing in Fig. 5(a).Significant amount of heat is produced as a result of oxidation within the nanoparticles in all three systems.The radial temperature distribution of the three systems at 0.3 ns is shown in Fig. 5(b).Temperature at the center of the core is higher than the shell temperature for the smallest system (D ¼ 26 nm) as shown in Fig. 5(b).However, the temperature at the center of the core is lower than the shell temperature at 0.3 ns for larger systems (D ¼ 36 and 46 nm).Comparison between radial temperature distribution and global temperature of all systems (at time less than 0.3 ns) shows that the shell temperature is higher than the average temperature of the system, which indicates that oxidation at core-shell interface dominates the heat release of the whole system at early stage.

C. Local aluminum atoms concentration
In order to study the distribution of aluminum atoms inside Al-NPs, we calculate local concentration of aluminum atoms from their two different origins, i.e., those from the alumina shell and from the aluminum core (see Fig. 6(a)).Local concentration at distance R from the center of Al-NP is calculated by averaging the concentration within a spherical shell of thickness dR ¼ 10 A ˚. We compute three concentrations as follows:  number of atoms, n Al is the total number of Al atoms, n shell Al is the number of Al atoms from alumina shell, n core Al is the number of Al atoms from aluminum core at distance R from the center of Al-NP within a spherical shell of thickness dR ¼ 10 A ˚. Figure 7 shows C Al as function of radius R for time t ¼ 0, 0.3 and 0.6 ns for D ¼ 26, 36 and 46 nm.In Figs.7(a)-7(c), at 0 ns, shown in blue lines, there are shoulders around R ¼ 10-13, 15-18, and 20-23 nm, respectively.These shoulders reflect the core-shell structure of Al-NPs.Namely, the center part, up to R ¼ 10, 15 and 20 nm, respectively, for the three systems, is pure Al.The shell part is Al 2 O 3 extends for 3 nm beyond the core.There are no Al atoms beyond the shell radius of the three systems.As time progresses, these shoulders become less pronounced indicating the intermixing of atoms between shell and core.
To distinguish origin of Al atoms, Figs. 8 and 9 show C shell Al and C core Al , respectively.In Fig. 8, the movement of shell Al atoms is observed from the radial distribution at different times.At 0 ns, blue lines, the local concentration is 0.4, reflecting the composition of Al 2 O 3 .At 0.3 ns, green lines, the alumina shells have moved outwards and become wider due to the expansion of the molten core.For D ¼ 26 nm, more Al atoms spread inwards, which is not observed in D ¼ 46 nm until 0.6 ns.At 0.6 ns, red lines, more shell Al atoms spread outwards than inwards for D ¼ 26 and 36 nm.At times beyond 0.6 ns, the same behavior is seen in D ¼ 46 nm.In Fig. 9, at 0 ps, blue lines, the center part is pure (metallic) Al, while there is no core aluminum in alumina shell or oxygen environment.At 0.6 ns, red lines, the core Al spreads outward.The core Al atoms in D ¼ 26 nm moves outwards earlier in time than other systems at the same instant.

D. Formation of Al x O y clusters
To understand the oxidation in the core-shell region for three different sizes of Al-NPs, we next investigate the chemical composition of oxide clusters using fragment analysis.9][50] Here, we take the fragment as the atoms that are covalently bonded.Except the large fragments (those having over a thousand atoms,

E. Evolution of the shell structure
During the oxidation, the shell structure of Al-NPs changes, in width due to core expansion and also in its chemical composition.Figure 13(a) shows the outer radius of the Al-NPs.Within 1 ns, D ¼ 26 and 36 nm systems show expansion followed by contraction, whereas 46 nm system only shows expansion.The contraction phase for the D ¼ 46 nm system is beyond 1 ns.The extent of expansion (the increment percentage of the radius from the initial value to the maximum value) of the D ¼ 26, 36 and 46 nm systems is 18.89%, 22.33% and 24.18%, respectively.The oxide shell radius is calculated by averaging distances from the center of the system to all atoms contained within the largest oxide fragment.Fig. 13(b) shows the ratio between oxygen atoms and aluminum atoms in the largest oxide fragments (i.e., the oxide shell).Initially, the O/Al ratio of 1.5 stands for the amorphous Al 2 O 3 shell.The formation of the shell structure is a dynamic process involving Al atoms coming into the shell from the core and Al atoms leaving the shell into the oxygen environment.The overall tendency of the ratio is declining.At early times, the O/Al ratio in the shell region decreases with time for all the three systems.For D ¼ 26 nm, oxygen incorporation from the environment increases the ratio until it peaks around 1.45 at 0.5 ns and then decreases.A similar behavior is observed for D ¼ 36 and 46 nm, but at later times and with different ratio values.By 1 ns, the ratio for all systems is under 1.4, which means that the largest fragment in the system changes from perfect amorphous oxide shell to an Al-rich shell.
In the process of oxidation, we observe some shell Al atoms eject out into the surrounding oxygen.In order to quantitatively determine the ejection process, after the oxide shell structure has been defined (the largest oxide fragment), we count the number of shell Al atoms that are located beyond the oxide shell radius.Figure 14(  shows that for all three Al-NP systems the shell temperature at the onset Al ejection time is the same, i.e., $2900 K.The dashed arrows specify the onset time of Al ejections as obtained from Fig. 14(a).The black dashed double arrows mark the shell temperature when Al ejections begin.Thus, we conclude that only after the shell temperature reaches 2900 K, the shell Al atoms will start their ejections into the oxygen environment.At 2900 K, the alumina shell is in molten state in all three systems.

F. Reaction delay
In order to explain the time dependence of reaction for the three systems, we plot temperature of Al-NPs including the surrounding oxygen at different time in Fig. 15(a) and the rate of change in temperature (dT system /dt) for the three systems in Fig. 15(b).Comparing Figs.15(a) and 15(c) (same as Fig. 14(b)), the temperature of the global systems and the oxide shells has the same increasing trend.This restates the fact that the temperature of the shell is higher than the average temperature of the whole system, which indicates that reaction in the core-shell region dominates the heat release as mentioned in Sec.IV A.
Figure 15(b) shows that the largest rate of change in temperature (dT system /dt) occurs at t 1 ¼ 0.2 ns, 0.32 ns and 0.51 ns for D ¼ 26, 36 and 46 nm, respectively.At these times, the reaction is most intensive.Once passing this peak (i.e., the largest rate of change in temperature (dT system /dt) max ), temperature still increases but at a slower rate.In Fig. 15(a), the temperature at the culmination point is given for each Al-NP (2500 K, 2700 K and 2800 K for D ¼ 26, 36 and 46 nm, respectively).
It is worth noting that the delay of the culmination points from the onset of shell-Al ejection is different for the three Al-NPs systems.The delay for D ¼ 26, 36 and 46 nm are (t 1 À t 0 ) ¼ 0.02, 0.04 and 0.09 ns, respectively.Namely, smaller system has shorter delay.Finally, we would like to point out the possible effects of oxygen density on the reaction dynamics.In this paper, we used a higher oxygen density so that the entire oxidation process can be observed during the simulation time.We expect that the earlier stage of the reaction is less sensitive to the oxygen density because the reaction is occurring inside the nanoparticle at the metal core and oxide shell.However, the reaction at the later stage involves the environmental oxygen, thus it should depend on the oxygen density. 51This can be seen in the D ¼ 26 nm Al-NP, where the depletion of environmental oxygen slows down the oxidation reaction.From Fig. 15(b), the oxygen from the environment for 26 nm system is almost consumed after 0.6 ns, which leads to a decrease of the reaction rate (dT system /dt), as shown in Fig. 15(c).

V. CONCLUSION
Reactive MD simulations results have provided quantitative information on the mechano-chemistry and the atomistic mechanism of oxidation for Al-NPs of different sizes, D ¼ 26, 36 and 46 nm in which the Al core diameter is varied from 20 to 40 nm and the alumina shell thickness is held constant at 3 nm.In the initial state of the system, the nanoparticles are uniformly heated to 1100 K in ambient oxygen environment.
We find that the oxidation dynamics starts within the Al-NP at the core-shell interface when molten aluminum starts reacting with alumina shell by taking some oxygen from expanded alumina shell to form aluminum rich Al x O y (x > y) clusters.This reaction at the aluminum/alumina interface generates heat first.By analyzing the structure of the oxide shell of all three Al-NPs in the process of this reaction, we found that the shell of Al-NPs does not break or shatter, but only deforms, where the deformations differ slightly depending on the size of the Al-NP.For smaller Al-NP, which has larger surface-to-volume ratio and lower ratio of core-radius-to-shell-thickness, the reaction is faster.In early stage of the reaction, a majority of Al x O y clusters formed are small, having less than six atoms.Most of the clusters in this stage are AlO and Al 2 O as shown in Figs.10-12.The induced heat is mainly as a result of the oxidation dynamics at the aluminum/alumina interface.When the temperature reaches $2900 K, the alumina shell melts.This shell melting is coincident with the start of aluminum ejections from the Al-NP into oxygen environment.The aluminum ejections start earlier in smaller Al-NP-at t 0 ¼ 0.18, 0.28 and 0.42 ns for D ¼ 26, 36 and 46 nm.We would like to emphasize here that in all the three systems the shell temperature at t 0 ¼ 0.18, 0.28, and 0.42 ns is 2900 K at which time the shell is in molten state.After the Al ejection starts from Al-NP, the oxidation dynamics also starts at the alumina/oxygen interface.The initial Al x O y clusters found in this outer region are oxygen rich (y > x).The start time of AlO 2 clusters production in the outer regions of the Al-NP is after the shell temperature has reached $2900 K, as can be inferred from Figs. 10(c)-12(c).We have also examined the rate of heat production in the three Al-NP systems.As the oxidation proceeds, the total system temperature (including the environmental oxygen) increases monotonically; however, the time derivative of the total system temperature, (dT system /dt), has a peak that occurs at t 1 ¼ 0.20, 0.32 and 0.51 ns for D ¼ 26, 36 and 46 nm.At this peak in the time derivative of the total temperature, in (dT system /dt), the shell temperature for the three Al-NPs is 3100 K, 3300 K, and 3500 K, respectively, as given in Table II.However, because all the oxygen is not uniformly heated due to continuing oxidation dynamics inside and outside the Al-NP, the full system temperature is lower than the shell temperature.At the peak in (dT system /dt), the system temperature for the three Al-NP systems is 2500 K, 2700 K, and 2800 K, respectively.The time lag between when the shell temperature reaches $2900 K and the peak in time derivative of the total system temperature, (dT system /dt) specifies maximum rate of heat production is t 1 À t 0 ¼ 0.02, 0.04 and 0.09 ns for D ¼ 26, 36 and 46 nm, clearly indicating the Al-NP size effect.

FIG. 1 .
FIG. 1. Radial temperature profile of D ¼ 26 nm Al-NP during (a) preparation of aluminum core and alumina shell system, and (b) preheating of the system to 1100 K.The dashed line indicates the interface between the Al-NP and environmental oxygen.The arrows indicate the time sequence of the preparation and preheating.

)CFIG. 7 .
FIG. 6.(a) Schematic of a quarter slice of the whole system.Dark blue, grey, and yellow circles are aluminum atoms from alumina shell, aluminum atoms from aluminum core, and oxygen atoms, respectively.(b) Local concentration is calculated at a distance R by averaging the concentration within a spherical shell of thickness dR ¼ 10 A åt radius R from the center of the nanoparticle.
FIG. 10.Number of (a) Al 2 O, (b) AlO, and (c) AlO 2 clusters as a function of time in D ¼ 26 nm.

FIG. 14
FIG. 14.(a) Semi-log plot of the number of ejected Al shell atoms as a function of time for D ¼ 26, 36, and 46 nm.(b) Oxide shell temperature during oxidation as a function of time for the three Al-NP systems.Ejection of Al atoms from the shell starts when the shell is in molten state at 2900 K.

TABLE I .
Reactive molecular dynamics simulation setup for aluminum (core) and amorphous alumina (shell) nanoparticles in oxygen environment.
Table II shows the detailed comparison of Figs.15(a) and 15(c).

TABLE II .
Relationship between the delay time and the nanoparticle size.