1. INTRODUCTION

In the relatively new technical field of nanoenergetics, a great deal of attention and research has been paid recently to investigating various types of nanoenergetic materials and nanostructured reactive systems, in efforts to develop novel applications. Most of the nanoenergetic materials and related system design concepts that have proven to be highly promising are metal based, which show extremely high energy density and reaction sensitivity [¹]. They could also be useful in improving the system performances of many weapons and spacecraft components and miniaturizing those pyrotechnic systems, compared to conventional energetic materials and systems.

Although metal based nanoenergetic materials and systems can be realized in several different forms, there are presently two categories being intensively investigated. The first are metallic nanoparticles and nanocomposites. As an example, highly reactive nano-sized aluminum powders have been used as metallic additives to solid propellants to expand rocket propulsion performance by increasing burning rate, higher specific impulse, and more complete combustion [²,³]. Also, composites of nanoscale metallic and metal based oxide particles commonly known as nano-thermite (such as aluminum nanopowders mixed with nano-sized molybdenum trioxide or copper oxide flakes), have demonstrated enhanced energetic performance, with extremely high exothermic heat release and considerably shorter ignition delay [⁴-⁶].

The other group of nanoenergetic concepts being extensively studied are based on nanoscale reactive multilayer foils (NRMF) or nano-laminates [⁷]. These are composed of multiple alternating nano-sized layers of two different metals, usually fabricated by physical vapor deposition methods of magnetron sputtering and electron beam evaporation, or sometimes by mechanical processing such as hot pressing and cold rolling of stacked foils [⁴,⁷]. The most popular intermetallic combinations of nanoscale multilayer foils are aluminum-nickel (Al-Ni) and borontitanium (B-Ti) nano-laminates. Several experimental reports have demonstrated the self-sustaining propagation of highly exothermic reaction waves, whose propagation speed easily exceeds ten meters per second. The ignition delay of such intermetallic reaction propagation is very small (much less than a microsecond) and efficiently controllable.

Because of its superior pyrotechnical characteristics, NRMF has been proposed for various applications, including the joining of metallic specimens or electronics die soldering, airbag initiators, and alternative heat sources for thermal batteries [⁷,⁸]. The NRMF is also uses potentially in a miniaturized ignition device to initiate a secondary reaction in nearby energetic materials in propellants, explosives, and other pyrotechnic materials [⁹-¹³].

There have been a significant amount of theoretical and computational studies attempting to understand the fundamental physical mechanism of the reaction ignition and self-propagating intermetallic reaction waves in bimetallic multilayer microsystems. These include studies with Al-Ni [¹⁴-¹⁷] and B-Ti nano-laminate foils [¹⁸] using the timedependent transient intermetallic reaction model, or the quasi-steady one-dimensional evolution approximation model. These previous numerical investigations were able to predict the reaction wave speed quite accurately compared to corresponding measurements, and reveal the reaction propagation mechanism. However, there are many characteristic details in NRMF studies to be explored, such as the periodic wave behavior in self-sustaining reaction wave propagation, or the effects of a pre-mixed zone at the bimetallic interface. Therefore, the present numerical investigation focuses on understanding those detailed characteristics of self-propagating reaction waves in Al-Ni based NRMF systems.

2. NUMERICAL MODELING

The ignition of the intermetallic reaction and the subsequent self-sustaining propagation of reaction waves in NRMF were numerically modeled and simulated in the present study. As graphically illustrated in Fig 1, the simple geometry of the NRMF consists of alternately placed thin metallic layers of two different metals whose layer thicknesses are in the nanoscale, and usually much less than 100 nm. In this study, the bimetallic multilayers were made of aluminum and nickel. This combination is expected to have an exothermic intermetallic reaction of Al + Ni → NiAl once ignited. The stoichiometric ratio is one to one on a molar basis. By considering the mass densities and molar mass of those metallic species, the bimetallic thickness ratio was found to be approximately Al:Ni = 1.5:1 or t_Al/t_Ni = 1.5.

