Simulation Study and Proper Orthogonal Decomposition Analysis of Buoyant Flame Dynamics and Heat Transfer of Wind-Aided Fires Spreading on Sloped Terrain

1. Introduction

The wildfire activity has increasingly been reported in recent years [1], causing the catastrophic damage biologically or ecologically [2] and posing unprecedented fire threat to the communities near the wildland-urban interface [3], the study of which with a focus on the fire spreading mechanism, therefore become imperative. Wildfire is a very complex phenomenon, involving highly nonlinear interactions between fire dynamics, vegetation chemistry and environmental conditions (wind, humidity, and rough terrain, etc.) [4,5]. To predict reliably the fire behaviour and spread for emergency response and mitigation of forest fire risks, a further understanding of the basic influential elements that are responsible for the extreme or erratic fire spreading is necessary.

Wind is commonly accepted to be important to the dynamic processes of wildfire initiation and spreading. During the fire-atmosphere interactions, the wind effect is usually determined by two components, i.e., the atmospheric wind and the fire-induced wind. A field experiment conducted by Clements and Seto [6] showed that the fire-induced circulation was formed because of the ambient shear and thermal instability related to the sensible heat flux from the fire. The strong fire-induced wind was responsible for the enhanced turbulence kinetic energy and a minimum in atmospheric pressure. The spectral analysis of the measured velocity field during the fire-front passage indicated that the increased energy in velocity and temperature spectra at high frequency was caused by the dynamics of small eddies near the fire plume [7]. Kochanski et al. [8] suggested that with a strong vertical shear in the ambient wind field, the fire may develop into an erratic behaviour [9]. The role of coupling between fire and fire-induced flow was also demonstrated in the numerical simulations by Sun et al. [10] They found that the head fire rate of spread (ROS) and area burnt were substantially enhanced as much as doubled when the fire-atmosphere coupling was considered. However, a more reliable simulation is required for further understanding the physics underlying the fire acceleration since in the work of Sun et al. [10], the fire sensible flux was coupled to the atmosphere in an empirical way. Additionally, since the extreme wildfire usually occurs in the complex terrains [1,11,12], those wind effects also need to be quantified together with slope effect to fully untangle the extreme fire dynamics under the real environment conditions [13].

For the slope effect, based on the laboratory-scale experiment, Dupuy et al.[14] reported that the increasing in slope can greatly increase the ROS, and it was hypothesized to be related to the stronger fire-induced wind at high slope condition. Silvani et al.[15] further observed the appearance of fire whirls rolling along the flank of the fire front under the steeper slope. The author then stressed the importance of fluid dynamics to explain the remarkable change of fire spread rate. Meanwhile, a combined radiative-convective heating mode was deemed to dominate the heat transfer and thus the fire spread when the slope was increased [16]. However, as indicated by the experimental results of Xie et al.[17], for a steeper slope, the increased likelihood of flame attachment and therefore the convective heating seems to be the mechanism explaining the fire acceleration.

The ambiguous conclusion on the main heating mode over the unburnt vegetation may be because the convective energy transfer is difficult to be measured in the experiment [18], but a more quantitative examination of the heat transport is necessary to accurately evaluate the wind and slope effects. For most fire spread models of operational use, for example, the Rothermel model [19], an analogy has often been invoked between wind and slope effects; an effective wind speed was explicitly incorporated to account for the slope terrains based on the assumption that the fire spread mainly result from the increase of incident radiation [20]. This however limits the applicability of most empirically derived fire model for predicting the extreme fire spread under the complex environment conditions [14,21].

In general, due to the limitations of experiment, some important mechanism of fire spreading associated with the heat transfer and fluid dynamics is still unclear. In recent years, with the advancement of numerical methods, the computational fluid dynamics (CFD)-based analysis of wildfire spread has gained a great attention. Mell et al.[22] conducted a pioneer work to study the grassland fires using the physics-based model of WFDS (wildland-urban interface fire dynamics simulator), which showed a favourable capability in predicting physically the fire properties compared to the semi-empirical or empirical models. Adopting WFDS model, Sánchez-Monroy et al.[23] investigated the upslope fire spread for fuel beds at various slopes. The simulations reproduced the experimental trend and identified the convection-dominated regime as the slope ranging from 30° to 45°. However, in those cases, the effect of wind speed was still unknown. Moinuddin et al.[24] discussed the wind effect alone and revealed a linear relationship between ROS and wind speed in the considered cases. Furthermore, for now, the fluid dynamics associated with the change of spread mode (e.g., flame spread acceleration) that were indicated in the previous experiments have not been fully investigated.

Another advantage of physics-based modeling that should be noticed is the large amount of meaningful flow data. Processing these data based on the data-driven strategies can be helpful in revealing the coherent motions that govern the erratic fire behaviour and extracting the meaningful insights. One of the most popular data-driven fluid analysis methods is the proper orthogonal decomposition (POD) [25]. Mathematically, the POD is a decomposition technique that determines a set of orthogonal basis or mode functions for a given input data (or the target field variables), $X(x,t)$; these basis functions are identified as the optimal candidates for collectively expressing the original data. The critical features of turbulent flow, for example, can then be efficiently captured through the combination of the basis functions that represent the most energetic flow structures. This method was, therefore, often applied to examine the dominant fluctuating and coherent structures within the various types of fluid flow[26].

