LES Simulation of Wind-Driven Wildfire Interaction with Idealized Structures in the Wildland-Urban Interface

Abstract

This paper presents a numerical investigation of the impact of a wind-driven surface fire, comparable to a large wildfire, on an obstacle located downstream of the fire source. The numerical modelling was conducted using FireFOAM, a coupled fire-atmosphere model underpinned by a large eddy simulation (LES) solver, which is based on the Eddy Dissipation Concept (EDC) combustion model and implemented in the OpenFOAM platform (an open source CFD tool). The numerical data were validated using the aerodynamic measurements of a full-scale building model in the absence of fire effects. The results highlighted the physical phenomena contributing to the fire spread pattern and its thermal impact on the building. In addition, frequency analysis of the surface temperature fluctuations ahead of the fire front showed that the presence of a building influences the growth and formation of buoyant instabilities, which directly affect the behaviour of the fire’s plume. The coupled fire-atmosphere modelling presented here constitutes a fundamental step towards better understanding the behaviour and potential impacts of large wind-driven wildland fires in wildland-urban interface (WUI) areas.

1. Introduction

Wildfires are a major natural hazard, which can have disastrous consequences from socio-economic and environmental points of view. Given the apparent increase in large wildfires [1] and continued urban expansion into wildland areas [2], there is a clear need for further study of the vulnerability of built assets to fires impacting the wildland-urban interface (WUI). Such studies underpin further development of comprehensive wildfire risk management strategies, including improved urban planning.

Modern Computational Fluid Dynamic codes have played an increasingly important role in understanding the physical processes driving fire behaviour and in improving our ability to predict it. In fact, computational tools are widely used to complement laboratory-scale studies and to compensate for the fact that the results of laboratory-scale studies are difficult to apply directly to large wildfire situations [3]. The first fully physical multiphase wildfire model was developed by Grishin [4]. This model has formed the basis for the development of more advanced, fully physical wildfire models capable of handling additional physicochemical phenomena (e.g., WFDS (Wildland Fire Dynamics Simulator) [5], and FIRETEC code [6]).

Indeed, over the last two decades, significant progress has been made in the development of physics-based wildfire models. Some examples include Meroney et al. [7], who presented a simplified computational model to investigate fire spread within the idealized generic porous models of 3D city structures and Filippi et al. [8,9,10], who proposed an approach to couple a fire area simulator to a mesoscale weather numerical model to simulate local fire/atmosphere interaction. The authors numerically compared the coupled approach between the meso-non-hydrostatic (NH) large eddy simulation (LES) mesoscale atmospheric model and the ForeFire wildland fire area simulator with experimental data to evaluate the performance of the proposed coupled approach in forecasting fine-scale properties of the dynamics of wildland fires.

ForeFire remains a front-tracking model similar to FARSITE, PROMETHEUS, ELMFIRE, etc., where the key principle is to consider a propagating front to represent how the active flame zones propagate over landscapes. This requires a parameterization of the rate of spread, and ForeFire can use either the Rothermel model [11], Finney et al. [12], or the Balbi model [13]. In fact, Rothermel’s model is considered as a semi-empirical rate-of-spread law, while Balbi’s model is considered as a quasi-physical rate-of-spread law. Still, both remain a simplified representation of the reality by using the fire front paradigm.

Generally, fire models range from tools based on the Rothermel model [11], fire spread rate formulas, such as BehavePlus [14] and FARSITE [15], suitable for operational prediction, to sophisticated 3D computational fluid dynamics and combustion simulations appropriate for research and reanalysis, such as FIRETEC and WFDS. PC-based solvers such as BehavePlus predict the fire spread rate at a single point of fuel and a set of environmental data; FARSITE employs the fire spread rate to provide a two-dimensional simulation on a PC; whereas FIRETEC and WFDS need a parallel supercomputer and run a lot slower than real time. Among coupled atmosphere-wildland fire solvers, WRF-Fire [16] consists of a fire-spread model, achieved by the level-set method, coupled with the Weather Research and Forecasting model (WRF). Kochanski et al. [17] evaluated and improved the performance of the coupled atmosphere-fire model WRF-SFIRE and showed WRF-SFIRE is capable of predicting head-fire rate of spread and fire plume vertical temperature up to 10 m above ground level.

Similarly, a series of papers [18,19] outlines the development and application of a numerical physics-based multiphase model called FIRESTAR3D. FIRESTAR3D is a 3D finite volume model based on a fully physical multiphase model to simulate the ignition and propagation of grassland fires, crown fire and wildfires at local scales (<500 m) [18]. The impact of wind velocity and terrain features on large scale fires was extensively studied using FireFOAM in a series of papers [20,21,22]. The authors successfully developed a numerical model to express the dynamic characteristics of the fire in terms of fire-induced pressure and viscous forces.