Depending on the input loading direction of the ignition source, the reaction waves may travel in a horizontal direction along the bimetallic layers, or in a transverse direction across the bimetallic multilayers. The case of reaction propagation in the transverse direction can be simplified into one-dimensional modeling [¹⁸,¹⁹]. However, this numerical study analyzes and elucidated the reaction propagation in the horizontal direction, as depicted in Fig 1. In this case the numerical model of the intermetallic reaction wave propagation should be two-dimensional. By employing the periodic geometry of the multilayers, the computational domain include only half of each metallic layer, and this is shown as the dotted rectangular area in the Fig 1.

The numerical model for the simulation of the ignition and self-propagation of the intermetallic reaction in Al-Ni based NRMF can be represented by the following governing equations for the atomic species concentration for aluminum and nickel, and the temperature distribution [⁷,¹⁵-¹⁸].

In the above equations of transient and two-dimensional mass and thermal diffusion, C(x,y,t) and T(x,y,t) represent the spatial distributions of the metallic species concentration and temperature, respectively, and t is the time. Note that x is the coordinate of the reaction propagation direction, while y represents the vertical coordinate across the multilayer. The last term in the thermal diffusion equation in Eq. (2) represents the local exothermic heat generation rate from the intermetallic reaction. Here, ΔH_f is the heat of reaction for the intermetallic reaction that occurs between Al and Ni, and its value is 330 cal/g or 1.38 MJ/kg [²⁰]. The subscript notations for mass density (ρ), specific heat (c), and thermal conductivity (k) in Eq. (2) specify the thermal properties in the corresponding metallic layers (i = Al, Ni). Since temperature varies significantly due to the high exothermic reaction in a very short width of reaction wave, the specific heat and thermal conductivity are regarded as functions of local temperature in order to accommodate the steep temperature gradient, as listed in Table 1. However, melting of the metallic layers is not considered here, assuming that its effect is negligible and mass density is given as a constant. The initial temperature distribution is uniformly given as 300 K at a time of t = 0.

For atomic species concentration, C = 1 is meant for pure aluminum and C = -1 is for pure nickel. The species concentration will go from 1 or -1 toward zero upon the progress and completion of Al-Ni intermetallic reaction. The initial and pre-defined distribution of atomic species concentration at the initial time before the reaction can be expressed in the form of a Meander function or square waves from the alternating arrangement of pure aluminum and pure nickel, as depicted in Fig 2(a). However, varying degrees of pre-mixing or inter-mixing of the two atomic species should exist at the interface between the two metallic layers even before the start of the reaction, due to unavoidable atomic diffusion. This would apparently decrease the ignition sensitivity and reaction intensity, as its effects could be significant especially in the present nano-sized system of the NRMF.

Thus, the numerical modeling in this investigation includes the effects of a pre-mixed zone on the formation and selfpropagation of the intermetallic reaction wave. This is illustrated as the thickness of 2w between the Al and Ni layers in Fig 1. Gavens et al. [¹⁴] gave a simple treatment by assigning C = 0 in the pre-mixed zone. In contrast, some previous studies on NRMF [¹⁵-¹⁷] have approximated the pre-mixed zone as linearly varying species concentration from 1 to -1, which is shown in Fig 2(b) and the following formulation.

where t_Al and t_Ni are the layer thicknesses of aluminum and nickel, respectively.

However, a more realistic variation within the pre-mixed zone would be a more gradual change in species concentration, due to the slow atomic diffusion at room temperature. Therefore, this study employs a sinusoidal form for the initial distribution of species concentration, as portrayed in Fig 2(c), and can be expressed as

The intermetallic reaction between the Al and Ni multilayers is assumed to follow a diffusion-limited reaction. Hence, the above atomic diffusion equation employs the following binary atomic diffusion coefficient based on the Arrhenius dependence with temperature [²⁰].

