- Email: [email protected]
Turbulence and fire-spotting effects into wild-land fire simulators
Accepted Manuscript
Turbulence and fire-spotting effects into wild-land fire simulators Inderpreet Kaur, Andrea Mentrelli, Fred ´ eric ´ Bosseur, Jean-Baptiste Filippi, Gianni Pagnini PII: DOI: Reference:
S1007-5704(16)30075-2 10.1016/j.cnsns.2016.03.003 CNSNS 3804
To appear in:
Communications in Nonlinear Science and Numerical Simulation
Received date: Revised date: Accepted date:
4 December 2015 29 February 2016 4 March 2016
Please cite this article as: Inderpreet Kaur, Andrea Mentrelli, Fred ´ eric ´ Bosseur, Jean-Baptiste Filippi, Gianni Pagnini, Turbulence and fire-spotting effects into wild-land fire simulators, Communications in Nonlinear Science and Numerical Simulation (2016), doi: 10.1016/j.cnsns.2016.03.003
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT
Highlights • A mathematical approach to model the random processes like turbupresented.
CR IP T
lence and fire-spotting into a DEVS based wild-land fire simulator is
• The performance of a Lagrangian and an Eulerian based wildfire sim-
ulators is compared in context of their application to simulate the ran-
AN US
dom behavior associated with the wildfires.
• The formulation is successful in reproducing the realistic situations associated with wildfires eg. flank fires and back fires, faster spread due to the effect of hot air and firebrands, fire crossing fuel break zones
M
and development of isolated secondary fires due to fire-spotting. • The results show the potential for implementation of the formulation
AC
CE
PT
Fire.
ED
into operational wildland fire simulators like WRF-SFIRE and Fore-
1
ACCEPTED MANUSCRIPT
Turbulence and fire-spotting effects into wild-land fire simulators
CR IP T
Inderpreet KAURa,∗, Andrea MENTRELLIb,a , Fr´ed´eric BOSSEURc , Jean-Baptiste FILIPPIc , Gianni PAGNINIa,d
AN US
a BCAM–Basque Center for Applied Mathematics Alameda de Mazarredo 14, E-48009 Bilbao, Basque Country – Spain b Department of Mathematics & 2 (AM) –Alma Mater Research Center on Applied Mathematics University of Bologna, via Saragozza 8, I-40123 Bologna, Italy c Laboratoire SPE–Sciences Pour l’Environnement CNRS - UMR6134 University of Corsica, BP 52, F-20250 Corte, Corsica – France d Ikerbasque–Basque Foundation for Science Calle de Mar´ a D´ıaz de Haro 3, E-48013 Bilbao, Basque Country – Spain
M
Abstract
This paper presents a mathematical approach to model the effects and the
ED
role of phenomena with random nature such as turbulence and fire-spotting into the existing wildfire simulators. The formulation proposes that the propagation of the fire-front is the sum of a drifting component (obtained from an
PT
existing wildfire simulator without turbulence and fire-spotting) and a random fluctuating component. The modelling of the random effects is embod-
CE
ied in a probability density function accounting for the fluctuations around the fire perimeter which is given by the drifting component. In past, this
AC
formulation has been applied to include these random effects into a wildfire simulator based on an Eulerian moving interface method, namely the Level ∗
Corresponding author Email address: [email protected] (Inderpreet KAUR)
Preprint submitted to Commun. Nonlinear Sci. Numer. Simulat.
March 12, 2016
ACCEPTED MANUSCRIPT
Set Method (LSM), but in this paper the same formulation is adapted for a wildfire simulator based on a Lagrangian front tracking technique, namely
CR IP T
the Discrete Event System Specification (DEVS). The main highlight of the present study is the comparison of the performance of a Lagrangian and an Eulerian moving interface method when applied to wild-land fire propagation.
Simple idealised numerical experiments are used to investigate the potential
applicability of the proposed formulation to DEVS and to compare its be-
AN US
haviour with respect to the LSM. The results show that DEVS based wildfire
propagation model qualitatively improves its performance (e.g., reproducing flank and back fire, increase in fire spread due to pre-heating of the fuel by hot air and firebrands, fire propagation across no fuel zones, secondary fire generation, . . . ) when random effects are included according to the present
M
formulation. The performance of DEVS and LSM based wildfire models is comparable and the only differences which arise among the two are due to
ED
the differences in the geometrical construction of the direction of propagation. Though the results presented here are devoid of any validation exercise
PT
and provide only a proof of concept, they show a strong inclination towards an intended operational use. The existing LSM or DEVS based operational
CE
simulators like WRF-SFIRE and ForeFire respectively can serve as an ideal basis for the same. Keywords: Wildland fire propagation, fire simulators, Level Set Method,
AC
Discrete Event System Specification, ForeFire, random phenomena, turbulence, fire-spotting
3
ACCEPTED MANUSCRIPT
1
1. Introduction and motivations Modelling of wild-land fire propagation serves as a crucial tool to combat
3
the high economic and environmental damage associated with the wildfires.
4
An effective and swift modelling can aid in building up an efficient plan for
5
fire suppression and restrict the life and property damage to a minimum.
6
Modelling the propagation of wild-land fire is a complex, multi-scale, multi-
7
physics and multi-discipline process, and the crux of this simulation lies in
8
delivering either a versatile or a specialised model which is easy to imple-
9
ment and is capable of providing timely information. With the availability
10
of higher computational power, various new and improved modelling ap-
11
proaches have been developed over the last decade. An extensive review by
12
Sullivan [1, 2, 3] provides a comprehensive overview of the wide spectrum of
13
physical, empirical, statistical and mathematical analogue models available
14
in literature.
ED
M
AN US
CR IP T
2
Most of the statistical and mathematical models available in literature
16
comprise of a Rate of Spread (ROS) formulation and a moving interface
17
technique. The analytical formulation of the ROS is developed independent
18
of the moving interface scheme and can be characterised in terms of the
19
wind speed, slope, fuel characterisation, combustion properties, along with
20
other experimental data. Various formulations for the ROS are available in
21
literature [4, 5, 6, 7, 8], and they are usually versatile enough to be used
AC
CE
PT
15
22
with most moving interface schemes. Among the different moving interface
23
methods available, Eulerian Level Set Method (LSM) is extensively used in
24
wildfires to represent the propagating fire-line [9, 6, 10]. LSM is a scheme
25
which represents a moving interface on a simple Cartesian grid. This method 4
ACCEPTED MANUSCRIPT
is particularly appropriate for wildfire simulations, as it permits an accurate
27
computation of the front normal vector which is used to propagate the front.
28
Other approaches available in literature which focus on wildfire propaga-
29
tion modelling include the Lagrangian Discrete Event System Specification
30
(DEVS) [11, 12, 13, 14] and Huygens’ principle [15, 16, 17, 18, 19, 20]. The
31
Lagrangian DEVS approach simulates the propagation of the wildfire with-
32
out defining any underlying mesh to represent the burning state. DEVS
33
being an event driven simulation considers time as a continuous variable and
34
permits faster simulations over higher resolution. On the other hand, wildfire
35
models based on Huygens’ principle utilise an elliptical spread at each point
36
of the fire-front and in some cases also benefit from some analytical results
37
[15, 16, 17, 18, 19, 20].
AN US
CR IP T
26
Wild-land fire propagation involves the heat transport to the surroundings
39
in accordance with the atmospheric wind, and other factors like fuel distribu-
40
tion, elevation, and orography. Different statistical and mathematical models
41
utilized for wild-fire propagation in operational simulators include a part of
42
the physical processes based on the experimental data or its mathematical
43
analogues, but are usually limited to the average value without taking into
44
account the random fluctuations. These models accommodate the variabil-
45
ity introduced by the changes in fuel type, wind conditions or orography but
46
fail in quantifying a sudden erratic output of the fire due to random effects,
AC
CE
PT
ED
M
38
47
e.g., due to turbulence or fire-spotting. The generation of heat introduces a
48
turbulent flow in the vicinity of the fire which in turn allows a wider contact
49
of the hot air with the unburned fuel [21, 22, 23, 24, 25, 10]. This convective
50
and radiative transfer of heat causes a rapid ignition of the unburned fuel 5
ACCEPTED MANUSCRIPT
bed and leads to a faster spread of the fire. Turbulence can also facilitate
52
the generation of spot fires when the firebrands are transported away from
53
the main fire due to advection [26, 27, 28, 29, 30]. The effect of such ran-
54
dom phenomena can range from subdued to catastrophic depending upon the
55
concurrent weather conditions and fuel characteristics. Existing operational
56
wildfire simulators are limited by their formulation to decipher these random
57
phenomena, and require implementation of additional schemes/variables to
58
represent such processes.
AN US
CR IP T
51
In past, few approaches have been proposed to model different random
60
processes within the wildfire propagation [31, 32, 33, 30] in order to provide
61
a more realistic fire spread. In Ref. [33], the authors incorporate the stochas-
62
tic fire-spotting phenomenon based on a continuous-time Markov chain on a
63
lattice. The state of the lattice site is controlled according to the local tran-
64
sition rate functions, where the physical processes like - fire spread, spotting
65
and burnout compete against each other. In Ref. [30], the authors integrate
66
the mathematical models with a cellular automata scheme to demonstrate
67
the variability in the landing patterns of firebrands. In Ref. [34], the authors
68
discuss a DEVS based model which incorporates the uncertainty by sampling
69
certain variables from arbitrary probability distributions.
ED
PT
The present authors in a number of their recent works [35, 36, 37, 38, 39,
40] have proposed a new approach to include the effects of random processes
AC
71
CE
70
M
59
72
in some operational wildfire propagation models. They derive a formulation
73
to modify the modelled position of the fire-front to include random processes.
74
The motion of the propagating fire-front is randomised by the addition of a
75
random displacement distributed according to a probability density function 6
ACCEPTED MANUSCRIPT
(PDF) corresponding to the turbulent transport of heat and to the fire-
77
spotting landing distance. The driving equation of the resulting averaged
78
process is analogous to an evolution equation of the reaction-diffusion type,
79
where the ROS controls the source term. If the motion is considered without
80
such random processes, the diffusive part disappears, and the spread of the
81
fire-line is identical to the one given by the chosen operational simulator
82
and driven only by the ROS. The proposed formulation to include random
83
processes holds true for any method of determination of the ROS and the
84
statistical spread is determined by the PDF of displacements of the random
85
contour points marked as active burning points. The ensemble averaging
86
of the displacements smoothens out of the noisy fluctuations and results
87
in a smooth fire-line contour. This formulation can be implemented into
88
existing wildfire models as a post-processing scheme at each time step without
89
calling for any major changes in the original framework. The efficacy of this
90
formulation by using a LSM based fire propagation model has already been
91
shown [35, 36, 38, 39, 40].
ED
M
AN US
CR IP T
76
In this study, the authors proceed with their preliminary attempt [41]
93
to include the proposed formalism for random processes (turbulence and
94
fire-spotting), already discussed and implemented into an Eulerian LSM
95
based wild-land fire simulator [35, 36, 37, 38, 39, 40], into a Lagrangian
96
DEVS based wild-land fire simulation model: ForeFire (https://github.
AC
CE
PT
92
97
com/forefireAPI/firefront) [13]. ForeFire is an open-source wildfire sim-
98
ulator developed by University of Corsica, France [42] to serve as a basis
99
for an operational simulator and provides flexibility to include new physics
100
options. ForeFire is based on a Lagrangian technique and integrates the so7
ACCEPTED MANUSCRIPT
called Balbi model for the ROS [5] with a DEVS [43] based front tracking
102
method to simulate the wildfire propagation. Coupling the front tracking
103
method with DEVS permits the utilisation of time as a continuous vari-
104
able; and hence removes the limitation of establishing a trade off between
105
the temporal resolution and the scale of the simulation according to the
106
Courant–Friedrichs–Lewy (CFL) condition. ForeFire model is developed de-
107
void of any description of random processes such as turbulent heat transfer
108
and fire-spotting and therefore hinders the evolution of the fire across fuel
109
break zones.
AN US
CR IP T
101
The aim of this article is the implementation of the mathematical for-
111
mulation for turbulence and fire-spotting in a DEVS based wildfire model
112
(ForeFire) to improve its efficiency. Operational utilisation of ForeFire for
113
active wildfire propagation is the main driving factor for improving its output
114
specially while dealing with realistic situations. This article tries to provide
115
a proof of concept through a series of idealised numerical simulations to
116
show the potential efficacy of the method in adapting the fire-line under the
117
effect of turbulence and fire-spotting. Besides this, the investigation also
118
provides a unique scope to compare the performance of Eulerian and La-
119
grangian approach based moving interface methods. Although a number of
120
studies on fluid flows focus on the comparison of performance of Eulerian and
121
Lagrangian methods for particle transport [44, 45, 46], none of the studies
AC
CE
PT
ED
M
110
122
concentrate on their comparison in context of wildfire propagation. In the
123
present article, various test cases are described to compare the performance
124
of these two schemes in application to wildfires.
125
The remaining paper is organised as follows. The salient features of DEVS 8
ACCEPTED MANUSCRIPT
and LSM based wildfire simulators are described in section 2. The descrip-
127
tion of the formulation for random effects is provided in section 3 and the
128
numerical set-up to include this formulation into the DEVS based simulator
129
is described in section 4. Section 5 describes the different simulation set-ups
130
for evaluation. A discussion of the numerical behaviour of the two methods
131
and a detailed comparison of their performance is provided in section 6. A
132
brief sensitivity analysis of the spatial distribution of the fire-brand landing
133
patterns is also presented in section 7. Concluding remarks are provided in
134
section 8.
135
2. Overview of wildfire simulators
AN US
CR IP T
126
Modelling of wildfire spread focusses on the development of a phenomeno-
137
logical and theoretical formula for a suitable representation of the ROS of
138
fire and the simulation of the propagation of the fire perimeter to recon-
139
struct the shape of fire over time. Various ROS models available in literature
140
[4, 5, 6, 7, 8] provide an analytical formulation of the ROS for given wind con-
141
ditions, slope and fuel characteristics. The available ROS models are usually
142
developed independent of the moving interface methods and user discretion
143
is advised in selecting one according to the requirements. Since this study
144
aims at studying the effect of the random processes using wildfire simulators
145
based over different moving interface methods, an analysis of different ROS
AC
CE
PT
ED
M
136
146
models is not presented in this article.
147
Wildfire simulators based on the DEVS and the LSM techniques are in-
148
herently different, but the concept of Lagrangian markers and the constant
149
perimeter resolution in the DEVS front tracking method can be considered 9
ACCEPTED MANUSCRIPT
analogous to the concept of active fuel points and constant action arc length
151
in LSM tracking method as discussed in Ref. [39]. This section provides a
152
brief overview of the salient features of the DEVS and LSM based wildfire
153
simulators.
154
2.1. DEVS based wildfire simulator
CR IP T
150
The temporal scheme used to simulate the fire perimeter in wild-land
156
fire simulator ForeFire is based on the DEVS [47]. DEVS handles the time
157
advancement in terms of the increment of physical quantities instead of a
158
discrete time step. The resulting front is a polygon whose marker points
159
have real (i.e.. non-discrete) coordinates instead of being located on nodes
160
of a regular mesh or grid. This computationally in-expensive Lagrangian
161
technique simulates the evolution of an interface without any underlying grid
162
representing the state of the system. Each fire-front is discretised into a series
163
of markers connected to the next marker by a piecewise linear segment. The
164
fire-break zones like roads and rivers are also represented as arbitrary shaped
165
polygons. Each marker is associated with a propagation speed and direction,
166
and each time a marker moves, the intersection with the neighbourhood is
167
checked to take care of collision and topological changes.
M
ED
PT
The propagation speed of the marker is defined by the ROS, while the
CE
168
AN US
155
propagation direction is defined by a front normal function. While spline
170
interpolation may be used to estimate this normal, the bisector angle made
AC
169
171
by the marker with its immediate left and right neighbours is used as a
172
computationally effective approximation of this outward normal.
173
In DEVS, time is treated as a continuous parameter and each marker
174
evolves according to its own independent time step. The time advancement 10
ACCEPTED MANUSCRIPT
of the markers is event based and all markers do not share the same time
176
step. Due to different time advances for each marker, the CFL condition
177
applies only locally and the markers move asynchronously. All events, trig-
178
gering marker movement, are time sorted in an event list and processed by
179
a scheduler. This self adjusting temporal resolution in Lagrangian formu-
180
lation is computationally efficient and provides an efficient way to simulate
181
the spatially in-homogeneous problem. The simulation in DEVS advances as
182
new events are generated; the generation of new events is managed by the
183
following two criteria:
184
Collision criterion: A collision happens when a marker moves into a differ-
185
ent area, e.g., from an unburned to a burned area or from an active fire
186
to a fuel break. Each collision generates a modification of the shape.
AN US
M
187
CR IP T
175
Quantum distance criterion: Quantum distance ∆q is defined as the maximum distance each marker is allowed to cover during advancement.
189
The actual resolution of the simulation is limited by this parameter and
190
details smaller than the quantum distance ∆q may not be accounted
191
for.
PT
Three type of events can be defined to control the front propagation:
CE
192
ED
188
decomposition, regeneration and coalescence. The decomposition function is
194
activated when a marker enters a different area (e.g., while approaching a fuel
AC
193
195
break zone). As soon as the marker enters a new area, two new markers are
196
created on the boundary of the new area. In regeneration, the markers are
197
redistributed to refine the shape; if two markers are separated by a distance
198
greater than the perimeter resolution ∆c, a new marker is generated between 11
ACCEPTED MANUSCRIPT
the two markers. The coalescence function reconstructs the fire-perimeter by
200
merging the markers. All markers with separation less than ∆r are merged
201
together. For stability and to avoid cross over of two markers, ∆r = ∆c/2
202
is assumed and the condition ∆c ≥ 2∆q should be respected. The precision
203
of the method is highly dependent on the choice of ∆q and ∆c. Quantum
204
distance ∆q, should be of a much higher resolution than the wind data for
205
minimal error. A detailed description of the DEVS front tracking method
206
can be found in Ref. [13].
AN US
CR IP T
199
ForeFire library is developed with C/C++/Java and Fortran bindings
208
for a UNIX compatible environment and compiled using a SWIG (http:
209
//swig.org/) build platform. NetCDF library with legacy C++ interface
210
is required to build ForeFire, and SCons (http://www.scons.org/) python
211
tool is used to make the library. A python based interface is used to launch
212
the simulations with software calls to the main code. A command line mode
213
interpreter is also available to run simulations.
214
2.2. LSM based wildfire simulator
ED
M
207
LSM is one of the widely used and successful tools for tracking fronts
216
in any dimension [48, 9]. It is particularly suitable for problems where the
217
speed of the evolving interface is dependent on the interface properties and
218
the boundary conditions at the interface. Wild-land fire propagation is one
219
of problems where LSM finds a frequent application to represent the evolving
AC
CE
PT
215
220
fire-fronts [6, 49, 10]. In particular, it is the basis for the operational software
221
WRF-SFIRE (http://github.com/jbeezley/wrf-fire/).
222
Let Γ be a simple closed curve representing the fire-front interface in
223
two-dimensions, and let Ω be the region bounded by the fire-front Γ. If 12
ACCEPTED MANUSCRIPT
the interface is made up of more than one closed surface, the domain Ω is
225
not simply connected and represents more than one independently evolving
226
bounded areas. Let γ : S ×[0, +∞[→ R be an implicit function defined on the
CR IP T
224
227
domain of interest S ⊆ R2 such that the level set γ(x, t) = γ∗ coincides with
228
the evolving front, i.e., Γ(t) = {x ∈ S | γ(x, t) = γ∗ }. If Γ is an ensemble
229
of n surfaces, the ensemble of the n interfaces is considered as an interface.
230
The evolution of the isocontour γ = γ∗ in time is governed by
232
233
(1)
If the motion of the interface is assumed to be directed towards the front b normal n then, Eq. (1) reduces to
b=− n
∇γ , k ∇γ k
M
231
AN US
Dγ ∂γ dx Dγ∗ = + · ∇γ = = 0. Dt ∂t dt Dt
(3)
ED
∂γ = V(x, t) k ∇γ k . ∂t
(2)
In application to wild-land fire modelling, an indicator function φ(x, t)
235
is introduced to allow a convenient identification of the burned area and the
236
fire perimeter
238
CE
φ(x, t) =
1, x ∈ Ω(t) ,
0, x 6∈ Ω(t) .
(4)
The boundary of Ω can be identified as the front-line contour of the wildfire
AC
237
PT
234
and V denotes the ROS of the fire-front.
239
The LSM code makes use of a general-purpose library which aims at pro-
240
viding a robust and efficient tool for studying the evolution of co-dimensional
241
fronts propagating in one-, two- and three-dimensional system. The library, 13
ACCEPTED MANUSCRIPT
written in Fortran2008/OpenMP, along with standard algorithms useful for
243
the calculation of the front evolution by means of the classical LSM, in-
244
cludes Fast Marching Method algorithms. All the results presented in this
245
paper are obtained by making use of the above-mentioned software; Matlab
246
(www.mathworks.com/products/matlab/) and open source softwares such
247
as SciPy [50] and Matplotlib [51] in the IPython framework [52] have been
248
used for visualisation purposes.
249
3. Random effects formulation
AN US
CR IP T
242
The proposed approach is based on the idea to split the motion of the
251
front position into a drifting part and a fluctuating part [35, 36, 37, 38,
252
39, 40, 53, 54]. This splitting allows specific numerical and physical choices
253
which can be useful to improve the algorithms and the models. In particular,
254
the drifting part can be related to existing methods for moving interfaces,
255
for example the Eulerian LSM [9] or the Lagrangian DEVS [11, 12]. The
256
fluctuating part is the result of a comprehensive statistical description of the
257
system which includes the random effects in agreement with the physical
258
properties of the system. As a consequence, the fluctuating part can have a
259
non-zero mean, implying that the drifting part does not correspond to the
260
average motion. Due to this fact, the present splitting is distinct from the
261
well-known Reynolds decomposition frequently adopted in turbulent flows.
AC
CE
PT
ED
M
250
262
The motion of the each burning point can be random due to the effect
263
of turbulence and/or fire-spotting. Let X ω (t, x0 ) be the ω realisation of
264
the random trajectory of an active burning point at time t with the initial
265
condition X ω (0, x0 ) = x0 . The trajectory of a single interface particle can be 14
ACCEPTED MANUSCRIPT
described by the one-particle PDF f ω (x; t) = δ(x − X ω (t, x0 )), where δ(x)
267
is the Dirac delta function. The random front contour can be derived by
268
using the sifting property of the Dirac delta function peaking at a random
269
position. This random position represents the stochastic trajectory which
270
describes the random motion of the front line and includes both drifting and
271
fluctuating part of the front position. In the ω realisation, the evolution in
272
time of the function γ ω (x, t), which embeds the random fire-line Γω is given
273
by [35, 39, 53, 54] γ (x, t) =
Z
S
274
γ(x, t)δ(x − X ω (t, x) dx .
(5)
The effective indicator of the burned area enclosed by the random front φe (x, t) : S × [0, +∞[→ [0, 1] may be defined as Z ω δ(x − X (t, x)) dx φe (x, t) = Ω(t) Z = hδ(x − X ω (t, x))i dx Ω(t) Z = f (x; t|x) dx Ω(t) Z = φ(x, t)f (x; t|x) dx ,
PT
ED
M
275
AN US
ω
CR IP T
266
(6) (7)
S
where f (x; t|x) = hδ(x − X ω (t, x))i is the PDF of the displacement of the
277
active burning points around the position x which is obtained from existing
278
wildfire simulators.
AC
CE
276
279
280
281
An arbitrary threshold value φth e is chosen to serve as the criterion to mark the region as burned or unburned, i.e., Ωe (x, t) = x ∈ S | φe (x, t) > φth e . The evolution of the effective indicator φe (x, t) can be obtained by using
15
ACCEPTED MANUSCRIPT
283
284
the Reynold transport theorem to Eq. (7) [55, 35] Z Z ∂f ∂φe = dx + ∇x · [V (x, t)f (x; t|x)] dx . ∂t Ω(t) ∂t Ω(t)
(8)
Since f (x; t|x) is a PDF, its evolution in time is here assumed to be governed by an equation of the type ∂f = εf , ∂t
CR IP T
282
(9)
where ε is an operator acting on space variable x, and not on x and t.
286
For example, the operator ε is the Laplacian when f (x; t|x) is the Gaussian
287
AN US
285
density function. Using (9), equation (8) can be written as Z ∂φe = εφe + ∇x · [V (x, t)f (x; t|x)] dx . ∂t Ω(t)
(10)
The velocity of the fire-line propagation is controlled by the ROS V and
289
defined through V = V(x, t)b n, while turbulence and fire-spotting phenomena
290
are modelled by modifying the PDF function. This method is compatible
291
with the estimation of the ROS by any technique.
ED
M
288
The modelling of the effects of the random processes is handled by the
293
PDF f (x; t|x), accounting for the two independent random variables repre-
294
senting turbulence and fire-spotting respectively. It should be noted that for
295
brevity, here fire-spotting is assumed to be an independent downwind phe-
296
nomenon; hence the effect of fire-spotting is accounted only for the leeward
297
part of the fire-line. Taking into account these assumptions the PDF for the
CE
random processes can be defined as Z ∞ b U ; t)q(`; t) d` , n b ·n bU ≥ 0 , G(x − x − ` n 0 f (x; t|x) = G(x − x; t), otherwise ,
AC
298
PT
292
16
(11)
ACCEPTED MANUSCRIPT
300
301
b U is the unit vector aligned with the mean wind direction. where n
Turbulent diffusion is assumed to be isotropic and modelled by a bi-
variate Gaussian PDF 1 G(x − x; t) = exp 4πDt
CR IP T
299
(x − x)2 + (y − y)2 4Dt
,
(12)
where D is the turbulent diffusion coefficient such that h(x − x)2 i = h(y − y)2 i =
303
2Dt. In the present model, the whole effect of the turbulent processes over
304
different scales is assumed to be parametrised only by the turbulent diffu-
305
sion coefficient. Also, only the turbulent fluctuations with respect to the
306
mean transport are considered. Hence, the characterisation of turbulence
307
is independent of the mean characteristics of the motion, and in particular
308
the estimation of the turbulent diffusion coefficient D is independent of the
309
mean wind. Moreover, since the simulations are performed with a flat terrain
310
and without boundaries in the spatial domain, horizontal isotropy is also as-
311
sumed. Within this framework, the mean wind does not break the isotropy,
312
because only the fluctuations are considered, but a non-trivial orography can
313
indeed break the isotropy, e.g., valleys or slopes. Even if an exact estima-
314
tion of D is out of the scope of the present study, D is one of the necessary
315
parameters for simulations and a quantitative estimation is required. The
316
desired diffusion coefficient D corresponds to the turbulent heat convection
317
generated by the fire. In literature, the ratio between the total heat trans-
318
fer and the molecular conduction of heat is known as the Nusselt number:
319
Nu = (D + χ)/χ. Where, χ is the thermal diffusivity of the air at ambi-
320
ent temperature and is well-known in literature: χ = 2 × 10−5 m2 s−1 . It is
321
also experimentally known that the relation between the Nusselt number Nu
322
and the Rayleigh number Ra is given by Nu ' 0.1 Ra1/3 [56]. The Rayleigh
AC
CE
PT
ED
M
AN US
302
17
ACCEPTED MANUSCRIPT
number is the ratio between the convection and the conduction of heat and
324
it is defined as Ra = γ ∆T g h3 /(νχ), where γ is the thermal expansion co-
325
efficient, ∆T the temperature difference between the bottom and the top of
326
the convective cell, h the dimension of the convection cell, g the acceleration
327
gravity and ν the kinematic viscosity. From above, D can be computed by
328
the formula:
CR IP T
323
D ' 0.1 χ [γ ∆T g h3 /(νχ)]1/3 − χ,
(13)
where values γ = 3.4 × 10−3 K−1 , g = 9.8 ms−2 and ν = 1.5 × 10−5 m2 s−1
330
are taken from literature. Since the heat transfer is considered in the hori-
331
zontal plane, perpendicular to the vertical “heating wall” embodied by the
332
fire, the length scale of the convective cell is assumed to be h = 100 m with
333
a temperature difference ∆T = 100 K. Considering these values, the scale of
334
the turbulent diffusion coefficient turns out to be approximately 104 times
335
the thermal diffusivity of air at ambient temperature (χ).
M
AN US
329
The jump length of the firebrands varies according to the meteorological
337
conditions, the intensity of the fire and the fuel characteristics [57, 58, 59].
338
A precise formulation of the landing distributions is difficult due to the lim-
339
ited small scale experimental results [60] and the complexities involved in
340
characterising the detailed aspect of firebrand generation and landing. The
341
firebrands land only in the positive direction and their landing distribu-
342
tion increases with distance to reach a maximum before dropping to zero
AC
CE
PT
ED
336
343
[61, 62]. In past, researchers have utilized different statistical distributions
344
to fit/reproduce the distribution pattern of particles deposited on ground
345
from a given source, e.g., the lognormal distribution is used to depict the dis-
346
tribution of the short-range firebrands (up to 4500 m) [27] while other studies 18
ACCEPTED MANUSCRIPT
347
also indicate the use of a Rayleigh distribution [63] or a Weibull distribution
348
[28] to describe the landing distributions.
350
In this study, a lognormal distribution is used to describe the downwind
CR IP T
349
distribution of the firebrands:
1 (ln `/`0 − µ(t))2 q(`; t) = √ exp − , 2 σ(t)2 2π σ(t) `
(14)
where µ(t) = hln `/`0 i and σ(t) = h(ln `/`0 − µ(t))2 i are the mean and the
352
standard deviation of ln `/`0 respectively, and `0 is a unit reference length.
AN US
351
Since the fuel ignition due to hot air and firebrands is not an instanta-
354
neous process, a suitable criterion related to an ignition delay is introduced
355
for application to wildfire modelling. This ignition delay was previously con-
356
sidered as a heating-before-burning mechanism due to the hot air [35, 36]
357
and generalised to include fire-spotting [39]. In particular, since the fuel can
358
burn because of two pathways, i.e., hot-air heating and firebrand landing,
359
the resistance analogy suggests that the resulting ignition delay can be ap-
360
proximately computed as resistances acting in parallel. Hence, representing
361
τh and τf as the ignition delay due to hot air and firebrands respectively, the
362
joint ignition delay τ is
ED
PT
CE
1 1 1 τh + τf = + = . τ τh τf τh τf
(15)
In this study, the delay effect is portrayed as a heating-before-burning
AC
363
M
353
364
mechanism and described by an accumulative process of the effective fire-
365
front over time, i.e., ψ(x, t) =
Z
t
φe (x, η)
0
19
dη , τ
(16)
ACCEPTED MANUSCRIPT
where, ψ(x, 0) = 0 corresponds to the initial unburned fuel. This accumula-
367
tion can be understood as an accumulation of heat that results in an increase
368
of the fuel temperature T (x, t). Finally, the function ψ(x, t) may be related
369
to the fuel temperature T (x, t) as follows: ψ(x, t) ∝
T (x, t) − Ta (x) , Tign − Ta (x)
CR IP T
366
(17)
where Tign is the ignition temperature, T ((x), 0) = Ta (x) and T (x, t) ≤ Tign .
371
Let ∆t be the ignition delay due to the effects of the random processes,
372
then after the elapse of the ignition delay time T (x, ∆t) = Tign , and with
373
the assumptions Ta Tign and the constant of proportionality equal to 1,
374
Eq. (17) reduces to
AN US
370
ψ(x, ∆t) = 1 .
(18)
Hence, as a cumulative effect of hot air and firebrands on the unburned fuel,
376
an ignition point x at time t is created when ψ(x, t) = 1, and the condition
377
φ(x, t) = 1 is also imposed.
378
4. Turbulence and fire-spotting effects into DEVS based wildfire
PT
ED
M
375
simulator
380
ForeFire is a DEVS based wild-land fire model developed bereft of any
CE
379
additional terms or parameters to model turbulence or fire-spotting. As
382
described earlier, DEVS is a Lagrangian method where each marker (active
383
burning point in the fire-line) advances by the quantum distance and the
384
advancement of each marker is managed by decomposition, regeneration or
385
coalescence functions. When one such marker arrives to a new region, the
386
marker decomposes to generate new markers at the interface between the two
AC
381
20
ACCEPTED MANUSCRIPT
regions. If this new region belongs to a fuel break zone, the new markers are
388
assigned null propagation velocity and their further evolution is inhibited.
389
But additional information about turbulence and fire spotting can allow the
390
fire-line to cross the fuel breaks. In this article, the DEVS based ForeFire
391
model is extended to include the effects of turbulence and fire-spotting by
392
post-processing the output obtained from the model at each time step. The
393
present method does not call for any major changes in the original framework
394
of DEVS, and can be versatilely used to reconstruct its output to include
395
these two processes. The detailed steps of the numerical procedure to re-
396
construct the output of DEVS based ForeFire are described as follows:
397
Step 1. Beginning with an initial fire-line, the DEVS front tracking method
398
is used to estimate the propagation of the markers to build up a new
399
fire perimeter for the next time step. This output is modified to
400
include the effects of turbulence and fire spotting by post-processing
401
numerical enrichment. This post-processing step is independent of
402
the definition of ROS.
ED
M
AN US
CR IP T
387
Step 2. The fire perimeter obtained from the active burning points is used to
404
construct the indicator function φ(x, t) (4). The spatial information
405
contained in φ(x, t) is sufficient to modify the fire-line with respect
to turbulence and fire-spotting and serves as the input to the post-
AC
407
CE
406
PT
403
408
processing step. Indicator function φ(x, t) has a value 0 over the region of the available fuel and 1 over the burned area.
409
Step 3. The spatial coordinates of the fire perimeter obtained from DEVS are
410
discrete values spread all over the domain and cannot be represented
411
by a specific mesh/grid. For brevity, the effective indicator function 21
ACCEPTED MANUSCRIPT
φe (x, t) (7) is generated over a Cartesian grid to facilitate the com-
413
putation of the function ψ(x, t) (16) over the same grid. Point in
414
polygon is used to generate the indicator function φ(x, t) from the lo-
415
cus of points representing the fire perimeter in DEVS. It is remarked
416
that φ(x, t), and φe (x, t) can also be computed from the real coor-
417
dinates representing the fire perimeter in DEVS, but the inherent
418
construction of the function ψ(x, t) mandates its computation over
419
a regular grid.
AN US
CR IP T
412
Step 4. The value of the effective indicator φe (x, t) is computed through the
421
numerical integration of the product of the indicator function φ(x, t)
422
and the PDF of fluctuations according to Eq. (7). The effect of
423
turbulence or fire-spotting is included by choosing the corresponding
424
PDF (11).
425
M
420
Step 5. The function ψ(x, t) is updated for each grid point by integration in time with the current value of φe (x, t) (16).
ED
426
Step 6. All points which satisfy the condition ψ(x, t) ≥ 1 are labelled as new
428
ignition points. The post-processing enrichment of the input fire-
429
front is completed at this step and the new markers describing the
430
enriched fire-front are used to define a new front valid at the next time step in DEVS framework. With the inclusion of heating-beforeburning mechanism, the new fire-front perimeter represents a larger
AC
432
CE
431
PT
427
433
area than the input fire perimeter.
434
Step 7. At the next time step, the new markers progress according to the
435
DEVS front tracking method, and the updated perimeter is again
436
subjected to the post-processing to enrich the fire front with the 22
ACCEPTED MANUSCRIPT
random fluctuations pertaining to turbulence and fire-spotting. The
438
sequence is repeated till the final “event time” step or till the fire
439
reaches the end of the domain.
CR IP T
437
Step 8. The DEVS models time as a continuous parameter and updates the
441
system in accordance with the events. The system is updated at user
442
defined “event times”. In between the consecutive events, no change
443
in the system occurs and a list of pending events is maintained in
444
chronological order. The event list is updated at the next “event
445
time”. Since, DEVS provides freedom to select an “event time”, it is
446
chosen to coincide with the maximum allowable time step according
447
to the CFL criterion. It is remarked that this choice of “event time”
448
limits the efficiency of DEVS computation (increase in the compu-
449
tation time due to excessive computations), but is indispensable in
450
order to introduce the heating-before-burning mechanism (16) in the
451
current formulation.
ED
M
AN US
440
It is remarked that the direction of the propagation vector for DEVS
453
is decided by the bisector of the angle formed between a marker and its
454
immediate neighbouring markers on left and right. This direction vector is a
455
weak approximation of the normal and can be very different from the actual
456
normal direction in the case of a non-circular profile.
AC
CE
PT
452
457
5. Simulation set-up
458
To estimate the performance of the DEVS front tracking method with the
459
inclusion of turbulence and fire-spotting, a series of numerical experiments
460
are performed. Similar set of experiments are also carried out with the LSM 23
ACCEPTED MANUSCRIPT
based simulator to enable a fair comparison between the two techniques. A
462
formulation developed in Ref. [6] and Ref. [10] is followed for LSM, while
463
for the DEVS method, ForeFire fire simulator [13] is used. To allow a direct
464
comparison between the two, both models are parametrised in an identical
465
set-up.
CR IP T
461
In the present study, for brevity, no particular type of vegetation is de-
467
fined and simulations are carried out with a pre-defined constant value of
468
ROS. It is assumed that the ROS remains constant for a particular domain.
469
The parametrisation of the numerical simulations is oversimplified and the
470
test cases are chosen with the purpose of highlighting the potential appli-
471
cability of the mathematical formulation. It is remarked that the chosen
472
parametrisations for the numerical simulations do not reproduce any of the
473
forest-fire or prescribed burn experiments but to emphasise the applicability
474
of the formulation, the values of various parameters are chosen to approxi-
475
mately lie in the valid range. The present scope of this work is to provide a
476
first look into the investigation of the performance of LSM and DEVS based
477
fire simulators in the context of inclusion of the effects due to turbulence and
478
fire-spotting.
PT
ED
M
AN US
466
A flat area of hypothetical homogeneous vegetation spread over a domain
480
size of 5000 m × 5000 m is selected for simulations. Different values of the
ROS are utilised for different test cases. The ROS is assumed to be 0.05 ms−1
AC
481
CE
479
482
483
484
b in no wind conditions, while in the presence of along the normal direction n wind, it is estimated by the 3% model [14], i.e.,
b, ROS = 0.03 U · n
(19)
where U is the mean wind velocity. Since, 3% model limits the propagation 24
ACCEPTED MANUSCRIPT
only towards the mean wind direction, a ROS formulation provided by Mallet
486
et al. [6] is followed to prescribe non-zero spread velocities for the flank and
487
rear fires: ROS(U, θ) =
√ εo + a U cosn θ,
CR IP T
485
if |θ| ≤
εo (α + (1 − α)| sin θ|), if |θ| >
π 2
,
π 2
,
(20)
where U 2 = U · U , εo is the flank velocity, (εo α) is the rear velocity with
489
α ∈ [0, 1], and θ is defined as the angle between the normal to the front
490
and the wind direction. In the present set-up, the parameters are defined as
491
α = 0.8, n = 3, a = 0.5 m1/2 s−1/2 , εo = 0.2 ms−1 , and U = 3 ms−1 [6].
AN US
488
In the case of the LSM based simulator, the domain is discretised with
493
a Cartesian grid of 20 m both in x and y directions, while for DEVS the
494
resolution of the simulation is established in the terms of quantum distance
495
∆q and perimeter resolution ∆c [13]. Quantum distance ∆q is defined as
496
the maximum allowable distance to be covered by a particle at each advance,
497
while a measure of ∆c is used to decompose/regenerate/coalesce two particles
498
on propagation. The choice of ∆q and ∆c is dependent on the type of problem,
499
and in the present study, two sets of values are utilized. The simulations
500
are performed with ∆q = 4 m, ∆c = 18 m for zero wind conditions, and
501
∆q = 2 m, ∆c = 8 m in the presence of wind. In DEVS based simulator, the
502
locus of the points constituting the initial front are wind up in the clockwise
503
direction to represent an expanding front (anticlockwise order represents a
504
contracting front). For each numerical simulation presented here, the initial
505
fire-front is represented by 200 markers. A brief study focussing on the
506
numerical behaviour of DEVS and LSM based techniques is also presented
507
in the next section.
AC
CE
PT
ED
M
492
25
ACCEPTED MANUSCRIPT
The mean wind, wherever used, is assumed to be constant in magnitude
509
and direction. The turbulent diffusion coefficient D, and ignition delays cor-
510
responding to hot air and firebrands heating are also assumed to be constant
511
throughout the simulations. In this study, the ignition delay correspond-
512
ing to the hot air τh and to the firebrand τf are assumed to be 10 min and
513
1 min respectively. Physically, a measure of the ignition delay can be inferred
514
from the fire intensity, fuel properties and the length scale at which the fuel
515
is exposed to the heat source [64]. According to Ref. [64], for the same
516
fuel characteristics and low fire intensities the ignition delay due to hot air is
517
higher than the firebrand landing, but in the cases of high fire intensities, the
518
ignition delay at different length scales can be comparable due to the func-
519
tional relationship between the ability of the fuel surface to absorb energy
520
and the ability of the fuel to conduct heat inward. Flammability classifica-
521
tion of different natural fuels available in literature [65, 66, 67] can aid in the
522
selection of an appropriate value for ignition delay.
ED
M
AN US
CR IP T
508
Appropriate values of D are chosen by using formula 13 to assess the
524
effect of increasing level of turbulence and to highlight its role in comparison
525
to the firebrand landing. In this study, three different values of D are utilised:
526
D = 0.30 m2 s−1 , D = 0.15 m2 s−1 , D = 0.075 m2 s−1 . Also, with reference to the lognormal distribution (14), the unit reference
length `0 in feet (1 ft = 0.30 m).
AC
528
CE
527
PT
523
529
Multiple idealised simulation tests are performed both in the presence and
530
the absence of wind by neglecting and considering the effects of the random
531
phenomena. The first case evaluates an isotropic growth of the fire-line in
532
the absence of wind by neglecting all the random processes. A fire-break 26
ACCEPTED MANUSCRIPT
zone 60 m wide is also considered in this case. In the second test, the spread
534
of fire-line for different ROS in different directions is studied. The third test
535
discusses the propagation of the fire-line with wind when the ROS is defined
536
according to the 3% of the mean wind, formula (19), and when it is defined
537
according to formula (20) . The random processes are neglected for the
538
first three test cases. The pure LSM and pure DEVS based wildfire models
539
fail while managing the realistic cases of fire propagation across the fire-
540
break zones, hence the fourth test evaluates the performance of DEVS when
541
turbulent processes are also considered in the model. The effect of turbulent
542
processes is studied both in the presence and the absence of the wind. The
543
fifth test evaluates the performance of the two simulators when fire-spotting
544
is also included along with turbulence. Fire-break zones are also introduced
545
in the last two tests (the fourth and the fifth) to observe the propagation of
546
the fire-line while encountering areas of null fuel. To simplify the simulation
547
and maintain equivalence between the two simulators, the region across and
548
behind the centre of the initial fire-line is demarcated as the leeward side and
549
the windward side respectively.
PT
ED
M
AN US
CR IP T
533
As an extension of the study, a brief analysis of the firebrand landing
551
distribution is also made. With reference to the lognormal distribution (14)
552
for firebrand landing, different values of µ and σ are chosen motivated by a
553
possible correlation with some physical aspects of wildfire like changing fire
AC
CE
550
554
intensity and meteorological conditions. For these simulations, only the LSM
555
based simulator is used and the wind speed is assumed to be 3 ms−1 , while
556
the ROS is given by the Rothermel model [4, 10].
557
All the simulations have been performed at BCAM–Basque Center for Ap27
ACCEPTED MANUSCRIPT
plied Mathematics (www.bcamath.org) in Bilbao, Basque Country – Spain.
559
6. Discussion
560
6.1. Numerical behaviour of DEVS and LSM
CR IP T
558
The stability of DEVS method relies on the choice of quantum distance
562
∆q and the perimeter resolution ∆c but the numerical simulations involving
563
a uniform ROS in all directions (e.g., a constant 0.05 ms−1 irrespective of
564
wind or direction) are unaffected by the choice of these two parameters. In a
565
homogenous domain, a uniform ROS forces all the markers to advance with
566
an identical speed at a same time step and leads to generation of synchronous
567
events. In this case, all the markers are driven by the same CFL condition,
568
and hence, the motion of the markers is analogous to an Eulerian technique.
569
On the contrary, markers with different speed generate different number of
570
discrete events and the resolution of ∆q and ∆c comes into picture. In con-
571
text of the idealised test cases involving wind, an evaluation of performance
572
of the method concerning different values of ∆q and ∆c is demonstrated in
573
Fig. 1. The numerical simulations show the evolution of the fire-line in the
574
presence of an east wind of 3 ms−1 with ROS as 3% of the normal wind, for
575
different pairs of ∆q and ∆c. Each contour represents the fire-front at 134
576
min. The perimeter resolution should be at least twice of the quantum dis-
577
tance (∆c ≥ 2∆q) to ensure that two markers do not cross-over [68]. Taking
AC
CE
PT
ED
M
AN US
561
578
into account this restriction, the upper plot of Fig. 1 highlights the differences
579
that arise in the evolution of the fire-front with respect to different values
580
of ∆c. DEVS employs a neighbourhood check for each marker to identify
581
its two immediate neighbours and utilises these neighbours to regulate the 28
ACCEPTED MANUSCRIPT
shape of the front by merging and generating new markers. A large value of
583
∆c can consider disconnected markers as neighbours and can lead to a unde-
584
sired generation or merging of markers. The loss of information by unsought
585
merging can cause the front to collapse. On the other hand, decreasing
586
the perimeter resolution ∆c without increasing the resolution of the marker
587
jumps causes the regeneration function to be evoked frequently and leads to
588
creation of a large number of markers with small separations. Over time,
589
too many markers at a coarse spatial scale leads to excessive unnecessary
590
computations and small numerical errors introduced by the computation can
591
escalate quickly over time. A finer scale of ∆q is desired with small values
592
of ∆c for a complete optimisation of the numerical simulation. The lower
593
plot in Fig. 1 shows the variation of the fire-front with respect to increasing
594
resolution of ∆q. Though a finer resolution of ∆q resolves the fire map at a
595
better resolution, it has no significant impact on the simulation results; on
596
the other hand, it increases the computation cost of the simulation. In the
597
present scenario, the homogeneous vegetation and topography of the present
598
domain causes the Lagrangian movement of markers to be driven predom-
599
inantly by regeneration and the fine scale information about the domain is
600
redundant. But on the other hand, in an inhomogeneous domain (involving
601
diverse vegetation, irregular topography, etc), where both decomposition and
602
coalescence play an important role, a finer marker jump length is necessary
CE
PT
ED
M
AN US
CR IP T
582
to account for the finer details of the domain. Usually, the value of ∆q is
604
taken to be around 2-3 m to ensure a correct representation of the vari-
605
able vegetation, changing topography and fire break zones like roads, rivers
606
etc. Concerning a pure quantitative analysis of the total error in DEVS, the
AC 603
29
ACCEPTED MANUSCRIPT
607
interested readers are referred to [68] and [69]. The numerical behaviour of the LSM is dependent on the time step (dt),
609
the grid resolution (dx, dy) and the maximum absolute speed of all points
610
on the grid. The iterative time step involved while solving the differential
611
equations needs to satisfy the CFL criteria for a converging solution. The
612
CFL condition defines the maximum limit of the time step in order to limit
613
the evolution of the contour to at most one grid cell at each time step. This
614
upper limit on the time step is important for the stability of the system and
615
a converging solution. In one-dimension, the CFL condition can be described
616
as:
AN US
CR IP T
608
dt ≤ k
min(dx, dy) vmax
(21)
where, vmax is the maximum absolute speed, and 0 < k < 1 is a constant. To
618
demonstrate the numerical behaviour of LSM, similar numerical simulations
619
with 3% model and an east wind of 3 ms−1 are repeated for different values
620
of dt and dx . The upper figure in Fig. 2 shows the fire-line contours at 380
621
min for various values of dt, and when dx = dy = 20 m. With values of dt
622
smaller than the threshold defined by CFL condition, the fire-line evolves
623
without any instabilities. Similar results are achieved with further decrease
624
in the time step but with an increase in the computation time. On the
625
other hand, when for the same speed function, dt is chosen to violate the
626
CFL criteria, the fire-line evolves with an undesirable behaviour and this
AC
CE
PT
ED
M
617
627
obscurity magnifies over time. In the same way, CFL condition can also
628
be interpreted as the constraint on the grid resolution. The lower figure in
629
Fig. 2 shows the evolution of the fire-line with respect to different values of
630
dx. The increase in the grid resolution over an identical time step introduces 30
ACCEPTED MANUSCRIPT
numerical instability and causes the fire-lines to evolve erroneously. For the
632
same velocity field, shorter time steps are required over small grid spacings
633
in order to respect the CFL condition, though it leads to an increase in the
634
computation costs. Hence, accounting for the numerical behaviour of the
635
LSM and to avoid the excessive computation, in the present study, dt is
636
defined according to (21), with k = 0.5.
CR IP T
631
It is remarked here that since the mathematical formulation for turbulence
638
and fire-spotting is implemented in DEVS and LSM based wildfire simulators
639
as a post-processing scheme, it preserves the original framework of the wildfire
640
simulators and does not affect the numerical stability and convergence of both
641
these methods.
642
6.2. Comparison of DEVS and LSM
M
AN US
637
In the first three test scenarios, the evolution of the fire-line without the
644
impact of random processes in both LSM and DEVS approaches is analysed.
645
Fig. 3 presents the propagation of the fire-line in no wind conditions for both
646
the simulators. The advancement of the initial circular fire-line follows an
647
isotropic growth and with identical initial conditions the two different moving
648
interface schemes provide a similar evolution of the fire-line. But the fire-line
649
fails to propagate across the fire-break zone in both the techniques. The
650
evolution is shown only up to 140 min, but an extended run up to 250 min
651
(not shown) indicates the limitation of the two approaches to reproduce the
AC
CE
PT
ED
643
652
653
fire jump across the no fuel zone. Fig. 4 presents the growth of an initial spot fire but with a non-homogeneous
654
ROS in absence of any wind. The different values of the ROS can be at-
655
tributed to different fuel types, but here no effort has been made to describe 31
ACCEPTED MANUSCRIPT
the fuel type. The fire-line propagates with different speed in the three dif-
657
ferent directions, but the isotropy is preserved. The results from these two
658
test cases indicate that in the absence of wind, the two advancing schemes:
659
LSM and DEVS, show an identical behaviour in simulating situations with
660
constant ROS and hence, this observation provides a crucial background and
661
scope for the comparison of situations with increasing variability and com-
662
plexity.
CR IP T
656
Fig. 5 shows the evolution of the fire-front of a circular initial profile with
664
radius 300 m in case of a weak wind of 3 ms−1 directed in the east direction.
665
The isochronous fronts are plotted at every 20 min and follow an oval shape
666
for both the simulators. DEVS fire contours diverge slightly from the mean
667
wind direction and develop into an increasing flanking fire over time. This
668
divergence in the evolution of fire-front occurs due to differences in the com-
669
putation of the normal between the two approaches. This fact can be very
670
well appreciated when an initial square profile is considered. Fig. 6 shows the
671
progress of the fire when the initial fire perimeter is considered to be a square
672
with side 600 m. Under the effect of a constant zonal wind and restricted 3%
673
model (19), the evolution in LSM strictly follows the initial square shape,
674
but in DEVS the active burning markers at the corners advance spuriously
675
to provide an additional flanking spread. The DEVS markers at the corners
676
are affected by the approximate normal and this leads to a larger flank fire.
AC
CE
PT
ED
M
AN US
663
677
Whereas, the normal direction for markers pointing towards the direction of
678
wind is identical for both the methods, hence the head fires arrive to same
679
location in both the simulations. The flank fires generated by the approxi-
680
mate construction of the normal seem to be in a more qualitative agreement 32
ACCEPTED MANUSCRIPT
681
with a realistic case, but within the stated analytical formulation of the ROS
682
(19), the lateral progress of the flank fires is an anomalous behaviour. The ROS formulation developed by Mallet et al. [6], here formula (20),
684
provides a simple way to study the evolution of flank fire by introducing
685
different ROS for the head, flank and rear directions. In the original paper,
686
θ is defined as the angle between the normal to the front and the wind
687
direction. Since, the normal computation in DEVS approach is approximate,
688
two separate tests are performed to evaluate the effect of such approximation
689
on the spread: firstly when θ is computed according to the original definition
690
[6], and secondly when θ is assumed to be the angle between the line joining
691
the front and the wind direction to enforce an identical definition of θ in
692
both the methods. Fig. 7 shows the evolution of fire when θ is computed
693
according to the original definition; the advancement of the head and rear
694
fires are identical, but the flank fire spread is spuriously larger in DEVS
695
because of the approximate computation of the front normal direction. On
696
the other hand, it is evident from Fig. 8 that in the second case, when the
697
computation of θ is not affected by the approximate computation of the
698
front normal in DEVS, the fire-line advancement is identical to LSM in all
699
directions.
AN US
M
ED
PT
As shown in Fig. 3, in case of fire-break zones, pure LSM and DEVS
based simulators are inherently unable to simulate the realistic situations of
AC
701
CE
700
CR IP T
683
702
fire overcoming a fire break. But Fig. 9 shows that with the introduction of
703
turbulence, the simulators can model the effect of hot air to overcome fire-
704
break zones. Here the value of turbulent diffusion coefficient D is assumed
705
to be 0.15 m2 s−1 . The evolution of the fire-line is almost similar for both 33
ACCEPTED MANUSCRIPT
the simulators, though a slight underestimation can be visibly observed in
707
the DEVS approach. A cross section of φ(x, t) and φe (x, t) at y = 2500 m
708
(Fig. 10) provides a more detailed illustration of the differences in the evolu-
709
tion of the fire-line between the two approaches. The introduction of marginal
710
differences with the inclusion of turbulence can be explained by the contrast
711
in the construction of the indicator function between the two simulators. In
712
LSM, the indicator function is defined over a Cartesian grid by definition,
713
while on the other hand, in DEVS, a gridded indicator function is constructed
714
from the original real coordinates. The point in polygon technique is used to
715
transform the real coordinates to a Cartesian grid. This technique provides
716
only an approximate measure of the same; small errors can be introduced
717
by the improper classification of the points lying very close to the boundary
718
[70], and can lead to a slight underestimation/overestimation of the original
719
fire-line.
M
AN US
CR IP T
706
Fig. 11 shows the results for the two simulators when the turbulent dif-
721
fusion coefficient is D = 0.30 m2 s−1 . Stronger turbulence causes a rapid
722
propagation of the fire-line and an earlier ignition across the fire-break zone.
723
A detailed analysis of the effect of varying turbulence over long-term propa-
724
gation with the LSM can be found in [35, 36].
PT
Fig. 12 presents the effect of inclusion of turbulence with a non-zero wind
profile and 3% model (19). The effect of turbulence is most pronounced in
AC
726
CE
725
ED
720
727
the direction of the wind and it can be clearly observed that with turbulence,
728
the spread over the head, flank and back fire-lines is faster and the head fire
729
is also able to overcome the fuel break. The quantum of increase in the
730
fire spread can also be appreciated in comparison to Fig. 5, which presents 34
ACCEPTED MANUSCRIPT
the ideal case with the same initial conditions. Both simulators show almost
732
similar characteristics in the spread of fire, though the flank fire has a slightly
733
larger spread in DEVS.
CR IP T
731
Another aspect contributing towards the increase in the fire spread due
735
to new fire ignitions generated by the firebrands is presented in Fig. 13. With
736
the inclusion of fire-spotting along with turbulence, the evolution of the fire-
737
front is faster in comparison to the effect of turbulence alone (Fig. 12). The
738
flank fire and the head fire are also well simulated in both the approaches,
739
and again a larger spread out the flanking fires is observed in the DEVS
740
scheme. The proposed mathematical formulation of spot-fires along with
741
turbulence is effective in mimicking the rapid advancement of the fire-line
742
due to the generation of new ignition points by the firebrands. It is also
743
effective in simulating the new ignition points developed across the no fuel
744
zones. Within the current parametrisation, the firebrand landing distances
745
lie close to the main fire and are not long enough to develop a new secondary
746
fire, but in the next section, a sensitivity study of the spatial pattern of the
747
firebrand landing distances is also presented.
PT
ED
M
AN US
734
Fig. 14 shows the time evolution of the fire-front in presence of a mean
749
wind velocity of 3 ms−1 and two fire-break zones. The results are obtained
750
when only the effects of turbulence are included in the formulation and the
751
diffusion coefficient D = 0.075 m2 s−1 . Similarly, Fig. 15 shows the time evo-
AC
CE
748
752
lution of fire-front but when both turbulence and fire-spotting are included
753
in the formulation. Both approaches follow an identical evolution pattern as
754
observed in other cases. As expected, the spread of the fire-front is faster and
755
is able to overcome both the fire-break zones when turbulence is included; 35
ACCEPTED MANUSCRIPT
but inclusion of fire-spotting along with turbulence causes the fire to propa-
757
gate at a much faster rate and results in a larger spread. Again, the major
758
differences which arise in the evolution of flank fires between the two meth-
759
ods are due to the dissimilarity in the construction of the direction of the
760
propagation of the active burning points.
761
7. Insight on firebrand distribution
CR IP T
756
In the simulations discussed previously, the phenomenon of fire-spotting
763
contributes towards an increase in the burning area and a faster propagation
764
of fire, but the landing distances are not sufficiently long enough to develop
765
new secondary fires detached from the primary fire.
AN US
762
Generally, the firebrands generated in low intensity fires fall close to the
767
primary fire and although they contribute towards a faster fire spread, the
768
separation from the primary fire is not large enough for the new ignitions
769
to develop a new secondary fire. But high intensity fires and extreme mete-
770
orological conditions such as strong wind/extremely high temperatures/low
771
humidity can facilitate the transport of embers over longer distances and the
772
area under firebrand attack can be sufficiently far away from the primary
773
fire [71, 72, 73, 74, 75]. In Ref. [76] and references therein, a brief account
774
of the impact of various historical forest fires is provided emphasising the
775
importance of strong winds, and extreme meteorological conditions on the
AC
CE
PT
ED
M
766
776
fire-spotting phenomenon. They describe different historical fires where low
777
intensity fires coupled with strong winds and low relative humidity caused
778
a large scale fire-spotting, while on the other hand, fires with high intensi-
779
ties in calm meteorological conditions had a negligible contribution from the 36
ACCEPTED MANUSCRIPT
780
firebrands. In lognormal distribution (14), a suitable choice of the parameters µ and
782
σ can be utilized to modify the firebrand landing distributions to represent
783
different fire intensities and wind conditions. In this study, the main interest
784
lies in the mathematical modelling of the firebrands landing distributions,
785
hence constant meteorological conditions are assumed for all the simulations
786
presented here. It is also assumed that each active ignition point is capable
787
of generating firebrands and all embers are controlled by identical transport
788
forces.
AN US
CR IP T
781
Fig. 16 and Fig. 17 show the variation of the area under the firebrand
790
attack with different values of µ and σ. A constant σ and an increasing µ
791
points towards a decrease in the effectiveness of the firebrands, and slows
792
down the fire propagation. A slight deviation from the pattern is observed
793
for µ = 15 and σ = 5, when the firebrands travel far enough to create a new
794
secondary fire.
ED
M
789
For an increasing value of σ (Fig. 17a-c) with constant µ, the fire evolu-
796
tion increases as larger number of firebrands land away from the source, but
797
the landing distances are still not long enough to delineate the new ignitions
798
from the main fire. With further increase in σ, a secondary fire emerges away
799
from the main source. Equal contribution to the firebrands by each flaming
800
point leads to a regular distribution pattern of the secondary fire. As time
AC
CE
PT
795
801
progresses, the primary fire also advances and covers up the gap with the
802
secondary fire. With further increase in σ, the probability of the firebrands
803
landing at farther distances increases, and this can lead to a potentially ex-
804
plosive firebrands attack. Fig. 17f shows the effect of fire-spotting at shorter 37
ACCEPTED MANUSCRIPT
scales contributing to a faster rate of spread, but the development of new
806
secondary fires is not observed. For this particular value of σ, the statistical
807
distribution of the long range landing distances lies outside the domain and
808
is impossible to resolve them in the limited time and spatial scales considered
809
in this simulation.
CR IP T
805
As a consequence of the idealised set-up for simulations, the new sec-
811
ondary fires evolve with the same intensity as the primary fire and also serve
812
as a source of generation of newer fire brands and henceforth the tertiary
813
fires. Fig. 18a-f demonstrates the generation of new secondary fires by the
814
primary fire, and tertiary fires by the secondary fires. As time progresses,
815
the each of the fire-line evolves independently along with a small contribution
816
from its parent fire-line, and eventually they merge among each other.
AN US
810
In a lognormal distribution, an increasing value of the µ causes a right
818
shift of the maximum, indicating towards a larger jump length, but the sim-
819
ulations suggest a decreasing trend in the landing distance. This can be
820
explained by the fact that the ignition of the unburned fuel is not an instan-
821
taneous process, but rather an accumulative process. The fuel starts burning
822
only when the temperature is high enough to allow spontaneous combustion.
823
Though a larger µ allows the firebrands to fall at larger distances but it limits
824
the frequency of landing. Insufficient number of firebrands landing at farther
825
distances diminishes the efficacy of their attack to cause combustion; hence
AC
CE
PT
ED
M
817
826
leads to a smaller contribution from fire-spotting. For an effective combus-
827
tion, sufficient number of firebrands should be available for an adequate heat
828
exchange with the unburned fuel. On the contrary, an increasing σ causes the
829
lognormal distribution to have a lower maximum and an elongated right tail 38
ACCEPTED MANUSCRIPT
of the distribution, but with a high probability of landing far away from the
831
main source. The increase in the number of firebrands landing away from
832
the maximum leads to a rapid ignition of the unburned fuel and a faster
833
propagation of the fire.
CR IP T
830
According to the increasing and decreasing trend of the landing distance
835
with respect to increasing σ and decreasing µ respectively, within the log-
836
normal characterisation assumed here for the firebrand landing distribution,
837
the reference length scale L of the landing distance can be stated as:
AN US
834
L∝ 838
σα , µβ
(22)
where α > 0 and β > 0 are exponents to be established.
From relation (22) it can be deduced that neither the maximum nor the
840
mean of the lognormal distribution can be considered as an estimation of
841
the landing distance length scale L. Physically, the landing distance can be
842
described in terms of the dimensions of the convective column, wind condi-
843
tions, fire intensity and fuel characteristics. Further analysis towards a phys-
844
ical parametrisation of fire-spotting is postponed for a future paper, where
845
firebrand distributions different from the lognormal will be also considered.
846
8. Conclusions
ED
PT
CE
847
M
839
This article describes a mathematical formulation to model the effects
of turbulence and fire-spotting in wild-land fire propagation models. Previ-
849
ously, this formulation has been applied to an Eulerian LSM based wildfire
850
propagation method, and this article describes the implementation of the
851
mathematical formulation into a Lagrangian DEVS based wildfire propaga-
852
tion method. This formulation splits the motion of the front into a drifting
AC
848
39
ACCEPTED MANUSCRIPT
part and a fluctuating part. The drifting part is independent of the fluctuat-
854
ing part and is given by the fire perimeter determined from one of the various
855
wildfire propagation methods that exist in literature (DEVS and LSM in the
856
present study); while the fluctuating part is generated by a comprehensive
857
statistical description of the system and includes the effects of random pro-
858
cesses in agreement with the physical properties of the process. This study
859
also provides a unique opportunity for the comparison of two simple wildfire
860
simulators based on the Eulerian LSM and on the Lagrangian DEVS scheme.
861
Numerical simulations based on identical set-ups are used to compare the
862
performance of the two approaches.
AN US
CR IP T
853
A series of idealised numerical simulations show that the proposed ap-
864
proach for random effects emerges to be suitable for both LSM and DEVS
865
based simulators to manage the real world situations related to random char-
866
acter of fire, e.g., increase in fire spread due to pre-heating of the fuel by hot
867
air, vertical lofting and transporting of firebrands, fire overcoming no fuel
868
zones, and development of new isolated secondary fires. The simulations
869
from both LSM and DEVS techniques follow nearly identical patterns, and
870
the only difference which distinguishes the outputs from the two techniques
871
is due to the distinction in the treatment of the direction of propagation. In
872
DEVS, the active burning points propagate in a direction given by the bisec-
873
tor of the angle made by each marker with its immediate left and right mark-
AC
CE
PT
ED
M
863
874
ers, in contrast to the normal direction of propagation used in LSM. Such
875
differences result in “spurious” flanking fires in DEVS, which are anomalous
876
with respect to the definition of the ROS, but qualitatively provide a more
877
“realistic” representation of the fire contour. Due to this fact, the two meth40
ACCEPTED MANUSCRIPT
ods can be considered complementary to each other for simple situations,
879
but with increasing complexity and the introduction of random processes a
880
validation exercise is necessary to single out the best representation of the
881
fire propagation.
CR IP T
878
A brief study on the sensitivity of the firebrand landing distribution to
883
the changes in the shape of the lognormal distribution indicates that the
884
shape parameters (µ and σ) can control the spread of the firebrand landing
885
and the lognormal distribution is capable of simulating new secondary fires
886
separated from the main fire perimeter. In future, the formulation for the
887
firebrand landing distance would be modified to include information about
888
the wind speed, fire intensity, radius of the firebrands and fuel characteristics.
889
Though without a validation test case, these simulations do provide an
890
insight into the flexibility of this formulation to incorporate effects of random
891
processes such as the turbulent heat diffusion and fire-spotting, specially
892
while characterising the landing distance of firebrands. It can be deduced
893
that the mathematical formulation is independent of the choice of the method
894
used to ascertain the drifting part and can serve as a versatile addition to
895
the existing fire spread simulators to include the effects of random spread.
896
The results from this study provide a support to the implementation of the
897
proposed approach for the effects of random phenomena into operational
898
softwares such as WRF-SFIRE and ForeFire. The source code utilized for all
AC
CE
PT
ED
M
AN US
882
899
the simulations in this article is hosted at https://github.com/ikaur17/
900
firefronts.
41
ACCEPTED MANUSCRIPT
901
Acknowledgement This research is supported by MINECO under Grant MTM2013-40824-P,
903
by Bizkaia Talent and European Commission through COFUND programme
904
under Grant AYD-000-226, and also by the Basque Government through the
905
BERC 2014-2017 program and by the Spanish Ministry of Economy and
906
Competitiveness MINECO: BCAM Severo Ochoa accreditation SEV-2013-
907
0323. The authors would also like to thank the two anonymous reviewers for
908
their insightful comments which helped in improving the manuscript both in
909
presentation and scientific content.
910
References
AN US
CR IP T
902
[1] A. L. Sullivan, Wildland surface fire spread modelling, 1990–2007. 1:
912
Physical and quasi-physical models, Int. J. Wildland Fire 18 (2009) 349–
913
368.
ED
M
911
[2] A. L. Sullivan, Wildland surface fire spread modelling, 1990–2007. 2:
915
Empirical and quasi-empirical models, Int. J. Wildland Fire 18 (2009)
916
369–386.
918
Simulation and mathematical analogue models, Int. J. Wildland Fire 18 (2009) 387–403.
AC
919
[3] A. L. Sullivan, Wildland surface fire spread modelling, 1990–2007. 3:
CE
917
PT
914
920
[4] R. C. Rothermel, A mathematical model for predicting fire spread in
921
wildland fires, Research Paper INT-115, USDA Forest Service, Inter-
922
mountain Forest and Range Experiment Station, Ogden, Utah 84401,
923
available at: http://www.treesearch.fs.fed.us/pubs/32533 (1972). 42
ACCEPTED MANUSCRIPT
924
925
[5] J. H. Balbi, F. Morandini, X. Silvani, J. B. Filippi, F. Rinieri, A physical model for wildland fires, Combust. Flame 156 (2009) 2217–2230. [6] V. Mallet, D. E. Keyes, F. E. Fendell, Modeling wildland fire propaga-
927
tion with level set methods, Comput. Math. Appl. 57 (2009) 1089–1101.
928
[7] M. Finney, Fire growth using minimum travel time methods, Can. J.
930
931
932
933
For. Res. 32 (2002) 420–1424.
[8] M. Finney, Calculation of fire spread rates across random landscapes,
AN US
929
CR IP T
926
Int. J. Wildland Fire 12 (2003) 167–174.
[9] J. A. Sethian, P. Smereka, Level set methods for fluid interfaces, Ann. Rev. Fluid Mech. 35 (2003) 341–372.
[10] J. Mandel, J. D. Beezley, A. K. Kochanski, Coupled atmosphere-
935
wildland fire modeling with WRF 3.3 and SFIRE 2011, Geosci. Model.
936
Dev. 4 (2011) 591–610.
ED
M
934
[11] H. Karimabadi, J. Driscoll, Y. A. Omelchenko, N. Omidi, A new asyn-
938
chronous methodology for modeling of physical systems: breaking the
939
curse of Courant condition, J. Comp. Phys. 205 (2005) 755–775.
CE
PT
937
[12] Y. A. Omelchenko, H. Karimabadi, HYPERS: A unidimensional asyn-
941
chronous framework for multiscale hybrid simulations, J. Comp. Phys.
AC
940
942
231 (2012) 1766–1780.
943
[13] J. B. Filippi, F. Morandini, J. H. Balbi, D. Hill, Discrete event front
944
tracking simulator of a physical fire spread model, Simulation 86 (2009)
945
629–646. 43
ACCEPTED MANUSCRIPT
[14] J. B. Filippi, V. Mallet, B. Nader, Evaluation of forest fire models on
947
a large observation database, Nat. Hazards Earth Syst. Sci. 14 (2014)
948
3077–3091.
950
951
952
[15] D. Anderson, E. Catchpole, N. de Mestre, T. Parkes, Modelling the spread of grass fires., J. Aus. Math. Soc. 23 (1982) 451–466.
[16] G. D. Richards, An elliptical growth model of forest fire fronts and its numerical solution, Int. J. Numer. Meth. Eng. 30 (1990) 1163–1179.
AN US
949
CR IP T
946
953
[17] G. D. Richards, The properties of elliptical wildfire growth for time
954
dependent fuel and meteorological conditions, Combust. Sci. Technol.
955
95 (1993) 357–383.
958
959
M
957
[18] G. D. Richards, A general mathematical framework for modelling twodimensional wildland fire spread, Int. J. Wildland Fire 5 (1995) 63–72. [19] J. Glasa, L. Halada, On elliptical model for forest fire spread modeling
ED
956
and simulation, Math. Comput. Simulat. 78 (2008) 76–88. [20] O. Rios, W. Jahn, G. Rein, Forecasting wind-driven wildfires using an
961
inverse modelling approach, Nat. Hazards Earth Syst. Sci. 14 (2014)
962
1491–1503.
CE
[21] T. L. Clark, M. A. Jenkins, J. Coen, D. Packham, A coupled
AC
963
PT
960
964
atmospheric-fire model: convective feedback on fire-line dynamics, J.
965
Appl. Meteor. 35 (1996) 875–901.
966
967
[22] B. E. Potter, A dynamics based view of atmosphere-fire interactions, Int. J. Wildland Fire 11 (2002) 247–255. 44
ACCEPTED MANUSCRIPT
[23] B. E. Potter, Atmospheric interactions with wildland fire behaviour - I.
969
Basic surface interactions, vertical profiles and synoptic structures, Int.
970
J. Wildland Fire 21 (2012) 779–801.
971
972
CR IP T
968
[24] B. E. Potter, Atmospheric interactions with wildland fire behaviour - II. Plume and vortex dynamics, Int. J. Wildland Fire 21 (2012) 802–817.
[25] C. B. Clements, S. Zhong, X. Bian, W. E. Heilman, D. W. Byun, First
974
observations of turbulence generated by grass fires, J. Geophys. Res. 113
975
(2008) D22102.
AN US
973
976
[26] N. Sardoy, J. L. Consalvi, B. Porterie, A. C. Fernandez-Pello, Modeling
977
transport and combustion of firebrands from burning trees, Combust.
978
Flame 150 (2007) 151–169.
[27] N. Sardoy, J. L. Consalvi, A. Kaiss, A. C. Fernandez-Pello, B. Porterie,
980
Numerical study of ground-level distribution of firebrands generated by
981
line fires, Combust. Flame 154 (2008) 478–488.
ED
M
979
[28] S. Kortas, P. Mindykowski, J. L. Consalvi, H. Mhiri, B. Porterie, Exper-
983
imental validation of a numerical model for the transport of firebrands,
984
Fire Safety J. 44 (2009) 1095–1102.
CE
PT
982
[29] H. A. Perryman, A mathematical model of spot fires and their manage-
986
ment implications, Master’s thesis, Humboldt State University, Arcata,
AC
985
987
CA (2009).
988
[30] H. A. Perryman, C. J. Dugaw, J. M. Varner, D. L. Johnson, A cellular
989
automata model to link surface fires to firebrand lift-off and dispersal,
990
Int. J. Wildland Fire 22 (2013) 428–439. 45
ACCEPTED MANUSCRIPT
992
993
994
[31] C. Favier, Percolation model of fire dynamic, Phys. Lett. A 330 (2004) 396–401. [32] H. Hunt, A new conceptual model for forest fires based on percolation
CR IP T
991
theory, Complexity 13 (2007) 12–17.
[33] D. Boychuk, W. J. Braun, R. J. Kulperger, Z. L. Krougly, D. A. Stan-
996
ford, A stochastic forest fire growth model, Environ. Ecol. Stat. 16 (2009)
997
133–151.
998
999
1000
AN US
995
[34] L. Ntaimo, B. P. Zeigler, M. J. Vasconcelos, B. Khargharia, Forest fire spread and supression in DEVS, Simulation 80 (2004) 479–500. [35] G.
Pagnini,
fire
Massidda,
method,
1003
dinia,
1004
http://publications.crs4.it/pubdocs/2012/PM12a/pagnini massidda-
1005
levelset.pdf and arXiv:1408.6129 (July 2012).
revised
version:
randomized
CRS4,
August
Pula 2014,
level-set
(CA),
Sar-
available
at:
PT
ED
Italy,
2012/PM12a,
the
effects
1002
Rep.
by
turbulence
in
Tech.
propagation
Modelling
1001
M
wildland
L.
[36] G. Pagnini, L. Massidda, The randomized level-set method to model
1007
turbulence effects in wildland fire propagation, in: D. Spano, V. Bacciu,
1008
M. Salis, C. Sirca (Eds.), Modelling Fire Behaviour and Risk. Proceedings of the International Conference on Fire Behaviour and Risk. ICFBR
AC
1009
CE
1006
1010
2011, Alghero, Italy, October 4–6 2011, 2012, pp. 126–131, ISBN 978-
1011
88-904409-7-7.
1012
[37] G. Pagnini, A model of wildland fire propagation including random ef-
1013
fects by turbulence and fire spotting, in: Proceedings of XXIII Congreso 46
ACCEPTED MANUSCRIPT
1014
de Ecuaciones Diferenciales y Aplicaciones XIII Congreso de Matem´atica
1015
Aplicada. Castell´o, Spain, 9–13 September 2013, 2013, pp. 395–403. [38] G. Pagnini, Fire spotting effects in wildland fire propagation, in:
1017
F. Casas, V. Mart´ınez (Eds.), Advances in Differential Equations and
1018
Applications, Vol. 4 of SEMA SIMAI Springer Series, Springer Interna-
1019
tional Publishing Switzerland, 2014, pp. 203–214, dOI 10.1007/978-3-
1020
319-06953-1 20. ISBN: 978-3-319-06952-4 (eBook: 978-3-319-06953-1).
AN US
CR IP T
1016
[39] G. Pagnini, A. Mentrelli, Modelling wildland fire propagation by track-
1022
ing random fronts, Nat. Hazards Earth Syst. Sci. 14 (2014) 2249–2263.
1023
[40] G. Pagnini, A. Mentrelli, The randomized level set method and an as-
1024
sociated reaction-diffusion equation to model wildland fire propagation,
1025
in: Proceedings of The 18th European Conference on Mathematics for
1026
Industry, ECMI 2014, Taormina, Italy, June 9–13, 2014, ECMI book
1027
subseries of Mathematics in Industry, Springer Heidelberg, accepted.
ED
M
1021
[41] I. Kaur, A. Mentrelli, F. Bosseur, J. B. Filippi, G. Pagnini, Wildland
1029
fire propagation modelling: A novel approach reconciling models based
1030
on moving interface methods and on reaction-diffusion equations, in:
1031
ˇıstek, T. Vejchodsk´ J. Brandts, S. Korotov, M. Kˇr´ıˇzek, K. Segeth, J. S´ y (Eds.), Proceedings of the International Conference on Applications of Mathematics 2015; Prague, Czech Republic, 18–21 November (2015),
AC
1033
CE
1032
PT
1028
1034
Institute of Mathematics – Academy of Sciences, Czech Academy of
1035
Sciences, Prague, Czech Republic, 2015, pp. 85–99.
1036
[42] J. B. Filippi, F. Bosseur, D. Grandi, ForeFire: open-source code for 47
ACCEPTED MANUSCRIPT
wildland fire spread models, in: D. X. Viegas (Ed.), Advances in forest
1038
fire research, Imprensa da Universidade de Coimbra, 2014, pp. 275–282.
1039
URL https://digitalis.uc.pt/handle/10316.2/34013
1040
1041
CR IP T
1037
[43] B. P. Zeigler, T. G. Kim, H. Praehofer, Theory of Modeling and Simulation, 2nd Edition, Academic Press, Inc., Orlando, FL, USA, 2000.
[44] J. C. Carvalho, M. T. Vilhena, D. M. Moreira, Comparison between Eu-
1043
lerian and Lagrangian semi-analytical models to simulate the pollutant
1044
dispersion in the PBL, Appl. Math. Model. 31 (2007) 120–129.
AN US
1042
[45] Z. Zhang, Q. Chen, Comparison of the Eulerian and Lagrangian methods
1046
for predicting particle transport in enclosed spaces, Atmos. Environ. 41
1047
(2007) 5236–5248.
M
1045
[46] M. Saidi, M. Rismanian, M. Monjezi, M. Zendehbad, S. Fatehiboroujeni,
1049
Comparison between Lagrangian and Eulerian approaches in predicting
1050
motion of micron-sized particles in laminar flows, Atmos. Environ. 89
1051
(2014) 199–206.
1054
object-oriented environment, Simulation 49 (1987) 219–230. [48] S. Osher, R. Fedkiw, Level Set Methods and Dynamic Implicit Surfaces, Springer–Verlag, New York, Berlin, Heidelberg, 2003.
AC
1055
PT
1053
[47] B. P. Zeigler, Hierarchical, modular discrete-event modelling in an
CE
1052
ED
1048
1056
1057
[49] R. G. Rehm, R. J. McDermott, Fire-front propagation using the level set method, Tech. Note 1611, Natl. Inst. Stand. Technol. (April 2009).
48
ACCEPTED MANUSCRIPT
1058
[50] E. Jones, T. Oliphant, P. Peterson, et al., SciPy: Open source scientific tools for Python (2001–).
1060
URL http://www.scipy.org/
1061
1062
1063
CR IP T
1059
[51] J. D. Hunter, Matplotlib: A 2D graphics environment, Comput. Sci. Eng. 9 (3) (2007) 90–95.
[52] F. P´erez, B. E. Granger, IPython: a system for interactive scientific computing, Comput. Sci. Eng. 9 (3) (2007) 21–29.
1065
URL http://ipython.org
1066
1067
AN US
1064
[53] A. Mentrelli, G. Pagnini, Random front propagation in fractional diffusive systems, Commun. Appl. Ind. Math. 6 (2) (2015) e–504. [54] A. Mentrelli, G. Pagnini, Front propagation in anomalous diffusive me-
1069
dia governed by time-fractional diffusion, J. Comp. Phys. 293 (2015)
1070
427–441.
1073
1074
ED
[56] J. J. Niemela, K. R. Sreenivasan, Turbulent convection at high Rayleigh numbers and aspect ratio 4, J. Fluid Mech. 557 (2006) 411–422.
[57] F. A. Albini, Potential spotting distance from wind-driven surface fires,
AC
1075
combustion, Phys. Rev. Lett. 107 (2011) 044503.
PT
1072
[55] G. Pagnini, E. Bonomi, Lagrangian formulation of turbulent premixed
CE
1071
M
1068
1076
Research Paper INT-309, U.S. Department of Agriculture Forest Service
1077
Intermountain Forest and Range Experiment Station (1983).
1078
[58] G. A. Morris, A simple method for computing spotting distances from
1079
wind-driven surface fires, Research Note INT-374, U.S. Department of 49
ACCEPTED MANUSCRIPT
1080
Agriculture Forest Service Intermountain Forest and Range Experiment
1081
Station (1987). [59] W. Catchpole, Fire properties and burn patterns in heterogeneous
1083
landscapes, in: R. A. Bradstock, J. E. Williams, M. A. Gill (Eds.),
1084
Flammable Australia: The Fire Regimes and Biodiversity of a Conti-
1085
nent, Cambridge University Press, 2002, Ch. 3, pp. 49–75.
CR IP T
1082
[60] M. E. Houssami, E. Mueller, J. C. Thomas, A. Simeoni, A. Filkov,
1087
N. Skowronski, M. R. Gallagher, K. Clark, R. Kremens, Experimental
1088
procedures characterising firebrand generation in wildfires, Fire Tech-
1089
nol.(Available online). doi:10.1007/s10694-015-0492-z.
1091
[61] K. D. Hage, On the dispersion of large particles from a 15-m source in the atmosphere, J. Appl. Meteor. 18 (1961) 534–539.
M
1090
AN US
1086
[62] S. L. Manzello, J. R. Shields, T. G. Cleary, A. Maranghides, W. E.
1093
Mell, J. C. Yang, Y. Hayashi, D. Nii, T. Kurita, On the development
1094
and characterization of a firebrand generator, Fire Safety J. 43 (2008)
1095
258–268.
1097
PT
[63] H. H. Wang, Analysis on downwind distribution of firebrands sourced
CE
1096
ED
1092
from a wildland fire, Fire Technol. 47 (2011) 321–340.
[64] H. E. Anderson, Forest fuel ignitibility, Fire Technol. 6 (1970) 312–319.
1099
[65] A. Dimitrakopoulos, K. K. Papaioannou, Flammability assessment of
AC
1098
1100
Mediterranean forest fuels, Fire Technol. 37 (2) (2001) 143–152.
50
ACCEPTED MANUSCRIPT
[66] K. J. M. Dickinson, J. B. Kirkpatrick, The flammability and energy
1102
content of some important plant species and fuel components in the
1103
forests of Southeastern Tasmania, J. Biogeogr. 12 (1985) 121–134.
1104
[67] D. M. J. S. Bowman, B. A. Wilson, Fuel characteristics of coastal mon-
1105
soon forests, Northern Territory, Australia, J. Biogeogr. 15 (1988) 807–
1106
817.
CR IP T
1101
[68] J. B. Filippi, F. Bosseur, C. Mari, C. Lac, P. L. Moigne, B. Cuenot,
1108
D. Veynante, D. Cariolle, J. H. Balbi, Coupled atmosphere-wildland fire
1109
modelling, J. Adv. Model. Earth Syst. 1 (2009) art. #11.
1112
1113
1114
tinuous systems, Simulation 78 (2002) 76–89.
M
1111
[69] E. Kofman, A second order approximation for DEVS simulation of con-
[70] K. Hormann, A. Agathos, The point in polygon problem for arbitrary polygons, Comp. Geom.-Theo. Appl. 20 (2001) 131–144.
ED
1110
AN US
1107
[71] R. Wilson, The Devil Wind and Wood Shingles: The Los Angeles Conflagration of 1961, National Fire Protection Assn., 1962.
1116
URL https://books.google.es/books?id=bWFCywAACAAJ
1118
sin Magazine of History 54 (1971) 246–272. URL http://www.jstor.org/stable/4634648
AC
1119
[72] P. Pernin, The great Peshtigo fire: An eyewitness account, The Wiscon-
CE
1117
PT
1115
1120
1121
[73] P. J. Pagni, Causes of the 20 October 1991 Oakland Hills conflagration, Fire Safety J. 21 (1993) 331–339.
51
ACCEPTED MANUSCRIPT
[74] J. Trelles, P. J. Pagni, Fire-induced winds in the 20 October 1991 Oak-
1123
land Hills fire, in: Y. Hasemi (Ed.), Proc. of Fifth International Sympo-
1124
sium on Fire Safety Science, Melbourne, Australia, International Asso-
1125
ciation of Fire Safety Science: Boston, MA, 1997, pp. 911–922.
1126
1127
CR IP T
1122
[75] J. G. Quintiere, Principles of Fire Behavior, Career Education Series, Delmar Publishers, Albany, N.Y., 1998.
[76] E. Koo, P. J. Pagni, D. R. Weise, J. P. Woycheese, Firebrands and
1129
spotting ignition in large-scale fires, Int. J. Wildland Fire 19 (2010)
1130
818–843.
AC
CE
PT
ED
M
AN US
1128
52
ACCEPTED MANUSCRIPT
3500
Y (m)
"c = 4 "c = 8 "c = 12 "c = 18
1500 "q = 1, "q = 2, "q = 3, "q = 4,
"c = 8 "c = 8 "c = 8 "c = 8
2500
AN US
500 0 2000
CR IP T
"q = 2, "q = 2, "q = 2, "q = 2,
2500
3000
3500
4000
4500
X (m)
Figure 1: The fire-lines showing the numerical behaviour of DEVS for different pairs of ∆c and ∆q. Each contour represents the fire-line at 134 min. The initial fire is a circle of
AC
CE
PT
ED
M
radius 300 m and represented by 200 markers in the model.
Figure 2: The fire-lines showing the numerical behaviour of LSM for different pairs of dt
and dx. Each contour represents the fire-line at 380 min and the initial fire is a circle of radius 300 m.
53
ACCEPTED MANUSCRIPT
(a)
120 200 280
Y (m)
80
40 0
1000 1000
2000
3000
0 16
X (m)
3000
2000
4000
1000 1000
2000
3000
4000
X (m)
AN US
2000
Time [min]
0 40 80 120 160 180 220 240
Y (m)
240
3000
(b)
4000
CR IP T
4000
Figure 3: Evolution in time of the fire-line contour without random processes in absence of wind for the a) LSM and b) DEVS based simulators. The initial fire-line is a circle of
M
radius 300 m. The fire-break zone is 60 m wide.
(a)
ED
8000
6000
CE
2000
4000 X (m)
4000
140
6000
2000
AC
20000
20 40
0 20 40 60 80 100 120 140 160 180
Y (m)
PT
Y (m)
60
160 120 80 180
Time [min]
6000
100
4000
(b)
8000
0
2000
4000 X (m)
6000
Figure 4: Evolution in time of the fire-line contour without random processes in absence of wind with a non-homogeneous ROS for a) LSM and b) DEVS based simulators. The initial fire-line is a circle of radius 30 m.
54
ACCEPTED MANUSCRIPT
(a)
3500
(b)
3000
60
80 120
20
0 40
Y (m)
160
2500 140
100
2000
2500
2000
2500
3000 X (m)
3500
4000
1500 2000
2500
3000 X (m)
3500
4000
AN US
1500 2000
0 20 40 60 80 100 120 140 160
CR IP T
3000 Y (m)
Time [min]
3500
Figure 5: Evolution in time of the fire-line contour without random processes with an east wind of 3 ms−1 for a) LSM and b) DEVS based simulators. The initial fire-line is a circle
ED
M
of radius 300 m.
(a)
PT
3500
(b)
3000
40 60
3000
AC
2500
Y (m)
2500
120 140
2000
1500 2000
0 20 40 60 80 100 120 140 160
160
0
CE
2500
100
80
Y (m)
20
Time [min]
3500
2000
3000 X (m)
3500
4000
1500 2000
2500
3000 X (m)
3500
4000
Figure 6: Same as Fig. 5, but when the initial fire-line is a square of side 600 m.
55
ACCEPTED MANUSCRIPT
(a)
6000
40 0 20 0 1 30 7050
60
Y (m)
Y (m)
6000 4000 2000
4000
2000
2000
4000 X (m)
6000
8000
0
0
Time [min]
2000
4000 X (m)
6000
0 10 20 30 40 50 60 70
8000
AN US
00
(b)
8000
CR IP T
8000
Figure 7: Evolution in time of the fire-line contour without random processes with ROS given by formula (20) and when θ is the angle between the outward normal in a contour point and the mean wind direction for a) LSM and b) DEVS based simulators. The mean
ED
M
wind velocity is 3 ms−1 in the positive x-direction.
(a)
8000
PT
0 10 20 30 40 50 60 70
6000
Y (m)
20
0
10
Time [min]
60
4000
30
50
Y (m)
6000
(b)
8000
4000
CE
40
AC
2000
00
2000
2000
70
4000 X (m)
6000
8000
0
0
2000
4000 X (m)
6000
8000
Figure 8: Same as Fig. 7, but when θ is the angle between the line joining a contour point
and the mean wind direction.
56
ACCEPTED MANUSCRIPT
Time [min]
4000
4000
3000
2000
2000
1000
1000
0
1000
2000 3000 X (m)
4000
5000
0
0
1000
2000 3000 X (m)
AN US
0
Time [min]
0 20 40 60 80 100 120 140
Y (m)
Y (m)
3000
(b)
5000
0 20 40 60 80 100 120 140
CR IP T
(a)
5000
4000
5000
Figure 9: Evolution in time of the fire-line contour with turbulence in absence of wind for a) LSM and b) DEVS based simulators. The initial fire-line is a circle of radius 300 m.
ED
M
The turbulent diffusion coefficient is D = 0.15 m2 s−1 .
1
AC
CE
?, ?e
PT
0.8
0.6
0.4 ? - DEVS ? e - DEVS ? - LSM ? e - LSM
0.2
0
0
1000
2000
3000
4000
5000
X (m)
Figure 10: Cross section of the indicator function φ(x, t) and φe (x, t) at y = 2500 m and
t = 140 min, corresponding to the front evolution showed in Fig. 9.
57
ACCEPTED MANUSCRIPT
Time [min]
4000
4000
3000
2000
2000
1000
1000
0
1000
2000 3000 X (m)
4000
5000
0
0
1000
2000 3000 X (m)
AN US
0
Time [min]
0 20 40 60 80 100 120 140
Y (m)
Y (m)
3000
(b)
5000
0 20 40 60 80 100 120 140
CR IP T
(a)
5000
4000
5000
0 20 40 60 80 100 120 140 160
ED
4000
3000
PT
Y (m)
Time [min]
CE
2000
0
1000
2000 3000 X (m)
Time [min]
0 20 40 60 80 100 120 140 160
4000
3000
2000
4000
1000
5000
AC
1000
(b)
5000
Y (m)
(a)
5000
M
Figure 11: Same as Fig. 9, but with turbulent diffusion coefficient D = 0.30 m2 s−1 .
0
1000
2000 3000 X (m)
4000
5000
Figure 12: Evolution in time of the fire-line contour with turbulence and a north wind of 3 ms−1 for a) LSM and b) DEVS based simulators. The turbulent diffusion coefficient is D = 0.15 m2 s−1 .
58
ACCEPTED MANUSCRIPT
Time [min]
Y (m)
4000
3000
3000
2000
1000
2000 3000 X (m)
4000
1000
5000
0
1000
2000 3000 X (m)
4000
5000
AN US
0
Time [min]
0 20 40 60 80 100 120 140 160
4000
2000
1000
(b)
5000
Y (m)
0 20 40 60 80 100 120 140 160
CR IP T
(a)
5000
Figure 13: Same as Fig. 9, but when both turbulence and fire-spotting are included with
M
D = 0.15 m2 s−1 , µ = 2.69 and σ = 1.25.
ED
PT
Y (m)
4000
3000
Time [min]
0 20 40 60 80 100 120 140 160 180 200 220
CE
2000
0
1000
2000 3000 X (m)
Time [min]
0 20 40 60 80 100 120 140 160 180 200 220
4000
3000
2000
4000
1000
5000
AC
1000
(b)
5000
Y (m)
(a)
5000
0
1000
2000 3000 X (m)
4000
5000
Figure 14: Evolution in time of the fire-line contour with turbulence, a north wind of 3 ms−1 and two fire-break zones for a) LSM and b) DEVS based simulators. The initial fire-line is a circle of radius 300 m. The turbulent diffusion coefficient is D = 0.075 m2 s−1 .
59
AN US
CR IP T
ACCEPTED MANUSCRIPT
(a)
Time [min]
0 20 40 60 80 100 120 140 160 180 200 220
3000
2000
0
1000
2000 3000 X (m)
4000
5000
ED
1000
M
Y (m)
4000
5000
(b)
Time [min]
0 20 40 60 80 100 120 140 160 180 200 220
4000
Y (m)
5000
3000
2000
1000
0
1000
2000 3000 X (m)
4000
5000
Figure 15: Same as Fig. 14, but when both turbulence and fire-spotting are included with
AC
CE
PT
D = 0.075 m2 s−1 , µ = 2.69 and σ = 1.25.
60
ACCEPTED MANUSCRIPT
10000
10000
8000
8000 150
00
2000
4000 6000 X (m)
00
8000 10000
(a)
2000
25
0
2000
4000 6000 X (m)
00
8000 10000
ED
00
2000
(c)
(d)
00
2000
100
0
4000 6000 X (m)
(e)
150
4000 2000
0
00
8000 10000
5
6000
17
Y (m)
175
125
150
PT
8000
0 10 25 50
2000
8000 10000
10000
75 25 50
CE
Y (m)
8000
4000 6000 X (m)
75 125
10000
175
0
2000
M
2000
125
4000
0 10 25 50
50
75
6000
1 75 50
75
125
Y (m)
100
4000
AC
8000 10000
8000
6000
4000
4000 6000 X (m)
AN US
10000
8000
6000
150
0
(b)
10000
Y (m)
50 100
2000
CR IP T
Y (m)
125
50 100 0
2000
4000
75 25
75 25
4000
6000 125
Y (m)
6000
2000
4000 6000 X (m)
8000 10000
(f)
Figure 16: Series of figures showing the variation in the fire propagation due to increasing value of µ. From a-f, µ varies from 13 to 18, in increments of 1, while σ = 5.
61
ACCEPTED MANUSCRIPT
10000
10000
8000
8000
50
25
2000
4000 6000 X (m)
0
2000 00
8000 10000
2000
(a) 8000
25
25
50
0
2000
4000 6000 X (m)
00
8000 10000
2000
4000 6000 X (m)
425
565
150
PT
Y (m)
4000
360
360 Y (m)
CE
5
100
X (m)
32
250 50
AC
8000
0
0
4000
475 400
0
50
530 50 4
35 0
12000
200
2000
50
16000
350
10050
250
00
300 200 1
5
32
4000
20000
350
16000
8000
(d)
ED
(c)
20000
8000 10000
300
00
0
2000
M
2000
5 12 71500
4000
50
75
6000
150
75
125
Y (m)
100
4000
12000
0
8000
15
AN US
10000
6000
8000 10000
(b)
10000
Y (m)
4000 6000 X (m)
375
00
25
50 100
4000
0
2000
100
5
75
25
75
CR IP T
4000
150 1
6000
Y (m)
200 150 12
175
Y (m)
6000
6000
8000
00
10000
(e)
4000
8000
12000 X (m)
16000
20000
(f)
Figure 17: Series of figures showing the variation in the fire propagation due to increasing value of σ. From a-f, σ varies from 4 to 6.5, in increments of 0.5, while µ = 15.
62
CR IP T
7000
Y (m)
5000
75
5000
75
5000
AN US
7000 Y (m)
7000
9000
100
9000
95
9000
75
Y (m)
ACCEPTED MANUSCRIPT
95
M
ED
1000 2000
4000 6000 X (m)
(b)
8000
1000 2000
25 75
25 75
25 75
(a)
8000
0
4000 6000 X (m)
50
1000 2000
3000
0
12
0
50
3000
50
3000
0
4000 6000 X (m)
8000
(c)
Figure 18: Series of figures showing the generation of secondary fires due to fire-spotting
AC
CE
PT
and their subsequent merging with the primary fire when µ = 15 and σ = 5.
63
Turbulence and fire-spotting effects into wild-land fire simulators Inderpreet Kaur, Andrea Mentrelli, Fred ´ eric ´ Bosseur, Jean-Baptiste Filippi, Gianni Pagnini PII: DOI: Reference:
S1007-5704(16)30075-2 10.1016/j.cnsns.2016.03.003 CNSNS 3804
To appear in:
Communications in Nonlinear Science and Numerical Simulation
Received date: Revised date: Accepted date:
4 December 2015 29 February 2016 4 March 2016
Please cite this article as: Inderpreet Kaur, Andrea Mentrelli, Fred ´ eric ´ Bosseur, Jean-Baptiste Filippi, Gianni Pagnini, Turbulence and fire-spotting effects into wild-land fire simulators, Communications in Nonlinear Science and Numerical Simulation (2016), doi: 10.1016/j.cnsns.2016.03.003
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT
Highlights • A mathematical approach to model the random processes like turbupresented.
CR IP T
lence and fire-spotting into a DEVS based wild-land fire simulator is
• The performance of a Lagrangian and an Eulerian based wildfire sim-
ulators is compared in context of their application to simulate the ran-
AN US
dom behavior associated with the wildfires.
• The formulation is successful in reproducing the realistic situations associated with wildfires eg. flank fires and back fires, faster spread due to the effect of hot air and firebrands, fire crossing fuel break zones
M
and development of isolated secondary fires due to fire-spotting. • The results show the potential for implementation of the formulation
AC
CE
PT
Fire.
ED
into operational wildland fire simulators like WRF-SFIRE and Fore-
1
ACCEPTED MANUSCRIPT
Turbulence and fire-spotting effects into wild-land fire simulators
CR IP T
Inderpreet KAURa,∗, Andrea MENTRELLIb,a , Fr´ed´eric BOSSEURc , Jean-Baptiste FILIPPIc , Gianni PAGNINIa,d
AN US
a BCAM–Basque Center for Applied Mathematics Alameda de Mazarredo 14, E-48009 Bilbao, Basque Country – Spain b Department of Mathematics & 2 (AM) –Alma Mater Research Center on Applied Mathematics University of Bologna, via Saragozza 8, I-40123 Bologna, Italy c Laboratoire SPE–Sciences Pour l’Environnement CNRS - UMR6134 University of Corsica, BP 52, F-20250 Corte, Corsica – France d Ikerbasque–Basque Foundation for Science Calle de Mar´ a D´ıaz de Haro 3, E-48013 Bilbao, Basque Country – Spain
M
Abstract
This paper presents a mathematical approach to model the effects and the
ED
role of phenomena with random nature such as turbulence and fire-spotting into the existing wildfire simulators. The formulation proposes that the propagation of the fire-front is the sum of a drifting component (obtained from an
PT
existing wildfire simulator without turbulence and fire-spotting) and a random fluctuating component. The modelling of the random effects is embod-
CE
ied in a probability density function accounting for the fluctuations around the fire perimeter which is given by the drifting component. In past, this
AC
formulation has been applied to include these random effects into a wildfire simulator based on an Eulerian moving interface method, namely the Level ∗
Corresponding author Email address: [email protected] (Inderpreet KAUR)
Preprint submitted to Commun. Nonlinear Sci. Numer. Simulat.
March 12, 2016
ACCEPTED MANUSCRIPT
Set Method (LSM), but in this paper the same formulation is adapted for a wildfire simulator based on a Lagrangian front tracking technique, namely
CR IP T
the Discrete Event System Specification (DEVS). The main highlight of the present study is the comparison of the performance of a Lagrangian and an Eulerian moving interface method when applied to wild-land fire propagation.
Simple idealised numerical experiments are used to investigate the potential
applicability of the proposed formulation to DEVS and to compare its be-
AN US
haviour with respect to the LSM. The results show that DEVS based wildfire
propagation model qualitatively improves its performance (e.g., reproducing flank and back fire, increase in fire spread due to pre-heating of the fuel by hot air and firebrands, fire propagation across no fuel zones, secondary fire generation, . . . ) when random effects are included according to the present
M
formulation. The performance of DEVS and LSM based wildfire models is comparable and the only differences which arise among the two are due to
ED
the differences in the geometrical construction of the direction of propagation. Though the results presented here are devoid of any validation exercise
PT
and provide only a proof of concept, they show a strong inclination towards an intended operational use. The existing LSM or DEVS based operational
CE
simulators like WRF-SFIRE and ForeFire respectively can serve as an ideal basis for the same. Keywords: Wildland fire propagation, fire simulators, Level Set Method,
AC
Discrete Event System Specification, ForeFire, random phenomena, turbulence, fire-spotting
3
ACCEPTED MANUSCRIPT
1
1. Introduction and motivations Modelling of wild-land fire propagation serves as a crucial tool to combat
3
the high economic and environmental damage associated with the wildfires.
4
An effective and swift modelling can aid in building up an efficient plan for
5
fire suppression and restrict the life and property damage to a minimum.
6
Modelling the propagation of wild-land fire is a complex, multi-scale, multi-
7
physics and multi-discipline process, and the crux of this simulation lies in
8
delivering either a versatile or a specialised model which is easy to imple-
9
ment and is capable of providing timely information. With the availability
10
of higher computational power, various new and improved modelling ap-
11
proaches have been developed over the last decade. An extensive review by
12
Sullivan [1, 2, 3] provides a comprehensive overview of the wide spectrum of
13
physical, empirical, statistical and mathematical analogue models available
14
in literature.
ED
M
AN US
CR IP T
2
Most of the statistical and mathematical models available in literature
16
comprise of a Rate of Spread (ROS) formulation and a moving interface
17
technique. The analytical formulation of the ROS is developed independent
18
of the moving interface scheme and can be characterised in terms of the
19
wind speed, slope, fuel characterisation, combustion properties, along with
20
other experimental data. Various formulations for the ROS are available in
21
literature [4, 5, 6, 7, 8], and they are usually versatile enough to be used
AC
CE
PT
15
22
with most moving interface schemes. Among the different moving interface
23
methods available, Eulerian Level Set Method (LSM) is extensively used in
24
wildfires to represent the propagating fire-line [9, 6, 10]. LSM is a scheme
25
which represents a moving interface on a simple Cartesian grid. This method 4
ACCEPTED MANUSCRIPT
is particularly appropriate for wildfire simulations, as it permits an accurate
27
computation of the front normal vector which is used to propagate the front.
28
Other approaches available in literature which focus on wildfire propaga-
29
tion modelling include the Lagrangian Discrete Event System Specification
30
(DEVS) [11, 12, 13, 14] and Huygens’ principle [15, 16, 17, 18, 19, 20]. The
31
Lagrangian DEVS approach simulates the propagation of the wildfire with-
32
out defining any underlying mesh to represent the burning state. DEVS
33
being an event driven simulation considers time as a continuous variable and
34
permits faster simulations over higher resolution. On the other hand, wildfire
35
models based on Huygens’ principle utilise an elliptical spread at each point
36
of the fire-front and in some cases also benefit from some analytical results
37
[15, 16, 17, 18, 19, 20].
AN US
CR IP T
26
Wild-land fire propagation involves the heat transport to the surroundings
39
in accordance with the atmospheric wind, and other factors like fuel distribu-
40
tion, elevation, and orography. Different statistical and mathematical models
41
utilized for wild-fire propagation in operational simulators include a part of
42
the physical processes based on the experimental data or its mathematical
43
analogues, but are usually limited to the average value without taking into
44
account the random fluctuations. These models accommodate the variabil-
45
ity introduced by the changes in fuel type, wind conditions or orography but
46
fail in quantifying a sudden erratic output of the fire due to random effects,
AC
CE
PT
ED
M
38
47
e.g., due to turbulence or fire-spotting. The generation of heat introduces a
48
turbulent flow in the vicinity of the fire which in turn allows a wider contact
49
of the hot air with the unburned fuel [21, 22, 23, 24, 25, 10]. This convective
50
and radiative transfer of heat causes a rapid ignition of the unburned fuel 5
ACCEPTED MANUSCRIPT
bed and leads to a faster spread of the fire. Turbulence can also facilitate
52
the generation of spot fires when the firebrands are transported away from
53
the main fire due to advection [26, 27, 28, 29, 30]. The effect of such ran-
54
dom phenomena can range from subdued to catastrophic depending upon the
55
concurrent weather conditions and fuel characteristics. Existing operational
56
wildfire simulators are limited by their formulation to decipher these random
57
phenomena, and require implementation of additional schemes/variables to
58
represent such processes.
AN US
CR IP T
51
In past, few approaches have been proposed to model different random
60
processes within the wildfire propagation [31, 32, 33, 30] in order to provide
61
a more realistic fire spread. In Ref. [33], the authors incorporate the stochas-
62
tic fire-spotting phenomenon based on a continuous-time Markov chain on a
63
lattice. The state of the lattice site is controlled according to the local tran-
64
sition rate functions, where the physical processes like - fire spread, spotting
65
and burnout compete against each other. In Ref. [30], the authors integrate
66
the mathematical models with a cellular automata scheme to demonstrate
67
the variability in the landing patterns of firebrands. In Ref. [34], the authors
68
discuss a DEVS based model which incorporates the uncertainty by sampling
69
certain variables from arbitrary probability distributions.
ED
PT
The present authors in a number of their recent works [35, 36, 37, 38, 39,
40] have proposed a new approach to include the effects of random processes
AC
71
CE
70
M
59
72
in some operational wildfire propagation models. They derive a formulation
73
to modify the modelled position of the fire-front to include random processes.
74
The motion of the propagating fire-front is randomised by the addition of a
75
random displacement distributed according to a probability density function 6
ACCEPTED MANUSCRIPT
(PDF) corresponding to the turbulent transport of heat and to the fire-
77
spotting landing distance. The driving equation of the resulting averaged
78
process is analogous to an evolution equation of the reaction-diffusion type,
79
where the ROS controls the source term. If the motion is considered without
80
such random processes, the diffusive part disappears, and the spread of the
81
fire-line is identical to the one given by the chosen operational simulator
82
and driven only by the ROS. The proposed formulation to include random
83
processes holds true for any method of determination of the ROS and the
84
statistical spread is determined by the PDF of displacements of the random
85
contour points marked as active burning points. The ensemble averaging
86
of the displacements smoothens out of the noisy fluctuations and results
87
in a smooth fire-line contour. This formulation can be implemented into
88
existing wildfire models as a post-processing scheme at each time step without
89
calling for any major changes in the original framework. The efficacy of this
90
formulation by using a LSM based fire propagation model has already been
91
shown [35, 36, 38, 39, 40].
ED
M
AN US
CR IP T
76
In this study, the authors proceed with their preliminary attempt [41]
93
to include the proposed formalism for random processes (turbulence and
94
fire-spotting), already discussed and implemented into an Eulerian LSM
95
based wild-land fire simulator [35, 36, 37, 38, 39, 40], into a Lagrangian
96
DEVS based wild-land fire simulation model: ForeFire (https://github.
AC
CE
PT
92
97
com/forefireAPI/firefront) [13]. ForeFire is an open-source wildfire sim-
98
ulator developed by University of Corsica, France [42] to serve as a basis
99
for an operational simulator and provides flexibility to include new physics
100
options. ForeFire is based on a Lagrangian technique and integrates the so7
ACCEPTED MANUSCRIPT
called Balbi model for the ROS [5] with a DEVS [43] based front tracking
102
method to simulate the wildfire propagation. Coupling the front tracking
103
method with DEVS permits the utilisation of time as a continuous vari-
104
able; and hence removes the limitation of establishing a trade off between
105
the temporal resolution and the scale of the simulation according to the
106
Courant–Friedrichs–Lewy (CFL) condition. ForeFire model is developed de-
107
void of any description of random processes such as turbulent heat transfer
108
and fire-spotting and therefore hinders the evolution of the fire across fuel
109
break zones.
AN US
CR IP T
101
The aim of this article is the implementation of the mathematical for-
111
mulation for turbulence and fire-spotting in a DEVS based wildfire model
112
(ForeFire) to improve its efficiency. Operational utilisation of ForeFire for
113
active wildfire propagation is the main driving factor for improving its output
114
specially while dealing with realistic situations. This article tries to provide
115
a proof of concept through a series of idealised numerical simulations to
116
show the potential efficacy of the method in adapting the fire-line under the
117
effect of turbulence and fire-spotting. Besides this, the investigation also
118
provides a unique scope to compare the performance of Eulerian and La-
119
grangian approach based moving interface methods. Although a number of
120
studies on fluid flows focus on the comparison of performance of Eulerian and
121
Lagrangian methods for particle transport [44, 45, 46], none of the studies
AC
CE
PT
ED
M
110
122
concentrate on their comparison in context of wildfire propagation. In the
123
present article, various test cases are described to compare the performance
124
of these two schemes in application to wildfires.
125
The remaining paper is organised as follows. The salient features of DEVS 8
ACCEPTED MANUSCRIPT
and LSM based wildfire simulators are described in section 2. The descrip-
127
tion of the formulation for random effects is provided in section 3 and the
128
numerical set-up to include this formulation into the DEVS based simulator
129
is described in section 4. Section 5 describes the different simulation set-ups
130
for evaluation. A discussion of the numerical behaviour of the two methods
131
and a detailed comparison of their performance is provided in section 6. A
132
brief sensitivity analysis of the spatial distribution of the fire-brand landing
133
patterns is also presented in section 7. Concluding remarks are provided in
134
section 8.
135
2. Overview of wildfire simulators
AN US
CR IP T
126
Modelling of wildfire spread focusses on the development of a phenomeno-
137
logical and theoretical formula for a suitable representation of the ROS of
138
fire and the simulation of the propagation of the fire perimeter to recon-
139
struct the shape of fire over time. Various ROS models available in literature
140
[4, 5, 6, 7, 8] provide an analytical formulation of the ROS for given wind con-
141
ditions, slope and fuel characteristics. The available ROS models are usually
142
developed independent of the moving interface methods and user discretion
143
is advised in selecting one according to the requirements. Since this study
144
aims at studying the effect of the random processes using wildfire simulators
145
based over different moving interface methods, an analysis of different ROS
AC
CE
PT
ED
M
136
146
models is not presented in this article.
147
Wildfire simulators based on the DEVS and the LSM techniques are in-
148
herently different, but the concept of Lagrangian markers and the constant
149
perimeter resolution in the DEVS front tracking method can be considered 9
ACCEPTED MANUSCRIPT
analogous to the concept of active fuel points and constant action arc length
151
in LSM tracking method as discussed in Ref. [39]. This section provides a
152
brief overview of the salient features of the DEVS and LSM based wildfire
153
simulators.
154
2.1. DEVS based wildfire simulator
CR IP T
150
The temporal scheme used to simulate the fire perimeter in wild-land
156
fire simulator ForeFire is based on the DEVS [47]. DEVS handles the time
157
advancement in terms of the increment of physical quantities instead of a
158
discrete time step. The resulting front is a polygon whose marker points
159
have real (i.e.. non-discrete) coordinates instead of being located on nodes
160
of a regular mesh or grid. This computationally in-expensive Lagrangian
161
technique simulates the evolution of an interface without any underlying grid
162
representing the state of the system. Each fire-front is discretised into a series
163
of markers connected to the next marker by a piecewise linear segment. The
164
fire-break zones like roads and rivers are also represented as arbitrary shaped
165
polygons. Each marker is associated with a propagation speed and direction,
166
and each time a marker moves, the intersection with the neighbourhood is
167
checked to take care of collision and topological changes.
M
ED
PT
The propagation speed of the marker is defined by the ROS, while the
CE
168
AN US
155
propagation direction is defined by a front normal function. While spline
170
interpolation may be used to estimate this normal, the bisector angle made
AC
169
171
by the marker with its immediate left and right neighbours is used as a
172
computationally effective approximation of this outward normal.
173
In DEVS, time is treated as a continuous parameter and each marker
174
evolves according to its own independent time step. The time advancement 10
ACCEPTED MANUSCRIPT
of the markers is event based and all markers do not share the same time
176
step. Due to different time advances for each marker, the CFL condition
177
applies only locally and the markers move asynchronously. All events, trig-
178
gering marker movement, are time sorted in an event list and processed by
179
a scheduler. This self adjusting temporal resolution in Lagrangian formu-
180
lation is computationally efficient and provides an efficient way to simulate
181
the spatially in-homogeneous problem. The simulation in DEVS advances as
182
new events are generated; the generation of new events is managed by the
183
following two criteria:
184
Collision criterion: A collision happens when a marker moves into a differ-
185
ent area, e.g., from an unburned to a burned area or from an active fire
186
to a fuel break. Each collision generates a modification of the shape.
AN US
M
187
CR IP T
175
Quantum distance criterion: Quantum distance ∆q is defined as the maximum distance each marker is allowed to cover during advancement.
189
The actual resolution of the simulation is limited by this parameter and
190
details smaller than the quantum distance ∆q may not be accounted
191
for.
PT
Three type of events can be defined to control the front propagation:
CE
192
ED
188
decomposition, regeneration and coalescence. The decomposition function is
194
activated when a marker enters a different area (e.g., while approaching a fuel
AC
193
195
break zone). As soon as the marker enters a new area, two new markers are
196
created on the boundary of the new area. In regeneration, the markers are
197
redistributed to refine the shape; if two markers are separated by a distance
198
greater than the perimeter resolution ∆c, a new marker is generated between 11
ACCEPTED MANUSCRIPT
the two markers. The coalescence function reconstructs the fire-perimeter by
200
merging the markers. All markers with separation less than ∆r are merged
201
together. For stability and to avoid cross over of two markers, ∆r = ∆c/2
202
is assumed and the condition ∆c ≥ 2∆q should be respected. The precision
203
of the method is highly dependent on the choice of ∆q and ∆c. Quantum
204
distance ∆q, should be of a much higher resolution than the wind data for
205
minimal error. A detailed description of the DEVS front tracking method
206
can be found in Ref. [13].
AN US
CR IP T
199
ForeFire library is developed with C/C++/Java and Fortran bindings
208
for a UNIX compatible environment and compiled using a SWIG (http:
209
//swig.org/) build platform. NetCDF library with legacy C++ interface
210
is required to build ForeFire, and SCons (http://www.scons.org/) python
211
tool is used to make the library. A python based interface is used to launch
212
the simulations with software calls to the main code. A command line mode
213
interpreter is also available to run simulations.
214
2.2. LSM based wildfire simulator
ED
M
207
LSM is one of the widely used and successful tools for tracking fronts
216
in any dimension [48, 9]. It is particularly suitable for problems where the
217
speed of the evolving interface is dependent on the interface properties and
218
the boundary conditions at the interface. Wild-land fire propagation is one
219
of problems where LSM finds a frequent application to represent the evolving
AC
CE
PT
215
220
fire-fronts [6, 49, 10]. In particular, it is the basis for the operational software
221
WRF-SFIRE (http://github.com/jbeezley/wrf-fire/).
222
Let Γ be a simple closed curve representing the fire-front interface in
223
two-dimensions, and let Ω be the region bounded by the fire-front Γ. If 12
ACCEPTED MANUSCRIPT
the interface is made up of more than one closed surface, the domain Ω is
225
not simply connected and represents more than one independently evolving
226
bounded areas. Let γ : S ×[0, +∞[→ R be an implicit function defined on the
CR IP T
224
227
domain of interest S ⊆ R2 such that the level set γ(x, t) = γ∗ coincides with
228
the evolving front, i.e., Γ(t) = {x ∈ S | γ(x, t) = γ∗ }. If Γ is an ensemble
229
of n surfaces, the ensemble of the n interfaces is considered as an interface.
230
The evolution of the isocontour γ = γ∗ in time is governed by
232
233
(1)
If the motion of the interface is assumed to be directed towards the front b normal n then, Eq. (1) reduces to
b=− n
∇γ , k ∇γ k
M
231
AN US
Dγ ∂γ dx Dγ∗ = + · ∇γ = = 0. Dt ∂t dt Dt
(3)
ED
∂γ = V(x, t) k ∇γ k . ∂t
(2)
In application to wild-land fire modelling, an indicator function φ(x, t)
235
is introduced to allow a convenient identification of the burned area and the
236
fire perimeter
238
CE
φ(x, t) =
1, x ∈ Ω(t) ,
0, x 6∈ Ω(t) .
(4)
The boundary of Ω can be identified as the front-line contour of the wildfire
AC
237
PT
234
and V denotes the ROS of the fire-front.
239
The LSM code makes use of a general-purpose library which aims at pro-
240
viding a robust and efficient tool for studying the evolution of co-dimensional
241
fronts propagating in one-, two- and three-dimensional system. The library, 13
ACCEPTED MANUSCRIPT
written in Fortran2008/OpenMP, along with standard algorithms useful for
243
the calculation of the front evolution by means of the classical LSM, in-
244
cludes Fast Marching Method algorithms. All the results presented in this
245
paper are obtained by making use of the above-mentioned software; Matlab
246
(www.mathworks.com/products/matlab/) and open source softwares such
247
as SciPy [50] and Matplotlib [51] in the IPython framework [52] have been
248
used for visualisation purposes.
249
3. Random effects formulation
AN US
CR IP T
242
The proposed approach is based on the idea to split the motion of the
251
front position into a drifting part and a fluctuating part [35, 36, 37, 38,
252
39, 40, 53, 54]. This splitting allows specific numerical and physical choices
253
which can be useful to improve the algorithms and the models. In particular,
254
the drifting part can be related to existing methods for moving interfaces,
255
for example the Eulerian LSM [9] or the Lagrangian DEVS [11, 12]. The
256
fluctuating part is the result of a comprehensive statistical description of the
257
system which includes the random effects in agreement with the physical
258
properties of the system. As a consequence, the fluctuating part can have a
259
non-zero mean, implying that the drifting part does not correspond to the
260
average motion. Due to this fact, the present splitting is distinct from the
261
well-known Reynolds decomposition frequently adopted in turbulent flows.
AC
CE
PT
ED
M
250
262
The motion of the each burning point can be random due to the effect
263
of turbulence and/or fire-spotting. Let X ω (t, x0 ) be the ω realisation of
264
the random trajectory of an active burning point at time t with the initial
265
condition X ω (0, x0 ) = x0 . The trajectory of a single interface particle can be 14
ACCEPTED MANUSCRIPT
described by the one-particle PDF f ω (x; t) = δ(x − X ω (t, x0 )), where δ(x)
267
is the Dirac delta function. The random front contour can be derived by
268
using the sifting property of the Dirac delta function peaking at a random
269
position. This random position represents the stochastic trajectory which
270
describes the random motion of the front line and includes both drifting and
271
fluctuating part of the front position. In the ω realisation, the evolution in
272
time of the function γ ω (x, t), which embeds the random fire-line Γω is given
273
by [35, 39, 53, 54] γ (x, t) =
Z
S
274
γ(x, t)δ(x − X ω (t, x) dx .
(5)
The effective indicator of the burned area enclosed by the random front φe (x, t) : S × [0, +∞[→ [0, 1] may be defined as Z ω δ(x − X (t, x)) dx φe (x, t) = Ω(t) Z = hδ(x − X ω (t, x))i dx Ω(t) Z = f (x; t|x) dx Ω(t) Z = φ(x, t)f (x; t|x) dx ,
PT
ED
M
275
AN US
ω
CR IP T
266
(6) (7)
S
where f (x; t|x) = hδ(x − X ω (t, x))i is the PDF of the displacement of the
277
active burning points around the position x which is obtained from existing
278
wildfire simulators.
AC
CE
276
279
280
281
An arbitrary threshold value φth e is chosen to serve as the criterion to mark the region as burned or unburned, i.e., Ωe (x, t) = x ∈ S | φe (x, t) > φth e . The evolution of the effective indicator φe (x, t) can be obtained by using
15
ACCEPTED MANUSCRIPT
283
284
the Reynold transport theorem to Eq. (7) [55, 35] Z Z ∂f ∂φe = dx + ∇x · [V (x, t)f (x; t|x)] dx . ∂t Ω(t) ∂t Ω(t)
(8)
Since f (x; t|x) is a PDF, its evolution in time is here assumed to be governed by an equation of the type ∂f = εf , ∂t
CR IP T
282
(9)
where ε is an operator acting on space variable x, and not on x and t.
286
For example, the operator ε is the Laplacian when f (x; t|x) is the Gaussian
287
AN US
285
density function. Using (9), equation (8) can be written as Z ∂φe = εφe + ∇x · [V (x, t)f (x; t|x)] dx . ∂t Ω(t)
(10)
The velocity of the fire-line propagation is controlled by the ROS V and
289
defined through V = V(x, t)b n, while turbulence and fire-spotting phenomena
290
are modelled by modifying the PDF function. This method is compatible
291
with the estimation of the ROS by any technique.
ED
M
288
The modelling of the effects of the random processes is handled by the
293
PDF f (x; t|x), accounting for the two independent random variables repre-
294
senting turbulence and fire-spotting respectively. It should be noted that for
295
brevity, here fire-spotting is assumed to be an independent downwind phe-
296
nomenon; hence the effect of fire-spotting is accounted only for the leeward
297
part of the fire-line. Taking into account these assumptions the PDF for the
CE
random processes can be defined as Z ∞ b U ; t)q(`; t) d` , n b ·n bU ≥ 0 , G(x − x − ` n 0 f (x; t|x) = G(x − x; t), otherwise ,
AC
298
PT
292
16
(11)
ACCEPTED MANUSCRIPT
300
301
b U is the unit vector aligned with the mean wind direction. where n
Turbulent diffusion is assumed to be isotropic and modelled by a bi-
variate Gaussian PDF 1 G(x − x; t) = exp 4πDt
CR IP T
299
(x − x)2 + (y − y)2 4Dt
,
(12)
where D is the turbulent diffusion coefficient such that h(x − x)2 i = h(y − y)2 i =
303
2Dt. In the present model, the whole effect of the turbulent processes over
304
different scales is assumed to be parametrised only by the turbulent diffu-
305
sion coefficient. Also, only the turbulent fluctuations with respect to the
306
mean transport are considered. Hence, the characterisation of turbulence
307
is independent of the mean characteristics of the motion, and in particular
308
the estimation of the turbulent diffusion coefficient D is independent of the
309
mean wind. Moreover, since the simulations are performed with a flat terrain
310
and without boundaries in the spatial domain, horizontal isotropy is also as-
311
sumed. Within this framework, the mean wind does not break the isotropy,
312
because only the fluctuations are considered, but a non-trivial orography can
313
indeed break the isotropy, e.g., valleys or slopes. Even if an exact estima-
314
tion of D is out of the scope of the present study, D is one of the necessary
315
parameters for simulations and a quantitative estimation is required. The
316
desired diffusion coefficient D corresponds to the turbulent heat convection
317
generated by the fire. In literature, the ratio between the total heat trans-
318
fer and the molecular conduction of heat is known as the Nusselt number:
319
Nu = (D + χ)/χ. Where, χ is the thermal diffusivity of the air at ambi-
320
ent temperature and is well-known in literature: χ = 2 × 10−5 m2 s−1 . It is
321
also experimentally known that the relation between the Nusselt number Nu
322
and the Rayleigh number Ra is given by Nu ' 0.1 Ra1/3 [56]. The Rayleigh
AC
CE
PT
ED
M
AN US
302
17
ACCEPTED MANUSCRIPT
number is the ratio between the convection and the conduction of heat and
324
it is defined as Ra = γ ∆T g h3 /(νχ), where γ is the thermal expansion co-
325
efficient, ∆T the temperature difference between the bottom and the top of
326
the convective cell, h the dimension of the convection cell, g the acceleration
327
gravity and ν the kinematic viscosity. From above, D can be computed by
328
the formula:
CR IP T
323
D ' 0.1 χ [γ ∆T g h3 /(νχ)]1/3 − χ,
(13)
where values γ = 3.4 × 10−3 K−1 , g = 9.8 ms−2 and ν = 1.5 × 10−5 m2 s−1
330
are taken from literature. Since the heat transfer is considered in the hori-
331
zontal plane, perpendicular to the vertical “heating wall” embodied by the
332
fire, the length scale of the convective cell is assumed to be h = 100 m with
333
a temperature difference ∆T = 100 K. Considering these values, the scale of
334
the turbulent diffusion coefficient turns out to be approximately 104 times
335
the thermal diffusivity of air at ambient temperature (χ).
M
AN US
329
The jump length of the firebrands varies according to the meteorological
337
conditions, the intensity of the fire and the fuel characteristics [57, 58, 59].
338
A precise formulation of the landing distributions is difficult due to the lim-
339
ited small scale experimental results [60] and the complexities involved in
340
characterising the detailed aspect of firebrand generation and landing. The
341
firebrands land only in the positive direction and their landing distribu-
342
tion increases with distance to reach a maximum before dropping to zero
AC
CE
PT
ED
336
343
[61, 62]. In past, researchers have utilized different statistical distributions
344
to fit/reproduce the distribution pattern of particles deposited on ground
345
from a given source, e.g., the lognormal distribution is used to depict the dis-
346
tribution of the short-range firebrands (up to 4500 m) [27] while other studies 18
ACCEPTED MANUSCRIPT
347
also indicate the use of a Rayleigh distribution [63] or a Weibull distribution
348
[28] to describe the landing distributions.
350
In this study, a lognormal distribution is used to describe the downwind
CR IP T
349
distribution of the firebrands:
1 (ln `/`0 − µ(t))2 q(`; t) = √ exp − , 2 σ(t)2 2π σ(t) `
(14)
where µ(t) = hln `/`0 i and σ(t) = h(ln `/`0 − µ(t))2 i are the mean and the
352
standard deviation of ln `/`0 respectively, and `0 is a unit reference length.
AN US
351
Since the fuel ignition due to hot air and firebrands is not an instanta-
354
neous process, a suitable criterion related to an ignition delay is introduced
355
for application to wildfire modelling. This ignition delay was previously con-
356
sidered as a heating-before-burning mechanism due to the hot air [35, 36]
357
and generalised to include fire-spotting [39]. In particular, since the fuel can
358
burn because of two pathways, i.e., hot-air heating and firebrand landing,
359
the resistance analogy suggests that the resulting ignition delay can be ap-
360
proximately computed as resistances acting in parallel. Hence, representing
361
τh and τf as the ignition delay due to hot air and firebrands respectively, the
362
joint ignition delay τ is
ED
PT
CE
1 1 1 τh + τf = + = . τ τh τf τh τf
(15)
In this study, the delay effect is portrayed as a heating-before-burning
AC
363
M
353
364
mechanism and described by an accumulative process of the effective fire-
365
front over time, i.e., ψ(x, t) =
Z
t
φe (x, η)
0
19
dη , τ
(16)
ACCEPTED MANUSCRIPT
where, ψ(x, 0) = 0 corresponds to the initial unburned fuel. This accumula-
367
tion can be understood as an accumulation of heat that results in an increase
368
of the fuel temperature T (x, t). Finally, the function ψ(x, t) may be related
369
to the fuel temperature T (x, t) as follows: ψ(x, t) ∝
T (x, t) − Ta (x) , Tign − Ta (x)
CR IP T
366
(17)
where Tign is the ignition temperature, T ((x), 0) = Ta (x) and T (x, t) ≤ Tign .
371
Let ∆t be the ignition delay due to the effects of the random processes,
372
then after the elapse of the ignition delay time T (x, ∆t) = Tign , and with
373
the assumptions Ta Tign and the constant of proportionality equal to 1,
374
Eq. (17) reduces to
AN US
370
ψ(x, ∆t) = 1 .
(18)
Hence, as a cumulative effect of hot air and firebrands on the unburned fuel,
376
an ignition point x at time t is created when ψ(x, t) = 1, and the condition
377
φ(x, t) = 1 is also imposed.
378
4. Turbulence and fire-spotting effects into DEVS based wildfire
PT
ED
M
375
simulator
380
ForeFire is a DEVS based wild-land fire model developed bereft of any
CE
379
additional terms or parameters to model turbulence or fire-spotting. As
382
described earlier, DEVS is a Lagrangian method where each marker (active
383
burning point in the fire-line) advances by the quantum distance and the
384
advancement of each marker is managed by decomposition, regeneration or
385
coalescence functions. When one such marker arrives to a new region, the
386
marker decomposes to generate new markers at the interface between the two
AC
381
20
ACCEPTED MANUSCRIPT
regions. If this new region belongs to a fuel break zone, the new markers are
388
assigned null propagation velocity and their further evolution is inhibited.
389
But additional information about turbulence and fire spotting can allow the
390
fire-line to cross the fuel breaks. In this article, the DEVS based ForeFire
391
model is extended to include the effects of turbulence and fire-spotting by
392
post-processing the output obtained from the model at each time step. The
393
present method does not call for any major changes in the original framework
394
of DEVS, and can be versatilely used to reconstruct its output to include
395
these two processes. The detailed steps of the numerical procedure to re-
396
construct the output of DEVS based ForeFire are described as follows:
397
Step 1. Beginning with an initial fire-line, the DEVS front tracking method
398
is used to estimate the propagation of the markers to build up a new
399
fire perimeter for the next time step. This output is modified to
400
include the effects of turbulence and fire spotting by post-processing
401
numerical enrichment. This post-processing step is independent of
402
the definition of ROS.
ED
M
AN US
CR IP T
387
Step 2. The fire perimeter obtained from the active burning points is used to
404
construct the indicator function φ(x, t) (4). The spatial information
405
contained in φ(x, t) is sufficient to modify the fire-line with respect
to turbulence and fire-spotting and serves as the input to the post-
AC
407
CE
406
PT
403
408
processing step. Indicator function φ(x, t) has a value 0 over the region of the available fuel and 1 over the burned area.
409
Step 3. The spatial coordinates of the fire perimeter obtained from DEVS are
410
discrete values spread all over the domain and cannot be represented
411
by a specific mesh/grid. For brevity, the effective indicator function 21
ACCEPTED MANUSCRIPT
φe (x, t) (7) is generated over a Cartesian grid to facilitate the com-
413
putation of the function ψ(x, t) (16) over the same grid. Point in
414
polygon is used to generate the indicator function φ(x, t) from the lo-
415
cus of points representing the fire perimeter in DEVS. It is remarked
416
that φ(x, t), and φe (x, t) can also be computed from the real coor-
417
dinates representing the fire perimeter in DEVS, but the inherent
418
construction of the function ψ(x, t) mandates its computation over
419
a regular grid.
AN US
CR IP T
412
Step 4. The value of the effective indicator φe (x, t) is computed through the
421
numerical integration of the product of the indicator function φ(x, t)
422
and the PDF of fluctuations according to Eq. (7). The effect of
423
turbulence or fire-spotting is included by choosing the corresponding
424
PDF (11).
425
M
420
Step 5. The function ψ(x, t) is updated for each grid point by integration in time with the current value of φe (x, t) (16).
ED
426
Step 6. All points which satisfy the condition ψ(x, t) ≥ 1 are labelled as new
428
ignition points. The post-processing enrichment of the input fire-
429
front is completed at this step and the new markers describing the
430
enriched fire-front are used to define a new front valid at the next time step in DEVS framework. With the inclusion of heating-beforeburning mechanism, the new fire-front perimeter represents a larger
AC
432
CE
431
PT
427
433
area than the input fire perimeter.
434
Step 7. At the next time step, the new markers progress according to the
435
DEVS front tracking method, and the updated perimeter is again
436
subjected to the post-processing to enrich the fire front with the 22
ACCEPTED MANUSCRIPT
random fluctuations pertaining to turbulence and fire-spotting. The
438
sequence is repeated till the final “event time” step or till the fire
439
reaches the end of the domain.
CR IP T
437
Step 8. The DEVS models time as a continuous parameter and updates the
441
system in accordance with the events. The system is updated at user
442
defined “event times”. In between the consecutive events, no change
443
in the system occurs and a list of pending events is maintained in
444
chronological order. The event list is updated at the next “event
445
time”. Since, DEVS provides freedom to select an “event time”, it is
446
chosen to coincide with the maximum allowable time step according
447
to the CFL criterion. It is remarked that this choice of “event time”
448
limits the efficiency of DEVS computation (increase in the compu-
449
tation time due to excessive computations), but is indispensable in
450
order to introduce the heating-before-burning mechanism (16) in the
451
current formulation.
ED
M
AN US
440
It is remarked that the direction of the propagation vector for DEVS
453
is decided by the bisector of the angle formed between a marker and its
454
immediate neighbouring markers on left and right. This direction vector is a
455
weak approximation of the normal and can be very different from the actual
456
normal direction in the case of a non-circular profile.
AC
CE
PT
452
457
5. Simulation set-up
458
To estimate the performance of the DEVS front tracking method with the
459
inclusion of turbulence and fire-spotting, a series of numerical experiments
460
are performed. Similar set of experiments are also carried out with the LSM 23
ACCEPTED MANUSCRIPT
based simulator to enable a fair comparison between the two techniques. A
462
formulation developed in Ref. [6] and Ref. [10] is followed for LSM, while
463
for the DEVS method, ForeFire fire simulator [13] is used. To allow a direct
464
comparison between the two, both models are parametrised in an identical
465
set-up.
CR IP T
461
In the present study, for brevity, no particular type of vegetation is de-
467
fined and simulations are carried out with a pre-defined constant value of
468
ROS. It is assumed that the ROS remains constant for a particular domain.
469
The parametrisation of the numerical simulations is oversimplified and the
470
test cases are chosen with the purpose of highlighting the potential appli-
471
cability of the mathematical formulation. It is remarked that the chosen
472
parametrisations for the numerical simulations do not reproduce any of the
473
forest-fire or prescribed burn experiments but to emphasise the applicability
474
of the formulation, the values of various parameters are chosen to approxi-
475
mately lie in the valid range. The present scope of this work is to provide a
476
first look into the investigation of the performance of LSM and DEVS based
477
fire simulators in the context of inclusion of the effects due to turbulence and
478
fire-spotting.
PT
ED
M
AN US
466
A flat area of hypothetical homogeneous vegetation spread over a domain
480
size of 5000 m × 5000 m is selected for simulations. Different values of the
ROS are utilised for different test cases. The ROS is assumed to be 0.05 ms−1
AC
481
CE
479
482
483
484
b in no wind conditions, while in the presence of along the normal direction n wind, it is estimated by the 3% model [14], i.e.,
b, ROS = 0.03 U · n
(19)
where U is the mean wind velocity. Since, 3% model limits the propagation 24
ACCEPTED MANUSCRIPT
only towards the mean wind direction, a ROS formulation provided by Mallet
486
et al. [6] is followed to prescribe non-zero spread velocities for the flank and
487
rear fires: ROS(U, θ) =
√ εo + a U cosn θ,
CR IP T
485
if |θ| ≤
εo (α + (1 − α)| sin θ|), if |θ| >
π 2
,
π 2
,
(20)
where U 2 = U · U , εo is the flank velocity, (εo α) is the rear velocity with
489
α ∈ [0, 1], and θ is defined as the angle between the normal to the front
490
and the wind direction. In the present set-up, the parameters are defined as
491
α = 0.8, n = 3, a = 0.5 m1/2 s−1/2 , εo = 0.2 ms−1 , and U = 3 ms−1 [6].
AN US
488
In the case of the LSM based simulator, the domain is discretised with
493
a Cartesian grid of 20 m both in x and y directions, while for DEVS the
494
resolution of the simulation is established in the terms of quantum distance
495
∆q and perimeter resolution ∆c [13]. Quantum distance ∆q is defined as
496
the maximum allowable distance to be covered by a particle at each advance,
497
while a measure of ∆c is used to decompose/regenerate/coalesce two particles
498
on propagation. The choice of ∆q and ∆c is dependent on the type of problem,
499
and in the present study, two sets of values are utilized. The simulations
500
are performed with ∆q = 4 m, ∆c = 18 m for zero wind conditions, and
501
∆q = 2 m, ∆c = 8 m in the presence of wind. In DEVS based simulator, the
502
locus of the points constituting the initial front are wind up in the clockwise
503
direction to represent an expanding front (anticlockwise order represents a
504
contracting front). For each numerical simulation presented here, the initial
505
fire-front is represented by 200 markers. A brief study focussing on the
506
numerical behaviour of DEVS and LSM based techniques is also presented
507
in the next section.
AC
CE
PT
ED
M
492
25
ACCEPTED MANUSCRIPT
The mean wind, wherever used, is assumed to be constant in magnitude
509
and direction. The turbulent diffusion coefficient D, and ignition delays cor-
510
responding to hot air and firebrands heating are also assumed to be constant
511
throughout the simulations. In this study, the ignition delay correspond-
512
ing to the hot air τh and to the firebrand τf are assumed to be 10 min and
513
1 min respectively. Physically, a measure of the ignition delay can be inferred
514
from the fire intensity, fuel properties and the length scale at which the fuel
515
is exposed to the heat source [64]. According to Ref. [64], for the same
516
fuel characteristics and low fire intensities the ignition delay due to hot air is
517
higher than the firebrand landing, but in the cases of high fire intensities, the
518
ignition delay at different length scales can be comparable due to the func-
519
tional relationship between the ability of the fuel surface to absorb energy
520
and the ability of the fuel to conduct heat inward. Flammability classifica-
521
tion of different natural fuels available in literature [65, 66, 67] can aid in the
522
selection of an appropriate value for ignition delay.
ED
M
AN US
CR IP T
508
Appropriate values of D are chosen by using formula 13 to assess the
524
effect of increasing level of turbulence and to highlight its role in comparison
525
to the firebrand landing. In this study, three different values of D are utilised:
526
D = 0.30 m2 s−1 , D = 0.15 m2 s−1 , D = 0.075 m2 s−1 . Also, with reference to the lognormal distribution (14), the unit reference
length `0 in feet (1 ft = 0.30 m).
AC
528
CE
527
PT
523
529
Multiple idealised simulation tests are performed both in the presence and
530
the absence of wind by neglecting and considering the effects of the random
531
phenomena. The first case evaluates an isotropic growth of the fire-line in
532
the absence of wind by neglecting all the random processes. A fire-break 26
ACCEPTED MANUSCRIPT
zone 60 m wide is also considered in this case. In the second test, the spread
534
of fire-line for different ROS in different directions is studied. The third test
535
discusses the propagation of the fire-line with wind when the ROS is defined
536
according to the 3% of the mean wind, formula (19), and when it is defined
537
according to formula (20) . The random processes are neglected for the
538
first three test cases. The pure LSM and pure DEVS based wildfire models
539
fail while managing the realistic cases of fire propagation across the fire-
540
break zones, hence the fourth test evaluates the performance of DEVS when
541
turbulent processes are also considered in the model. The effect of turbulent
542
processes is studied both in the presence and the absence of the wind. The
543
fifth test evaluates the performance of the two simulators when fire-spotting
544
is also included along with turbulence. Fire-break zones are also introduced
545
in the last two tests (the fourth and the fifth) to observe the propagation of
546
the fire-line while encountering areas of null fuel. To simplify the simulation
547
and maintain equivalence between the two simulators, the region across and
548
behind the centre of the initial fire-line is demarcated as the leeward side and
549
the windward side respectively.
PT
ED
M
AN US
CR IP T
533
As an extension of the study, a brief analysis of the firebrand landing
551
distribution is also made. With reference to the lognormal distribution (14)
552
for firebrand landing, different values of µ and σ are chosen motivated by a
553
possible correlation with some physical aspects of wildfire like changing fire
AC
CE
550
554
intensity and meteorological conditions. For these simulations, only the LSM
555
based simulator is used and the wind speed is assumed to be 3 ms−1 , while
556
the ROS is given by the Rothermel model [4, 10].
557
All the simulations have been performed at BCAM–Basque Center for Ap27
ACCEPTED MANUSCRIPT
plied Mathematics (www.bcamath.org) in Bilbao, Basque Country – Spain.
559
6. Discussion
560
6.1. Numerical behaviour of DEVS and LSM
CR IP T
558
The stability of DEVS method relies on the choice of quantum distance
562
∆q and the perimeter resolution ∆c but the numerical simulations involving
563
a uniform ROS in all directions (e.g., a constant 0.05 ms−1 irrespective of
564
wind or direction) are unaffected by the choice of these two parameters. In a
565
homogenous domain, a uniform ROS forces all the markers to advance with
566
an identical speed at a same time step and leads to generation of synchronous
567
events. In this case, all the markers are driven by the same CFL condition,
568
and hence, the motion of the markers is analogous to an Eulerian technique.
569
On the contrary, markers with different speed generate different number of
570
discrete events and the resolution of ∆q and ∆c comes into picture. In con-
571
text of the idealised test cases involving wind, an evaluation of performance
572
of the method concerning different values of ∆q and ∆c is demonstrated in
573
Fig. 1. The numerical simulations show the evolution of the fire-line in the
574
presence of an east wind of 3 ms−1 with ROS as 3% of the normal wind, for
575
different pairs of ∆q and ∆c. Each contour represents the fire-front at 134
576
min. The perimeter resolution should be at least twice of the quantum dis-
577
tance (∆c ≥ 2∆q) to ensure that two markers do not cross-over [68]. Taking
AC
CE
PT
ED
M
AN US
561
578
into account this restriction, the upper plot of Fig. 1 highlights the differences
579
that arise in the evolution of the fire-front with respect to different values
580
of ∆c. DEVS employs a neighbourhood check for each marker to identify
581
its two immediate neighbours and utilises these neighbours to regulate the 28
ACCEPTED MANUSCRIPT
shape of the front by merging and generating new markers. A large value of
583
∆c can consider disconnected markers as neighbours and can lead to a unde-
584
sired generation or merging of markers. The loss of information by unsought
585
merging can cause the front to collapse. On the other hand, decreasing
586
the perimeter resolution ∆c without increasing the resolution of the marker
587
jumps causes the regeneration function to be evoked frequently and leads to
588
creation of a large number of markers with small separations. Over time,
589
too many markers at a coarse spatial scale leads to excessive unnecessary
590
computations and small numerical errors introduced by the computation can
591
escalate quickly over time. A finer scale of ∆q is desired with small values
592
of ∆c for a complete optimisation of the numerical simulation. The lower
593
plot in Fig. 1 shows the variation of the fire-front with respect to increasing
594
resolution of ∆q. Though a finer resolution of ∆q resolves the fire map at a
595
better resolution, it has no significant impact on the simulation results; on
596
the other hand, it increases the computation cost of the simulation. In the
597
present scenario, the homogeneous vegetation and topography of the present
598
domain causes the Lagrangian movement of markers to be driven predom-
599
inantly by regeneration and the fine scale information about the domain is
600
redundant. But on the other hand, in an inhomogeneous domain (involving
601
diverse vegetation, irregular topography, etc), where both decomposition and
602
coalescence play an important role, a finer marker jump length is necessary
CE
PT
ED
M
AN US
CR IP T
582
to account for the finer details of the domain. Usually, the value of ∆q is
604
taken to be around 2-3 m to ensure a correct representation of the vari-
605
able vegetation, changing topography and fire break zones like roads, rivers
606
etc. Concerning a pure quantitative analysis of the total error in DEVS, the
AC 603
29
ACCEPTED MANUSCRIPT
607
interested readers are referred to [68] and [69]. The numerical behaviour of the LSM is dependent on the time step (dt),
609
the grid resolution (dx, dy) and the maximum absolute speed of all points
610
on the grid. The iterative time step involved while solving the differential
611
equations needs to satisfy the CFL criteria for a converging solution. The
612
CFL condition defines the maximum limit of the time step in order to limit
613
the evolution of the contour to at most one grid cell at each time step. This
614
upper limit on the time step is important for the stability of the system and
615
a converging solution. In one-dimension, the CFL condition can be described
616
as:
AN US
CR IP T
608
dt ≤ k
min(dx, dy) vmax
(21)
where, vmax is the maximum absolute speed, and 0 < k < 1 is a constant. To
618
demonstrate the numerical behaviour of LSM, similar numerical simulations
619
with 3% model and an east wind of 3 ms−1 are repeated for different values
620
of dt and dx . The upper figure in Fig. 2 shows the fire-line contours at 380
621
min for various values of dt, and when dx = dy = 20 m. With values of dt
622
smaller than the threshold defined by CFL condition, the fire-line evolves
623
without any instabilities. Similar results are achieved with further decrease
624
in the time step but with an increase in the computation time. On the
625
other hand, when for the same speed function, dt is chosen to violate the
626
CFL criteria, the fire-line evolves with an undesirable behaviour and this
AC
CE
PT
ED
M
617
627
obscurity magnifies over time. In the same way, CFL condition can also
628
be interpreted as the constraint on the grid resolution. The lower figure in
629
Fig. 2 shows the evolution of the fire-line with respect to different values of
630
dx. The increase in the grid resolution over an identical time step introduces 30
ACCEPTED MANUSCRIPT
numerical instability and causes the fire-lines to evolve erroneously. For the
632
same velocity field, shorter time steps are required over small grid spacings
633
in order to respect the CFL condition, though it leads to an increase in the
634
computation costs. Hence, accounting for the numerical behaviour of the
635
LSM and to avoid the excessive computation, in the present study, dt is
636
defined according to (21), with k = 0.5.
CR IP T
631
It is remarked here that since the mathematical formulation for turbulence
638
and fire-spotting is implemented in DEVS and LSM based wildfire simulators
639
as a post-processing scheme, it preserves the original framework of the wildfire
640
simulators and does not affect the numerical stability and convergence of both
641
these methods.
642
6.2. Comparison of DEVS and LSM
M
AN US
637
In the first three test scenarios, the evolution of the fire-line without the
644
impact of random processes in both LSM and DEVS approaches is analysed.
645
Fig. 3 presents the propagation of the fire-line in no wind conditions for both
646
the simulators. The advancement of the initial circular fire-line follows an
647
isotropic growth and with identical initial conditions the two different moving
648
interface schemes provide a similar evolution of the fire-line. But the fire-line
649
fails to propagate across the fire-break zone in both the techniques. The
650
evolution is shown only up to 140 min, but an extended run up to 250 min
651
(not shown) indicates the limitation of the two approaches to reproduce the
AC
CE
PT
ED
643
652
653
fire jump across the no fuel zone. Fig. 4 presents the growth of an initial spot fire but with a non-homogeneous
654
ROS in absence of any wind. The different values of the ROS can be at-
655
tributed to different fuel types, but here no effort has been made to describe 31
ACCEPTED MANUSCRIPT
the fuel type. The fire-line propagates with different speed in the three dif-
657
ferent directions, but the isotropy is preserved. The results from these two
658
test cases indicate that in the absence of wind, the two advancing schemes:
659
LSM and DEVS, show an identical behaviour in simulating situations with
660
constant ROS and hence, this observation provides a crucial background and
661
scope for the comparison of situations with increasing variability and com-
662
plexity.
CR IP T
656
Fig. 5 shows the evolution of the fire-front of a circular initial profile with
664
radius 300 m in case of a weak wind of 3 ms−1 directed in the east direction.
665
The isochronous fronts are plotted at every 20 min and follow an oval shape
666
for both the simulators. DEVS fire contours diverge slightly from the mean
667
wind direction and develop into an increasing flanking fire over time. This
668
divergence in the evolution of fire-front occurs due to differences in the com-
669
putation of the normal between the two approaches. This fact can be very
670
well appreciated when an initial square profile is considered. Fig. 6 shows the
671
progress of the fire when the initial fire perimeter is considered to be a square
672
with side 600 m. Under the effect of a constant zonal wind and restricted 3%
673
model (19), the evolution in LSM strictly follows the initial square shape,
674
but in DEVS the active burning markers at the corners advance spuriously
675
to provide an additional flanking spread. The DEVS markers at the corners
676
are affected by the approximate normal and this leads to a larger flank fire.
AC
CE
PT
ED
M
AN US
663
677
Whereas, the normal direction for markers pointing towards the direction of
678
wind is identical for both the methods, hence the head fires arrive to same
679
location in both the simulations. The flank fires generated by the approxi-
680
mate construction of the normal seem to be in a more qualitative agreement 32
ACCEPTED MANUSCRIPT
681
with a realistic case, but within the stated analytical formulation of the ROS
682
(19), the lateral progress of the flank fires is an anomalous behaviour. The ROS formulation developed by Mallet et al. [6], here formula (20),
684
provides a simple way to study the evolution of flank fire by introducing
685
different ROS for the head, flank and rear directions. In the original paper,
686
θ is defined as the angle between the normal to the front and the wind
687
direction. Since, the normal computation in DEVS approach is approximate,
688
two separate tests are performed to evaluate the effect of such approximation
689
on the spread: firstly when θ is computed according to the original definition
690
[6], and secondly when θ is assumed to be the angle between the line joining
691
the front and the wind direction to enforce an identical definition of θ in
692
both the methods. Fig. 7 shows the evolution of fire when θ is computed
693
according to the original definition; the advancement of the head and rear
694
fires are identical, but the flank fire spread is spuriously larger in DEVS
695
because of the approximate computation of the front normal direction. On
696
the other hand, it is evident from Fig. 8 that in the second case, when the
697
computation of θ is not affected by the approximate computation of the
698
front normal in DEVS, the fire-line advancement is identical to LSM in all
699
directions.
AN US
M
ED
PT
As shown in Fig. 3, in case of fire-break zones, pure LSM and DEVS
based simulators are inherently unable to simulate the realistic situations of
AC
701
CE
700
CR IP T
683
702
fire overcoming a fire break. But Fig. 9 shows that with the introduction of
703
turbulence, the simulators can model the effect of hot air to overcome fire-
704
break zones. Here the value of turbulent diffusion coefficient D is assumed
705
to be 0.15 m2 s−1 . The evolution of the fire-line is almost similar for both 33
ACCEPTED MANUSCRIPT
the simulators, though a slight underestimation can be visibly observed in
707
the DEVS approach. A cross section of φ(x, t) and φe (x, t) at y = 2500 m
708
(Fig. 10) provides a more detailed illustration of the differences in the evolu-
709
tion of the fire-line between the two approaches. The introduction of marginal
710
differences with the inclusion of turbulence can be explained by the contrast
711
in the construction of the indicator function between the two simulators. In
712
LSM, the indicator function is defined over a Cartesian grid by definition,
713
while on the other hand, in DEVS, a gridded indicator function is constructed
714
from the original real coordinates. The point in polygon technique is used to
715
transform the real coordinates to a Cartesian grid. This technique provides
716
only an approximate measure of the same; small errors can be introduced
717
by the improper classification of the points lying very close to the boundary
718
[70], and can lead to a slight underestimation/overestimation of the original
719
fire-line.
M
AN US
CR IP T
706
Fig. 11 shows the results for the two simulators when the turbulent dif-
721
fusion coefficient is D = 0.30 m2 s−1 . Stronger turbulence causes a rapid
722
propagation of the fire-line and an earlier ignition across the fire-break zone.
723
A detailed analysis of the effect of varying turbulence over long-term propa-
724
gation with the LSM can be found in [35, 36].
PT
Fig. 12 presents the effect of inclusion of turbulence with a non-zero wind
profile and 3% model (19). The effect of turbulence is most pronounced in
AC
726
CE
725
ED
720
727
the direction of the wind and it can be clearly observed that with turbulence,
728
the spread over the head, flank and back fire-lines is faster and the head fire
729
is also able to overcome the fuel break. The quantum of increase in the
730
fire spread can also be appreciated in comparison to Fig. 5, which presents 34
ACCEPTED MANUSCRIPT
the ideal case with the same initial conditions. Both simulators show almost
732
similar characteristics in the spread of fire, though the flank fire has a slightly
733
larger spread in DEVS.
CR IP T
731
Another aspect contributing towards the increase in the fire spread due
735
to new fire ignitions generated by the firebrands is presented in Fig. 13. With
736
the inclusion of fire-spotting along with turbulence, the evolution of the fire-
737
front is faster in comparison to the effect of turbulence alone (Fig. 12). The
738
flank fire and the head fire are also well simulated in both the approaches,
739
and again a larger spread out the flanking fires is observed in the DEVS
740
scheme. The proposed mathematical formulation of spot-fires along with
741
turbulence is effective in mimicking the rapid advancement of the fire-line
742
due to the generation of new ignition points by the firebrands. It is also
743
effective in simulating the new ignition points developed across the no fuel
744
zones. Within the current parametrisation, the firebrand landing distances
745
lie close to the main fire and are not long enough to develop a new secondary
746
fire, but in the next section, a sensitivity study of the spatial pattern of the
747
firebrand landing distances is also presented.
PT
ED
M
AN US
734
Fig. 14 shows the time evolution of the fire-front in presence of a mean
749
wind velocity of 3 ms−1 and two fire-break zones. The results are obtained
750
when only the effects of turbulence are included in the formulation and the
751
diffusion coefficient D = 0.075 m2 s−1 . Similarly, Fig. 15 shows the time evo-
AC
CE
748
752
lution of fire-front but when both turbulence and fire-spotting are included
753
in the formulation. Both approaches follow an identical evolution pattern as
754
observed in other cases. As expected, the spread of the fire-front is faster and
755
is able to overcome both the fire-break zones when turbulence is included; 35
ACCEPTED MANUSCRIPT
but inclusion of fire-spotting along with turbulence causes the fire to propa-
757
gate at a much faster rate and results in a larger spread. Again, the major
758
differences which arise in the evolution of flank fires between the two meth-
759
ods are due to the dissimilarity in the construction of the direction of the
760
propagation of the active burning points.
761
7. Insight on firebrand distribution
CR IP T
756
In the simulations discussed previously, the phenomenon of fire-spotting
763
contributes towards an increase in the burning area and a faster propagation
764
of fire, but the landing distances are not sufficiently long enough to develop
765
new secondary fires detached from the primary fire.
AN US
762
Generally, the firebrands generated in low intensity fires fall close to the
767
primary fire and although they contribute towards a faster fire spread, the
768
separation from the primary fire is not large enough for the new ignitions
769
to develop a new secondary fire. But high intensity fires and extreme mete-
770
orological conditions such as strong wind/extremely high temperatures/low
771
humidity can facilitate the transport of embers over longer distances and the
772
area under firebrand attack can be sufficiently far away from the primary
773
fire [71, 72, 73, 74, 75]. In Ref. [76] and references therein, a brief account
774
of the impact of various historical forest fires is provided emphasising the
775
importance of strong winds, and extreme meteorological conditions on the
AC
CE
PT
ED
M
766
776
fire-spotting phenomenon. They describe different historical fires where low
777
intensity fires coupled with strong winds and low relative humidity caused
778
a large scale fire-spotting, while on the other hand, fires with high intensi-
779
ties in calm meteorological conditions had a negligible contribution from the 36
ACCEPTED MANUSCRIPT
780
firebrands. In lognormal distribution (14), a suitable choice of the parameters µ and
782
σ can be utilized to modify the firebrand landing distributions to represent
783
different fire intensities and wind conditions. In this study, the main interest
784
lies in the mathematical modelling of the firebrands landing distributions,
785
hence constant meteorological conditions are assumed for all the simulations
786
presented here. It is also assumed that each active ignition point is capable
787
of generating firebrands and all embers are controlled by identical transport
788
forces.
AN US
CR IP T
781
Fig. 16 and Fig. 17 show the variation of the area under the firebrand
790
attack with different values of µ and σ. A constant σ and an increasing µ
791
points towards a decrease in the effectiveness of the firebrands, and slows
792
down the fire propagation. A slight deviation from the pattern is observed
793
for µ = 15 and σ = 5, when the firebrands travel far enough to create a new
794
secondary fire.
ED
M
789
For an increasing value of σ (Fig. 17a-c) with constant µ, the fire evolu-
796
tion increases as larger number of firebrands land away from the source, but
797
the landing distances are still not long enough to delineate the new ignitions
798
from the main fire. With further increase in σ, a secondary fire emerges away
799
from the main source. Equal contribution to the firebrands by each flaming
800
point leads to a regular distribution pattern of the secondary fire. As time
AC
CE
PT
795
801
progresses, the primary fire also advances and covers up the gap with the
802
secondary fire. With further increase in σ, the probability of the firebrands
803
landing at farther distances increases, and this can lead to a potentially ex-
804
plosive firebrands attack. Fig. 17f shows the effect of fire-spotting at shorter 37
ACCEPTED MANUSCRIPT
scales contributing to a faster rate of spread, but the development of new
806
secondary fires is not observed. For this particular value of σ, the statistical
807
distribution of the long range landing distances lies outside the domain and
808
is impossible to resolve them in the limited time and spatial scales considered
809
in this simulation.
CR IP T
805
As a consequence of the idealised set-up for simulations, the new sec-
811
ondary fires evolve with the same intensity as the primary fire and also serve
812
as a source of generation of newer fire brands and henceforth the tertiary
813
fires. Fig. 18a-f demonstrates the generation of new secondary fires by the
814
primary fire, and tertiary fires by the secondary fires. As time progresses,
815
the each of the fire-line evolves independently along with a small contribution
816
from its parent fire-line, and eventually they merge among each other.
AN US
810
In a lognormal distribution, an increasing value of the µ causes a right
818
shift of the maximum, indicating towards a larger jump length, but the sim-
819
ulations suggest a decreasing trend in the landing distance. This can be
820
explained by the fact that the ignition of the unburned fuel is not an instan-
821
taneous process, but rather an accumulative process. The fuel starts burning
822
only when the temperature is high enough to allow spontaneous combustion.
823
Though a larger µ allows the firebrands to fall at larger distances but it limits
824
the frequency of landing. Insufficient number of firebrands landing at farther
825
distances diminishes the efficacy of their attack to cause combustion; hence
AC
CE
PT
ED
M
817
826
leads to a smaller contribution from fire-spotting. For an effective combus-
827
tion, sufficient number of firebrands should be available for an adequate heat
828
exchange with the unburned fuel. On the contrary, an increasing σ causes the
829
lognormal distribution to have a lower maximum and an elongated right tail 38
ACCEPTED MANUSCRIPT
of the distribution, but with a high probability of landing far away from the
831
main source. The increase in the number of firebrands landing away from
832
the maximum leads to a rapid ignition of the unburned fuel and a faster
833
propagation of the fire.
CR IP T
830
According to the increasing and decreasing trend of the landing distance
835
with respect to increasing σ and decreasing µ respectively, within the log-
836
normal characterisation assumed here for the firebrand landing distribution,
837
the reference length scale L of the landing distance can be stated as:
AN US
834
L∝ 838
σα , µβ
(22)
where α > 0 and β > 0 are exponents to be established.
From relation (22) it can be deduced that neither the maximum nor the
840
mean of the lognormal distribution can be considered as an estimation of
841
the landing distance length scale L. Physically, the landing distance can be
842
described in terms of the dimensions of the convective column, wind condi-
843
tions, fire intensity and fuel characteristics. Further analysis towards a phys-
844
ical parametrisation of fire-spotting is postponed for a future paper, where
845
firebrand distributions different from the lognormal will be also considered.
846
8. Conclusions
ED
PT
CE
847
M
839
This article describes a mathematical formulation to model the effects
of turbulence and fire-spotting in wild-land fire propagation models. Previ-
849
ously, this formulation has been applied to an Eulerian LSM based wildfire
850
propagation method, and this article describes the implementation of the
851
mathematical formulation into a Lagrangian DEVS based wildfire propaga-
852
tion method. This formulation splits the motion of the front into a drifting
AC
848
39
ACCEPTED MANUSCRIPT
part and a fluctuating part. The drifting part is independent of the fluctuat-
854
ing part and is given by the fire perimeter determined from one of the various
855
wildfire propagation methods that exist in literature (DEVS and LSM in the
856
present study); while the fluctuating part is generated by a comprehensive
857
statistical description of the system and includes the effects of random pro-
858
cesses in agreement with the physical properties of the process. This study
859
also provides a unique opportunity for the comparison of two simple wildfire
860
simulators based on the Eulerian LSM and on the Lagrangian DEVS scheme.
861
Numerical simulations based on identical set-ups are used to compare the
862
performance of the two approaches.
AN US
CR IP T
853
A series of idealised numerical simulations show that the proposed ap-
864
proach for random effects emerges to be suitable for both LSM and DEVS
865
based simulators to manage the real world situations related to random char-
866
acter of fire, e.g., increase in fire spread due to pre-heating of the fuel by hot
867
air, vertical lofting and transporting of firebrands, fire overcoming no fuel
868
zones, and development of new isolated secondary fires. The simulations
869
from both LSM and DEVS techniques follow nearly identical patterns, and
870
the only difference which distinguishes the outputs from the two techniques
871
is due to the distinction in the treatment of the direction of propagation. In
872
DEVS, the active burning points propagate in a direction given by the bisec-
873
tor of the angle made by each marker with its immediate left and right mark-
AC
CE
PT
ED
M
863
874
ers, in contrast to the normal direction of propagation used in LSM. Such
875
differences result in “spurious” flanking fires in DEVS, which are anomalous
876
with respect to the definition of the ROS, but qualitatively provide a more
877
“realistic” representation of the fire contour. Due to this fact, the two meth40
ACCEPTED MANUSCRIPT
ods can be considered complementary to each other for simple situations,
879
but with increasing complexity and the introduction of random processes a
880
validation exercise is necessary to single out the best representation of the
881
fire propagation.
CR IP T
878
A brief study on the sensitivity of the firebrand landing distribution to
883
the changes in the shape of the lognormal distribution indicates that the
884
shape parameters (µ and σ) can control the spread of the firebrand landing
885
and the lognormal distribution is capable of simulating new secondary fires
886
separated from the main fire perimeter. In future, the formulation for the
887
firebrand landing distance would be modified to include information about
888
the wind speed, fire intensity, radius of the firebrands and fuel characteristics.
889
Though without a validation test case, these simulations do provide an
890
insight into the flexibility of this formulation to incorporate effects of random
891
processes such as the turbulent heat diffusion and fire-spotting, specially
892
while characterising the landing distance of firebrands. It can be deduced
893
that the mathematical formulation is independent of the choice of the method
894
used to ascertain the drifting part and can serve as a versatile addition to
895
the existing fire spread simulators to include the effects of random spread.
896
The results from this study provide a support to the implementation of the
897
proposed approach for the effects of random phenomena into operational
898
softwares such as WRF-SFIRE and ForeFire. The source code utilized for all
AC
CE
PT
ED
M
AN US
882
899
the simulations in this article is hosted at https://github.com/ikaur17/
900
firefronts.
41
ACCEPTED MANUSCRIPT
901
Acknowledgement This research is supported by MINECO under Grant MTM2013-40824-P,
903
by Bizkaia Talent and European Commission through COFUND programme
904
under Grant AYD-000-226, and also by the Basque Government through the
905
BERC 2014-2017 program and by the Spanish Ministry of Economy and
906
Competitiveness MINECO: BCAM Severo Ochoa accreditation SEV-2013-
907
0323. The authors would also like to thank the two anonymous reviewers for
908
their insightful comments which helped in improving the manuscript both in
909
presentation and scientific content.
910
References
AN US
CR IP T
902
[1] A. L. Sullivan, Wildland surface fire spread modelling, 1990–2007. 1:
912
Physical and quasi-physical models, Int. J. Wildland Fire 18 (2009) 349–
913
368.
ED
M
911
[2] A. L. Sullivan, Wildland surface fire spread modelling, 1990–2007. 2:
915
Empirical and quasi-empirical models, Int. J. Wildland Fire 18 (2009)
916
369–386.
918
Simulation and mathematical analogue models, Int. J. Wildland Fire 18 (2009) 387–403.
AC
919
[3] A. L. Sullivan, Wildland surface fire spread modelling, 1990–2007. 3:
CE
917
PT
914
920
[4] R. C. Rothermel, A mathematical model for predicting fire spread in
921
wildland fires, Research Paper INT-115, USDA Forest Service, Inter-
922
mountain Forest and Range Experiment Station, Ogden, Utah 84401,
923
available at: http://www.treesearch.fs.fed.us/pubs/32533 (1972). 42
ACCEPTED MANUSCRIPT
924
925
[5] J. H. Balbi, F. Morandini, X. Silvani, J. B. Filippi, F. Rinieri, A physical model for wildland fires, Combust. Flame 156 (2009) 2217–2230. [6] V. Mallet, D. E. Keyes, F. E. Fendell, Modeling wildland fire propaga-
927
tion with level set methods, Comput. Math. Appl. 57 (2009) 1089–1101.
928
[7] M. Finney, Fire growth using minimum travel time methods, Can. J.
930
931
932
933
For. Res. 32 (2002) 420–1424.
[8] M. Finney, Calculation of fire spread rates across random landscapes,
AN US
929
CR IP T
926
Int. J. Wildland Fire 12 (2003) 167–174.
[9] J. A. Sethian, P. Smereka, Level set methods for fluid interfaces, Ann. Rev. Fluid Mech. 35 (2003) 341–372.
[10] J. Mandel, J. D. Beezley, A. K. Kochanski, Coupled atmosphere-
935
wildland fire modeling with WRF 3.3 and SFIRE 2011, Geosci. Model.
936
Dev. 4 (2011) 591–610.
ED
M
934
[11] H. Karimabadi, J. Driscoll, Y. A. Omelchenko, N. Omidi, A new asyn-
938
chronous methodology for modeling of physical systems: breaking the
939
curse of Courant condition, J. Comp. Phys. 205 (2005) 755–775.
CE
PT
937
[12] Y. A. Omelchenko, H. Karimabadi, HYPERS: A unidimensional asyn-
941
chronous framework for multiscale hybrid simulations, J. Comp. Phys.
AC
940
942
231 (2012) 1766–1780.
943
[13] J. B. Filippi, F. Morandini, J. H. Balbi, D. Hill, Discrete event front
944
tracking simulator of a physical fire spread model, Simulation 86 (2009)
945
629–646. 43
ACCEPTED MANUSCRIPT
[14] J. B. Filippi, V. Mallet, B. Nader, Evaluation of forest fire models on
947
a large observation database, Nat. Hazards Earth Syst. Sci. 14 (2014)
948
3077–3091.
950
951
952
[15] D. Anderson, E. Catchpole, N. de Mestre, T. Parkes, Modelling the spread of grass fires., J. Aus. Math. Soc. 23 (1982) 451–466.
[16] G. D. Richards, An elliptical growth model of forest fire fronts and its numerical solution, Int. J. Numer. Meth. Eng. 30 (1990) 1163–1179.
AN US
949
CR IP T
946
953
[17] G. D. Richards, The properties of elliptical wildfire growth for time
954
dependent fuel and meteorological conditions, Combust. Sci. Technol.
955
95 (1993) 357–383.
958
959
M
957
[18] G. D. Richards, A general mathematical framework for modelling twodimensional wildland fire spread, Int. J. Wildland Fire 5 (1995) 63–72. [19] J. Glasa, L. Halada, On elliptical model for forest fire spread modeling
ED
956
and simulation, Math. Comput. Simulat. 78 (2008) 76–88. [20] O. Rios, W. Jahn, G. Rein, Forecasting wind-driven wildfires using an
961
inverse modelling approach, Nat. Hazards Earth Syst. Sci. 14 (2014)
962
1491–1503.
CE
[21] T. L. Clark, M. A. Jenkins, J. Coen, D. Packham, A coupled
AC
963
PT
960
964
atmospheric-fire model: convective feedback on fire-line dynamics, J.
965
Appl. Meteor. 35 (1996) 875–901.
966
967
[22] B. E. Potter, A dynamics based view of atmosphere-fire interactions, Int. J. Wildland Fire 11 (2002) 247–255. 44
ACCEPTED MANUSCRIPT
[23] B. E. Potter, Atmospheric interactions with wildland fire behaviour - I.
969
Basic surface interactions, vertical profiles and synoptic structures, Int.
970
J. Wildland Fire 21 (2012) 779–801.
971
972
CR IP T
968
[24] B. E. Potter, Atmospheric interactions with wildland fire behaviour - II. Plume and vortex dynamics, Int. J. Wildland Fire 21 (2012) 802–817.
[25] C. B. Clements, S. Zhong, X. Bian, W. E. Heilman, D. W. Byun, First
974
observations of turbulence generated by grass fires, J. Geophys. Res. 113
975
(2008) D22102.
AN US
973
976
[26] N. Sardoy, J. L. Consalvi, B. Porterie, A. C. Fernandez-Pello, Modeling
977
transport and combustion of firebrands from burning trees, Combust.
978
Flame 150 (2007) 151–169.
[27] N. Sardoy, J. L. Consalvi, A. Kaiss, A. C. Fernandez-Pello, B. Porterie,
980
Numerical study of ground-level distribution of firebrands generated by
981
line fires, Combust. Flame 154 (2008) 478–488.
ED
M
979
[28] S. Kortas, P. Mindykowski, J. L. Consalvi, H. Mhiri, B. Porterie, Exper-
983
imental validation of a numerical model for the transport of firebrands,
984
Fire Safety J. 44 (2009) 1095–1102.
CE
PT
982
[29] H. A. Perryman, A mathematical model of spot fires and their manage-
986
ment implications, Master’s thesis, Humboldt State University, Arcata,
AC
985
987
CA (2009).
988
[30] H. A. Perryman, C. J. Dugaw, J. M. Varner, D. L. Johnson, A cellular
989
automata model to link surface fires to firebrand lift-off and dispersal,
990
Int. J. Wildland Fire 22 (2013) 428–439. 45
ACCEPTED MANUSCRIPT
992
993
994
[31] C. Favier, Percolation model of fire dynamic, Phys. Lett. A 330 (2004) 396–401. [32] H. Hunt, A new conceptual model for forest fires based on percolation
CR IP T
991
theory, Complexity 13 (2007) 12–17.
[33] D. Boychuk, W. J. Braun, R. J. Kulperger, Z. L. Krougly, D. A. Stan-
996
ford, A stochastic forest fire growth model, Environ. Ecol. Stat. 16 (2009)
997
133–151.
998
999
1000
AN US
995
[34] L. Ntaimo, B. P. Zeigler, M. J. Vasconcelos, B. Khargharia, Forest fire spread and supression in DEVS, Simulation 80 (2004) 479–500. [35] G.
Pagnini,
fire
Massidda,
method,
1003
dinia,
1004
http://publications.crs4.it/pubdocs/2012/PM12a/pagnini massidda-
1005
levelset.pdf and arXiv:1408.6129 (July 2012).
revised
version:
randomized
CRS4,
August
Pula 2014,
level-set
(CA),
Sar-
available
at:
PT
ED
Italy,
2012/PM12a,
the
effects
1002
Rep.
by
turbulence
in
Tech.
propagation
Modelling
1001
M
wildland
L.
[36] G. Pagnini, L. Massidda, The randomized level-set method to model
1007
turbulence effects in wildland fire propagation, in: D. Spano, V. Bacciu,
1008
M. Salis, C. Sirca (Eds.), Modelling Fire Behaviour and Risk. Proceedings of the International Conference on Fire Behaviour and Risk. ICFBR
AC
1009
CE
1006
1010
2011, Alghero, Italy, October 4–6 2011, 2012, pp. 126–131, ISBN 978-
1011
88-904409-7-7.
1012
[37] G. Pagnini, A model of wildland fire propagation including random ef-
1013
fects by turbulence and fire spotting, in: Proceedings of XXIII Congreso 46
ACCEPTED MANUSCRIPT
1014
de Ecuaciones Diferenciales y Aplicaciones XIII Congreso de Matem´atica
1015
Aplicada. Castell´o, Spain, 9–13 September 2013, 2013, pp. 395–403. [38] G. Pagnini, Fire spotting effects in wildland fire propagation, in:
1017
F. Casas, V. Mart´ınez (Eds.), Advances in Differential Equations and
1018
Applications, Vol. 4 of SEMA SIMAI Springer Series, Springer Interna-
1019
tional Publishing Switzerland, 2014, pp. 203–214, dOI 10.1007/978-3-
1020
319-06953-1 20. ISBN: 978-3-319-06952-4 (eBook: 978-3-319-06953-1).
AN US
CR IP T
1016
[39] G. Pagnini, A. Mentrelli, Modelling wildland fire propagation by track-
1022
ing random fronts, Nat. Hazards Earth Syst. Sci. 14 (2014) 2249–2263.
1023
[40] G. Pagnini, A. Mentrelli, The randomized level set method and an as-
1024
sociated reaction-diffusion equation to model wildland fire propagation,
1025
in: Proceedings of The 18th European Conference on Mathematics for
1026
Industry, ECMI 2014, Taormina, Italy, June 9–13, 2014, ECMI book
1027
subseries of Mathematics in Industry, Springer Heidelberg, accepted.
ED
M
1021
[41] I. Kaur, A. Mentrelli, F. Bosseur, J. B. Filippi, G. Pagnini, Wildland
1029
fire propagation modelling: A novel approach reconciling models based
1030
on moving interface methods and on reaction-diffusion equations, in:
1031
ˇıstek, T. Vejchodsk´ J. Brandts, S. Korotov, M. Kˇr´ıˇzek, K. Segeth, J. S´ y (Eds.), Proceedings of the International Conference on Applications of Mathematics 2015; Prague, Czech Republic, 18–21 November (2015),
AC
1033
CE
1032
PT
1028
1034
Institute of Mathematics – Academy of Sciences, Czech Academy of
1035
Sciences, Prague, Czech Republic, 2015, pp. 85–99.
1036
[42] J. B. Filippi, F. Bosseur, D. Grandi, ForeFire: open-source code for 47
ACCEPTED MANUSCRIPT
wildland fire spread models, in: D. X. Viegas (Ed.), Advances in forest
1038
fire research, Imprensa da Universidade de Coimbra, 2014, pp. 275–282.
1039
URL https://digitalis.uc.pt/handle/10316.2/34013
1040
1041
CR IP T
1037
[43] B. P. Zeigler, T. G. Kim, H. Praehofer, Theory of Modeling and Simulation, 2nd Edition, Academic Press, Inc., Orlando, FL, USA, 2000.
[44] J. C. Carvalho, M. T. Vilhena, D. M. Moreira, Comparison between Eu-
1043
lerian and Lagrangian semi-analytical models to simulate the pollutant
1044
dispersion in the PBL, Appl. Math. Model. 31 (2007) 120–129.
AN US
1042
[45] Z. Zhang, Q. Chen, Comparison of the Eulerian and Lagrangian methods
1046
for predicting particle transport in enclosed spaces, Atmos. Environ. 41
1047
(2007) 5236–5248.
M
1045
[46] M. Saidi, M. Rismanian, M. Monjezi, M. Zendehbad, S. Fatehiboroujeni,
1049
Comparison between Lagrangian and Eulerian approaches in predicting
1050
motion of micron-sized particles in laminar flows, Atmos. Environ. 89
1051
(2014) 199–206.
1054
object-oriented environment, Simulation 49 (1987) 219–230. [48] S. Osher, R. Fedkiw, Level Set Methods and Dynamic Implicit Surfaces, Springer–Verlag, New York, Berlin, Heidelberg, 2003.
AC
1055
PT
1053
[47] B. P. Zeigler, Hierarchical, modular discrete-event modelling in an
CE
1052
ED
1048
1056
1057
[49] R. G. Rehm, R. J. McDermott, Fire-front propagation using the level set method, Tech. Note 1611, Natl. Inst. Stand. Technol. (April 2009).
48
ACCEPTED MANUSCRIPT
1058
[50] E. Jones, T. Oliphant, P. Peterson, et al., SciPy: Open source scientific tools for Python (2001–).
1060
URL http://www.scipy.org/
1061
1062
1063
CR IP T
1059
[51] J. D. Hunter, Matplotlib: A 2D graphics environment, Comput. Sci. Eng. 9 (3) (2007) 90–95.
[52] F. P´erez, B. E. Granger, IPython: a system for interactive scientific computing, Comput. Sci. Eng. 9 (3) (2007) 21–29.
1065
URL http://ipython.org
1066
1067
AN US
1064
[53] A. Mentrelli, G. Pagnini, Random front propagation in fractional diffusive systems, Commun. Appl. Ind. Math. 6 (2) (2015) e–504. [54] A. Mentrelli, G. Pagnini, Front propagation in anomalous diffusive me-
1069
dia governed by time-fractional diffusion, J. Comp. Phys. 293 (2015)
1070
427–441.
1073
1074
ED
[56] J. J. Niemela, K. R. Sreenivasan, Turbulent convection at high Rayleigh numbers and aspect ratio 4, J. Fluid Mech. 557 (2006) 411–422.
[57] F. A. Albini, Potential spotting distance from wind-driven surface fires,
AC
1075
combustion, Phys. Rev. Lett. 107 (2011) 044503.
PT
1072
[55] G. Pagnini, E. Bonomi, Lagrangian formulation of turbulent premixed
CE
1071
M
1068
1076
Research Paper INT-309, U.S. Department of Agriculture Forest Service
1077
Intermountain Forest and Range Experiment Station (1983).
1078
[58] G. A. Morris, A simple method for computing spotting distances from
1079
wind-driven surface fires, Research Note INT-374, U.S. Department of 49
ACCEPTED MANUSCRIPT
1080
Agriculture Forest Service Intermountain Forest and Range Experiment
1081
Station (1987). [59] W. Catchpole, Fire properties and burn patterns in heterogeneous
1083
landscapes, in: R. A. Bradstock, J. E. Williams, M. A. Gill (Eds.),
1084
Flammable Australia: The Fire Regimes and Biodiversity of a Conti-
1085
nent, Cambridge University Press, 2002, Ch. 3, pp. 49–75.
CR IP T
1082
[60] M. E. Houssami, E. Mueller, J. C. Thomas, A. Simeoni, A. Filkov,
1087
N. Skowronski, M. R. Gallagher, K. Clark, R. Kremens, Experimental
1088
procedures characterising firebrand generation in wildfires, Fire Tech-
1089
nol.(Available online). doi:10.1007/s10694-015-0492-z.
1091
[61] K. D. Hage, On the dispersion of large particles from a 15-m source in the atmosphere, J. Appl. Meteor. 18 (1961) 534–539.
M
1090
AN US
1086
[62] S. L. Manzello, J. R. Shields, T. G. Cleary, A. Maranghides, W. E.
1093
Mell, J. C. Yang, Y. Hayashi, D. Nii, T. Kurita, On the development
1094
and characterization of a firebrand generator, Fire Safety J. 43 (2008)
1095
258–268.
1097
PT
[63] H. H. Wang, Analysis on downwind distribution of firebrands sourced
CE
1096
ED
1092
from a wildland fire, Fire Technol. 47 (2011) 321–340.
[64] H. E. Anderson, Forest fuel ignitibility, Fire Technol. 6 (1970) 312–319.
1099
[65] A. Dimitrakopoulos, K. K. Papaioannou, Flammability assessment of
AC
1098
1100
Mediterranean forest fuels, Fire Technol. 37 (2) (2001) 143–152.
50
ACCEPTED MANUSCRIPT
[66] K. J. M. Dickinson, J. B. Kirkpatrick, The flammability and energy
1102
content of some important plant species and fuel components in the
1103
forests of Southeastern Tasmania, J. Biogeogr. 12 (1985) 121–134.
1104
[67] D. M. J. S. Bowman, B. A. Wilson, Fuel characteristics of coastal mon-
1105
soon forests, Northern Territory, Australia, J. Biogeogr. 15 (1988) 807–
1106
817.
CR IP T
1101
[68] J. B. Filippi, F. Bosseur, C. Mari, C. Lac, P. L. Moigne, B. Cuenot,
1108
D. Veynante, D. Cariolle, J. H. Balbi, Coupled atmosphere-wildland fire
1109
modelling, J. Adv. Model. Earth Syst. 1 (2009) art. #11.
1112
1113
1114
tinuous systems, Simulation 78 (2002) 76–89.
M
1111
[69] E. Kofman, A second order approximation for DEVS simulation of con-
[70] K. Hormann, A. Agathos, The point in polygon problem for arbitrary polygons, Comp. Geom.-Theo. Appl. 20 (2001) 131–144.
ED
1110
AN US
1107
[71] R. Wilson, The Devil Wind and Wood Shingles: The Los Angeles Conflagration of 1961, National Fire Protection Assn., 1962.
1116
URL https://books.google.es/books?id=bWFCywAACAAJ
1118
sin Magazine of History 54 (1971) 246–272. URL http://www.jstor.org/stable/4634648
AC
1119
[72] P. Pernin, The great Peshtigo fire: An eyewitness account, The Wiscon-
CE
1117
PT
1115
1120
1121
[73] P. J. Pagni, Causes of the 20 October 1991 Oakland Hills conflagration, Fire Safety J. 21 (1993) 331–339.
51
ACCEPTED MANUSCRIPT
[74] J. Trelles, P. J. Pagni, Fire-induced winds in the 20 October 1991 Oak-
1123
land Hills fire, in: Y. Hasemi (Ed.), Proc. of Fifth International Sympo-
1124
sium on Fire Safety Science, Melbourne, Australia, International Asso-
1125
ciation of Fire Safety Science: Boston, MA, 1997, pp. 911–922.
1126
1127
CR IP T
1122
[75] J. G. Quintiere, Principles of Fire Behavior, Career Education Series, Delmar Publishers, Albany, N.Y., 1998.
[76] E. Koo, P. J. Pagni, D. R. Weise, J. P. Woycheese, Firebrands and
1129
spotting ignition in large-scale fires, Int. J. Wildland Fire 19 (2010)
1130
818–843.
AC
CE
PT
ED
M
AN US
1128
52
ACCEPTED MANUSCRIPT
3500
Y (m)
"c = 4 "c = 8 "c = 12 "c = 18
1500 "q = 1, "q = 2, "q = 3, "q = 4,
"c = 8 "c = 8 "c = 8 "c = 8
2500
AN US
500 0 2000
CR IP T
"q = 2, "q = 2, "q = 2, "q = 2,
2500
3000
3500
4000
4500
X (m)
Figure 1: The fire-lines showing the numerical behaviour of DEVS for different pairs of ∆c and ∆q. Each contour represents the fire-line at 134 min. The initial fire is a circle of
AC
CE
PT
ED
M
radius 300 m and represented by 200 markers in the model.
Figure 2: The fire-lines showing the numerical behaviour of LSM for different pairs of dt
and dx. Each contour represents the fire-line at 380 min and the initial fire is a circle of radius 300 m.
53
ACCEPTED MANUSCRIPT
(a)
120 200 280
Y (m)
80
40 0
1000 1000
2000
3000
0 16
X (m)
3000
2000
4000
1000 1000
2000
3000
4000
X (m)
AN US
2000
Time [min]
0 40 80 120 160 180 220 240
Y (m)
240
3000
(b)
4000
CR IP T
4000
Figure 3: Evolution in time of the fire-line contour without random processes in absence of wind for the a) LSM and b) DEVS based simulators. The initial fire-line is a circle of
M
radius 300 m. The fire-break zone is 60 m wide.
(a)
ED
8000
6000
CE
2000
4000 X (m)
4000
140
6000
2000
AC
20000
20 40
0 20 40 60 80 100 120 140 160 180
Y (m)
PT
Y (m)
60
160 120 80 180
Time [min]
6000
100
4000
(b)
8000
0
2000
4000 X (m)
6000
Figure 4: Evolution in time of the fire-line contour without random processes in absence of wind with a non-homogeneous ROS for a) LSM and b) DEVS based simulators. The initial fire-line is a circle of radius 30 m.
54
ACCEPTED MANUSCRIPT
(a)
3500
(b)
3000
60
80 120
20
0 40
Y (m)
160
2500 140
100
2000
2500
2000
2500
3000 X (m)
3500
4000
1500 2000
2500
3000 X (m)
3500
4000
AN US
1500 2000
0 20 40 60 80 100 120 140 160
CR IP T
3000 Y (m)
Time [min]
3500
Figure 5: Evolution in time of the fire-line contour without random processes with an east wind of 3 ms−1 for a) LSM and b) DEVS based simulators. The initial fire-line is a circle
ED
M
of radius 300 m.
(a)
PT
3500
(b)
3000
40 60
3000
AC
2500
Y (m)
2500
120 140
2000
1500 2000
0 20 40 60 80 100 120 140 160
160
0
CE
2500
100
80
Y (m)
20
Time [min]
3500
2000
3000 X (m)
3500
4000
1500 2000
2500
3000 X (m)
3500
4000
Figure 6: Same as Fig. 5, but when the initial fire-line is a square of side 600 m.
55
ACCEPTED MANUSCRIPT
(a)
6000
40 0 20 0 1 30 7050
60
Y (m)
Y (m)
6000 4000 2000
4000
2000
2000
4000 X (m)
6000
8000
0
0
Time [min]
2000
4000 X (m)
6000
0 10 20 30 40 50 60 70
8000
AN US
00
(b)
8000
CR IP T
8000
Figure 7: Evolution in time of the fire-line contour without random processes with ROS given by formula (20) and when θ is the angle between the outward normal in a contour point and the mean wind direction for a) LSM and b) DEVS based simulators. The mean
ED
M
wind velocity is 3 ms−1 in the positive x-direction.
(a)
8000
PT
0 10 20 30 40 50 60 70
6000
Y (m)
20
0
10
Time [min]
60
4000
30
50
Y (m)
6000
(b)
8000
4000
CE
40
AC
2000
00
2000
2000
70
4000 X (m)
6000
8000
0
0
2000
4000 X (m)
6000
8000
Figure 8: Same as Fig. 7, but when θ is the angle between the line joining a contour point
and the mean wind direction.
56
ACCEPTED MANUSCRIPT
Time [min]
4000
4000
3000
2000
2000
1000
1000
0
1000
2000 3000 X (m)
4000
5000
0
0
1000
2000 3000 X (m)
AN US
0
Time [min]
0 20 40 60 80 100 120 140
Y (m)
Y (m)
3000
(b)
5000
0 20 40 60 80 100 120 140
CR IP T
(a)
5000
4000
5000
Figure 9: Evolution in time of the fire-line contour with turbulence in absence of wind for a) LSM and b) DEVS based simulators. The initial fire-line is a circle of radius 300 m.
ED
M
The turbulent diffusion coefficient is D = 0.15 m2 s−1 .
1
AC
CE
?, ?e
PT
0.8
0.6
0.4 ? - DEVS ? e - DEVS ? - LSM ? e - LSM
0.2
0
0
1000
2000
3000
4000
5000
X (m)
Figure 10: Cross section of the indicator function φ(x, t) and φe (x, t) at y = 2500 m and
t = 140 min, corresponding to the front evolution showed in Fig. 9.
57
ACCEPTED MANUSCRIPT
Time [min]
4000
4000
3000
2000
2000
1000
1000
0
1000
2000 3000 X (m)
4000
5000
0
0
1000
2000 3000 X (m)
AN US
0
Time [min]
0 20 40 60 80 100 120 140
Y (m)
Y (m)
3000
(b)
5000
0 20 40 60 80 100 120 140
CR IP T
(a)
5000
4000
5000
0 20 40 60 80 100 120 140 160
ED
4000
3000
PT
Y (m)
Time [min]
CE
2000
0
1000
2000 3000 X (m)
Time [min]
0 20 40 60 80 100 120 140 160
4000
3000
2000
4000
1000
5000
AC
1000
(b)
5000
Y (m)
(a)
5000
M
Figure 11: Same as Fig. 9, but with turbulent diffusion coefficient D = 0.30 m2 s−1 .
0
1000
2000 3000 X (m)
4000
5000
Figure 12: Evolution in time of the fire-line contour with turbulence and a north wind of 3 ms−1 for a) LSM and b) DEVS based simulators. The turbulent diffusion coefficient is D = 0.15 m2 s−1 .
58
ACCEPTED MANUSCRIPT
Time [min]
Y (m)
4000
3000
3000
2000
1000
2000 3000 X (m)
4000
1000
5000
0
1000
2000 3000 X (m)
4000
5000
AN US
0
Time [min]
0 20 40 60 80 100 120 140 160
4000
2000
1000
(b)
5000
Y (m)
0 20 40 60 80 100 120 140 160
CR IP T
(a)
5000
Figure 13: Same as Fig. 9, but when both turbulence and fire-spotting are included with
M
D = 0.15 m2 s−1 , µ = 2.69 and σ = 1.25.
ED
PT
Y (m)
4000
3000
Time [min]
0 20 40 60 80 100 120 140 160 180 200 220
CE
2000
0
1000
2000 3000 X (m)
Time [min]
0 20 40 60 80 100 120 140 160 180 200 220
4000
3000
2000
4000
1000
5000
AC
1000
(b)
5000
Y (m)
(a)
5000
0
1000
2000 3000 X (m)
4000
5000
Figure 14: Evolution in time of the fire-line contour with turbulence, a north wind of 3 ms−1 and two fire-break zones for a) LSM and b) DEVS based simulators. The initial fire-line is a circle of radius 300 m. The turbulent diffusion coefficient is D = 0.075 m2 s−1 .
59
AN US
CR IP T
ACCEPTED MANUSCRIPT
(a)
Time [min]
0 20 40 60 80 100 120 140 160 180 200 220
3000
2000
0
1000
2000 3000 X (m)
4000
5000
ED
1000
M
Y (m)
4000
5000
(b)
Time [min]
0 20 40 60 80 100 120 140 160 180 200 220
4000
Y (m)
5000
3000
2000
1000
0
1000
2000 3000 X (m)
4000
5000
Figure 15: Same as Fig. 14, but when both turbulence and fire-spotting are included with
AC
CE
PT
D = 0.075 m2 s−1 , µ = 2.69 and σ = 1.25.
60
ACCEPTED MANUSCRIPT
10000
10000
8000
8000 150
00
2000
4000 6000 X (m)
00
8000 10000
(a)
2000
25
0
2000
4000 6000 X (m)
00
8000 10000
ED
00
2000
(c)
(d)
00
2000
100
0
4000 6000 X (m)
(e)
150
4000 2000
0
00
8000 10000
5
6000
17
Y (m)
175
125
150
PT
8000
0 10 25 50
2000
8000 10000
10000
75 25 50
CE
Y (m)
8000
4000 6000 X (m)
75 125
10000
175
0
2000
M
2000
125
4000
0 10 25 50
50
75
6000
1 75 50
75
125
Y (m)
100
4000
AC
8000 10000
8000
6000
4000
4000 6000 X (m)
AN US
10000
8000
6000
150
0
(b)
10000
Y (m)
50 100
2000
CR IP T
Y (m)
125
50 100 0
2000
4000
75 25
75 25
4000
6000 125
Y (m)
6000
2000
4000 6000 X (m)
8000 10000
(f)
Figure 16: Series of figures showing the variation in the fire propagation due to increasing value of µ. From a-f, µ varies from 13 to 18, in increments of 1, while σ = 5.
61
ACCEPTED MANUSCRIPT
10000
10000
8000
8000
50
25
2000
4000 6000 X (m)
0
2000 00
8000 10000
2000
(a) 8000
25
25
50
0
2000
4000 6000 X (m)
00
8000 10000
2000
4000 6000 X (m)
425
565
150
PT
Y (m)
4000
360
360 Y (m)
CE
5
100
X (m)
32
250 50
AC
8000
0
0
4000
475 400
0
50
530 50 4
35 0
12000
200
2000
50
16000
350
10050
250
00
300 200 1
5
32
4000
20000
350
16000
8000
(d)
ED
(c)
20000
8000 10000
300
00
0
2000
M
2000
5 12 71500
4000
50
75
6000
150
75
125
Y (m)
100
4000
12000
0
8000
15
AN US
10000
6000
8000 10000
(b)
10000
Y (m)
4000 6000 X (m)
375
00
25
50 100
4000
0
2000
100
5
75
25
75
CR IP T
4000
150 1
6000
Y (m)
200 150 12
175
Y (m)
6000
6000
8000
00
10000
(e)
4000
8000
12000 X (m)
16000
20000
(f)
Figure 17: Series of figures showing the variation in the fire propagation due to increasing value of σ. From a-f, σ varies from 4 to 6.5, in increments of 0.5, while µ = 15.
62
CR IP T
7000
Y (m)
5000
75
5000
75
5000
AN US
7000 Y (m)
7000
9000
100
9000
95
9000
75
Y (m)
ACCEPTED MANUSCRIPT
95
M
ED
1000 2000
4000 6000 X (m)
(b)
8000
1000 2000
25 75
25 75
25 75
(a)
8000
0
4000 6000 X (m)
50
1000 2000
3000
0
12
0
50
3000
50
3000
0
4000 6000 X (m)
8000
(c)
Figure 18: Series of figures showing the generation of secondary fires due to fire-spotting
AC
CE
PT
and their subsequent merging with the primary fire when µ = 15 and σ = 5.
63