Kostas et al.[27] utilized POD to analyze the instantaneous velocity field and vorticity field of flow over a backward-facing step, while Podvin et al.[28] established a connection between the first few high-energy POD modes and the shear layer instability in the cavity flow. The POD analysis was successfully applied also to understand the flow-flame interactions, such as in Shen et al.[29] where POD identified the different interaction modes between the precessing vortex core and the swirling flame under different inlet fuel conditions. A POD analysis of instantaneous OH field and velocity-temperature field also revealed the locations of ignition kernels and flame growth for the flames stabilized in the hot-coflow configuration[30]. For the wild fire, although Guelpa et al.[31] once used POD for model order reduction in order to improve the computational efficiency of physical models solving the wildfire behaviour. However, to the best of the our knowledge, there is no attempt yet with POD to understand the coherent structures that govern the fire-induced flow and heat transport.

Therefore, this study aims to investigate systematically the characteristics of fire spread under the various wind and slope effects, and a focus is given to the understanding of the mechanism underlying the fire erratic behaviour. This is achieved by performing a detailed parametric study and a quantitative analysis of the fluid mechanics and heat transfer based on the physics-based numerical simulations [32]. The modal analysis facilitated by Proper Orthogonal Decomposition is also employed for the first time to examine the coherent structures that govern the propagating flame front in wildfire.

2. Computational methods

A multiphase formulation of mass, momentum and energy equations is considered for the present wildfire simulation to account for phase coupling during the fire spreading on the solid vegetation. The thermally-driven gas flow is solved by a low-Mach number LES approach based on the open source code of Fire Dynamics Simulator (FDS 6.7.7) [32]. Combustion in gas flow is modelled using the Eddy Dissipation Concept (EDC), and the radiation, as the main driving mechanism for fire spreading, is solved by finite-volume method [32].

2.1. Gas-Phase Governing Equations

Considering the gas-solid interactions, the fluid conservation equations for mass ($\rho$), momentum ($\mathbf{u}$), species ($Z_{\alpha }$) and sensible energy ($h_s$) are given as,

$${\displaystyle \frac{\partial \rho }{\partial t}}+\nabla ·\left(\rho \mathbf{u}\right)={\dot{m}}_b^{\prime \prime }$$

(1)

$${\displaystyle \frac{\partial \mathbf{u}}{\partial t}}-\mathbf{u}\times \omega +\nabla H-\tilde{p}\nabla \left({\displaystyle \frac{1}{\rho }}\right)={\displaystyle \frac{1}{\rho }}\left[(\rho -{\rho }_0)\mathbf{g}+{\mathbf{f}}_{\mathrm{b}}+\nabla ·\tau \right]$$

(2)

$${\displaystyle \frac{\partial }{\partial t}}\left(\rho Z_{\alpha }\right)+\nabla ·\left(\rho Z_{\alpha }\mathbf{u}\right)=\nabla ·(\rho D_{\alpha }\nabla Z_{\alpha })+{\dot{m}}_{\alpha }^{\prime \prime \prime }+{\dot{m}}_{b,\alpha }^{\prime \prime \prime }$$

(3)

$${\displaystyle \frac{\partial }{\partial t}}\left(\rho h_s\right)+\nabla ·\left(\rho h_s\mathbf{u}\right)={\displaystyle \frac{D\overline{p}}{Dt}}+{\dot{q}}^{\prime \prime \prime }+{\dot{q}}_b^{\prime \prime \prime }+\nabla ·{\dot{\mathbf{q}}}^{\prime \prime }$$

(4)

where $\rho$ is the gas density, $\mathbf{u}$ the velocity vector. $\tau$ is the viscous stress, and the term ${\mathbf{f}}_{\mathrm{b}}$ represents the drag force exerted by the fuel particles, $\omega$ is vorticity vector, $\tilde{p}$ pressure perturbation, $D_{\alpha }$ the species diffusion coefficients. ${\dot{m}}_b^{\prime \prime }$ and ${\dot{m}}_{b,\alpha }^{\prime \prime \prime }$ are the mass increment due to particle degradation, and ${\dot{m}}_{\alpha }^{\prime \prime \prime }$ stems from the chemical reactions. ${\dot{q}}^{\prime \prime \prime }$ is heat release rate, ${\dot{q}}_b^{\prime \prime \prime }$ denotes the energy transferred to particles, and ${\dot{\mathbf{q}}}^{\prime \prime }$ accounts for subgrid heat fluxes of conduction and radiation.

2.2. Vegetation Fuel Model

In the wildfire, the fire spread is controlled by the close interplay between the combined slope & wind actions and the solid fuel (vegetation) combustion[33]. Therefore, considering the vegetative fuel that is usually composed of fine particles, the fuel element (FE) model is employed, which represents the vegetation as thermally thin fuel particles that are inserted within each numerical grid with assigned physiochemical properties. Compared to the alternative model available in FDS, i.e., boundary fuel (BF) model, theoretically, FE model can solve more reliably the fuel thermal-degradation process affected by both the convective and radiant heat transfers[32] that are one of the primary focus in this study; contrarily, BF model assumes the fuel to burn like a porous solid and in a radiation dominated mode.