While full physical (and other) models have been used extensively to gain a better understanding of the fundamental behaviour of wildfires [11,19,23,24,25], mechanisms of fire spread in WUI fires have not received as much research attention and so are not as well-understood. Recently, however, there have been some efforts to study fire configurations that are more relevant to WUI fires. Linn et al. [26] developed a fast-running numerical tool based on integration of a phenomenologically-based fire spread model with a fast-running atmospheric solver initially designed for application to urban environments. This tool allowed for more computationally efficient accounting of fire-atmospheric feedbacks and provides fire planners with the opportunity to compare, assess, and design prescribed burn plans to assist with WUI fire management more rapidly. Despite its relative simplicity, this model showed a similar capability to capture the underlying physics and basic trends in fire behaviour to more sophisticated models such as FIRETEC.

Hilton et al. [27] presented a novel low-fidelity numerical model for the calculation of radiative heat flux loadings on buildings in the WUI and demonstrated its utility in risk assessment of infrastructure in areas impacted by wildfires. The evolution of the fire front in this approach was simulated using a quasi-steady model for the fire’s rate of spread, which means that it does not properly account for dynamic fire behaviour and fire-atmosphere interactions. Similarly, Cohen [28] proposed a theoretical model based on ideal heat transfer characteristics to investigate direct flame heating leading to building ignition during wildland fires. This study also determined the critical flame-to-structure distance, at which wooden walls would readily ignite. Mell et al. [29] conducted a comprehensive review of approaches commonly used for addressing problems concerning WUI fires. They concluded that physics-based modelling approaches can be reliably used to determine structure exposure conditions, and to develop new building test methods and standards.

A detailed numerical study based on a fully physical model (Fire Dynamic Simulator, or FDS) focusing on the effect of wind-driven fires on a building was first carried out by He et al. [30]. The numerical model was validated in the absence of fire using mean pressure measurements made on the full-scale, cube-shaped Silsoe building [31]. It is worth mentioning that the Silsoe experiment was carried out in idealized situations (e.g., cubic building, flat terrain) and provides a set of important benchmark data about the fundamental aerodynamic behaviour of low-rise buildings under wind loads. Several attempts have also been made to predict high-risk areas in the WUI using high-fidelity physics-based CFD models. Fryanova et al. [32] used WFDS to study the effects of fire intensity and wind speed on the likelihood of structure ignition. A crown fire was represented as an equivalent heat and mass source at an appropriate height, and the resulting hydrodynamic and thermal interactions between the plume and wind flow were analysed to assess the possibility of ignition of buildings. Similarly, Pimont et al. [33] used FIRETEC to investigate the minimal “clearing” distance required between buildings and forest in WUI areas to maintain safe working zones for firefighters and to minimize the chance of structure ignition. They found that a clearing distance of 50 m is an acceptable threshold for both materials and firefighters for both thermal radiation and gas temperature.

Understanding flow pattern due to freestream velocity variation associated with the interaction of fire and cross-wind flow around the building, and also variation of the structure’s surface temperature due to radiative heat flux, are of great importance [34,35] because the consequences may have major implications in building design against wildland fire attacks [36,37,38]. However, a fundamental understanding of how the interaction of fire and wind can alter the time-dependent flow patterns of velocity variation and the surface temperature of the building has remained elusive.

In the CFD approach, the fire behaviour could be modelled through the computation of algebraic equations arising from fundamental physical laws and constitutive equations to accurately model fire and its impact on structure, properties and so on.

This study explores the fundamental mechanisms of how the interaction of horizontal cross wind flow with a vertical plume alters the velocity profile around the building and the surface temperature of the building downstream of the fire source. FireFOAM was used as a CFD solver in this study. This solver is a derivative of the OpenFOAM [39] platform, specifically designed for fire dynamic simulations. FireFOAM is a large eddy simulation (LES) fire dynamics transient solver for large-scale fires and turbulent diffusion flames with reacting particle clouds, surface film and pyrolysis modelling. Despite its capability in large-scale fire dynamic modelling, it has not yet been employed in wind-driven wildfire-urban interface simulations, hence this study aimed to provide further fundamental insights into the interaction between a wind-driven fire and a building, thereby improving our understanding of the mechanisms involved in fire propagation and structure loss in WUI areas using FireFOAM code. The code is used to investigate the thermal loading on a building exposed to a wildfire that is encroaching on the wildland-urban interface. Dynamic fire behaviour is also discussed, and observations of the flame structure are described. The present study provides valuable information relevant to the identification of sources for structure ignition, enhanced turbulence (important for ember generation, transport and accumulation) and strong winds/pressure gradients that could affect the integrity of structures, allowing for greater ember incursion and other damage. Indeed, it is worth noting that the effect of strong pyrogenic winds is not generally considered as part of existing risk management doctrine (including building standards, asset protection zones (APZs), etc.).

