RandomFront 2.3: a physical parameterisation of fire spotting for operational fire spread models – implementation in WRF-SFIRE and response analysis with LSFire+

Abstract

Fire spotting is often responsible for dangerous flare-ups in wildfires and causes secondary ignitions isolated from the primary fire zone, which lead to perilous situations. The main aim of the present research is to provide a versatile probabilistic model for fire spotting that is suitable for implementation as a post-processing scheme at each time step in any of the existing operational large-scale wildfire propagation models, without calling for any major changes in the original framework. In particular, a complete physical parameterisation of fire spotting is presented and the corresponding updated model RandomFront 2.3 is implemented in a coupled fire–atmosphere model: WRF-SFIRE. A test case is simulated and discussed. Moreover, the results from different simulations with a simple model based on the level set method, namely LSFire+, highlight the response of the parameterisation to varying fire intensities, wind conditions and different firebrand radii. The contribution of the firebrands to increasing the fire perimeter varies according to different concurrent conditions, and the simulations show results in agreement with the physical processes. Among the many rigorous approaches available in the literature to model firebrand transport and distribution, the approach presented here proves to be simple yet versatile for application to operational large-scale fire spread models.

Introduction

Fire spotting is an important phenomenon associated with wildfires . It is documented as a dominant aspect that has contributed to the rampant spread of fire in many devastating historical fires . Spot fires occur when fragments of the fuel tear off from the main fuel source and horizontal wind transports the burning embers beyond the zone of direct ignition. The burning embers/firebrands can develop new secondary ignition spots and lead to a perilous increase in the effective rate of spread (ROS) of the fire.

The main aim of the present research is to provide a versatile probabilistic model for fire spotting that is suitable for implementation as a post-processing scheme at each time step in any of the existing large-scale operational codes for simulating wildfire propagation, without calling for any major changes in the original framework.

Researchers have tried to understand the phenomenology of fire spotting through both experimental and theoretical aspects to update the existing wildfire management decision support systems. Most of the experimental procedures for studying the fire-spotting phenomenon focus on the characterisation of the generation of firebrands , the shape and size of the firebrands , drag forces and ignition processes . The short temporal and spatial scales of the experiments limit a detailed description of the landing distributions and flight paths of the firebrands. Conversely, firebrand transport models provide an estimate of the maximum landing distance and flight paths of firebrands through a simplified overview of the physical dynamics of fire behaviour, plume characteristics and the atmospheric conditions around the fire. and were the first to develop simplified plume models for the estimation of firebrand lifetimes, flight paths and potential fire-spotting distance. Beginning with their studies, there has been a paradigm shift in the development of firebrand transport models, with the latest models benefiting from advanced computational techniques and resources.

provide a model for the lofting of spherical and cylindrical firebrands using the plume model proposed by . They suggest analytical functions for the maximum loftable diameter and the maximum loftable height in terms of the fire intensity, atmospheric wind and the fuel characteristics. Numerical experiments by Sardoy and co-workers also analyse the effect of atmospheric conditions, fire properties and fuel properties on firebrand behaviour and provide a statistical estimate of the ground level distributions of disk shaped firebrands. Their results highlight that firebrands landing at short distances (up to $\mathrm{1000}$ m) from the source follow a lognormal distribution. A study by also provides a mathematical model to quantify the distribution and the mass of firebrands through a Rayleigh distribution function. In an another study, present a physics based multiphase transport model for wildfires (FIRETEC) to study firebrand transport. In a recent study, also discuss the underlying physical processes for firebrands in detail and derive a landing distribution based on these physical processes. Besides these statistical approaches, a few numerical models based on large eddy simulation (LES) or computational fluid dynamics (CFD) , small world networks and cellular automata models also exist in the literature. present one such study based on a coupled fire–atmosphere LES for predicting short-range fire spotting. They simulate multiple firebrand trajectories for analysing the sensitivity of the flight path to different particle sizes, release heights and wind conditions but also mention the limited applicability of such models to operational use due to the computational demands.

Despite the presence of multiple studies focusing on the detailed aspects of firebrand landing distributions, none of them are able to provide a comprehensive yet versatile approach for an application to operational fire spread models. The continuing demand for operational management tools is to provide a quick and efficient output with simple inputs but at the same time to take the most important parameters into consideration. A few operational fire spread models like FARSITE , BEHAVEPLUS and Prometheus incorporate the phenomenon of fire spotting through the Albini's model . However, the Albini's model only provides an estimate of the maximum distance for a spot fire and does not include any function for the ignition probability to model the spread of spot fires. The Australian wildfire simulator PHOENIX RapidFire is designed to model large fast moving fires and also includes a fire-spotting module; however, the formulations for fire spread in PHOENIX are calibrated for eucalyptus forests and a generic application to other types of fuels requires a recalibration . The new operational models like WRF-SFIRE and FOREFIRE are fast and allow coupling with atmospheric models for a better representation of the initial and concurrent atmospheric conditions; but they lack any specific module to tackle fire-spotting behaviour.

In this article, the authors proceed with the RandomFront statistical formulation for including the effects of random processes in wildfire simulators, namely turbulence and fire-spotting phenomena. The chronology of this approach refers to the following papers: v1.0 includes only turbulence, with no parameterisation ; v2.0 includes turbulence and fire spotting with literature parameterisation for fire spotting ; v2.1 includes turbulence and fire spotting with parameterisation for turbulence ; and v2.2 includes turbulence and fire spotting with a first physical parameterisation of fire spotting . Finally, in the present version v2.3 the parameterisation of fire spotting has been modified and corrected (also in view of a remark by one of the referees) with respect to the previous version.

The physical parameterisation of the probabilistic model is developed to incorporate the fire-spotting behaviour in terms of the fire intensity, wind conditions and fuel characteristics. This formulation is independent of the method used for the fire-line propagation and the definition of the ROS, and it is versatile enough to be utilised with any of the existing operational fire spread models. In their previous work , the authors demonstrate the applicability of the formulation for two wildfire models based on different fire-line propagation methods, i.e. a Eulerian moving interface method based on the level set method (LSM) that is the basis of the WRF-SFIRE model and a Lagrangian front tracking technique based on the discrete event system specification (DEVS) that is the basis of the FOREFIRE model. The aim of the present study is to provide a simple yet complete addition to operational fire spread models for representing the random behaviour of fire spotting through simple inputs related to wildfires. This probabilistic model is devised to provide a physical meaning to the spread of fire by virtue of firebrands. The proposed parameterisation has been implemented into WRF-SFIRE and a paradigmatic test case has been simulated. Moreover, the proposed formulation for including random process and the corresponding physical parameterisations are implemented into a much simpler fire spread simulator, also based on the LSM, namely LSFire+ . Results from different test cases with LSFire+ are presented to highlight the sensitivity of the simple parameterisation in simulating the generation of secondary fires by fire spotting under different wind conditions, fire intensities and firebrand radii.

The article is organised as follows. In Sect. 2 a very brief description of the mathematical model is presented, while the physical parameterisation of fire spotting within the framework of a lognormal distribution of the landing distance is described in Sect. 3. The implementation of RandomFront 2.3 in the computational environments WRF-SFIRE in addition to a test case and the corresponding discussion are reported in Sect. 4, and the implementation in LSFire+ and the corresponding response analysis are discussed in Sect. 5. Conclusions are drawn in Sect. 6.

Model formulation

A mathematical model to represent the random effects associated with the wildland fires has been developed by Pagnini and co-authors . This formulation describes the motion of the fire line as a composition of the drifting part and the fluctuating part. The drifting part represents the fire perimeter obtained through the definition of the ROS based on fuel characteristics and the averaged fire properties. The output from most of the existing operational fire spread models can be considered as the drifting part. Conversely, the fluctuating part is independent of the drifting part and represents the additional contribution to the fire perimeter as an effect of random processes like turbulence and fire spotting. This model can be implemented as a crucial addition to operational fire spread models through a post-processing application at each time step. The drifting component obtained from the output of any wildfire model can be updated with the fluctuating component at each time step to include the effects of turbulence and fire spotting. A brief overview of the mathematical details is provided in the following; for a detailed description readers are referred to Pagnini and Mentrelli (2014) and Kaur et al. (2016).

