## 2 Materials and Methods
The thermocouples are known to typically provide underestimated values of temperature due to their susceptibility to remove heat and interfere with the combustion process [22]. These simulations can offer some insight into the removal of heat, but the heat released will be unaffected by their presence. The domain comprises both fluid and solid mesh structures, with the respective thermal properties shown in Table 1. These properties are drawn from the ANSYS Fluent material databases and the GreenTape(tm) DuPont DP 951 material data sheet. The silicon, reactive foil, ceramic, platinum and silver structures have solid states assigned, whereas the air and solder have fluid status and utilise fluid models in the ANSYS Fluent solver. The melting/solidification model is activated in ANSYS Fluent for the solder, in order that the transitions between the molten and solid behaviour, and the associated latent heat transfer throughout, are modelled. The pure solvent melting heat for the solder is 58.5 kJ/kg. A liquidus temperature of 217 \({}^{\circ}\)C and solidus temperature of 220 \({}^{\circ}\)C were used for the solder. These values correspond closely to typical thermal properties of a SAC solder. The heat released by the reactive multilayer system was modelled in ANSYS Fluent through a heat source profile in the form of Equation (1). This heat source profile is a probability density function (PDF) propagating in the positive x direction, starting at the leading edge of the reactive multilayer. The reaction is assumed to initiate at time = 0 s and continually propagates as a function of time. \[Ae^{B(x+C-Dr)^{2}} \tag{1}\]  
The coefficient \(A\) corresponds to the amplitude of the reaction wave, 6 \(\times\) 10\({}^{12}\), \(B\) the width of the PDF, \(-\)8 \(\times\) 10\({}^{6}\), \(C\) is an offset to align with the leading edge, 2 mm in the \(x\) direction, and \(D\) is the velocity of the reaction. For the results presented here, a reaction speed, \(D\), of 1 m/s has been modelled. The reaction speed is known to vary with respect to several factors, for example the ignition potential [21]. In [20], the general range of these multilayer reactions was quoted as being between 0.1 and 100 m/s, therefore a reaction speed of 1 m/s is comfortably within this range. The entire domain was initialised to a temperature of 300 K, in other words no preheating was used. Figure 7: Detailed dimensions in x–y plane. In ANSYS Fluent, the pressure-coupled transient solver is used in combination with the implicitly formulated Volume-of-Fluid (VOF) model complete with 2 Eulerian phases, namely air and solder. The 2-equation k-\(\Omega\) SST model is used for turbulence modelling with the energy equation also solved numerically. The melting/solidification model, which makes use of the enthalpy-porosity approach [36], was used with a mushy zone parameter of 100,000. This method does not track the melt front explicitly but rather computes a liquid fraction, which indicates the fraction of the cell volume that is in the molten state, i.e., a liquid fraction of 1 corresponds to molten solder, whereas 0 is solder in the solid state. An inlet boundary condition was specified on the y-max boundary with a velocity of 0.1 m/s, with a corresponding outlet pressure boundary condition at 0 Pa (gauge pressure) on the opposing y min extremity. The turbulence intensity and viscosity ratios were specified to be the same, 0.1% and 10, on each boundary. These boundaries enabled the domain to have incoming and outgoing air, and this was found to be beneficial from a stability standpoint. Given that this crossflow inlet velocity was low, it was considered that artificial crossflow cooling effects would be negligible. The PISO pressure-velocity coupling method was used with 2 nd-order upwind schemes for spatial discretisation. Likewise, for the transient formulation, 2 nd-order schemes were used, only this time they were bounded implicit schemes.