2. Model Description

In the present paper, a 6 × 6 × 6 m cube is used to represent an idealized building impacted by a wind-driven wildfire. The geometric dimensions are the same as those used in the full-scale Silsoe cube experiment [31]. The cubic building is surrounded by a 50 × 30 × 25 m computational domain as illustrated in Figure 1. According to the recommendations made by Richards and Norris [40], the domain boundaries in this configuration are far enough from the building to minimize adverse boundary condition effects. A fire source running across the entire domain was placed 20 m upstream of the building to mimic a line fire configuration. Methane with a fixed mass flow rate corresponding to a specified heat release rate (HRR) of 180 MW was selected as the fuel source. In the wildfire literature, it is more common to use fire line intensity instead of HRR [41]. As such, in this study we consider a fire with a fixed fire line intensity of 6 MW/m, equivalent to a typical wildfire scenario involving a fire propagating with a rate of spread of 0.75 m/s through a fuel load of 0.4 kg/m², which represents “average” grassland fuels.

A power-law inflow velocity profile based on Equation (1) was imposed at the domain inlet. Random noise was also superimposed with the mean velocity to reconstruct the turbulent fluctuations.

To lower the cost linked with simulation of the transition process, simulation of spatially evolving turbulence begins with an inflow boundary that is shifted in a shorter distance upstream of the domain of interest [42]. In an ideal situation, having adequately precise inflow conditions at this boundary that result in a more realistic turbulent boundary layer with accurate skin friction is obtained within a short distance downstream. However, practically this is not usually feasible, and the inflow boundary might have to be moved further upstream to let error relaxation assist in estimating the inflow conditions. Having such a “development section” increases the overall cost of the simulations [42]. One simple and yet effective method to specify turbulent inflow conditions and to minimize the cost associated with generating the inflow data is to superpose random fluctuations on a desired mean velocity profile [43]. Despite its simplicity, the random fluctuation approach to introduce perturbations to the flow has been used with varying degrees of success [44,45].

$${\mathrm{u}}^{*}=\frac{\mathrm{U}\left(\mathrm{Z}\right)}{{\mathrm{U}}_{\mathrm{ref}}}={\left(\frac{\mathrm{Z}}{{\mathrm{Z}}_{\mathrm{ref}}}\right)}^{\mathsf{\gamma }}$$

(1)

In Equation (1), Zref is the reference height, which is taken as equal to the building’s height (6 m) and Uref is the reference velocity equal to 6 m/s. The exponent $\mathsf{\gamma }$ is a function of terrain characteristics and is assumed to be 0.16 in the present study based on the terrain category of the experiment site [46]. A typical atmospheric pressure condition (free inflow/outflow at atmospheric pressure) was applied for the top and outlet boundaries, while the side boundaries were treated as free slip boundaries. A wall model as suggested by Launder and Spalding [47] was applied for the near-wall treatment of turbulent flow. This self-adapting wall function allows direct resolution of the near-wall eddies for a fine wall grid corresponding to the dimensionless wall distance parameter of less than 5 (i.e., in the case of y+ < 5 m, the turbulent boundary layer is fully resolved up to viscous sublayer). It should be noted that y+ is a dimensionless parameter representing the first cell size near the wall. Further information regarding the y+ criterion can be found in [48]. For more information regarding implementation of boundary conditions, readers are referred to [19,49].

A high-quality grid was generated using successive zonal refinement. Following this strategy, a subdomain of size 22 × 20 × 12 m was defined around the building to ensure that the complex vortical flow structures formed behind the building were captured accurately (see Figure 1). The second level of grid refinement was set to resolve the near-wall regions around the building. In addition, a spanwise refinement was defined throughout the entire domain providing a high-resolution mesh near the floor. Further information on the grid characteristics will be given in the following section.

2.1. Mathematical and Numerical Modelling