Here, E_a and D₀ are the activation energy and the Arrhenius pre-exponential factor and R is the universal gas constant, respectively. The characteristics of the reaction wave, such as reaction initiation and propagation speed, are expected to heavily depend on those factors. Therefore, it is important to provide accurate values for a specific set of intermetallic reaction. Those reaction factors are usually estimated by linear fits of measured data from a series of high precision reaction experiments [¹⁴], but it is quite difficult to obtain such data for various bimetallic combination of NRMF [¹⁸]. However, the estimation of reaction factors has been very well established for the Al-Ni combination, and the Arrhenius reaction factors for the Al-Ni reaction used here are given as E_a = 137 kJ/mol and D₀ = 2.18×10^-6 m²/s [¹⁵].

The present study assumes that reaction wave propagates through a semi-infinite multilayer structure of NRMF, but the actual computational domain should be finite. Thus, the computational domain is made sufficiently long in order to capture the characteristics of the self-propagating reaction waves unaffected by the finite size of the domain length in the x-direction. The boundary conditions for the governing equations are given such that the spatial gradients of temperature and species concentration are set to zero on each side of the computational domain, which indicates there is no diffusive loss out of the computational domain (in the x-direction) or geometric symmetry (in the y-direction).

For ignition of the intermetallic reaction in NRMF, several different approaches have been successfully applied in experimental tests, including electrostatic spark discharge, DC electrical heating, laser irradiation, and mechanical impact [⁷]. Most of those methods can be regarded as a thermal ignition, where the intense heating and sufficiently elevated temperature on the targeted NRMF surface initiates the reaction ignition.

Several previous investigations have modeled thermal ignition by imposing an initial profile of high temperature on a certain boundary for a prearranged period of time [¹⁴-¹⁷]. However, imposing heat flux would be a more realistic way of simulating reaction initiation, since the external heat sources usually employed in experimental tests consist of joule heating of a bridgewire or laser heating on the target surface of the NRMF specimen [¹⁸]. Accordingly, in this numerical study, a constant heat flux is prescribed on the boundary of x = 0 for a predetermined initial time period for all simulation cases. Note that this boundary becomes adiabatic for the entire remainder of the computational time after this time period.

3. RESULTS AND DISCUSSION

In order to investigate and better understand the propagation characteristics of the self-sustaining intermetallic reaction wave in Al-Ni based NRMF, a series of numerical simulations were conducted using the aforementioned transient reaction model on the two-dimensional computational domain of the Al-Ni multilayer system. Although this study tested the effects of various Al-Ni layer thicknesses, a bilayer spacing of L = 20 nm (t_Al + t_Ni) was set for the base case, and accordingly the individual layer thickness of Al and Ni were t_Al= 12 nm and t_Ni= 8 nm, respectively, since the stoichiometric thickness ratio of Al-Ni reaction is t_Al/t_Ni = 1.5. The computational domain consisted of each half of the Al and Ni layers, since the reaction wave travels along the layer direction. The length of the domain should be long enough to eliminate the effect of an arbitrary adiabatic condition, and here was set to be 200 μm. The width of the domain was half of the bilayer spacing, or 10 nm for the base case, so the aspect ratio was 2000 which should be sufficiently long.

The computational mesh system is a structured quadrilateral type with a size of 40×1000. While the 1000 meshes were placed uniformly in the x-direction, finer meshes were applied in the y-direction near the interface of the two metallic layers in order to adequately resolve the temporal change in the steep gradients for temperature and species concentration which occur during the reaction at the bimetallic interface. The mesh sensitivity was thoroughly tested to ensure that the present mesh system was acceptable, by checking the propagation speed of the reaction wave with several different mesh systems, of 20 × 500, 20 × 1000, 40 × 500, 40 × 2000, and 80 × 2000. All of the numerical computations were efficiently carried out by employing the multi-purpose FEM-based software package, COMSOL Multiphysics 3.5.