Therefore, a three-step thermal degradation is considered to model complex fuel pyrolysis and its interaction with the air flow. It includes three successive reactions, i.e., the water evaporation and then the pyrolysis of dry vegetation that are both endothermic, followed by the exothermic char oxidation:

Drying:

$$\mathrm{Vegetation}\Rightarrow {\gamma }_{{\mathrm{H}}_2\mathrm{O}}{\mathrm{H}}_2\mathrm{O}+(1-{\gamma }_{{\mathrm{H}}_2\mathrm{O}})\mathrm{Dry}\mathrm{Vegetation}$$

(5)

Pyrolysis:

$$\mathrm{Dry}\mathrm{Vegetation}\Rightarrow {\gamma }_{\mathrm{char}}\mathrm{Char}+(1-{\gamma }_{\mathrm{char}})\mathrm{Fuel}\mathrm{Gas}$$

(6)

Char oxidation:

$$\mathrm{Char}+{\gamma }_{\mathrm{char},{\mathrm{O}}_2}{\mathrm{O}}_2\Rightarrow (1+{\gamma }_{\mathrm{char},{\mathrm{O}}_2}-{\gamma }_{\mathrm{ash}}){\mathrm{CO}}_2+{\gamma }_{\mathrm{ash}}\mathrm{Ash}$$

(7)

where the stoichiometric constant of ${\gamma }_{{\mathrm{H}}_2\mathrm{O}}=M/(1+M)$, and M is the fuel moisture content determined based on a dry weight. ${\gamma }_{\mathrm{char}}$ is the mass fraction of dry fuel decomposed into the char, and ${\gamma }_{\mathrm{ash}}$ is the mass fraction of char that is converted to the ash during oxidation process. ${\gamma }_{\mathrm{char},{\mathrm{O}}_2}=1.65$ is the oxygen consumption coefficient [34]. The above three-step degradation of solid is solved like the element reactions in gas-phase, and then to calculate the mass generation rate of the product in each reactions, the Arrhenius type kinetics are used.

Associated with the mass change, the variation of solid vegetation temperature ($T_{\mathrm{s}}$) is modeled following,

$${\rho }_{\mathrm{b}}{c}_{\mathrm{p},\mathrm{b}}{\displaystyle \frac{d{T}_{\mathrm{s}}}{dt}}=-{\dot{Q}}_{\mathrm{s}.\mathrm{vap}}^{\prime \prime }-{\dot{Q}}_{\mathrm{s}.\mathrm{kin}}^{\prime \prime }-{\dot{q}}_{\mathrm{sc}}^{\prime \prime }-\nabla ·{\dot{q}}_{\mathrm{sr}}^{\prime \prime }$$

(8)

here, ${\rho }_{\mathrm{b}}$, $c_{\mathrm{p},\mathrm{b}}$ are the bulk density and specific heat for the solid phase. ${\dot{Q}}_{\mathrm{s}.\mathrm{vap}}^{\prime \prime }$ denotes the endothermic effect of water evaporation, and ${\dot{Q}}_{\mathrm{s}.\mathrm{kin}}^{\prime \prime }$ accounts the heats associated with pyrolysis and char oxidation. The convective heat transfer to the fuel surface is determined by ${\dot{q}}_{\mathrm{sc}}^{\prime \prime }={\sigma }_s{\beta }_{s}h(T_{\mathrm{s}}-T_{\mathrm{g}})$, where ${\sigma }_s$ is the surface area to the volume ratio, ${\beta }_s$ is the packing ratio defined as the volume of solid needles divided by the volume they occupy, and $T_{\mathrm{g}}$ is extracted from the first gas-phase grid point adjacent to the surface. h is the convective heat transfer coefficient [32].

The radiative heat flux is given by,

$$\nabla ·{\dot{q}}_{\mathrm{sr}}^{\prime \prime }=\kappa (4\sigma T_{\mathrm{s}}^4-U)$$

(9)

where the absorption coefficient is $\kappa =C_{\mathrm{s}}{\sigma }_s{\beta }_s$ and $C_{\mathrm{s}}$ is the shape factor. $\sigma$ the Stefan-Boltzman constant and U is the integrated radiation intensity in all directions.

In addition to the mass and heat transfers, the fuel elements also impose the drag force on the surrounding flow that is expressed as a force term (${\mathrm{f}}_{\mathrm{b}}$) in the gaseous momentum equation, given by

$${\mathrm{f}}_{\mathrm{b}}={\displaystyle \frac{\rho }{2}}{C}_d{C}_s{\beta }_s{\sigma }_s\mathbf{u}\left|\mathbf{u}\right|$$

(10)

where $\rho$ is the air density and $C_d$ the drag coefficient.

Due to the solution of fuel thermal degradation in a detailed way, a physics-based simulation and numerical analysis of fire spread is conducted. More detailed fuel thermo-properties adopted in this work can refer to the Table 1.

2.3. The Principle of Proper Orthogonal Decomposition