Numerical simulations were conducted using FireFOAM [50], which is a powerful large eddy simulation (LES)-based finite volume solver that is available as a C++ library in the OpenFOAM CFD package. FireFOAM incorporates a number of efficient CFD submodels to describe processes such as radiant heating, pyrolysis, combustion and turbulence. It has been successfully used in many practical applications such as solid fuel pyrolysis [51], fire suppression [52] and fire-wall interaction [53]. FireFOAM has also been demonstrated to be an effective tool in wildfire modelling [25]. The Favre-filtered formulation of the fully compressible Navier-Stokes equations representing the fire dynamics in the most common form is written as a set of conservation equations of mass, momentum, energy and chemical species mass fraction. A reduced single-step combustion composed of the ${\mathrm{CH}}_4,{\mathrm{O}}_2,{\mathrm{N}}_2,{\mathrm{CO}}_2,{\mathrm{H}}_2\mathrm{O}$ components has been assumed (CH₄ + 2O₂ + 7.5N₂ → CO₂ + 2H₂O + 7.5N₂ [54]).

$$\frac{\partial \overline{\mathsf{\rho }}}{\partial \mathrm{t}}+\nabla \cdot \left(\overline{\mathsf{\rho }}\tilde{\mathrm{u}}\right)=0$$

(2)

$$\frac{\partial \overline{\mathsf{\rho }}\tilde{\mathrm{u}}}{\partial \mathrm{t}}+\nabla \cdot \left(\overline{\mathsf{\rho }}\tilde{\mathrm{u}}\tilde{\mathrm{u}}\right)=-\nabla \overline{\mathrm{p}}+\nabla \cdot \left[\left(\mathsf{\mu }+{\mathsf{\mu }}_{\mathrm{sgs}}\right)\left(\nabla \tilde{\mathrm{u}}+{\left(\nabla \tilde{\mathrm{u}}\right)}^{\mathrm{T}}-\frac{2}{3}\left(\nabla \cdot \tilde{\mathrm{u}}\right)\mathrm{I}\right)\right]+\overline{\mathsf{\rho }}\mathrm{g}$$

(3)

$$\frac{\partial \overline{\mathsf{\rho }}\widetilde{{\mathrm{Y}}_{\mathrm{k}}}}{\partial \mathrm{t}}+\nabla \cdot \left(\overline{\mathsf{\rho }}(\tilde{\mathrm{u}}+{\tilde{\mathrm{u}}}_{\mathrm{c}}){\tilde{\mathrm{Y}}}_{\mathrm{k}}\right)=\nabla \cdot \left[\left(\frac{{\mathsf{\mu }}_{\mathrm{sgs}}}{\mathrm{Sct}}+\overline{{\mathrm{j}}_{\mathrm{k}}}\right)\nabla \widetilde{{\mathrm{Y}}_{\mathrm{k}}}\right]+\overline{{\dot{{\mathsf{\omega }}_{\mathrm{k}}}}^{‴}}\mathrm{k}=1,\dots ,{\mathrm{N}}_{\mathrm{s}}-1$$

(4)

$$\frac{\partial \overline{\mathsf{\rho }}\widetilde{{\mathrm{h}}_{\mathrm{s}}}}{\partial \mathrm{t}}+\nabla \cdot \left({\overline{\mathsf{\rho }}\tilde{\mathrm{u}}\tilde{\mathrm{h}}}_{\mathrm{s}}\right)=\frac{\mathrm{D}\overline{\mathsf{p}}}{\mathrm{Dt}}+\nabla \cdot \left[\left(\mathsf{\alpha }+{\mathsf{\alpha }}_{\mathrm{sgs}}\right)\nabla {\tilde{\mathrm{h}}}_{\mathrm{s}}\right]-\nabla \cdot \overline{{\dot{{\mathrm{q}}_{\mathrm{r}}}}^{″}}+\overline{{\dot{{\mathrm{q}}_{\mathrm{c}}}}^{‴}}+\nabla \cdot (\mathsf{\alpha }{\displaystyle \sum }_{\mathrm{k}=1}^{{\mathrm{N}}_{\mathrm{s}}}\left({\mathrm{h}}_{\mathrm{s},\mathrm{k}}\nabla {\tilde{\mathrm{Y}}}_{\mathrm{k}}\right))+\nabla \cdot ({\displaystyle \sum }_{\mathrm{k}=1}^{{\mathrm{N}}_{\mathrm{s}}}\left(\overline{{\mathrm{j}}_{\mathrm{k}}}{\mathrm{h}}_{\mathrm{s},\mathrm{k}}\right))$$

(5)