Knepper et al. [²²] carried out experiments on Al-Ni based NRMF made by DC magnetron sputtering with an Al-Ni thickness ratio of 1.5, as mentioned above. They prepared Al-Ni multilayered nanofoils whose bilayer spacing ranged from 10 to 200 nm. They used a DC spark as the ignition method and measured the speed of the self-propagating reactions by detecting the visible radiation emission. The experimental study also provided the thickness of the pre-mixed zone at the Al-Ni interface, which was 2.4 ± 0.3 nm of the total thickness or w = 1.2 nm in the present model.

This numerical study will use their experimental reaction wave propagation speed results to the validate computational model. To simulate the ignition of the intermetallic reaction by thermal stimuli sources such as pulsed laser ignition, the incoming heat flux of q_in = 10 MW/cm² was given at the left boundary of x = 0. The time period of spatially uniform heating was set to be 20 ns for the base case. This is equivalent to a laser fluence of 200 mJ/cm², which is in the same order of laser ignition intensity for the reaction initiation tests of the reactive multilayer system by Picard et al. [²³].

The computational results for a base case of bilayer spacing L = 20 nm and a pre-mixed zone thickness w = 1.2 nm are shown in Fig 3, in a time progression form of temperature contours on the computational domain at 12 time steps. The intermetallic reaction was initiated quite immediately and the ignition delay was very small. The reaction wave can be easily identified by the elevated temperature behind the thin reaction zone, close to the Al-Ni reaction temperature of 1912 K. It then propagates at a constant speed of 12.7 m/s into the unreacted domain. Since the thickness of the Al-Ni bilayer is relatively small with high thermal conductivity, the temperature is uniform across the layer, and the individual metallic layers are not distinguishable in the temperature contours.

The remnants of the prescribed heat flux for the reaction ignition can be seen at a time of 1 μs, but it disappears quickly as time progresses. The ascribed heat flux input for ignition was carefully minimized to a level just enough to initiate the reaction but not large enough to affect the establishment of self-sustaining reaction waves in the present computational domain of finite size. For a heat flux input lower than a 12 ns time period, the reaction wave propagation does not occur, and the heat load of 100 mJ/cm² dissipates without triggering the reaction ignition.

Figure 4 shows the spatial changes in the atomic species concentration distribution at 11 lateral locations (x = 93 to 103 μm) around the reaction zone at the time that the reaction wave passes by (t = 7.6 μs). It can be seen from the Fig 4 that, the step function shape of atomic species concentration evolves into a flat distribution in the reaction zone, and the thickness of the reaction zone is estimated to be approximately 5 μm.

The effects of the pre-mixed zone thickness were investigated for a bilayer spacing of L = 20 nm, both with linear and sinusoidal models of atomic pre-mixing, as shown in Fig 5. The experimental study [²²] reported a measured speed of approximately 13 m/s for the reaction wave with a pre-mixed zone thickness of w = 1.2 nm. The sinusoidal model produced a reaction wave speed of 12.7 m/s, which is in excellent agreement with the measured value, while the linear model yielded a significant lower speed of approximately 9.9 m/s. As the thickness of the pre-mixed zone increased, reaction propagation speed dropped quite significantly. When there was no or a marginally narrow premixed zone (less than w = 0.02 nm), the reaction wave speed approached approximately 24.4 m/s in the simulation, and it was twice that speed for the case of w = 1.2 nm. When the pre-mixed zone thickness increased to over w = 3 nm, the reaction wave speed dropped below 5 m/s and went to zero or no ignition at all.

The experiments by Knepper et al. [²²] also tested different sizes of bilayer spacing in Al-Ni NRMF, from L = 10 to 200 nm. In Fig 6, the reaction wave speeds measured in those experiments were compared with the numerical results from several computations with different bilayer spacing, but with the identical pre-mixed zone thickness of w = 1.2 nm. As expected, the increase in bilayer spacing resulted in a substantial decrease in reaction wave speed.