In this work, the Singular Value Decomposition (SVD) POD[25] is applied for a direct decomposition of the selected variable X of buoyant flame. As shown in the diagram of SVD method (Fig. Figure A1), the input matrix $\overline{X}$ is the instantaneous physical field and can be divided into average field X and fluctuation part $X^{\prime }$. Decomposed by the SVD, $X^{\prime }$ can be transformed into the product of three matrices:

$$X^{\prime }=U\times S\times V^T$$

(11)

where U is a left singular vector matrix, $V^T$ is a right singular vector matrix, and S is a diagonal matrix containing eigenvalue[35], in which the eigenvalue $\lambda$ is arranged in a descending order. The value of ${\lambda }_i$ on the diagonal is highly correlated with the energy contained in each mode and is considered as the amplification factors[36]. Energy content of each mode is calculated by ${ϵ}_k={\displaystyle \frac{{\lambda }_k}{{\sum }_{i=1}^N{\lambda }_i}}$.

$A=U\times S$ is the time coefficient matrix of POD, representing the change of each mode over time[35]. Usually, the energy associated with the first several modes is much higher. As such, the $X^{\prime }$ can then be reconstructed effectively by using a limited number of modes through

$$X^{\prime }\approx \sum _{i=1}^N{a}_i{\varphi }_i(k<N)$$

(12)

where $A=[a_1,a_2,\cdots ,a_N],V=[{\varphi }_1,{\varphi }_2,\cdots ,{\varphi }_N,\cdots ]$. In addition, those chosen modes often correspond to the critical characteristics of the original physical field.

2.4. Experiments and Numerical Setup

Referring to the experiment of Silvani et al.[15] and the DESIRE bench used in the INRA laboratory[14], a large-scale bench with 10 m long and 4 m wide is considered for modelling. Therefore, as depicted in Figure 1, the chosen computational domain for this study is a block region (13 m×6 m×5 m), and all borders are set as open boundaries except the bottom ground. Within this domain, there is a 8 cm thick plate of 10 m×4 m to hold the fuel bed. The leading edge of the plate is 1 m away from the open boundary, and the same is for its height from the ground.

The wind is introduced from the left side of the domain with the turbulence imposed using the Synthetic Eddy Methodology[22], where the turbulence intensity is set as 20% of the mean velocity. The wind is set as parallel to the surface of fuel bed even when the slope is changed. The fuel bed measuring 3 m×7 m is ignited as a line fire along the left edge by an ignition line of 3 m×0.1 m, which delivers a power of 500 kW for a duration of 5 s. To ensure a fully developed wind field across the domain, the ignition is delayed with 25 s after the start of the simulation. The calculated profile of wind speed at the downstream location of 6 m is shown in the Figure 2. The radiation fraction is determined from the smoke points ($L_{sp}$) as:

$${\chi }_{rad}=0.41-0.85L_{sp}$$

(13)

where $L_{sp}$ is set based on the measurement in the Fire Propagation Apparatus [37](i.e., $L_{sp}$ = 0.08 m).

The rate of fire spread is the important parameter to characterize how fast the fire propagates into the unburned region, for which the fire front needs to be first properly defined. However, there is still no rigorous definition due to the complex fuel and violent fire, although some prior works attempted to measure it based on the pyrolysis temperature (400 K) [33] or ignition temperature (623 K) [38]. In this work, due to the use of a detailed pyrolysis model, the pyrolysis characteristics of the fuel can be derived (as shown in Figure 3). We can observe that above the 560 K, the drying process of virgin fuel is initialized, and shortly afterwards, the intense decomposition reaction is induced, which indicates the usefulness of this temperature as the threshold to differentiate the thermal state of fuel bed. Thus, this temperature value has been adopted in this work to localize the instantaneous fireline.

To determine the appropriate grid resolution, a thorough analysis of the simulation results based on different sizes of mesh and domain is first performed: three sets of domain and grid sizes are examined, respectively, with the smallest mesh size being 0.02 m and the largest domain extending to 14 m in streamwise direction. It has been observed that (shown in Figure 4) the computed results are more sensitive to the mesh size compared to the dimension of the domain. Therefore, after a parametric study, a computational domain of 13×6×5 m and the grid size of 0.05×0.05×0.04 m have been chosen for the later discussion. It has been suggested that to resolve the fire spread on fuel bed, the grid size needs to be smaller than one-third of extinction length [24], which is 0.0547 m in the present study. The final mesh used is smaller than this value. Figure 5 depicts the simulated fire front at different instants under the ${30}^{\circ }$ slope and windless condition. In general, the simulation reproduces the time evolution of wildfire perimeter observed in experiments. Finally, to have a systematic study of the combined slope-wind effects on the fire dynamics, simulation cases with five different slopes (ranging from $0^{\circ }$ to ${30}^{\circ }$) and three wind speeds (0.5, 1.0 and 1.5 m/s) are considered (see Table 2). We also compared the computed firelines using the newest version of FDS (6.9.1), which shows a close match with present computations. Overall, the present fire modeling is adequate for a thorough analysis of fire dynamics under different wind & slope conditions.

3. Results and Discussion

3.1. Fire Perimeter and Rate of Spread (ROS)