In these equations, $\mathsf{\rho }$, p, u, ${\mathrm{h}}_{\mathrm{s}}$, Y and Sct are the density, pressure, velocity, sensible enthalpy, mass fraction and turbulent Schmidt number, respectively. The symbols − and ~ denote spatial and Favre filtering [55], respectively. $\overline{{\dot{{\mathrm{q}}_{\mathrm{r}}}}^{″}}$ is the total radiative heat flux (W/m²) of the gas mixture. $\overline{{\dot{{\mathrm{q}}_{\mathrm{c}}}}^{‴}}$ is the heat generated by combustion and is expressed as the reaction rate of fuel multiplied by the heat of combustion:

$$\overline{{\dot{\mathrm{q}}}_{\mathrm{c}}^{‴}}=\overline{{\dot{\mathsf{\omega }}}_{\mathrm{F}}^{‴}}.\mathsf{\Delta }{\mathrm{H}}_{\mathrm{c}}$$

(6)

Finally, J_k in Equation (4) denotes the mass-diffusive flux accounting for the Soret effect (thermodiffusion) and is calculated using Hirschfelder-Curtiss correlation [56], in which both advective and diffusive mass fluxes are taken into consideration. This avoids possible underestimation error in the mass flow rate calculation at the inlet fuel boundary [50]. The mass diffusion coefficient for each component is written as

$${\mathrm{D}}_{\mathrm{k}}=\frac{\mathsf{\mu }}{{\overline{\mathsf{\rho }}\mathrm{Sc}}_{\mathrm{k}}}$$

(7)

It should be noted that the mass diffusion coefficients, D_k, were assumed to be constant and written as a function of the corresponding Schmidt number, Sc_k. The values of the Schmidt number for all components in the methane combustion process are reported in [57]. In the present LES approach, the wall-adapting local eddy (WALE)-viscosity method was employed for subgrid scale (SGS) modelling. This method has shown remarkable results, particularly in wall-bounded flows [58]. According to this model, the SGS kinetic energy, SGS viscosity and corresponding turbulent mixing time are defined as follows:

$${\mathrm{k}}_{\mathrm{sgs}}=(\frac{{\mathrm{C}}_{\mathrm{w}}{}^2{\mathsf{\Delta }}^2}{{\mathrm{C}}_{\mathrm{k}}})\frac{{({\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}}{\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}})}^3}{{({({\tilde{\mathrm{S}}}_{\mathrm{ij}}{\tilde{\mathrm{S}}}_{\mathrm{ij}})}^{5/2}+{({\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}}{\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}})}^{5/4})}^2}$$

(8)

$${\mathsf{\nu }}_{\mathrm{sgs}}=\left({\mathrm{C}}_{\mathrm{w}}{\mathsf{\Delta }}^2\right)\frac{{\left({\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}}{\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}}\right)}^{3/2}}{{\left({\tilde{\mathrm{S}}}_{\mathrm{ij}}{\tilde{\mathrm{S}}}_{\mathrm{ij}}\right)}^{5/2}+{\left({\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}}{\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}}\right)}^{5/4}}$$

(9)

$${\mathsf{\tau }}_{\mathrm{sgs}}=\frac{\left({\mathrm{C}}_{\mathrm{w}}{\mathsf{\Delta }}^2\right)}{{\mathsf{\nu }}_{\mathrm{sgs}}}=\frac{{\left({\tilde{\mathrm{S}}}_{\mathrm{ij}}{\tilde{\mathrm{S}}}_{\mathrm{ij}}\right)}^{5/2}+{\left({\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}}{\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}}\right)}^{5/4}}{{\left({\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}}{\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}}\right)}^{3/2}}$$

(10)

where ${\overline{\mathrm{S}}}_{\mathrm{ij}}$ is the rate of strain tensor, ${\mathrm{S}}_{\mathrm{ij}}{}^{\mathrm{d}}$ is a special tensor defined in [58], ${\mathrm{C}}_{\mathrm{w}}$ and ${\mathrm{C}}_{\mathrm{k}}$ are the model constants, ∆ is the filter width and ${\mathrm{S}}_{\mathrm{ij}}^{\mathrm{d}}$ is a specific form of the rotation rate tensor [58]. The well-established eddy dissipation concept (EDC) [59] was applied for combustion modelling. The EDC model has shown great superiority in modelling well-ventilated diffusion flames [50]. In the EDC, the combustion is controlled by turbulent mixing assuming infinitely fast reaction kinetics. Thus, the fuel mass reaction rate, $\overline{{\dot{\mathsf{\omega }}}_{\mathrm{F}}^{‴}}$, is calculated as