However, when pre-mixing exists, a critical bilayer spacing for maximum obtainable reaction wave speed may exist. Both the measurements and present modeling with a sinusoidal model of atomic pre-mixing indicated the reaction wave speed reaches a maximum of approximately 13.4 m/s at the critical bilayer spacing between 15 and 20 nm. Once the Al and Ni layers become thinner than this range, the reaction wave speed sharply decreases toward zero or no ignition. In contrast, the linear model does not yield reasonable predictions when bilayer spacing is less than 40 nm. This comparison confirms that the sinusoidal model of interface pre-mixing is more suitable than the linear model in previous modeling studies, and for that reason the following computational results will only use the sinusoidal model for treating the pre-mixed zone.

In Fig 7, the effects of various pre-mixed zone thicknesses are shown for bilayer spacings up to 180 nm. If the premixed zone is nonexistent or relatively small (w = 0 or 0.5 nm in the figure), a much faster propagation of the reaction wave can be expected. However, in reality, the maximum reaction wave speed would be severely limited by the premixing zone thickness. For example, if the pre-mixed zone thickness was as large as w = 3 nm, the maximum reaction wave speed would occur at L = 40 nm, and be slightly greater than 5 m/s.

In order to explore the detailed characteristics of the reaction wave propagation in NRMF, the time-evolving location of the reaction wave was estimated by assessing the positions where the atomic species concentration first reached values of +0.5 for the Al layer, and -0.5 for the Ni layer, counting from the reaction starting point of x = 0. It is revealed from the Fig 8, for a bilayer spacing of L = 20 nm with varying pre-mixed zone thicknesses, the difference in reaction wave locations based on C = +0.5 in the Al layer or C = -0.5 in the Ni layer is negligible. In fact, the reaction wave speeds listed in Figs 5Fig. 6.-7 were evaluated by applying linear fitting on the time-progressing reaction wave locations in such a manner. As discussed before, the presence of the pre-mixed zone significantly slows the propagation of the reaction wave.

In case of no pre-mixed zone (w = 0), the temporal change in the reaction wave location is very linear. Fig 8 also shows the distinguishing characteristic of periodic wave unsteadiness in the reaction propagation in some cases. This periodic nature of repeated acceleration and deceleration in the reaction wave propagation is more clearly visible, especially when the premixed zone thickness is relatively larger.

To examine the details of periodic reaction wave behavior, time-space contours of atomic species concentration (along the symmetrical centerline of Al layer) and temperature (along the Al-Ni interface line) were plotted for three cases of pre-mixed zone thickness (w = 0, 1.2, and 2 nm) in Figs 9(a) to (f). Figs 9(a), (c), and (e) show that the thickness of the reaction zone becomes larger with higher pre-mixing at the interface. When the periodic reaction wave unsteadiness occurred for cases of w = 1.2 and 2 nm, small hot spots behind the reaction wave could be observed, as shown in Figs 9(d) and (f), and the size of the hot spots grows with higher degrees of pre-mixing.

At this time, an even higher pre-mixed zone thickness of w = 3 nm was numerically tested in the computations with an identical bilayer spacing of L = 20 nm. As discussed in Fig 5, the maximum allowable pre-mixed zone thickness for the given bilayer spacing would be roughly 3.5 nm. The ignition of the self-sustaining reaction wave would fail with a premixed zone thickness higher than this value, even if ignition inducing heat flux increases further. w = 3 nm was quite close to this allowable level of pre-mixing, and in fact the reaction wave propagation slowed down significantly, as shown in Fig 5. In this case, the periodicity of the reaction wave is expected to be prominent, and one cycle of periodic reaction unsteadiness occurs in a mu ch wider lateral range of the computational domain. This makes it necessary to expand the length of the computational domain from 200 μm (as in previous computational cases) to 400 μm to accurately capture the much widened reaction wave periodicity.