Figure 6 first presents the time evolution of the fire perimeter extracted from the representative slope - wind conditions and the topmost cross-section of fuel bed (A complete results are provided in the Appendix of Figure A2), which are also overlaid with the iso-lines of convective and radiative heat fluxes received by the fuel bed. A dramatical change of fire topology is clearly seen for different ambient conditions. In most cases, a common pattern of U-shaped fire line [14] is observed, in particular, for the intermediate value of slope from ${10}^{\circ }$ to ${20}^{\circ }$. With the increase of wind speed (Figure 6b,c), the flame front is apparently broadened in space and meanwhile more stripe burning zones appear behind the front line. This indicates a faster propagation but incomplete combustion under the strong wind condition. As shown by Porterie et al. [39], the side wind can cause the continuous oxygenation of the fuel bed, promoting the combustion, but at the same time, the strong wind may results in the cooling effect behind the flame front.

A deep change of fire line shape is found for the slope effect (see Figure 6b,d). Apart from the transformation of fire shape from a U to a V, the flame front appears to be no longer compact but more distributed within the central region, indicating an acceleration of fire spread along the axis on the sloped terrain. As discussed in the earlier works[14,15], this fire behaviour could be ascribed to the tilted and elongated flame volume as well as a coupled effect from the fire-induced turbulent motions. The latter one may enhance the heat transport through the convection. As shown in plots, the convective isolines are more fluctuated under this condition. Furthermore, its distribution is close to the fire line, confirming the dominated role of turbulent heat convection in the vicinity of fire line. This is an important supplement to the existing studies since it is difficult to be measured in experiment. The radiation has often been assumed to govern the fire spread, but in most cases, the temperature rise of fuel particles and their subsequent ignition is convection controlled because the heat transfer coefficients for free and forced convection depend on the inverse function of characteristic surface length [40].

The larger spatial distribution of radiative heat flux (red isoline) in a comparison to the convection part (black isoline) beyond the fire front also implies that in this study, the radiation may dominates the preheating of the combustible bed. At the same time, it is interesting to note that the enlargement of the dominated area of radiation heat flux is more sensitive to the increase of wind speed than that of the slope.

This study also discovers another intriguing phenomenon: under low wind speed and sloped environments, the fire front deforms into a W-shaped structure with two flame heads (Figure 6a). This will be elaborated on in detail later.

Subsequently, to quantitatively evaluate the influence of slope-wind on the fire spread, the rate of spread (ROS) is calculated by post-processing the time-averaged fire lines under different ambient slope-wind conditions. Figure 7 illustrates the ROS computed for various slope and wind cases. Additionally, results calculated using the empirical models[22], such as the Rothermel model, CSIRO, MarkIII and MarkV[24] are also presented. The experimental data are only available for the windless condition. In general, ROS is increasing with the increase in slope, and this trend is almost same for three considered wind velocities. Specifically, the change of ROS shows a linear function at lower slope and when the slope is higher than ${10}^{\circ }$, a non-linear increase is observed for all side winds, which indicates the existence of a universal critical slope regardless of wind speed to enhance the conjunct effects of slope and wind and hence induce the eruptive fire[17].

The predictions from the empirical models generally depart from the present physics-based simulations. However, in comparison, except for the apparent deviation after ${10}^{\circ }$, MarkIII shows a closer match at small slope conditions ($0^{\circ },{10}^{\circ }$). This could be attributed to the simplified assumption and the specific experimental scenarios used for deriving these empirical formulas.

Figure 8 shows the variation of ROS as increasing the wind speed conditioned on different slopes. Compared to the observation shown in Figure 7, the fireline spreads more linearly with the increased wind. However, the acceleration of spreading fire seems to be less sensitive to the wind increment at the higher slope (e.g., ${30}^{\circ }$), which indicates the prevailing influence of slope even the wind speed is high. It also implies that the mechanism controlling the fire spreading by the increase of slope and wind, respectively, is different, which would be further explored in the later discussion.

3.2. Flow Field and Perturbation Pressure

To better understand the individual wind and slope effects on the fire spread process, Figure 9 illustrates the velocity field and perturbation pressure contour in x-y cross-section computed under three wind speeds without slope and three different slopes with 1 m/s wind speed, respectively. It can be observed that without slope, under only the effect of ambient wind, the fireline tends to exhibit a concave shape at the center region, forming a W-shaped structure with two heading flames. This can be attributed to the pair of vortices and a negative pressure zone formed in front of leading edge. With the increase of wind speed, however, the apparent big eddy dissipates because of the strengthened wind speed, and at the same time, the more pronounced enlargement of the negative pressure zone shows that the wind effect can be understood from the induced turbulent flow field, which develops a strong fluctuating buoyant plume with the fresh air largely entrained from the front side of flame. The low pressure zone is also observed experimentally by Clements and Seto [6]. They pointed out that the pressure minimum may be responsible for the observed acceleration of horizontal flows into the fire and the strong updraft, which will be shown in Figure 10.