In a domain $S$, let $\mathrm{\Omega }\subseteq S$ represent the burnt area and let $X^{\mathit{\omega }}=X+{\mathit{\eta }}^{\mathit{\omega }}$ represent the trajectory of each active fire point as the sum of a drifting part $X$ and a fluctuating part ${\mathit{\eta }}^{\mathit{\omega }}$. The drifting part $X$ is obtained from the output of a wildfire propagation model, while the fluctuations in the fire-line are included through a probability density function (PDF) corresponding to the type of random process under consideration. Let the area enclosed by the drifting part be described through an indicator function $I_{\mathrm{\Omega }}(\mathit{x},t)=\mathrm{1}$ when $\mathit{x}$ is inside the domain $\mathrm{\Omega }$, and $I_{\mathrm{\Omega }}(\mathit{x},t)=\mathrm{0}$ when $\mathit{x}$ is outside. Considering the ensemble average of the active burning points, a new effective indicator function is defined as follows: $${\mathit{\varphi }}_{\text{e}}(\mathit{x},t)=\underset{S}{\int }{I}_{\mathrm{\Omega }}(\stackrel{\mathrm{‾}}{\mathit{x}},t)f(\mathit{x};t|\stackrel{\mathrm{‾}}{\mathit{x}})\text{d}\stackrel{\mathrm{‾}}{\mathit{x}},$$ where $f(\mathit{x};t|\stackrel{\mathrm{‾}}{\mathit{x}})$ represents the PDF which accounts for the fluctuations of the random effects. The effective indicator ${\mathit{\varphi }}_{\text{e}}\in [\mathrm{0},\mathrm{1}]$, and an arbitrary threshold is fixed to mark points as burned, i.e. ${\mathrm{\Omega }}_e(\mathit{x},t)=\mathit{\{}\mathit{x}\in S\mid {\mathit{\varphi }}_{\text{e}}(\mathit{x},t)>{\mathit{\varphi }}_{\text{e}}^{\text{th}}\mathit{\}}$. The ignition of fuel by firebrands involves heat exchange over a sufficient period of time; therefore, a sufficient delay is also incorporated in the model through another function: $\mathit{\psi }$. The function $\mathit{\psi }$ simulates the ignition of fuel by hot air and burning embers as an accumulative process over time (the heating-before-burning mechanism): $$\mathit{\psi }(\mathit{x},t)=\underset{\mathrm{0}}{\overset{t}{\int }}{\mathit{\varphi }}_{\text{e}}(\mathit{x},\mathit{\eta }){\displaystyle \frac{\text{d}\mathit{\eta }}{\mathit{\tau }}},$$ where $\mathit{\tau }$ is the ignition delay. At each time $t$, all points $\mathit{x}$ that satisfy the conditions $\mathit{\psi }(\mathit{x},t)>\mathrm{1}$ or ${\mathit{\varphi }}_{\text{e}}(\mathit{x},t)>{\mathit{\varphi }}_{\text{e}}^{\text{th}}$ are labelled as burned.