$$\overline{{\dot{\mathsf{\omega }}}_{\mathrm{F}}^{‴}}={\mathrm{C}}_{\mathrm{EDC}}\frac{\overline{\mathsf{\rho }}}{{\mathsf{\tau }}_{\mathrm{EDC}}}\min \left({\tilde{\mathrm{Y}}}_{\mathrm{F}},\frac{{\tilde{\mathrm{Y}}}_{{\mathrm{O}}_2}}{\mathrm{rs}}\right)$$

(11)

where ${\mathrm{C}}_{\mathrm{EDC}}$ is the model coefficient with a default value of 4 and ${\mathsf{\tau }}_{\mathrm{EDC}}$ is the mixing time-scale. ${\tilde{\mathrm{Y}}}_{\mathrm{F}},{\tilde{\mathrm{Y}}}_{{\mathrm{O}}_2}$ and rs are the fuel mass fraction, the oxygen mass fraction, and the stoichiometric oxygen-to-fuel mass ratio, respectively.

Radiation modelling was conducted using the finite volume discrete ordinate method (fvDOM) [50]. To solve the radiative transfer equation (RTE) accurately, 64 solid angles were used for angular discretization. In addition, the weighted-sum-of-grey-gases model (WSGGM) was used to calculate absorption/emission coefficients [60].

A fully implicit second-order Gauss linear scheme was used in space while a first-order Euler method was used as the temporal scheme. FireFOAM takes advantage of the PIMPLE algorithm for pressure-velocity coupling [61]. Making use of a SIMPLE iteration loop embedded in a PISO (Pressure-Implicit with Splitting of Operators) algorithm, the PIMPLE algorithm ensures numerical stability with large Courant numbers (i.e., CFL = U∆t/∆x). This is found to be of significant importance in LES calculations where the grid resolution is high and local Courant numbers might exceed the stability limit (CFL = 1). In the current study, 5 iterations in the outer SIMPLE loop plus 2 iterations over the inner PISO loop (with pressure correction) are considered to ensure the required level of numerical convergence. In addition, the convergence criterion is met when the L2-norms of all transport equation residuals (in normalized form) reach 10⁻⁵.

2.2. Computational Model Verification

Two sets of experimental data were employed to validate the numerical model of the current study. The first experiment involved measuring the pressure on the vertical and horizontal centrelines of the Silsoe cube of Richards and Hoxey [35]. In this work [35] the full-scale data were sorted out in non-overlapping blocks of cube surface pressure along with the reference upstream flow method measured at the cube height. Because in a full-scale system natural wind conditions are considered, each block of data is exclusive and there is minor or zero chance to accurately repeat a test multiple times, which is not the case in the wind tunnel. As a result, these researchers could process the data in a way that constructs a large number of blocks recorded. The second experiment used the experimental investigation of Castro and Robins [62], who examined the flow around a surface-mounted cube in uniform, irrotational and sheared turbulent flows. The third benchmarking case involved the results of a numerical study on bushfires, wind and structure interactions by He et al. [30], who predicted pressure coefficient distribution around a building under no-fire conditions and compared this with the experimental results measured by others [30].

As mentioned, the Richards and Hoxey experimental data pertaining to the full-scale Silsoe cubic building [31] were used to validate the numerical model in the absence of fire effects. As stated earlier, since the EDC is mostly controlled by turbulent mixing, the fire dynamic characteristics are closely tied with the building aerodynamics. Thus, the pressure distribution over the building in a no-fire scenario is an appropriate indicator showing the accuracy and validity of the numerical results. To minimize numerical uncertainties, a grid sensitivity analysis was conducted. To this end, three sets of different grid sizes were generated. A comparison in terms of numerical settings and computing time on a 24-core node for different grid resolutions is presented in

Table 1. Computational grid specifications and numerical settings.

ID	Number of Cells (Millions)	y+	Timestep (s)	Computing Time (Hours)
grid A	4.6	37	0.004	16
grid B	7.8	1.65	0.0015	62
grid C	9.5	0.95	0.0015	71

. All simulations were carried out for a total time of 20 s after fuel injection. It is worth mentioning that a 20 s preliminary run was performed to achieve quasi-fully developed turbulence before the fire source was activated.