In contrast, for the sloped terrain, even with ambient wind present, the region ahead of the fire front tends to be under the influence of positive pressure instead of negative one. Moreover, the strength of positive pressure increases with a steeper terrain. Therefore, it can be inferred that the concavity of the fireline at $0^{\circ }$ is mainly caused by the entrainment of lateral airflow and the backflow induced by the turbulent fire plume; the corresponding pressure drop suppresses the central fire spreading and alternatively, promotes the movement of the lateral heading fire. However, with an increase in slope angle, due to the inclination of flame, flame region expands more closer to the unburned area ahead, resulting in the higher possibility of high-temperature gas flowing over the unburned vegetation. This leads to the suppression of the backflow and a formation of high pressure zone. A faster spreading of the central fire heads is then expected, forcing the fireline to develop into a pointed V shape. As indicated by Hassan et al. [41], for the V shaped fire line, the interaction between the two junction arms will be another critical factor, enhancing the fire dynamics.

The flow field in x-z cross-section shown in Figure 10 clearly illustrates the complex change of flow pattern above the unburned fuel bed under the different ambient conditions. The dominate influence of wind (left-side plots) creates a complex turbulent flow, interacting with the flame-zone buoyancy. The slope effect, on the other hand, can be seen as the suppression of the indraft flow ahead of the flame due to the strong buoyant flow, and thus the promotion of the hot-temperature gases convected from the fire-front region to the inert vegetation fuel. In the field of establishing the operational model for fire spread prediction, the slope effect is usually modeled assuming that it acts like the ambient wind [42]. However, the present analysis implies that although both wind and slope could lead to the flame inclination, the fluid dynamics induced near the fire front would be totally different, especially the difference in flow of frontal area could impact the influential region of convective heating or cooling. Proper modeling of this difference may be the key to capture the erratic fire spread at higher slopes.

In addition, it is interesting to note the convergence of incoming wind field from the profiles of streamline in x-y cross-section (especially the case of ${30}^{\circ }$-1.0 m/s in Figure 9). It is associated with the stream-wise streaks of flame behind the flame front as reported by Finney et al. [40]. A saw-toothed flame geometry is then observed which is related to the complex flow dynamics of Taylor-Görtler instability and the vegetation smoldering combustion due to the fast fire spreading.

3.3. Flame Morphology

To understand the abrupt change in the speed of spreading fire shown in Figure 7, it is useful to analyse the flame topology that is key to determine the heat transfers to the vegetative fuel [15,17]. As revealed above, the ambient condition shows the determining effect on the fire propagation through the marked change of flame dynamics. Therefore, the specific parameters pertaining to the flame morphology and to easy understanding of the slope - wind effect are discussed. First, there are various methods to define the flame region, and in the experiments, the averaged luminous region captured by the camera is most often adopted to extract the flame shape[12].

However, in simulation, numerous parameters can serve as the indicator of flame zone, including the temperature, heat release rate, combustible gases and oxygen concentration[43]. Given that the combustion model is an Eddy Dissipation Concept (EDC) combustion model based on mixed fraction, where the chemical reactions are primarily depending on the mixedness of vapor fuel and the oxygen, we choose, therefore, the fuel gases as the variable to define the flame region. Based on the analysis of the flame profiles, it was found that the exothermic reaction zone is primarily located within the region where the fuel mass fraction $Y_f$ is less than 0.001. Therefore, this critical value is adopted as the threshold to define the flame shape.

By employing $Y_f$ = 0.001 as the threshold to create the binary image and extract the instantaneous flame front, and adopting the right-corner of flame base in x-z plane as the reference point, the contours of time-averaged flame probability for different conditions are obtained, as shown in Figure 11. It is seen that in general, apart from the case of $0^{\circ }$ - 0.5 m/s, the flames at other conditions tilt further towards the unburned area and meantime, the flame is elongated significantly.

The flame angle (FA) and length (FL) (illustrated in Figure 11) are two critical geometry parameters that can help illustrate the heat transfer processes, the driving force for the fire spreading. Following the earlier works [14,44], they are calculated in this study based on the same process as to define flame probability contours. Figure 12 shows the profiles of computed flame angle and length in a comparison to the experimental results of Guo et al.[44], Morandini et al.[45] and Mendes-Lopeser et al.[46] It can be observed that the increase in slope and wind speed leads to an increase in both flame length and angle, and particularly, after an almost linear change between $0^{\circ }$ and ${20}^{\circ }$, a sharp increase is found after ${20}^{\circ }$, indicating the existence of a critical slope angle for the transition of normal fire spread to the abrupt acceleration[17].

Comparatively, the present simulations reproduce the general trend reported in earlier measurements. This is particularly notable given that errors may arise from the differences in fuel bed parameters (e.g., fuel depth and bench size) and the method used to determine the flame shape. However, due to a more close match between the simulation and the experiments for FA, the result implies that the flame angle might has a smaller sensitivity to the fuel bed properties as also indicated in the work of Sánchez-Monroy et al.[23].

Overall, together with the variations observed on the flame probability contours, the remarkable changes of flame length and angle indicate that the higher wind and slope allows the flame to approach further the unburned fuel bed, increasing significantly the contact area for radiative heat transfer. This is manifested visually in Figure 6 for a larger impact area of radiative heating under the high slope & wind conditions. Meanwhile, the expanding flame volume and the plume-induced circulation flow would induce a substantial amount of high-temperature gas to flow over the unburned fuel ahead of fireline. This significantly enhances the convective heating. Thus, the enhancement of both convective and radiative heat transfer modes is often regarded as the driving force for the nonlinear behavior of fire spread under different ambient conditions[15,17,23]. In addition, at higher slope (30°), the flame inclines almost parallel (FA >60°) with respect to the fuel bed, which promotes the probability of flame attachment on the unburned fuel surface.