The shape of the PDF $f(\mathit{x};t|\stackrel{\mathrm{‾}}{\mathit{x}})$ is established by analysing the random processes under consideration. The diversity in the shapes of the PDF provides the model a multifaceted outlook. Assuming fire spotting to be a downwind phenomenon occurring in turbulent atmosphere, the shape of the PDF is defined as follows: $$f(\mathit{x};t|\stackrel{\mathrm{‾}}{\mathit{x}})=\left\{\begin{array}{ll}{\displaystyle }{\int }_{\mathrm{0}}^{\mathrm{\infty }}G(\mathit{x}-\stackrel{\mathrm{‾}}{\mathit{x}}-l{\widehat{\mathit{n}}}_U;t)q\left(l\right)\text{d}l,& \text{downwind},\\ G(\mathit{x}-\stackrel{\mathrm{‾}}{\mathit{x}};t),& \text{upwind}.\end{array}\right.$$ The distribution function $G(\mathit{x}-\stackrel{\mathrm{‾}}{\mathit{x}};t)$ is an isotropic bivariate Gaussian and provides for the effect of the turbulent heat fluxes in fire propagation, while the distribution function $q\left(l\right)$ represents the firebrand landing distribution. The strength of the turbulence around the fire is parameterised through a turbulent diffusion coefficient $\mathfrak{D}$. A short description of the physical characterisation of $\mathfrak{D}$ is presented in the next section. A precise description of the landing distributions through experimental observations is difficult due to temporal and spatial constraints. But the experimental results analysing the flight paths, shape and landing distributions of firebrands have shown that the frequency of firebrands landing in the positive direction from the source increases with distance to a maximum value and then gradually decays to zero . The landing distributions of firebrands have also been studied though the numerical solution of the energy balance equations . Among the different transport models proposed in the literature, both and describe the lognormal density function as an approximate fit to the landing distribution of firebrands. , in comparison, proposes a Rayleigh distribution for the same. In this article, the shape of $q\left(l\right)$ is defined by a lognormal distribution to describe the frequency profile of the fallen firebrands: $$q\left(l\right)={\displaystyle \frac{\mathrm{1}}{\sqrt{\mathrm{2}\mathit{\pi }}\mathit{\sigma }l}}\exp {\displaystyle \frac{-(\ln l/\mathit{\mu }{)}^{\mathrm{2}}}{\mathrm{2}{\mathit{\sigma }}^{\mathrm{2}}}},$$ where $\mathit{\mu }$ is the ratio between the square of the mean of the landing distance $l$ and its standard deviation, while $\mathit{\sigma }$ is the standard deviation of $\ln l/\mathit{\mu }$.

Physical parameterisation of fire spotting

The firebrands generated from vegetation face strong buoyant forces and those with a size less than the maximum loftable size are vertically uplifted in the convective column. These firebrands rise to a maximum height until the buoyant and the gravitational forces counterbalance one another. Once the firebrands are expelled from the column, they are steered by atmospheric wind and fly at their terminal velocity of fall. Simplified models for the landing distance assume that the ejection of firebrands from the vertical convective column is a random process affected by the turbulence in the environment around the fire. Among other factors, the strength of the convective column, the atmospheric conditions and the dimensions of firebrands play a vital role in governing the trajectory of the firebrands. In this section, the landing distribution of firebrands based on a lognormal probability function is combined with the physical characterisation of firebrand transport. The parameterisation presented here is simplified and only includes the vital ingredients necessary to describe firebrand transport. Each firebrand is assumed to be spherical and for a particular set of concurrent atmospheric conditions and fuel characteristics, and the size is assumed to be constant. Any modification in the flight of the firebrand due to rotation of the firebrand or collision with other firebrands is also neglected. Preliminary results were discussed in .

Literature studies identify the maximum spotting distance as a numerical measure to assess the severity of the fire-spotting danger under different circumstances . Recognising the importance of maximum-spotting distance, we select to parameterise the mathematical model in terms of the $p$th percentile of a lognormal distribution as a measure of the maximum landing distance. The $p$th percentile for a lognormal distribution is described by its location and shape parameters $\mathit{\mu }$ and $\mathit{\sigma }$, respectively: $$\mathcal{L}=\mathit{\mu }\exp \left(z_p\mathit{\sigma }\right),$$ where the value of $z_p$ corresponding to the $p$th percentile can be estimated from the $Z$-table (see, for example, http://www.itl.nist.gov/div898/handbook/eda/section3/eda3671.htm, last access: 19 December 2018). We hypothesise that the $p$th percentile represents the maximum landing distance for firebrands under different situations and that no ignition is possible beyond this cut-off. To ascertain the value of the “cut-off” percentile, it is assumed that the effective contribution of the firebrands ceases to be meaningful when its probability reduces to 20 times its peak value. Thereafter, the ability of the firebrands to cause an ignition is assumed to be negligible. The cut-off criteria is chosen empirically, but a sufficiently large number (like 20) ensures that we do not miss any considerable fire-spotting behaviour that exists outside this range.

For this particular distribution, the cut-off for the 50th percentile lies way beyond the point denoting $\mathrm{1}/\mathrm{20}$th of the maximum probability. In order to define a generalized value of the cut-off percentile for all of the simulation cases presented in this paper, the value of $z_p$ is set to 0.45, which corresponds to the $\mathrm{67}$th percentile point.

The process of fire spotting can be roughly segregated into generation, lofting and transport of firebrands. The generation of firebrands from a burning canopy is a random and dynamic process, while the lofting and transport of firebrands is regulated by the firebrand geometry, fuel combustion rates, plume dynamics and ambient wind conditions. The different sub-processes involved in firebrand lofting and transport interact with each other and affect the maximum spotting distances. An explicit modelling of the coupled processes is difficult and different approximations and assumptions are often used to simplify the physical processes. An example of important works on fire-spotting distributions and maximum spotting distances are those by Tarifa and co-workers . In their different publications they describe spotting distributions and maximum spotting distances by combining a series of experimental and theoretical approaches. We follow these existing approaches and formulate the physical parameterization for our model. Below we provide a brief discussion of the different processes which are considered in the physical parameterization:

• Firebrand lofting

  • Vertical gas flow: in a convective column, the updraught introduced by the fire lifts the firebrands in the convective column. The strength of vertical gas flow $U_{\text{gas}}$ increases with the fire intensity $I$ and is empirically expressed as in Muraszew (1974):$$U_{\text{gas}}=\mathrm{9.35}{\left({\displaystyle \frac{I}{H_{\text{c}}}}\right)}^{\mathrm{1}/\mathrm{3}},$$where $H_{\text{c}}$ is the heat of combustion of wildland fuels.

  • Size of firebrands: the convective activity inside the plume regulates the maximum size of the firebrand that can be lofted. The terminal velocities of the loftable firebrands can not exceed the vertical gas flow rate. As the vertical gas-flow increases with increasing fire rate, heavier firebrands can be uplifted into the plume. In the literature, the maximum loftable radius for spherical firebrands is expressed as follows :$$r_{\text{max}}={\displaystyle \frac{\mathrm{3}}{\mathrm{2}}}{C}_{\text{d}}{\displaystyle \frac{{\mathit{\rho }}_{\text{a}}}{{\mathit{\rho }}_{\text{f}}}}{\displaystyle \frac{U_{\text{gas}}^{\mathrm{2}}}{g}},$$where ${\mathit{\rho }}_{\text{a}}$ and ${\mathit{\rho }}_{\text{f}}$ represent the density of the ambient air and wild-land fuels, respectively, $C_{\text{d}}$ is the drag coefficient and $g$ is the acceleration due to gravity.

  • Maximum loftable height: and parameterise the maximum loftable height for spherical firebrands in terms of the radius of firebrand $r$, constrained by the maximum radius of the firebrands $r_{\text{max}}$:$$\mathcal{H}=\mathrm{1.46}\left({\displaystyle \frac{{\mathit{\rho }}_{\text{f}}}{{\mathit{\rho }}_{\text{a}}{C}_{\text{d}}}}\right){\displaystyle \frac{r_{\text{max}}^{\mathrm{5}/\mathrm{2}}}{r^{\mathrm{3}/\mathrm{2}}}}.$$

• Horizontal transport:

  • Maximum landing distance: assuming the shape of the firebrands to be spherical, the study by combines both experimental and theoretical approaches to characterise the maximum landing distances of firebrands. Based on these results, provides an approximation of the maximum travel distance for spherical firebrands from a vertical convective column in terms of the following: maximum loftable height $\mathcal{H}$; the meteorological mean wind $\mathcal{U}=\left|{\mathit{U}}_{\text{h}}\right|$, where ${\mathit{U}}_{\text{h}}$ is the horizontal wind vector field at fire height; and the radius of the firebrands $r$. This results in the following equation:$$\mathcal{L}=\mathcal{H}{\displaystyle \frac{\mathcal{U}}{U_{\text{gas}}}}{\left({\displaystyle \frac{r_{\text{max}}}{r}}\right)}^{\mathrm{1}/\mathrm{2}}.$$

• Ignition probability: as described in the previous section, the probability that fuel is ignited by burning embers is modelled using the function $\mathit{\psi }$ (Eq. ). Here we assume that the fuel conditions are homogeneous and the ignition probability only depends on an ignition delay $\mathit{\tau }$. No other local variables are taken into account.

• Secondary fire-lines: the secondary emissions generated during fire-spotting modelling are assumed as new sources of fire with a proper fire intensity. These new fires act as additional input along with the primary fire towards the generation of other secondary fires; it is assumed that these new sources are capable of generating firebrands of the same size as the primary source.

Small-scale processes, such as the mass loss of a firebrand due to combustion, affect the fire-spotting phenomenon by generating random fluctuations in the firebrand trajectory. These fluctuations are embodied by the use of a distribution for the landing distance.

Finally, the above large-scale sub-processes under lofting and transport mechanisms are linked through Eq. () to obtain the physical parameterisation of $\mathit{\mu }$ and $\mathit{\sigma }$ of the lognormal distribution. Using Eqs. () and (), we express the shape and location parameters as follows: $$\mathit{\sigma }={\displaystyle \frac{\mathrm{1}}{\mathrm{2}{z}_p}}\ln \left({\displaystyle \frac{{\mathcal{U}}^{\mathrm{2}}}{rg}}\right),$$ $$\mathit{\mu }=\mathcal{H}{\left({\displaystyle \frac{\mathrm{3}}{\mathrm{2}}}{\displaystyle \frac{{\mathit{\rho }}_{\text{a}}}{{\mathit{\rho }}_{\text{f}}}}{C}_{\text{d}}\right)}^{\mathrm{1}/\mathrm{2}}.$$

List of the symbols used to describe the model quantities and the parameters used in this study.

Model quantities		Units
${\mathit{\varphi }}_{\text{e}}$	Effective indicator	–
$\mathit{\psi }$	Ignition function	–
$\mathit{x}=(x,y)$	Horizontal space variable	m
$t$	Time	s
$f$	Probability density function of the random processes	m${}^{-\mathrm{2}}$
$G(\mathit{x};t)$)	Isotropic bivariate Gaussian probability density	m${}^{-\mathrm{2}}$
$q\left(l\right)$	Lognormal distribution of firebrand landing	m${}^{-\mathrm{1}}$
Parameters
Static parameters		Value
$\mathit{\mu }$	Parameter of $q\left(l\right)$	–
$\mathit{\sigma }$	Parameter of $q\left(l\right)$	–
$\mathfrak{D}$	Turbulent diffusion coefficient of $G$	0.15 m${}^{\mathrm{2}}$ s${}^{-\mathrm{1}}$
${\mathit{\rho }}_{\text{a}}$	Density of air	1.2 kg m${}^{-\mathrm{3}}$
${\mathit{\rho }}_{\text{f}}$	Density of wildland fuel (Pinus ponderosa)	542 kg m${}^{-\mathrm{3}}$
$C_{\text{d}}$	Drag coefficient	0.45
$z_p$	$p$th percentile	0.45
$g$	Acceleration due to gravity	9.8 m s${}^{-\mathrm{1}}$
$H_{\text{c}}$	Heat of combustion of wildland fuels	18 620 kJ kg${}^{-\mathrm{1}}$
${\mathit{\omega }}_{\mathrm{0}}$	Oven-dry mass of fuel	2.243 kg m${}^{-\mathrm{2}}$
$H$	Fuel low heat of combustion	22 000 kJ kg${}^{-\mathrm{1}}$
Dynamic parameters		Units
$\mathit{U}=(U,V,W)$	Wind vector at fire height	m s${}^{-\mathrm{1}}$
${\mathit{U}}_{\text{h}}=(U,V)$	Horizontal wind vector field at fire height	m s${}^{-\mathrm{1}}$
$\mathit{\tau }$	Ignition delay of firebrands	s
$I$	Fire-line intensity	MW m${}^{-\mathrm{1}}$
$U_{\text{gas}}$	Vertical gas flow	m s${}^{-\mathrm{1}}$
$r$	Radius of spherical firebrand	m
$r_{\text{max}}$	Maximum loftable radius for spherical firebrand	m
$\mathcal{H}$	Maximum loftable height for spherical firebrands	m
${\mathcal{V}}_{\text{ros}}$	Rate of spread	m s${}^{-\mathrm{1}}$

We chose such a parameterization of $\mathit{\mu }$ and $\mathit{\sigma }$ in order to delineate the governing parameters for lofting and transport mechanisms, respectively. We hypothesise that the definition of $\mathit{\mu }$ covers the essential input parameters needed to describe the lofting mechanism of firebrands inside the convective column. The relative density ${\mathit{\rho }}_{\text{a}}/{\mathit{\rho }}_{\text{f}}$ and atmospheric drag quantify the buoyant forces experienced by the firebrand; hence, it is appropriate to include these quantities in the definition of $\mathit{\mu }$. Substituting maximum loftable height from Eq. () in $\mathit{\mu }$ gives $$\mathit{\mu }=\mathrm{3.52}\times {\mathrm{10}}^{\mathrm{5}}{\left({\displaystyle \frac{{\mathit{\rho }}_{\text{a}}}{{\mathit{\rho }}_{\text{f}}}}{C}_{\text{d}}\right)}^{\mathrm{2}}{\left({\displaystyle \frac{I}{H_{\text{c}}}}\right)}^{\mathrm{5}/\mathrm{3}}{r}^{-\mathrm{3}/\mathrm{2}}{g}^{-\mathrm{5}/\mathrm{2}}.$$ The radius of the firebrand $r$ and the fuel density are important factors regarding the determination of the height of the lofted firebrands. Conversely, $\mathit{\sigma }$ is hypothesised to define the transport of firebrands under the effect of horizontal wind after ejection from the convective column. The definition of $\mathit{\sigma }$ includes a dimensionless ratio $\mathcal{F}={\mathcal{U}}^{\mathrm{2}}/\left(rg\right)$ which is analogous to the Froude number, which quantifies the balance between inertial and gravitational forces experienced by firebrands. All firebrands with $r\le {\mathcal{U}}^{\mathrm{2}}/g$ can be transported by horizontal wind.

In this model, the phenomenon of fire spotting is assumed to occur together with the turbulent heat flux around the fire, and the turbulent diffusion coefficient $\mathfrak{D}$ is utilised as a measure of the turbulent heat transfer generated by the fire. It is parameterised in terms of the Nusselt number Nu. The Nusselt number defines the ratio between the convective and conductive heat transfer in fluids and is defined as $\mathit{Nu}=(\mathfrak{D}+\mathit{\chi })/\mathit{\chi }$, where $\mathit{\chi }=\mathrm{2}\times {\mathrm{10}}^{-\mathrm{5}}$ m${}^{\mathrm{2}}$ s${}^{-\mathrm{1}}$ is the thermal diffusivity of air at ambient temperature. Experimentally, it is shown that the Nusselt number is related to the Rayleigh number as $\mathit{Nu}\simeq \mathrm{0.1}{\mathit{Ra}}^{\mathrm{1}/\mathrm{3}}$ . The Rayleigh number is defined as $\mathit{Ra}=\mathit{\gamma }\mathrm{\Delta }Tg{h}^{\mathrm{3}}/\left(\mathit{\nu }\mathit{\chi }\right)$, where $\mathit{\gamma }=\mathrm{3.4}\times {\mathrm{10}}^{-\mathrm{3}}$ K${}^{-\mathrm{1}}$ is the thermal expansion coefficient, $h$ is the dimension of the convective cell, $\mathit{\nu }=\mathrm{1.5}\times {\mathrm{10}}^{-\mathrm{5}}$ m${}^{\mathrm{2}}$ s${}^{-\mathrm{1}}$ is the kinematic viscosity and $\mathrm{\Delta }T$ is the temperature gradient between the top and bottom faces of the convective cell. Finally, we estimate the turbulent diffusion coefficient using the following formula: $$\mathfrak{D}\simeq \mathrm{0.1}{\left({\displaystyle \frac{\mathit{\gamma }g}{\mathit{\nu }}}\right)}^{\mathrm{1}/\mathrm{3}}{\mathit{\chi }}^{\mathrm{2}/\mathrm{3}}\mathrm{\Delta }{T}^{\mathrm{1}/\mathrm{3}}h-\mathit{\chi },$$ and assuming $\mathrm{\Delta }T\simeq \mathrm{1000}$ K and $h\simeq \mathrm{100}$ m, we have $\mathfrak{D}\sim {\mathrm{10}}^{-\mathrm{1}}$ m${}^{\mathrm{2}}$ s${}^{-\mathrm{1}}$. For all of the simulations presented in this article, the value of the turbulent diffusion coefficient $\mathfrak{D}$ is assumed to be $\mathrm{0.15}$ m${}^{\mathrm{2}}$ s${}^{-\mathrm{1}}$.

Values of the main parameters for numerical simulations performed with LSFire+.

	Unit of	First	Second
Quantity	measurement	test case	test case
$\mathfrak{D}$	m s${}^{-\mathrm{1}}$	$\mathrm{0.15}$	$\mathrm{0.15}$
$U$	m s${}^{-\mathrm{1}}$	$\mathrm{10}$	$\mathrm{2}/\mathrm{26}$
$I$	MW m${}^{-\mathrm{1}}$	$\mathrm{5}/\mathrm{100}$	$\mathrm{50}$
$r$	m	$\mathrm{0.015}$	$\mathrm{0.015}/\mathrm{0.03}$
$\mathit{\tau }$	s	$\mathrm{1}$	$\mathrm{1}$

The simple design of the physical parameterisation makes the model computationally less expensive, and the requirement of defining only few vital parameters to execute any simulation also serves as an added advantage to the operational users. Static and dynamic input parameters of the model are reported in Table .

Numerical simulations

The detailed steps of the numerical procedure are as follows :

• Beginning with an initial fire line, an operational code, i.e. WRF-SFIRE in this section and LSFire+ in the next section, is used to estimate the propagation of the front and to generate a new fire perimeter for the next time step. This output is modified to include the effects of turbulence and fire spotting using a “post-processing numerical procedure”. This post-processing step is independent of the definition of the ROS.

• The fire perimeter obtained from the chosen operational code is used to construct the indicator function $I_{\mathrm{\Omega }}(\mathit{x},t)$, i.e. the indicator function $I_{\mathrm{\Omega }}(\mathit{x},t)$ has a value of $\mathrm{1}$ inside the domain surrounded by the fire-line and $\mathrm{0}$ outside of the domain. The spatial information contained in $I_{\mathrm{\Omega }}(\mathit{x},t)$ is necessary to modify the fire line with respect to turbulence and fire spotting and serves as an input to the post-processing step.

• The effective indicator function ${\mathit{\varphi }}_{\text{e}}(\mathit{x},t)$ (Eq. ) is generated over a Cartesian grid to facilitate the computation of the function $\mathit{\psi }(\mathit{x},t)$ (Eq. ) over the same grid.

• The value of the effective indicator ${\mathit{\varphi }}_{\text{e}}(\mathit{x},t)$ is computed through the numerical integration of the product of the indicator function $I_{\mathrm{\Omega }}(\mathit{x},t)$ and the PDF of fluctuations according to Eq. (). The effect of turbulence or fire spotting is included by choosing the corresponding PDF (Eq. ).

• The function $\mathit{\psi }(\mathit{x},t)$ is updated for each grid point by integration in time with the current value of ${\mathit{\varphi }}_{\text{e}}(\mathit{x},t)$.

• All points that satisfy the condition $\mathit{\psi }(\mathit{x},t)\ge \mathrm{1}$ are labelled as new ignition spots. The post-processing procedure is completed at this step.

• At the next time step, the new fire perimeter evolves according to the chosen operational code, and the updated perimeter is again subjected to the post-processing procedure to enrich the fire front with the random fluctuations pertaining to turbulence and fire spotting. The sequence is repeated until the final “event time” step or until the fire reaches the boundaries of the domain.

Implementation of RandomFront 2.3 in WRF-SFIRE

In order to prove the viability of the proposed formulation within an operational code, we implemented RandomFront 2.3 in the framework of the WRF-SFIRE simulator . WRF-SFIRE is a coupled fire–atmosphere model, which operates in the computational environment of a well-known public domain numerical weather prediction model: WRF v3.4 (Weather Research and Forecasting; ; https://www.mmm.ucar.edu/weather-research-and-forecasting-model , last access: 19 December 2018).

The fire module embedded into WRF simulates a surface fire and takes a two-way coupling with the atmospheric model into account. The near-surface winds from the atmospheric model are interpolated on the fire grid and are used along with fuel properties and local terrain gradients to compute both the ROS and the outward front-direction, which are further used as an input to the front propagation routines through a LSM scheme. Fuel consumption is responsible for the release of sensible and latent heat into the lowest layers of the atmosphere, and this has a role in the computation of the boundary-layer meteorology. Recently, the model has been equipped with a fuel-moisture sub-model and a chemistry sub-model (WRF Chem), which contribute to reproducing and investigating the effects of the fire–atmosphere coupling.

points out that fires generally start from a horizontal extent much smaller than the size of the fire mesh-cell. The same may be argued for the secondary ignitions related to fire-spotting phenomenon. In this respect, propose (and explain in detail) an algorithm for a punctual or line ignition that actually runs on WRF-SFIRE. The purpose of this algorithm is two-fold: (i) it guarantees from a physical point of view that the ignition starts at sub-grid scale without generating unrealistically large initial heat flux and an accelerated ignition; (ii) this procedure is numerically robust because it is fully integrated into the representation in terms of a “signed distance function” of the LSM .

In the proposed formulation, a punctual ignition occurs whenever the condition $\mathit{\psi }(\mathit{x},t)\ge \mathrm{1}$ holds true. This procedure is not computationally viable so we set a threshold distance $R_{\text{th}}=\mathrm{200}$ m for separating each pair of punctual ignitions. In particular, let $\mathcal{P}$ be the set of point-wise fire-spotting ignitions, the actual algorithm performed at each time step within WRF-SFIRE model is reported in Algorithm .

algorithm for point-wise ignition due to fire spotting

[Figure omitted. See PDF]

Discussion of the test case with WRF-SFIRE

In this paper we consider a slight modification of the $\mathtt{hill}$ test case (https://github.com/openwfm/wrf-fire/blob/master/wrfv2_fire/test/em_fire/hill/namelist.input.hill, last access: 19 December 2018). In order to simplify the underlying dynamics, but maintain the fire–atmosphere coupling, the hill is suppressed (fire_mountain_type = 0) and we have a square-grid simulation over a flat domain with a side length of $\mathrm{2.5}$ km. The horizontal atmospheric grid-spacing at terrain-height is $\mathrm{60}$ m, while the fire spread grid-spacing is $\mathrm{15}$ m. The simulation starts at $t=\mathrm{0}$ min and ends at $t=\mathrm{20}$ min. The fire line is initially located along the segment joining the points (1900 m, 1500 m) and (1900 m, 1800 m), and the initial wind field at the fire height is ($U(\mathit{x}$; $t=\mathrm{0}$), $V(\mathit{x};t=\mathrm{0})$) $=$ ($-\mathrm{6.4}$ m s${}^{-\mathrm{1}}$, $-\mathrm{3.6}$ m s${}^{-\mathrm{1}}$) at all points in the simulation domain.

The fuel has been set equal to fire “Type 9”, i.e. “FM 9 hardwood litter” according to Anderson classification . This fuel type may represent a terrain covered by P. ponderosa trees. The radius of spherical embers was set equal to $r=\mathrm{12.5}$ mm following a size considered by with respect to the same vegetation.

Concerning the proposed formulation, the fire-line intensity $I$ and the wind field are computed using the WRF-SFIRE model; this allows for a varying field of both parameters $\mathit{\sigma }$ and $\mathit{\mu }$ according to Eqs. () and (), respectively. In particular, following the latest advancements of the SFIRE model environment in , the spatial representation of the potential-fire characteristics are available from which a field extension of the fire-line intensity $I$ is available in order to have a space and time varying field of $\mathit{\mu }$ for the fire-spotting routines. Parameters $\mathfrak{D}$ and $\mathit{\tau }$ are set as $\mathfrak{D}=\mathrm{0.15}$ m${}^{\mathrm{2}}$ s${}^{-\mathrm{1}}$ and $\mathit{\tau }=\mathrm{8}$ s, without using estimations from WRF-SFIRE.

Wind vector components ($U$, $V$ $W$) performed with WRF-SFIRE at times $t=\mathrm{6}$, 10 and 20 min. The fire front is represented by a dashed line.

[Figure omitted. See PDF]

Fire intensity $I$ and PDF shape parameter $\mathit{\mu }$ performed with WRF-SFIRE at times $t=\mathrm{6}$, 10 and 20 min. The fire front is represented by a dashed line.

[Figure omitted. See PDF]

The PDF shape parameter $\mathit{\sigma }$ and horizontal wind squared magnitude ($|{\mathit{U}}_{\text{h}}{|}^{\mathrm{2}}$) performed with WRF-SFIRE at times $t=\mathrm{6}$, 10 and 20 min. The fire front is represented by a dashed line.

[Figure omitted. See PDF]

Figures – display the simulation results. In each figure the fire front is represented by a dashed line at the following instants $t=\mathrm{6}$ min, $t=\mathrm{10}$ min and $t=\mathrm{20}$ min. The selected instants allow for the observation of the propagation of the main fire alone, the generation of a secondary fire and the multi-generation of secondary fires. In particular, in Fig.  the evolution of the fire line is shown in relation to the three components of the wind field; Fig.  shows the relationship with the parameter $\mathit{\mu }$ and the fire intensity field; Fig.  shows the relationship with the parameter $\mathit{\sigma }$ and the squared norm of the horizontal wind.

Overall, we observe that the fire-line propagation is “pulled” in the direction of the maximum value of the squared norm of the horizontal wind (see right column in Fig. ); this direction is induced by the fire itself as a feedback on the weather as is shown by the patterns of the atmospheric observables when secondary fires are generated. The geometrical profile of the fire perimeter always plays an important role in determining fire-spotting behaviour. The asymmetry in the fire perimeter at 20 min along the prominent direction of propagation, causes the first secondary fire to appear in the top-left part of the domain. With time, as the fire-activity increases, the differences between the maximum value of the squared norm of the horizontal wind and its surroundings increase and the fire-line becomes symmetrical with respect to the main direction of propagation. This has a direct influence on the fire-spotting action, and the new secondary fires appear increasingly aligned towards the main direction of propagation.

The secondary fires are equally as important as the primary fire with respect to influencing the weather around the fire. The plots clearly show the influence of fire–atmosphere coupling, and feedback dynamics from secondary fires to primary fire can be also observed. The secondary fires affect the wind (see Fig. ) and the parameter $\mathit{\sigma }$ (Fig. ), which implies a refinement in fire-spotting characteristics for further ignitions. A point worth noting is that the first secondary fire occurs at a distance of almost $\mathrm{1500}$ m from the main fire. This observation supports the viability of the proposed formulation to simulate the fire-spotting mechanism in studies of large-scale fires. However, if the implementation in WRF-SFIRE allows for a comprehensive picture including the physical features of a multi-scale and multi-physics process, the complexity of the model, the number of parameters and the numerical cost increases.

In the next section, we perform a response analysis of our parameterisation, by using the simple finite difference code LSFire+, which allows for extensive simulations in terms of the spatial domain and simulation time. A sensitivity analysis on the inputs and an uncertainty quantification on the outputs of RandomFront implemented in LSFire+ will be considered in a separate paper.

Response analysis

Implementation of RandomFront 2.3 in LSFire+

A few idealised simulations are carried out to highlight the potential applicability of the formulation. For all the simulations, a flat domain with a homogeneous coverage of a P. ponderosa ecosystem is selected, as seen in the WRF-SFIRE implementation. The simulations are run using a basic set-up of the wildfire model LSFire+ which involves a moving interface method based on the LSM . The Byram formula is used to estimate the ROS of the fire line: $${\mathcal{V}}_{\text{ros}}(\mathit{x},t)={\displaystyle \frac{I(\mathrm{1}+f_{\text{w}})}{H\mathit{\alpha }{\mathit{\omega }}_{\mathrm{0}}}},$$ where $I$ is the fire intensity, $H=$ 22 000 kJ kg${}^{-\mathrm{1}}$ is the fuel low heat of combustion, ${\mathit{\omega }}_{\mathrm{0}}=\mathrm{2.243}$ kg m${}^{-\mathrm{2}}$ is the oven-dry mass of the fuel and the functional dependence on the wind is included through the factor $f_{\text{w}}$. The user has the flexibility to introduce a different ecosystem into the simulations by modifying the parameters $H$ and ${\mathit{\omega }}_{\mathrm{0}}$. The parameter $\mathit{\alpha }$ is chosen to guarantee that the maximum ROS is always equal to the ROS prescribed by the Byram formulation.

The response of the formulation with respect to depicting the different firebrand landing distributions is highlighted through two sets of test cases. In the first test case, the wind conditions and the size of the firebrands are assumed to be constant as the fire intensity changes. In the second test case, the fire intensity is assumed to be constant and the simulations for different wind conditions are carried out. The second test case is also repeated for different firebrand radii.

(a) Line contours showing the fire perimeter at different time steps with the presence of a fire barrier. The wind velocity is $\mathrm{10}$ m s ${}^{-\mathrm{1}}$, the fire intensity is $\mathrm{25}$ MW m${}^{-\mathrm{1}}$ and the diffusion coefficient is $\mathrm{0.15}$ m${}^{\mathrm{2}}$ s ${}^{-\mathrm{1}}$. The $x$ and $y$ axes of the plot are scaled by a factor of 4. The same plot is proposed in (b), but with $U=\mathrm{20}$ m s ${}^{-\mathrm{1}}$ and no barrier. (c) A comparison of the total burned area at different time steps when only turbulence is considered (black) and when both turbulence and fire spotting are included (red). The total burned area is simply the number of burned grid points at any each instant. For both line plots, the wind velocity is $\mathrm{10}$ m s ${}^{-\mathrm{1}}$, the fire intensity is $\mathrm{25}$ MW m${}^{-\mathrm{1}}$ and the diffusion coefficient is $\mathrm{0.15}$ m${}^{\mathrm{2}}$ s ${}^{-\mathrm{1}}$.

[Figure omitted. See PDF]

For speeding-up all simulations presented in this section, the domain has been scaled by a factor of $\mathrm{4}$ to reduce the computation time. In the scaled mode each grid cell represents $\mathrm{4}\times \mathrm{4}$ times the area of each grid cell of the original domain. This scaling also affects the ROS and the turbulent diffusion coefficient: their value is reduced by a factor of $\mathrm{4}$. Fire intensity and the wind speed remain unaffected by the rescaling. It is noted that such scaling has no effect on the outputs of the simulations but helps to reduce the computation time.

It is remarked that, in the simulations presented in this paper, the firebrands are considered to be a sphere of constant radius for each simulation; however, in real situations all shapes and sizes of firebrands are produced from the fuel. It is also emphasized that the selection of the domain and other parameters do not correspond to any real fire, but an effort is made to chose values of the different parameters that lie within a valid range (Table ).

Discussion of the test cases with LSFire+

The mathematical formulation of the random effects presented in this article considers the combined impacts of turbulence and fire spotting. The main highlight of this formulation is its ability to incorporate the generation of secondary fires. Figure 4a–b show the evolution of the fire perimeter under the effect of turbulence and fire spotting. Figure 4a shows the effect of fire spotting in the presence of a barrier. This barrier is a fuel free zone with zero probability of ignition. The fire break zone stops the spread of the fire, but at 50 min, a new secondary fire appears beyond the barrier. This new fire line is completely detached from the main parent line, although it originates from the fire-spotting effects of the main fire line. The parent fire shows negligible growth in the head and cross wind directions, but the secondary fires increase quickly under the effect of wind. The wavy pattern of the fire perimeter results from the merging of multiple secondary fires in an idealised set-up of constant conditions. Figure 4b shows the scenario where secondary fires appear without the presence of a barrier. The effect of firebrands first appears around 130 min, and the new secondary fire behaves as a separate fire and evolves accordingly. In LSFire+ the effects of fire spotting occur in conjunction with turbulence, and both processes contribute towards the fire propagation. It is difficult to separate the effects of both processes individually, but a comparison of the increase in burned area due to turbulence and turbulence + fire spotting is presented in Fig. c. The total number of burned points is plotted at different times for two simulations: only turbulence, and turbulence $+$ fire spotting. All of the simulation parameters remain the same in both of the simulations. As the fire starts evolving, the line plots for both of the simulations overlap signifying that fire spotting has no visible contribution, but after $\mathrm{50}$ min the fire-spotting effect picks up and the burned area rapidly increases. At $\mathrm{140}$ min the increase in the burned area under the combined effect of the two random processes is almost three times the effect of turbulence alone.

(a) Line plot showing the sensitivity of the formulation to different values of fire intensity over constant wind conditions ($\mathrm{10}$ m s ${}^{-\mathrm{1}}$) and a constant firebrand radius ($\mathrm{0.015}$ m). The sensitivity is measured in terms of the total increase in the burned area when both fire spotting and turbulence are included over the case when no random effects are considered. The parameter ${\mathit{\beta }}_{\text{e}}$ is defined as ${\mathit{\beta }}_{\text{e}}=({\mathit{x}}_{\text{random}}-{\mathit{x}}_{\text{no-random}})/{\mathit{x}}_{\text{no-random}}$. The diffusion coefficient $\mathfrak{D}$ is $\mathrm{0.15}$ m${}^{\mathrm{2}}$ s ${}^{-\mathrm{1}}$. (b) The line plots show the lognormal distribution for selected values of fire intensity $I$ but constant values of wind speed $U=\mathrm{10}$ m s ${}^{-\mathrm{1}}$ and firebrand radius $r=\mathrm{0.03}$ m. According to the physical parameterisation, the plots can also be interpreted as the behaviour of the lognormal distribution for varying values of the parameter $\mathit{\mu }$ and a constant value of the parameter $\mathit{\sigma }$.

[Figure omitted. See PDF]

For the response analysis, a separate set of simulations are carried out, and the response is evaluated using the parameter ${\mathit{\beta }}_{\text{e}}$, which describes the effective increase in the burned area as follows: $${\mathit{\beta }}_{\text{e}}=({\mathit{x}}_{\text{random}}-{\mathit{x}}_{\text{no-random}})/{\mathit{x}}_{\text{ no-random}}.$$ ${\mathit{\beta }}_{\text{e}}$ is simply the increase in the number of burned grid points with respect to the simulation when no random effects are considered.

The simulated domain for the response analysis is chosen as a rectangle with the following dimensions (0 m, 6000 m )$\times$ (0 m, 6000 m). The simulations start at time $t=\mathrm{0}$ min and end at time $t=\mathrm{140}$ min. The grid spacing is $\mathrm{\Delta }x=\mathrm{\Delta }y=\mathrm{20}$ m. At time $t=\mathrm{0}$ min the initial fire line is a circle with a radius of $\mathrm{180}$ m centred at ${\mathit{x}}_c=$ (720 m, 3000 m). The horizontal wind is assumed as a constant field parallel to the vector $\mathit{j}=(\mathrm{1},\mathrm{0})$ and with modulus $\left|{\mathit{U}}_{\text{h}}\right|=\left|\right(U,V\left)\right|$ in this simulation set-up.

Response analysis to fire intensity

An increase in the fire intensity causes an increase in the burned area (see Eq.  for the definition of the ROS); simultaneously, the fire-spotting behaviour is also affected by any change in $I$. The parameter ${\mathit{\beta }}_{\text{e}}$ allows us to identify the contribution of fire spotting to the fire-propagation. Figure 5a shows the change in the burned area under the combined effect of turbulence and fire spotting with an increase in the fire-intensity. The two line plots correspond to a constant wind speed ($\mathrm{10}$ m s${}^{-\mathrm{1}}$) but two different firebrand radii, i.e. 0.015 and 0.030 m. According to the physical parameterisation of the lognormal shape parameters $\mathit{\mu }$ and $\mathit{\sigma }$, for this set of simulations (increase in fire intensity $I$), the parameter $\mathit{\mu }$ varies while the parameter $\mathit{\sigma }$ remains constant. For both sets of radii, an increase in the fire intensity shows a sharp increase in the burnt area for low fire intensities. A zoom-in of this sharp rise is also shown in the top right of Fig. 5. For the smaller firebrand radius, the fire-spotting effect shows a slight saturation between 15 and 25 MW m${}^{-\mathrm{1}}$, but with a further increase in the fire intensity, the contribution of fire spotting remains positive until 60 MW m${}^{-\mathrm{1}}$ and then it saturates again. Any further increase in the fire intensity causes a decrease in the importance of fire spotting. For the larger firebrand radius, a similar behaviour is observed for fire intensities lower than 15 MW m${}^{-\mathrm{1}}$, but with a further increase in $I$ the contribution from fire spotting declines before it starts to increase again. A zoom-in of the response analysis is important as the literature shows that for fire intensities around 8 MW m${}^{-\mathrm{1}}$ the vegetation type P. ponderosa is classified as a high “severity class” . For fire intensities ranging up to 10 MW m${}^{-\mathrm{1}}$ we can observe a rapid increase in the fire-spotting behaviour up to 4 MW m${}^{-\mathrm{1}}$. Any further increase in the fire intensity has a positive effect on the ROS and on the propagation of the main fire, such that less impact is observed on the fire line due to fire spotting. It is also interesting to note that for weak fires (less than 1 MW m${}^{-\mathrm{1}}$), the fire-spotting mechanism is independent of the firebrand radius.

In order to explain these observations, we plot the lognormal distribution for selected values of $I$ (Fig. b). These lognormal distribution plots show a general trend in the distribution with varying $\mathit{\mu }$ but constant $\mathit{\sigma }$. With an increase in $I$ (or $\mathit{\mu }$), the maximum probability increases but the distribution becomes increasingly skewed. For this particular set-up of parameters $\mathit{\mu }$ and $\mathit{\sigma }$ the skewness is more pronounced for fire intensities greater than 20 MW m${}^{-\mathrm{1}}$ (bottom right). For lower values of $I$ (less than 20 MW m${}^{-\mathrm{1}})$, the lognormal distribution tapers off slowly and the probability of a “long-range” ignition increases. This explains the large initial increase in the contribution of the fire-spotting behaviour for a low range of fire intensities. As the magnitude of the corresponding peak value also decreases with an increase in $\mathit{\mu }$ or $I$, this “long-range” ignition probability can only have a positive influence up to a certain threshold. This threshold is the range where parameter ${\mathit{\beta }}_{\text{e}}$ decreases or shows a saturation. Beyond this point, the effective contribution of the “long-range” probability diminishes, and the contribution of the “short-range” ignitions becomes increasingly important. The gradual increase in the effective burned area for both firebrand sizes (at fire intensities greater than 30 MW m${}^{-\mathrm{1}}$) can be attributed to this reason. For large values of $\mathit{\mu }$ the lognormal distributions tend to be similar though they retain their behaviour of becoming increasingly skewed (Fig. 5b right panel). Ideally, with increasing skewness, the “short-range” probability will lie much closer to the main fire line and the effective contribution of fire spotting should decrease. However, we only observe such behaviour for the smaller firebrands; this behaviour may exist for the large firebrand outside of the current range of simulations used in this study.

The effect of the fire intensity on fire-spotting behaviour can also be explained by the physical parameterisation provided in this paper. According to the physical parameterisation proposed here, an increase in the fire intensity increases the maximum loftable height; hence, the firebrands are ejected from elevated heights. A higher release height contributes to an increase in the firebrand activity at longer distances, and the initial increase in the fire perimeter follows this observation. Simultaneously, the increase in the firebrand ejection height over constant wind conditions causes the firebrands to travel further in the atmosphere before hitting the ground. Increased travel time for a firebrand promotes its combustion and the firebrand reaches the ground at a lower temperature (less “long-range” ignition probability) than its counterparts that are ejected at lower heights. The lower temperature of the firebrands leads to an inadequate heat exchange with the unburned fuel for a successful ignition; therefore, after reaching an area of maximum activity, the effective contribution of the “long-range” firebrands under same atmospheric conditions diminishes with increasing fire activity. This explains the initial decrease/saturation in firebrand activity. By the same token, the “short-range” firebrands have larger energy and become the dominant cause of the fire spreading. This range can be considered as the transition time when the “long-range” firebrands become less important and the “short-range” activity picks up. For the heavier firebrands a similar behaviour is expected, although the maximum loftable height under identical wind and fire conditions is lower. A lower loftable height decreases the maximum landing distance and the magnitude of the burned area is less, which is also evident from the lower magnitude of ${\mathit{\beta }}_{\text{e}}$ in the results.

(a) Line plot showing the sensitivity of the formulation to different wind conditions and radii when the fire intensity is constant ($\mathrm{50}$ MW m${}^{-\mathrm{1}}$). The measure of the effective increase in area (${\mathit{\beta }}_{\text{e}}$) and other simulation parameters are the same as defined in Fig. . (b) The line plots show the lognormal distribution for selected values of wind speed $U$ but constant values of fire intensity $I=\mathrm{5000}$ kW m${}^{-\mathrm{1}}$ and firebrand radius $r=\mathrm{0.03}$ m. According to the physical parameterisation, the plots can also be interpreted as the behaviour of the lognormal distribution for varying values of the parameter $\mathit{\sigma }$ and a constant value of the parameter $\mathit{\mu }$.

[Figure omitted. See PDF]

Response analysis to wind speed

Figure 6b highlights the simulation results with an increasing value of the wind velocity over constant fire intensity (50 MW m${}^{-\mathrm{1}}$). The results for the two different radii (0.015 and $\mathrm{0.030}$ m) are presented. The fire-spotting mechanism over varying wind speeds shows similar behaviour for the different sized firebrands. For both radii, the effective burnt area increases with increasing wind speeds, but after a certain threshold an increase in the wind speed leads to a decline in the effective burnt area. The line plot for $r=\mathrm{0.015}$ m shows that after a value of around $\mathrm{10}$ m s ${}^{-\mathrm{1}}$, the contribution of the firebrands decreases. Similarly, in the blue line plot for the radius $r=\mathrm{0.030}$ m, the effective increase in the burned area follows an identical pattern but the total increase in the burned area is of a lesser magnitude and shows a saturation around the maxima before it starts to decrease (around $\mathrm{22}$ m s ${}^{-\mathrm{1}}$). For bigger firebrands the onset of the maximum occurs at higher wind speeds and it sustains longer.

The lognormal distributions for selected wind speed values (but fixed $r$ and $I$) are plotted in Fig. b. These two plots show two different aspects of the response behaviour of the lognormal distribution when parameter $\mathit{\sigma }$ is varied but parameter $\mathit{\mu }$ is constant. Firstly, from Fig. 6b (left) it is evident that with an increase in the wind speed (increasing $\mathit{\sigma }$, constant $\mathit{\mu }$) the lognormal distribution shifts towards the left but the tails taper off slowly, making the distribution wider around the maximum. The increase in the width of the lognormal distribution leads to a larger area of “long-range” and “short-range” probability; hence, this explains the initial increase in the burned area in Fig. 6a. At the same time, the increasing wind speed also causes a decrease in the magnitude of the probability; therefore, beyond a certain threshold, the overall contribution from fire spotting starts to fall. The saturation in fire-spotting behaviour can be explained by the second aspect of the lognormal response. Figure 6b (right) shows that after a certain threshold of parameter $\mathit{\sigma }$ or wind speed, the lognormal distributions become increasingly similar. This threshold depends upon the value of parameter $\mathit{\mu }$, and for smaller values of $\mathit{\mu }$ (smaller $r$ or $I$) we have an early onset of this threshold. As the lognormal distributions tend to have similar probability distributions with increasing wind speed (Fig. 6b, right panel), the contribution from fire spotting also becomes similar and explains the saturation in the fire-spotting behaviour of the bigger firebrand radius.

This response of the model over different wind velocities can also be explained through physical parameterisation of $\mathit{\sigma }$. In terms of the physical quantities used in the parameterisation, it can be argued that strong winds can carry the firebrands longer distances from the main source and result in a larger fire perimeter (increasing “long-range probability”). Historically it has been reported that strong winds coupled with extremely dry conditions formed the perfect recipe for long-range fire spotting. Strong wind speeds can loft the smaller firebrands longer distances, but with an increasing wind speed the combustion process quickens and the firebrands reach the ground at a lower temperature. This fact explains the reduced effect of fire spotting on the burned area at high wind conditions. Conversely, a larger firebrand size can sustain longer in the atmosphere which means that their “long-range” probability is relatively higher than for smaller firebrands. This explains the onset of maximum burned area at $\mathrm{15}$ m s ${}^{-\mathrm{1}}$ instead of 12 m s${}^{-\mathrm{1}}$ for the 0.015 m radius. The heavier mass of bigger firebrands restricts their flight to shorter distances in comparison to the lighter firebrands and a lower magnitude of the burned area is therefore observed. The longer saturation in firebrand activity for the larger firebrands is due to their ability to maintain longer in air without burning out.

Conclusions

A mathematical formulation complete with its physical parameterisation RandomFront 2.3 to reproduce/mimic fire-spotting behaviour is presented in this article. In particular, we provide a versatile probabilistic model for fire spotting that is suitable for implementation as a post-processing scheme at each time step in any of the existing operational large-scale wildfire propagation models, without calling for any major changes in the original framework. This simple physical parameterisation is also an added advantage for real-time application. In this respect, RandomFront 2.3 is implemented in the coupled fire–atmosphere model WRF-SFIRE, and a simple test case is discussed. Moreover, the proposed formulation and parameterisation can be extended to include further variables such as moisture, spatial distribution of combustible, orography in addition to atmospheric variables such as wind and pressure that are available by running the simulation with WRF-SFIRE. This allows for an examination of the interaction between the concurrent factors in wildfires.

Furthermore, simulations with a simple propagation model based on the level set method, namely LSFire+, are performed to highlight the different responses of the model to varying fire intensities and wind conditions, and constant climatic conditions are assumed throughout the entire set-up. Results using different firebrand radii are also shown.

In both implementations, i.e. WRF-SFIRE and LSFire+, the simulations are simplified to highlight the physical applicability of the model.

The parameterisation RandomFront 2.3 provides a simple yet versatile addition to operational fire spread models reproducing the different physical aspects of the firebrand landing behaviour and simulating the occurrence of secondary fires as a result. The new secondary fires are also capable of modulating the weather around the primary fire and clear interactions between the two are observed. The wind conditions and fire-perimeter play an important role in determining the occurrence of the secondary fires.

The results highlight that the parameterisation is successful in reproducing the different physical aspects of firebrand landing behaviour. In this model, the complexities related to the shape and density of the firebrands are not considered and for brevity they are assumed to be spherical with a diameter in the order of the “collapse diameter”. The model also does not include an explicit computation of the time taken to reach the charred oxidation state, although a heating-before-burning mechanism is introduced in the mathematical formulation to serve a similar purpose. The inferences made from the simulations clearly fit within the physical aspects of the fire-spotting process. The increase in the wind speed causes an initial flare-up in the fire perimeter, but in really high wind conditions the wind enhances the propagation of the main fire and this reduces the effective contribution of fire spotting. Similarly, with increasing fire intensities, the fire intensity enhances the ROS of the main fire such that new ignitions of the unburned fuel ahead of the main fire by fire spotting are reduced.

Whilst many other studies focus on the long range landing distributions of the firebrands, most of them include rigorous computational aspects like LES which limit their applicability to the operational models for wildfire propagation. The simple yet powerful probabilistic formulation presented in this paper obeys the physical aspects of the fire-spotting process and provides scope for its applicability to operational fire spread models.

The code LSFire+ is developed in C and Fortran where the model RandomFront 2.3 acts as a post-processing routine at each time step in a LSM code for the front propagation implemented using LSMLIB ; the ROS is computed using the FireLib library . The numerical library LSMLIB is written in Fortran2008/OpenMP and propagates the fire line through standard algorithms for the LSM, including the fast marching method algorithms. Furthermore, the RandomFront 2.3 routine has also been implemented in the latest released version of WRF-SFIRE (https://github.com/openwfm/wrf-fire/, last access: 19 December 2018) by introducing new ad hoc routines. Both implementations of RandomFront 2.3 in LSFire+ and WRF-SFIRE are freely available at the official git repository of BCAM, Bilbao (https://gitlab.bcamath.org/atrucchia/randomfront-wrfsfire-lsfire, last access: 19 December 2018).

Simulations using WRF-SFIRE were performed on an Intel(R) Core(TM) i5-4310M 2.70 GHz CPU laptop with 8 GB of RAM. Each simulation that spanned 20 min of physical time took about 100 min of computational time.

Simulations using LSFire+ were performed over the HYPATIA cluster of BCAM, Bilbao, using OpenMP shared memory parallelism, running over 24 cores inside an Intel(R) Xeon(R) CPU E5-2680 v3 2.50GHz node with 128GB RAM. The computational time for each simulation, which spanned 140 min of physical time, was about 45 min.

Eighty percent of the computational cost in both cases, i.e. WRF-SFIRE and LSFire+, was due to the post-processing routine RandomFront 2.3. Computational time can be reduced in the future through further code optimisation.

The research was planned and coordinated by GP in collaboration with IK. GP and IK formulated the physical parameterisation. IK performed the simulations of the test cases using LSFire+ and was the main contributor with respect to writing the text. AT implemented the RandomFront 2.3 routines into WRF-SFIRE and performed the corresponding simulations. AB contributed to running the LSFire+ simulations, and VE contributed to running the WRF-SFIRE simulations.

The authors declare that they have no conflict of interest.

Acknowledgements

This research was supported by the Basque Government through the BERC 2014–2017 and BERC 2018–2021 programs. It was also funded by the Spanish Ministry of Economy and Competitiveness MINECO via the BCAM Severo Ochoa SEV-2013-0323 and SEV-2017-0718 accreditations, the MTM2013-40824-P “ASGAL” and MTM2016-76016-R “MIP” projects, and the PhD grant “La Caixa 2014”. The research began and was primarily developed at BCAM, Bilbao, during the postdoc fellowship of Inderpreet Kaur and the internship period of Anton Butenko. Edited by: Tomomichi Kato Reviewed by: two anonymous referees