The numerical predictions of the mean pressure coefficient along the centrelines over the front, top and rear surfaces of the building (shown by 0–1–2–3 solid lines) are plotted in Figure 2. For a more conclusive comparison, the results of previous numerical studies are also plotted in Figure 2 with those calculated from the three different grids considered in this study. While there is a satisfactory agreement between the experimental and numerical data on both windward and leeward faces (0–1 and 2–3 lines), noticeable discrepancies still exist across the roof face (1–2 line). These discrepancies, which were also present in the other studies included in the figure, can be attributed to the high scatter of wind velocity and direction data, which considerably affects the corner separation region and pressure distribution over the roof. In addition, since the fluctuations of wind pressure on the building are mainly due to large-scale motions [60], the computational domain size could limit the size of the largest eddy and adversely affect the results. While the results for all grid densities give almost the same level of accuracy, it should be noted that smaller scales should also be resolved to simulate turbulent energy exchange in reaction zones accurately [61]. Thus, grid B with a high resolution near the walls was finally selected for further study while grid A was sufficient for aerodynamic performance prediction.

3. Results and Discussion

The present study investigated the simulated effect of a dynamically varying wind-field on a stationary heat source, representing a WUI fire configuration. The simulations permitted investigation of the transient fire behaviour that could be expected. Figure 3 shows the time-dependent flow patterns of velocity variation in the fire domain. As shown, the simulated fire front is tilted towards the ground downstream of the fire bed. This is mainly due to the fact that the inertia forces are dominant in this region. In addition, air entrainment into the turbulent surface fire creates a low-pressure region in the flame-ground spacing ahead of the fire. As a result, the fire plume is accelerated and forms a jet-like structure in a relatively small region almost attached to the floor immediately downwind of the fire source. This explanation is consistent with the calculated variation of the velocity profile shown in Figure 4. According to this figure, an intense velocity gradient is observed near the ground (z < 2m) for a specific streamwise length in which a dominant wind-driven fire regime exists. For this fire regime, a large portion of the incoming air passes over the fire plume and is mildly affected by the fire.

Moving further downstream of the fire location, the turbulent boundary flow enters a region where there is no clear balance between buoyant and inertial forcing and as a result, the flame shows random time-dependent motions. The pulsating behaviour of the fire plume propagation leading to the fire intermittency can be reflected in the extent of the attachment region. While the plume attachment region ahead of the building oscillates with a low amplitude, sharp fluctuations (in 27 m < x < 35 m) are observed in the flame-ground contact interface at the adjacent flat areas (see Figure 3b,d,f). Finally, the fire plume encounters a buoyancy-induced reverse flow (coloured in blue in Figure 3) and is driven upward. This phenomenon is often referred to as buoyant instability [63]. As a result, a strong recirculation zone appears, and the low-momentum fire plume eventually separates from the surface. Indeed, the buoyancy-driven fire regime becomes predominant in this region. This effect is further enhanced by the presence of the building blocking part of the incident flow. As shown in Figure 4, beside the building (y = 9 m), where the plume is unimpeded, the plume remains attached to the surface for as far as 35 m downwind from the fire front. In contrast, the plume in the midplane (y = 0 m), where the building acts to block the flow, has already separated from the surface at a distance of only about 17 m downwind from the fire front.

In order to visualize the three-dimensional flow features inside the domain, the Q-criterion [64] is employed. Figure 5 illustrates the vortical structures generated in the presence and absence of fire using Q-Criterion isocontours coloured by the spanwise velocity magnitude in the z-direction to better illustrate buoyancy effects. In the absence of fire, a horse-shoe vortex is formed upstream of the building. This is a common phenomenon as a bluff-body encounters a turbulent boundary layer. Other toroidal and hair-pin vortex cores shedding from the sides of the building are relatively small and are mostly aligned with the streamwise direction. However, as can be seen in Figure 5b, the flow pattern substantially changes in the presence of the fire. In the attachment region, the breakdown of the turbulent boundary layer into small streamwise vortices is highly controlled by the high-momentum flow near the surface. Similar patterns have been reported in experimental observations of laboratory-scale fires, for example [65]. Moving further downstream, as the buoyancy-induced forces become dominant, the flow structures are transformed into larger crosswise vortices creating uplift.

In order to better describe the intermittent behaviour of the unsteady flow in the vicinity of the building, the surface temperature was monitored at three different locations (P1, P2 and P3) in front of the building, as indicated in Figure 1. Figure 6a shows the time-dependent temperature variation obtained at the three monitoring points over a 10 s timeframe (from t = 10 s to 20 s). The temperature fluctuations are mainly caused by the stochastic movement of the flames in the attachment region, which directly affects the fire-surface contact frequency.