It was noted that the flame attachment could be the main reason for the development of the extreme or erratic fire spreading [17]. To have a better illustration, the flame attachment length is defined here by measuring the distance between the flame front position on the surface of combustible fuel particles and the pyrolysis front position, as shown in Figure 13 under different ambient conditions. In the abscissa, the $t_{\mathrm{ign}}$ denotes the duration of initial ignition stage. It can be seen that as expected, this attachment length is increasing with the increase in wind speed, however, a more pronounced change is observed for the slope effect especially when the slope angle increases from ${20}^{\circ }$ to ${30}^{\circ }$, which may explain the abrupt rise of ROS by factor 2−2.5 observed both in Figure 7 and in laboratory studies of Dupuy and Maréchal [14].

In addition, it is also worth noting that the FA and FL both seem to show high dependence on the slope of terrain while the FL is only more sensitive to the wind speed. This is consistent with the observations in Figure 9 and Figure 10; the wind effect can cause the backflow ahead of flame front, which however, is restrained by the high pressure zone at the sloped conditions and the slant flame or the tilted angle is therefore more pronounced for the slope effect.

3.4. Radiative and Convective Heat Transfers

A quantitative evaluation of thermal budget (i.e., convective and radiative heats) is conducted to further understanding how the heat transport is controlled by the slope and wind effects and its feedback on the fire-line dynamics. Figure 14 first presents the instantaneous profiles of the convective and radiative heat fluxes absorbed by the unburnt fuel near the fire front along the centerline. It is seen that in the initial stage of ignition, the absorbed heat budgets are large, after which the magnitude of received energies become smaller and their distributed area tends to be stable. Close to the fire front (denoted by the vertical dash-dot line), the influential domain of positive convection and radiant heating is almost overlaid with each other except the area far behind and ahead of the fire, where heat transport is dominated by radiation. The fire-induced turbulent convection enhances the intermittent presence of hot gases impinging on top of the inert fuel, which is even the overwhelming heat budget close to the fire front. This indicates the invalidation of the simplified assumption, i.e., the neglected convective heating, adopted in many empirical fire spread model [45]. In addition, the negative regime of both convection and radiation heating behind the fireline corresponds to the locations of exothermic char oxidation with high temperature, and the incoming wind flow introduces mainly the cooling effect.

To focus on the heating fluxes that are responsible for fire spreading, the time-average over a domain, ranging from the fire line to the isoline of zero convective heat flux ahead, is employed to obtain the averaged heat components (convection and radiation). Figure 15 presents the averaged heating fluxes of convection and radiation under different conditions. It is seen that the radiative heating shows a more dramatical increase with the rise of wind speed and slope. For case of 1.5 m/s wind velocity, however, its increase becomes more mildly with the change of slope, indicating the less dominated effect of slope at high wind speed condition. On the other hand, the convective component shows the considerable change only after the critical slope angle of 20°, which resembles more the trend indicated by flame angle in Figure 12.

This may be explained by that the different mechanisms could exist for the radiative and convective heat transports, respectively: radiation heat flux is controlled by the combined effects of flame length and view angle that determine together the amount of heat transport between the radiant surface of flame and the fuel bed. Conversely, the convective part of heating flux is dominated by the action of flame attachment and therefore, it is only more dependent on the tilted angle.

Therefore, we can see that the heat transport through the convection and radiation shows a different pathway during the change of wind speed and slope, within which the flame morphology plays a critical role. The ratio of their magnitude shown in Figure 15(c) further illustrates that the strength of convection is still smaller, but it is demonstrated that the enhancement of convection part of heating would be the main factor for accelerating the fire spread under the steeper slopes (especially above 25°)[17]; in case of 30°, the convective heating accounts for almost half of the heat gain. In addition, the decreased role of radiation part when slope angle rises above 20° is consistent with the observation in experiments of Dupuy and Maréchal [14].

To further elucidate the mechanism governing the heat transport for fire propagation, the Byram convective number ($N_c$) [47], denoting the ratio between flame buoyancy and inertial force, is calculated as follows:

$$N_c={\displaystyle \frac{2g{I}_B}{\rho C_{p}T(u-ROS)}}$$

(14)

where g is the gravity, and $\rho$, $C_p$, T denote the density, specific heat capacity and temperature of ambient air, respectively. u the wind speed and $I_B$ is the fire line intensity. For the scenarios considered in this work, the calculated $N_c$ values are presented in Table 3. For $N_c$ <2, the fire is wind-driven, and $N_c$ >10 denotes a buoyancy-dominated fire. A mixed regime for fire propagation is deemed to occur at 2 <$N_c$ <10.