Figure 10 shows the temporal changes of reaction wave speed locations in Al-based NRMF. In Fig 10, the numerical results of the above case (w = 3 nm) are shown as the time trace of the reaction wave location, along with the previous cases of no or smaller pre-mixed zone thicknesses (w = 0 to 2 nm). While the periodic wave unsteadiness shows up prominently here, the time-averaged reaction wave speed from the linear fitting is approximately 3.3 m/s, and the reaction propagation slows down very much compared to the previous cases. After the reaction ignition at an early time stage, it takes approximately three to four cycles of periodic unsteadiness to establish the regular periodic manner. One cycle of reaction wave periodicity spans approximately 30 μs in time or 75 μm in length.

Figures 11(a), (b), and (c) present time versus space contours w = 3.0 nm with three kinds of time-space contours: atomic species concentration, temperature, and reaction heat generation rate. The periodic occurrence of hot spots, discussed previously, can be seen more clearly here. Some of them are marked A and C in these figures. The hot spots are associated with a sudden burst in the intermetallic reaction and a subsequent much higher level of heat generation, as shown in Fig 11(c). This leads to the abrupt propulsion or acceleration of reaction wave propagation to great speed (from A to B). The reaction heat diffuses backwards from the hot spots to the reacted zone, which can be observed in Fig 11(b). After the reaction wave acceleration wears down, the reaction moves forward very little with much decreased reaction and heat generation rates (from B to C) until the heat builds up around the slow moving reaction wave front and the other hot spot abruptly shows up again.

The periodic manner of reaction wave propagation can be appreciated in Fig 12 as time progressing temperature contours at 14 time steps, which cover approximately one cycle of periodic wave unsteadiness (from A to C in Fig 11). The hot spot starts to occur faintly at t = 40 μs and becomes clearly noticeable at t = 41-43 μs. It propels the reaction wave to high speed, up to t = 50 μs. This high speed movement of the reaction wave due to hot spot burst may lead to an incomplete reaction in the area swept by the fast moving reaction wave. This is followed by a relatively inactive period with little moving reaction wave in the time period of 50 to 68 μs, until the time of another hot spot at t = 70 μs.

The present computational simulations on the NRMF system of Al-Ni nano-laminates positively confirm that the numerical predictions are consistent with the corresponding measurements. It also demonstrates an extremely high level of reactivity, with the propagation speed of the reaction wave reaching 10 to 20 m/s.

The existence of atomic pre-mixing affects the reaction wave speed significantly, and induces complicated patterns in the reaction wave, so it should be carefully considered in application developments, such as microscale reactive bridge initiators and others.

4. CONCLUSIONS

This numerical study focused on computational simulations and predictions of self-sustaining reaction wave propagation in micron-sized NRMF systems. A computational model involving intermetallic atomic species diffusion and thermal diffusion was established, along with exothermic heat generation which assumes Arrhenius reaction kinetics. The atomic pre-mixing at the bimetallic interface was also modeled while assuming linear or sinusoidal changes in the pre-mixed zone thickness.

The numerical results of the reaction wave speed in the Al-Ni NRMF system with an Al-Ni thickness ratio of 1.5 were compared with the measured values, for bilayer spacing from 10 to 180 nm. There was good agreement between computations and experiments. The sinusoidal treatment of atomic pre-mixing provided more accurate predictions than the other pre-mixing models. For a pre-mixed zone thickness of ±1.2 nm, there is a critical bilayer spacing, which resulted in a maximum reaction wave speed of 13.4 m/s, consistent with measurements. For a bilayer spacing of 20 nm, the effects of pre-mixing were studied with a series of computations. With small or without pre-mixing, the reaction wave propagated linearly. However, as the pre-mixed zone thickness increased up to ±3 nm, the reaction wave slowed down significantly, and the numerical results revealed periodic unsteadiness in the reaction wave propagation. This periodic cycle features superfast reaction wave acceleration, due to a sudden burst of hot spots behind the reaction front, followed by a period of reaction wave deceleration, possibly due to incomplete reaction during the acceleration period.