The temperature-time series were subject to fast Fourier transform (FFT) analysis, and the dominant frequencies of the temperature signals were calculated. For points P2 and P3, where the fire is spreading over a flat surface (no building), the peak frequencies are nearly the same and occur at f₂ = 0.42 Hz and f₃ = 0.46 Hz. However, the peak frequency at P1, located in front of the building, is approximately twice (f₁ = 0.93 Hz) the frequency at P1 and P2, which suggests that the presence of the building downwind of the fire front results in more frequent pulsations in the fire plume. This can be attributed to premature formation of buoyant instabilities, which manifest in the form of rising crosswise vortices (see Figure 5b) in front of the building. The transient plume behaviour with a higher frequency of pulsation results in more local heating (both radiative and convective heating) around the building, which increases the potential of new ignitions in and around the building; for example, structural extremities (e.g., decking) or sources of fuel that may be present due to lack of household maintenance.

The higher frequency pulsations also mean that the building will be exposed to more frequent convective heating, with less opportunity for convective cooling [65]. This has implications for the likelihood of ignition of finer fuel particles in the vicinity of the structure, including building extremities. It also has implications for the validity of the way flames are represented in risk mitigation frameworks, such as building standards. For example, the Australian Standard for building in bushfire prone areas (AS3959 [66]) bases its calculation of building exposure on a “design fire” concept, in which flames are represented essentially as a static feature rather than the pulsating flames evident in the simulations.

Figure 7 compares the incident radiative heat flux on the building at different times. As expected, the frontal façade of the building receives most of the radiant heat flux over the whole simulation time. At t = 5 s, when the convective plume has not yet fully developed, the heat flux gradually increases with the building height, which is mainly due to the fact that the view factors from the heated ground increase at larger heights. After t = 5 s, a consistent fire plume structure is established, which approximates the effect of a sustained fire front near the building. As a result of the fire’s intermittent pulsing upwind of the building, which occurs after t = 10 s, the pattern of incident radiation becomes highly irregular. The maximum local radiation intensity is observed near the points where the flame extensions and pockets of hot gases come in close contact with the building. On the contrary, the lateral walls and the roof are not exposed to direct fire impingement and therefore experience low exposure to radiation. It is worth mentioning that the local incident heat flux is very close to the critical heat flux for wood pyrolysis [67] (i.e., approximately 10 kW/m²). This suggests that the frontal façade of the building is exposed to a very high risk of thermal pyrolysis compared to the other sides, which are less affected by the fire’s radiant heat.

4. Conclusions

A numerical investigation of the dynamic behaviour of wind-driven surface fires has been presented, with an emphasis on building-fire interaction relevant to consideration of WUI fire impacts. The study was based on LES modelling of a wind-driven fire using the FireFOAM code. The numerical results in the absence of fire were shown to compare well with the data obtained from the full-scale Silsoe building experiments. A comprehensive analysis of the dynamics of the fire’s plume revealed that the building can substantially change the plume attachment length. In other words, the building causes premature buoyant instabilities and the fire plume lift-off distance decreases in front of the building. In addition, an FFT analysis was performed to examine the characteristics of plume pulsations near the building. It was found that the building increases the plume’s intermittency and local heat exchange at some locations upstream of the building. Thus, potentially risky zones susceptible to the generation of new fire sources are formed near the building. In addition, the capability of the present model to predict the thermal response of the building to the fire was studied. These findings provide insights into the likely impacts of wind-driven wildfires applicable to WUI fire scenarios.

Moreover, the results presented here raise some questions about the validity of risk management measures, such as AS3959 [66], in which radiant heat exposure is calculated based on a static “design fire”. In AS3959, flames are essentially considered as a hot surface tilted at a constant angle determined by the balance between wind strength and fire intensity, the latter of which is estimated using quasi-steady fire behaviour models. In contrast, the simulations presented in this study highlight the transient (dynamic) nature of fire propagation and indicate that a fire immediately upwind of a structure can produce flames that alternate between being attached to the surface and more vertically oriented. Flames that are attached to the surface are more likely to bathe a structure in convective heat, and the higher frequency pulsations (between attached/separated regimes) that occur upwind of a structure suggest that convective heating can have a significant influence on the heat exposure of building façades directly in the path of an approaching fire front.

The present study offers a useful framework for evaluating the fire exposure of WUI communities. While it would incur a high computational cost to fully implement the proposed framework at larger scales, this approach could offer a reliable solution to reduce risks from future WUI fires. These findings can be extended by conducting similar studies with fires exhibiting different characteristics (e.g., more or less intense), and realistic building geometries or multiple buildings. In addition, further research is needed to explore fire dynamic behaviour and its impact on smoke dispersion and ember generation in complex urban environments. Such studies are required to inform comprehensive wildfire risk management strategies, including improved urban planning, and WUI building construction design.