It is seen from the Table III that the effects of wind and slope differ from each other: with the increase of wind speed, the propagation mode tends to be wind-driven because of the decrease in $N_c$, and conversely, the fire is more strongly buoyant with the increase of $N_c$ when the slope is larger. In the meantime, the transition of the fire mode appears to be more sensitive to the change of wind at low slope condition (<${20}^{\circ }$), in which a more dramatical change of $N_c$ is observed. In general, for most cases, the fire considered in this work is plume-dominated with $N_c$ greater than 10. The exception is indicated for cases with wind speed 1.5 m/s, where the $N_c$ is relatively small and the flame propagation is driven more likely by the wind. This may explain the mild change of thermal budget for 1.5 m/s wind speed shown in Figure 15(a)(b) at different slope conditions. In addition, the abrupt change of $N_c$ at slope ${30}^{\circ }$ indicates a stronger fire plume that would increase the flame attachment and accelerate the flame spread.

3.5. Proper Orthogonal Decomposition (POD) Analysis

In this paper, POD analysis is carried out for this works, as shown by the energy content of each mode (k) in Figure 16(b), the energy associated with the first several modes is usually much higher.

Figure 16(a) illustrates the reconstructed temperature field by mode 1 of propagating fire front at slope-wind condition (${20}^{\circ }$-0.5 m/s). It reveals the dominant pulsating features of gas temperature (in red and blue) at the flame base, interacting with the coherent vortex structures (in yellow).

To understand the coherent flame structures under different environmental conditions and their interaction with the turbulent flow fields, the temperature-velocity fields at the wind speed of 1 m/s and angles of ${10}^{\circ }$, ${20}^{\circ }$, and ${30}^{\circ }$ and three wind speeds without slope are selected as the input data for the POD analysis. The energy distributions of the first 15 modes and their cumulative energy obtained for three different conditions are shown in Figure 16(b). The mode 1 occupies an absolutely dominant position, with the second and third modes contributing to a partial share of total energy (approaching to 10%). Therefore, in this study, the POD analysis mainly relies on these first three modes.

The reconstructed temperature and velocity fields are then depicted in Figure 17 for various wind speed conditions. The dominant mode 1 clearly reveals the primary structural feature of the fluctuating flame front that is in a back and forth mode (in blue and red colours, respectively). The lower energetic modes 2 and 3 show a similar structure, but with three distinct regions and smaller spatial distribution. In general, the spatial distribution of different modes increases and the slant of flame becomes apparent when the side-wind becomes stronger. In the 1.5 m/s case, it is interesting to note that the positive part of pulsating temperature locates in front. This means that the influx of fresh air and their cooling effect caused by backflow ahead of the flame front is suppressed at higher wind condition, which is also indicated in Figure 10.

For different sloped terrains, it is shown in Figure 18 that the positive and negative coherent structures under different slope angle are similar to the stronger wind condition in Figure 17, and in the meantime, for the same wind speed, the mode pattern of ${10}^{\circ }$, 1.0 m/s is distinct from its counterpart with flat terrain in Figure 17. This indicates under the slope effect, the intense burning of fuel bed exists near the fire front and thus the buoyancy flow and the incoming wind will enhance the temperature rise in the direction of fire plume. As the slope increases, it is seen that the flame zones tilt further with an enhanced forward pulsation.

For the modes 2 and 3, however, the dominant structure is not observed and it seems that for those modes, the corresponding flame dynamics are controlled mainly by the small scale of vortex that generates the fluctuation in the direction of buoyancy flow, which is different from that observed in Fig. Figure 17. This further demonstrates that the mechanisms of slope effect controlling fire spread is different from that of wind speed. For the steeper slope condition, those structures of turbulent fluctuation tend to enlarge in space, indicating a higher turbulent intensity, which will then strength the effects of turbulent convection and heating in front of fire zones. This explains why there exists a critical slope for the acceleration behavior of propagating fire.

Under the condition of 30° and 1 m/s (Figure 18(d)), the velocity field confirms the similar phenomenon with fluctuating structures resembling those of the temperature field. The same phenomenon can be seen at different wind speeds.

4. Conclusion

This paper presented a detailed numerical simulation of wildfire, investigating comprehensively the main characteristics of the fires spreading under the combined effects of wind and slope. More specifically, the rate of spread, flame morphology, convective and radiative heat transfers as well as the fluid dynamics were studied quantitatively in detail, focusing on the fire spreading mechanisms controlled by slopes and wind speeds, for which the physics-based modeling of interactions between vegetal fuel and flame as well as the POD analysis were employed. The main conclusions are as follows.

1.

A power-law relationship was indicated existing between ROS and the slope angle. It was revealed that at high slope conditions, the convergence of incoming wind and the weakened indraft air from the frontal area made a significant contribution to the abrupt rise of ROS and the eruptive spread of head fire.

2.

The enlarged volume of fire plume was deemed to enhance the radiation heat transfer, and in contrast, the higher possibility of flame attachment at higher slopes (especially >20°) led to the prominent role of convective heating.

3.

The investigation into the joint temperature-velocity field utilizing POD approach revealed an increased forward pulsation of the flame front with the escalating slope, leading to a higher energy density in the pre-combustion zone ahead of the fireline that further explained the mechanism underlying the accelerated flame propagation.

4.

For wildfire modeling, the more decent model should distinguish the respective roles of wind & slope, where the slope has a more profound effect in terms of determining flame structures and convective heat; the unsteady feature of flame puffing could be incorporated, considering the dominated mode pattern of back-forward pulsation.