# Parameter Estimation of Fire Propagation Models
Using Level Set Methods

Angelo Alessandri  [ Patrizia Bagnerini  [ Mauro Gaggero  [ Luca Mantelli  [
Department of Mechanical Engineering, University of Genoa, I-16145 Genoa,
Italy  Institute of Marine Engineering, National Research Council of Italy,
I-16149 Genoa, Italy

###### Abstract

The availability of wildland fire propagation models with parameters estimated
in an accurate way starting from measurements of fire fronts is crucial to
predict the evolution of fire and allocate resources for firefighting. Thus,
we propose an approach to estimate the parameters of a wildland fire
propagation model combining an empirical rate of spread and level set methods
to describe the evolution of the fire front over time and space. The
estimation of parameters in the rate of spread is performed by using fire
front shapes measured at different time instants as well as wind velocity and
direction, landscape elevation, and vegetation distribution. Parameter
estimation is performed by solving an optimization problem, where the
objective function to be minimized is the symmetric difference between
predicted and measured fronts at different time instants. Numerical results
coming from by the application of the proposed method are reported in two
simulated scenarios and in a case study based on data originated by the 2002
Troy fire in Southern California. The obtained results showcase the
effectiveness of the proposed approach both from qualitative and quantitative
viewpoints.

###### keywords:

Wildland fire propagation model, level set methods, parameter estimation,
optimization.

mytitlenotemytitlenotefootnotetext: This work was supported by the Air Force
Office of Scientific Research under Grant FA9550-15-1-0530.

<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
##  1 Introduction

Wildland fires are a serious threat for human beings and protected natural
areas. Hence, the capability to forecast the evolution of fires over time and
space and evaluate the hazard level of a certain site is crucial. Many
different modeling approaches have been developed for this purpose, including
full-physics models, which provide the most accurate results by taking into
account many of the phenomena involved in fire propagation together with their
mutual interactions [1]. Such models consider regional and local weather as
well as terrain characteristics to forecast fire behavior, and also predict
production and dispersion of smoke. However, full-physics models may be too
computationally demanding for real-time applications, short-term forecasting,
and parameter estimation. As an alternative, empirical models have been
developed in the literature, such as the one proposed by Rothermel [2], where
the rate of spread of the fire is described through correlations based on
experimental results. The complexity of empirical models may be very
different, ranging from models where wind is constant [3] to models that
account for local interactions between heat sources and wind [4]. Indeed, the
level of detail adopted to describe topological features and vegetation (or
fuel) distribution may vary so much over space to undermine simulation
accuracy [5].

Based on the knowledge of the rate of spread, the evolution of the fire front
can be described by means of different techniques, which can be split into two
main categories depending on the implementation, i.e., raster and vector.
Examples of raster techniques are those based on cellular automata [6, 7, 8,
9, 10, 11], while the most popular vector approaches are level set methods
[12, 13, 14]. In raster implementations, the fire evolution only depends on
the interaction among contiguous cells, which can be denoted as burnt,
burning, or not burning. The fire propagation from one cell to its neighbors
is defined by a set of rules. Such methods are computationally efficient, but
the fire front can be subject to significant distortions. Many solutions have
been proposed to mitigate this issue, such as the increase of the number of
possible spread directions [7], particular kinds of discretization of the
spatial domain (e.g., using hexagonal grids) [9], and the adoption of suitable
correction factors [8]. In vector implementations, the fire front is defined
explicitly by a given amount of points or implicitly by the solution of a
partial differential equation (PDE). In both cases, the fire shape is updated
at different time steps based on the rate of spread. These methods are
accurate, but are more computationally demanding than raster implementations.
This issue worsens in the presence of merging fire fronts or unburnt areas
when using explicit methods (sometimes referred to as Lagrangian), but it is
less significant for implicit methods (also called Eulerian) such as level set
methods. The reader interested in a deeper comparison between Lagrangian and
Eulerian methods in the context of fire propagation is referred to [15] and
[16]. Owing to their properties, level set methods are widely employed also in
many other application fields, like image processing, computational fluid
dynamics, and material science [17, 18, 19, 20]. Instead, their use in the
context of estimation and control has received much less attention from the
research community (among the available results, see [21, 22, 23, 24, 25]).

In this paper, we focus on a wildland fire propagation model based on the
combination of an empirical description of the rate of spread and the use of
level set methods to account for fire front evolution. For this model, we
propose an approach for the optimal estimation of the main parameters that
results from the minimization of a least-squares cost fitting the available
measures of the fire front at different time instants. In our case, the cost
to be minimized is not continuously differentiable, and therefore we have to
resort to non-derivative methods such as the generalized pattern search (GPS)
algorithm [26]. After obtaining a fire propagation model with accurately
estimated parameters, the evolution of the fire front can be predicted in real
time, thus allowing the use of the model in a decision support system for fire
extinguishing.

The estimation of parameters involved in wildland fire propagation models
starting from available measures is a topic that has attracted the interest of
researchers in the last decades. To this purpose, various terminology has been
used, such as “model calibration” or “parameter identification.” As pointed
out in [27], wildland fire model calibration is typically accomplished
manually, and only very recently optimization algorithms have been adopted to
automate the estimation process, either for Lagrangian [28] or Eulerian [27]
models. In more detail, [29] focuses on parametric uncertainty quantification
of the rate of spread. Moreover, references are available that perform
parameter estimation using Kalman filtering approaches [30, 31, 32]. Other
studies combine parameter and state estimation based on a Luenberger observer
framework [33] or adopt other data-driven techniques based on fire propagation
simulators, also exploiting historical data [34, 35]. Such methods perform
dynamic estimation over time and are quite complex. Further, they are based on
specific fire propagation models that are different from the one considered in
this work. Instead, the method presented in this paper proposes a simple yet
effective method for calibration of the parameters that relies only on
observations of the fire front over time and wind measurements to perform
estimation of several parameters (some of them depending on vegetation) that
are involved in the rate of spread. Toward this end, the minimization of a
cost function measuring the symmetric difference between the observed and
predicted fronts is performed111This idea was exploited by the authors also in
[23, 24, 25] with the goal of finding optimal control policies for propagating
fronts (in different settings than wildland fires). Here the context is
completely different since we focus on parameter estimation rather than
control.. The considered fire propagation model combines the rate of spread
proposed in [36, 14] and level set methods to model front evolution. To the
best of the authors’ knowledge, no previous attempts concerning parameter
estimation for the model adopted in this paper are available in the
literature. This model was used also in [37] to predict fire evolution in a
real-scale shrubland fire scenario using machine learning techniques, but no
parameter estimation was performed. However, the focus on the rate of spread
reported in [36, 14] is not a limitation for the proposed approach, which
could be extended also to more complex expressions with no conceptual
difficulties. The effectiveness of the proposed method is tested from
qualitative and quantitative viewpoints in different scenarios using both
synthetic and real data. Another approach focused on parameter estimation
based on level set methods is proposed in [27], but it relies on different
equations and parameters as compared to our work. In more detail, the main
contribution of [27] is the development of ELMFire, a geospatial model for
simulating wildland fires using level set methods. Several equations are
inserted into the model to account for surface fire acceleration, crown fire
initiation, dead fuel moisture content, spot fire formation, as well as
adjustment factors due to wind. Then, a simple genetic optimization approach
is presented to automatically calibrate baseline model parameters, by adopting
a fitness function based on a scoring policy that compares burned and not
burned areas in a pixel-by-pixel fashion, differently from our approach where
we evaluate the symmetric difference between simulated and measured front.
However, unlike our work, in [27] no quantitative evaluation of the
effectiveness of proposed approach is given by comparing several intermediate
fronts during fire evolution.

The rest of this paper is structured as follows. Section 2 describes the
considered model for fire propagation. The proposed technique for parameter
estimation is presented in Section 3. Section 4 reports numerical results on
model testing and parameter estimation in two simulated case studies and in a
real one. Conclusions are drawn in Section 5.

##  2 Wildland Fire Front Propagation Modeling

In this section, we investigate the propagation of a wildland fire with
different environmental conditions by combining an empirical approach to
compute the rate of spread and level set methods to represent the moving fire
front.

###  2.1 Fire Rate of Spread

We consider the model of wildland fire rate of spread proposed by Lo in [36],
which is derived by previous works of Rothermel [2] and others [14]. The main
assumption of such model is that, in the presence of wind, if we consider a
flat landscape with homogeneous fuel, the fire spread velocity is different at
the head, rear, and flanks of the front. As showcased in Fig. 1, by “head” we
mean the region where fire propagation is in the wind direction, by “rear” we
denote the region where it is in the opposite direction of wind, and by
“flanks” we indicate the region where it is across of wind.

Figure 1: Propagation of a wildland fire front from time $t_{0}$ to $t_{f}$
under the influence of wind.

(a) (b)

Figure 2: Sketches of the definitions of the angles $\gamma$ (a) and $\chi$ of
the rate of spread (b). A red color of the terrain elevation isolines in (a)
denotes a greater elevation with respect to a blue one.

The head is the area with the highest velocity since wind blows high-
temperature combustion products towards unburnt regions. Instead, at the rear
such products are driven over already-burnt vegetation, and therefore velocity
is lower. At the flanks, the component of the wind speed that is perpendicular
to the front is zero, but from experimental tests it has been observed that
front propagation is faster than in the absence of wind. However, according to
[36], the effect of wind on the spread rate at the flanks is negligible since
its impact on the front is small222The interested reader can find a deeper
discussion on the modeling of flank fire using level set methods in [38, 13,
39]. A detailed investigation of this aspect is outside the scope of this
paper.. If we assume that a wildland fire propagates only in the normal
direction to the front, the front rate of spread $F(U,\theta,\gamma,\chi)$ is
given by

$F(U,\theta,\gamma,\chi):=\begin{cases}\\!\epsilon+a\sqrt{U\cos^{n}(\theta)}+\epsilon\psi(\chi)\cos(\gamma),\quad\mbox{if}\>\>\theta\in\Big{[}0,\frac{\pi}{2}\Big{]},\\\
\epsilon(\alpha+(1-\alpha)|\sin(\theta)|)+\epsilon\psi(\chi)\cos(\gamma),\quad\mbox{if}\>\>\theta\in\Big{(}\frac{\pi}{2},\pi\Big{]},\end{cases}$
(1)

with

$\psi(\chi):=5.275\,\beta^{-0.3}\,(\tan(\chi))^{2},$ (2)

where $U\geq 0$ is the module of the wind speed, $\theta\in[0,\pi]$ is the
angle between the wind and the normal direction to the front head (see Fig.
1), and $\gamma\in[0,\pi]$ is the angle between the elevation gradient
direction and the normal to the fire front (see Fig. 2(a)). The model also
depends on parameters $a\geq 0$ taking into account fuel, $\epsilon\geq 0$
accounting for the velocity of fire at the flanks, $\alpha\geq 0$ representing
the ratio between the velocity of fire at the rear and $\epsilon$, and $n\geq
0$ influencing the way the front propagates. The term $\psi(\chi)$ takes into
account the terrain effect and includes the fuel bed packing ratio parameter
$\beta\geq 0$ and the landscape slope $\tan(\chi)$, where
$\chi\in[-\pi/2,+\pi/2]$ is the angle between the terrain and the horizontal
line (see Fig. 2(b)). Since fire naturally tends to move upward, velocity is
higher when propagation is uphill, while it is lower when it spreads downhill
[14, 36]. Toward this end, the term $\psi(\chi)$ in (1) is multiplied by the
cosine of the angle $\gamma$. The coefficients $5.275$ and $-0.3$ in (2) have
been fixed according to the reference literature [2, 36].

The accuracy of the model (1)-(2) strictly depends on the chosen values for
the parameters $n$, $\epsilon$, $a$, $\alpha$, and $\beta$. Thus, in Section 3
we will propose an approach for their optimal estimation based on observations
of the fire front. Instead, the quantities $U$, $\theta$, $\gamma$, and $\chi$
in (1), (2) will be supposed to be known since they depend on wind and terrain
elevation.

###  2.2 Level Set Methods

We account for the evolution of the fire front by means of an Eulerian
approach based on level set methods. Toward this end, let
$\Omega\subset{\mathop{\mathbb{R}}}^{2}$ and $t\in[t_{0},t_{f}]$ be a compact
two-dimensional space domain and the time, respectively, where $t_{0}\geq 0$
and $t_{f}>t_{0}$ are given initial and final instants, respectively. Level
set methods represent a moving front at each time $t$, i.e., a curve in two
dimensions separating two regions, as the zero level set of a function
$\phi:\Omega\times[t_{0},t_{f}]\rightarrow\mathbb{R}$ (see Fig. 3). The front
$\mbox{\boldmath$x$}(t,s)$ is given at time $t$ by the points such that
$\phi(\mbox{\boldmath$x$}(t,s),t)=0$, where $s$ is the arc-length parameter of
the initial curve $\mbox{\boldmath$x$}(t_{0},s)$. If we differentiate with
respect to $t$, we obtain

$\phi_{t}(\mbox{\boldmath$x$},t)+v(\mbox{\boldmath$x$},t)\cdot\nabla\phi(\mbox{\boldmath$x$},t)=0.$
(3)

Equation (3) is a Hamilton-Jacobi PDE, where
$v(\mbox{\boldmath$x$},t):=\frac{d}{dt}\mbox{\boldmath$x$}(t,s)$ is the
Lagrangian particle velocity giving the direction of propagation of the front
at the point $\mbox{\boldmath$x$}(t,s)$ and $\nabla$ denotes the gradient with
respect to space. In the case of wildland fires, we focus on the normal flow
equation since we assume that propagation only occurs in the normal direction
to the fire front (see Section 2.1). This corresponds to choose
$v(\mbox{\boldmath$x$},t)$ proportional to the normal $n$ to the front, i.e.,

$v(\mbox{\boldmath$x$},t):=S\,\mbox{\boldmath$n$}=S\,\frac{\nabla\phi(\mbox{\boldmath$x$},t)}{|\nabla\phi(\mbox{\boldmath$x$},t)|},\,$
(4)

where $S$ is the speed of front propagation and $n:=\nabla\phi/|\nabla\phi|$
is the normal to the front. In our case, $S$ is given by the rate of spread
$F(U,\theta,\gamma,\chi)$ defined in (1). If we replace (4) in (3) and $S$
with $F(U,\theta,\gamma,\chi)$, we get the following normal flow equation:

$\phi_{t}(\mbox{\boldmath$x$},t)+F(U,\theta,\gamma,\chi)\,|\nabla\phi(\mbox{\boldmath$x$},t)|=0.$
(5)

Since $F(U,\theta,\gamma,\chi)$ depends on the angle $\theta$ between the
front (i.e., the zero level set of $\phi$) and the wind, it can be considered
as a function of $\phi$. Thus, with a little abuse of notation, we can re-
write (5) as follows:

$\phi_{t}(\mbox{\boldmath$x$},t)+F(\phi)\,|\nabla\phi(\mbox{\boldmath$x$},t)|=0.$
(6)

Equation (6) is associated with initial conditions
$\phi_{0}:\Omega\rightarrow{\mathop{\mathbb{R}}}$, i.e.,
$\phi(\mbox{\boldmath$x$},t_{0})=\phi_{0}(\mbox{\boldmath$x$})$ for all
$\mbox{\boldmath$x$}\in\Omega$. Usually, $\phi_{0}$ is the signed distance to
the initial front, and it is conventionally positive inside the fire front
(the burnt area) and negative outside (the unburnt area).

At each time $t$, the fire front is represented by the zero level set of the
function $\phi$ that solves (6), i.e., it is a set-valued mapping
$\Gamma:[t_{0},t_{f}]\rightrightarrows\mathcal{C}\subset{\mathop{\mathbb{R}}}^{2}$,
where
$\Gamma(t):=\left\\{\mbox{\boldmath$x$}\in\Omega:\phi(\mbox{\boldmath$x$},t)=0\right\\}$.

Figure 3: Example of function $\phi$ and its zero level set. Figure 4:
Schematic diagram of the parameter estimation procedure.

Since the rate of spread $F$ defined in (1) depends on the angle between the
fire front and the wind, it is not defined everywhere in the domain $\Omega$,
but only for the points $x$ such that $\phi(\mbox{\boldmath$x$},t)=0$.
However, it is possible to extend the expression of $F$ in all $\Omega$ by
considering the angle formed by all level sets of $\phi$ with the wind. In
this way, $F$ is defined for all points $\mbox{\boldmath$x$}\in\Omega$. The
value of $\theta$ in (1) is computed at each time $t$ as the angle between the
normal to the front and the wind direction $U$, i.e.,

$\theta=\arccos\left(\frac{\nabla\phi(\mbox{\boldmath$x$},t)\cdot\mbox{\boldmath$U$}}{|\nabla\phi(\mbox{\boldmath$x$},t)|\,U}\right),$

where $\cdot$ denotes the usual scalar product and $U$ is the module of the
vector $U$.

Level set methods are well-suited to describing wildland fire evolution since
no front reconstruction algorithms are required to determine the front shape
at a given time. In fact, the fire front can be obtained implicitly from
$\phi$, which is defined continuously on the spatial domain at any time
instant. For this reason, level set methods are able to easily manage changes
of topology, such as merging of fire fronts, or large deformations.

##  3 Parameter Estimation

The model of fire propagation described in Section 2 depends on several
parameters describing the characteristics of the region where fire propagates,
such as the type of vegetation or terrain elevation. In this section, we
investigate how to compute optimal estimates of quantities related to the kind
of fuel, while the quantities depending on terrain elevation will be assumed
known. In more detail, we focus on estimation of the parameters $n$,
$\epsilon$, $a$, $\alpha$, and $\beta$ involved in the rate of spread. A
precise knowledge of the values of such parameters is crucial to accurately
forecast the evolution of wildland fires and therefore plan effective
firefighting procedures.

The estimation of parameters in the rate of spread equation (1) is performed
by using measurements of the fire front shape at different time instants,
velocity and direction of wind, landscape elevation, and fuel distribution.
Without loss of generality, the wind is assumed to be constant (as regards
both velocity and direction) in the domain $\Omega$ over the whole simulation
time. Moreover, we assume that two different kinds of fuels are present in
$\Omega$, denoted as “fuel A” and “fuel B.” We suppose to know the fuel
distribution over $\Omega$, i.e., whether fuel A or fuel B is present in the
various regions of the domain (such information can be obtained, for instance,
from aerial observations [40]), but its properties are supposed to be unknown.
In other words, the values of the parameters $\epsilon$, $a$, and $\beta$ for
the two fuels, denoted by $\epsilon^{A}$, $a^{A}$, and $\beta^{A}$ for fuel A
and $\epsilon^{B}$, $a^{B}$, and $\beta^{B}$ for fuel B, are assumed to be
unknown. The parameters $n$ and $\alpha$ are unknown as well, while the
quantities $\theta$, $\gamma$, and $\tan(\chi)$ can be easily computed by
combining the available information on the considered region (i.e., landscape
elevation) and observations of the fire front. Hence, there is no need to
estimate them. Summarizing, we collect all unknown parameters in a
$8$-dimensional vector $p$, as follows:

$\mbox{\boldmath$p$}:=\left(n,\epsilon^{A},\epsilon^{B},a^{A},a^{B},\alpha,\beta^{A},\beta^{B}\right)\in{\mathop{\mathbb{R}}}^{8}.$

Let $\Gamma^{\rm meas}(t)$ be the fire front measured from observations at
time $t$ and let $\phi^{\rm meas}(x,t)$ be the corresponding level set
function obtained as the signed distance to the front (conventionally,
positive inside the front and negative outside). Moreover, let
$\Gamma(t,\mbox{\boldmath$p$})$ be the fire front predicted by the model
described in Section 2 corresponding to a given value for the parameter vector
$p$, and let $\phi(x,t,\mbox{\boldmath$p$})$ be the corresponding level set
function obtained as the signed distance to the front using the same
convention described above. The best estimate $\hat{\mbox{\boldmath$p$}}$ of
$p$ can be obtained by minimizing a least-squares fitting cost function as
follows:

$\hat{\mbox{\boldmath$p$}}:=\arg\>\min_{\mbox{\boldmath$p$}\in[\mbox{\boldmath$p$}_{\rm
min},\mbox{\boldmath$p$}_{\rm max}]}J(\mbox{\boldmath$p$})$

where $\mbox{\boldmath$p$}_{\rm min}\in{\mathop{\mathbb{R}}}^{8}$ and
$\mbox{\boldmath$p$}_{\rm max}\in{\mathop{\mathbb{R}}}^{8}$ are given lower
and upper bounds for the vector $p$, respectively, and

$J(\mbox{\boldmath$p$}):=\int_{t_{0}}^{t_{f}}\eta\left(\Gamma(t,\mbox{\boldmath$p$})\>\,\Delta\>\,\Gamma^{\rm
meas}(t)\right)\,dt.$ (7)

The operator $\Delta$ in (7) is the symmetric difference, i.e.,
$\Gamma(t,\mbox{\boldmath$p$})\>\>\Delta\>\>\Gamma^{\rm
meas}(t):=(\Gamma(t,\mbox{\boldmath$p$})\cup\Gamma^{\rm
meas}(t))\setminus(\Gamma(t,\mbox{\boldmath$p$})\cap\Gamma^{\rm meas}(t))$,
and $\eta$ corresponds to an outer measure on ${\mathop{\mathbb{R}}}^{2}$.

(a) (b)

(c) (d)

Figure 5: Fire front propagation without wind on a flat landscape (a), with
wind on a flat landscape (b), without wind in the presence of a hill (c), and
without wind in the presence of a depression (d) (fire fronts are depicted in
magenta and elevation isolines are represented in black).

Owing to the definition of the symmetric difference and the convention used
for the sign of the function $\phi$, the cost (7) can be rewritten as follows
to highlight its dependence on $\phi$:

$J(\mbox{\boldmath$p$}):=\int_{t_{0}}^{t_{f}}\\!\\!\int_{\Omega}\left(H(\phi(\mbox{\boldmath$x$},t,\mbox{\boldmath$p$}))\\!-\\!H(\phi^{\rm
meas}(\mbox{\boldmath$x$},t))\right)^{2}d\mbox{\boldmath$x$}\,dt,$ (8)

where $H(\cdot)$ is the Heaviside step function.

At each iteration of the optimization algorithm chosen to minimize the cost
$J$ in (8), a certain value for the parameter vector $p$ is used to run a fire
propagation simulation in the time interval $[t_{0},t_{f}]$ and the domain
$\Omega$. Hence, the cost $J$ is computed by comparing the predicted fronts
and the measured ones. Then, an updated value for the parameter vector $p$ is
chosen by the optimization algorithm, and a new simulation of fire propagation
is performed to compute a new value for the cost. The procedure is iterated
until convergence of the optimization algorithm is obtained and stopping
criteria are satisfied. In more detail, they are satisfied if a “small” norm
of the difference between the estimated parameters or the values of the cost
function in two consecutive iterations is found, as well as if a maximum
number of iterations is exceeded. Fig. 4 sketches a schematic diagram of the
considered parameter estimation procedure, where $k$ denotes a generic
iteration of the optimization algorithm used to minimize (8).

Unfortunately, the cost function associated with the front evolution may be
affected by several local minima. Thus, optimization may be subject to local
minima trapping, which can provide poor fitting and inaccurate parameter
estimates.

##  4 Numerical Results

In this section, we first evaluate the effectiveness of the model considered
in Section 2 to predict fire fronts in simple case studies taken from the
reference literature [27]. Then, we focus on parameter identification using
the approach presented in Section 3. In particular, to evaluate the
effectiveness of the proposed technique for estimating the parameter vector
$p$, we consider three case studies denoted by “Valley,” “Hill,” and “Troy.”
The first and second ones are based on simulated data, while the third case
involves measurements from a real fire event.

Either for testing or parameter estimation, a suitable discretization method
is required for the numerical solution of the model of fire front propagation
described in Section 2. Toward this end, we discretize the spatial domain on a
squared regular grid, with cells characterized by given fuel properties and
elevation. We compute slope and elevation gradients along the North-South and
East-West directions from the elevation matrix using a finite-differences
central scheme for inner points and single-sided finite differences for the
boundary.

To numerically solve the Hamilton-Jacobi equation (6), we adopt the Matlab
toolbox of level set methods [41]. In more detail, an upwind third-order
essentially non-oscillatory scheme [20] is applied to approximate space
derivatives, while a three-steps second-order Runge-Kutta scheme is used to
compute time derivatives. From the values of $\phi$, it is straightforward to
obtain the fire front shapes by computing the zero level sets.

As a consequence of discretization, the cost $J(\mbox{\boldmath$p$})$ in (8)
does not vary in a continuous way, and therefore it is not possible to use
gradient-based optimization methods requiring the computation of derivatives
of the cost. Without loss of generality, among the various alternatives, in
this paper we focus on the mesh adaptive GPS algorithm owing to its
computational efficiency and accuracy of the solution [26]. The optimization
routine implementing this method is available in Matlab through the
patternsearch function of the Global Optimization Toolbox. The tolerance for
the stopping criteria is set to $10^{-4}$, while the maximum number of
iterations is equal to $2000$. The time integral in the cost function (8) is
discretized over time with a sampling interval $\Delta t$ that varies on a
case-by-case basis, as detailed in the following sections. All the tests are
carried out on a computer equipped with a 2.5 GHz Intel Xeon CPU with 32 GB of
RAM.

(a) (b)

Figure 6: Valley case study: 3D elevation map (a) and initial condition (b).
The fire front at $t_{0}$ is in red, the wind direction is in blue, the
elevation level curves are in black, fuel A is in light-yellow, and fuel B is
in light-green.

(a) (b)

Figure 7: Valley case study: measured fire fronts obtained with
$\mbox{\boldmath$p$}=\mbox{\boldmath$p$}^{*}$ (a) and predicted fire fronts
computed with $\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}$ (b).

(a) (b)

Figure 8: Valley case study: trends of the cost function (a) and of the
parameters (b) during the estimation procedure.

###  4.1 Model Testing

Along the lines of [27] and others, the model described in Section 2 is first
tested in four simple cases taken from the reference literature. In all cases,
a circle is considered as starting shape of the fire front, and a squared
$[-1,1]\times[-1,1]$ domain $\Omega$ is chosen and discretized by means of a
$101\times 101$ regular grid; the same fuel is considered in all the domain.
The following parameters are adopted: $n=3$, $\epsilon=0.4$, $a=0.5$,
$\alpha=0.5$, and $\beta=0.02$.

In the first case (see Fig. 5(a)), no wind is considered and a final time
$t_{f}$ equal to $0.5$ is fixed starting from $t_{0}=0$. Owing to the absence
of wind, the fire front maintains its circular shape while propagating, as
expected [36]. In the second case, depicted in Fig. 5(b), a wind with speed
components $(-0.5,-0.5)$ is considered over a time span from $t_{0}=0$ to
$t_{f}=0.5$. The fire front evolves predominantly in the wind direction with
an elliptical shape, as shown also in [14, 36, 42].

In the third and fourth cases (reported in Fig. 5(c) and Fig. 5(d),
respectively), no wind is considered, the landscape elevation is non-constant,
and the initial and final times are set equal to $0$ and $1.5$, respectively.
In more detail, a hill is present in the third case and, when the fire goes
uphill, its spread rate increases (the zero level set curves of $\phi$
considered at a constant time step are more distant one to the others). When
the fire is close to the top of the hill, its velocity decreases and it
spreads downhill very slowly (the zero level set curves are close one to the
others). In the fourth case, a depression is present and the fire behavior is
similar, with a lower spread rate when the front moves downhill the depression
and a higher one when it moves uphill. Again, such behaviors are in accordance
with [36, 14].

Hence, we conclude that the model described in Section 2 behaves as expected
and is able to predict future evolution of fire fronts. Thus, we can proceed
with the optimal estimation of its parameter to adapt the predicted fronts to
the available measures.

###  4.2 Parameter Estimation: Valley Case Study

In this case study, we consider a region with a terrain characterized by a
depression, as shown in Fig. 6(a). The initial shape of the fire front is a
circle (see Fig. 6(b)), and the wind speed components are equal to
$(-0.5,-2)$. The domain is given by $\Omega=[-1,1]\times[-1,1]$, and it is
discretized by means of a grid made up of $101\times 101$ cells. Two kinds of
fuels are present: the fuel A, represented by the light-yellow region in Fig.
6(b), is characterized by $\epsilon^{A}=0.8$, $a^{A}=0.7$, and
$\beta^{A}=0.03$, while the fuel B, given by the light-green region, is
characterized by $\epsilon^{B}=0.4$, $a^{B}=0.4$, and $\beta^{B}=0.08$. The
parameters $n$ and $\alpha$ are equal to 3 and 0.5, respectively, on the whole
domain. Hence, we define
$\mbox{\boldmath$p$}^{*}:=(3,0.8,0.4,0.7,0.4,0.5,0.03,0.08)\in{\mathop{\mathbb{R}}}^{8}$
as the vector of true parameters. Such true parameters are first used to run a
simulation with sampling, initial, and final times equal to $0.01$, $0$, and
$0.1$, respectively. As shown in Fig. 7(a), the wildland fire propagates
faster in the region with fuel A and slows down going downhill of the
depression.

Then, the obtained fire front shapes are saved with a sampling time $\Delta
t=0.01$ and used as measurement fronts $\Gamma^{\rm meas}(t)$ to perform
estimation of the parameter vector $p$ using the procedure described in
Section 3. Toward this end, the initial values of the identified parameters
vector $p$ are randomly chosen between $\mbox{\boldmath$p$}_{\rm
min}=(2,0.1,0.1,0.1,0.1,0.01,0.12,0.12)$ and $\mbox{\boldmath$p$}_{\rm
max}=(4,3,3,3,3,1,0.001,$ $0.001)$. The trends of the components of $p$ and of
the cost function $J$ during the optimization procedure are shown in Fig. 8.
Optimization was performed in $5409$ s.

It turns out that, after a transient behavior, the cost and the components of
the vector $p$ converge to stationary values, which represent the optimal
estimates. The final estimate of the parameter vector $p$ obtained by the
optimization algorithm is
$\hat{\mbox{\boldmath$p$}}=(3.44,0.91,0.44,0.78,0.53,0.5,0.035,0.066)$, while
the optimal cost $J(\hat{\mbox{\boldmath$p$}})$ is equal to 13.

(a) (b)

Figure 9: Hill case study: 3D elevation map (a) and initial condition (b). The
fire front at $t_{0}$ is in red, the wind direction is in blue, the elevation
level curves are in black, fuel A is in light-yellow, and fuel B is in light-
green.

(a) (b)

Figure 10: Hill case study: measured fire fronts obtained with
$\mbox{\boldmath$p$}=\mbox{\boldmath$p$}^{*}$ (a) and predicted fire fronts
computed with $\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}$ (b).

(a) (b)

Figure 11: Hill case study: trends of the cost function (a) and of the
parameters (b) during the estimation procedure.

The accuracy of the estimation can be quantified by computing the following
two performance indicators.

  * 1)

The relative estimation error $e$, defined as the norm of the difference
between $\hat{\mbox{\boldmath$p$}}$ and $\mbox{\boldmath$p$}^{*}$ with respect
to the norm of $\mbox{\boldmath$p$}^{*}$, i.e.,

$e:=\frac{\|\hat{\mbox{\boldmath$p$}}-\mbox{\boldmath$p$}^{*}\|}{\|\mbox{\boldmath$p$}^{*}\|}.$
(9)

  * 2)

The quantity $r$, defined as the ratio between the value of the optimal cost
$J(\hat{\mbox{\boldmath$p$}})$ and the overall number of cells $N$ used to
discretize the domain, i.e.,

$r:=\frac{J(\hat{\mbox{\boldmath$p$}})}{N}.$ (10)

It provides the number of cells of the discretization grid where the sign of
the functions $\phi^{\rm meas}$ and $\phi$ of the measured and predicted
fronts, respectively, differ with respect to the total number of cells.

Clearly, the lower are $e$ and $r$, the better is the parameter estimation.

In the valley case study, we have $e=0.14$ and $r=0.0013$. A simulation of
fire propagation with $\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}$ provides
the fire front snapshot predictions depicted in Fig. 7(b), which are very
close to the corresponding “measured” ones reported in Fig. 7(a). This
confirms the effectiveness of the proposed approach for parameter estimation
in this case study.

###  4.3 Parameter Estimation: Hill Case Study

(a) 13:30 (b) 13:46

(c) 14:05 (d) 14:25

(e) 14:32 (f) 14:39

(g) 14:46 (h) 14:52

(i) 15:01 (j) 15:14

(k) 15:21 (l) 15:28

(m) 15:35 (n) 15:42

(o) 15:48 (p) 15:55

Figure 12: Evolution of the Troy fire between 13:30 and 15:55. The red area
represents the burnt region, while the black contour is the fire front.

This case study is characterized by the presence of a hill instead of a valley
in the domain $\Omega$, as showcased in Fig. 6(a). The initial shape of the
fire front is still a circle (see Fig. 6(b)), and the wind speed components
are fixed to $(-1,-3)$. The domain is given again by
$\Omega=[-1,1]\times[-1,1]$, and it is discretized with a $101\times 101$
grid. Two fuels are present, referred to as A and B, as in the previous case
study. The fuel A is in the light-yellow region in Fig. 9(b), and it is
characterized by $\epsilon^{A}=0.8$, $a^{A}=0.6$, and $\beta^{A}=0.02$, while
the fuel B is in the light-green region, and we have $\epsilon^{B}=0.5$,
$a^{B}=0.4$, and $\beta^{B}=0.06$. On the whole domain, $n$ and $\alpha$ are
set equal to $3$ and $0.5$, respectively. Thus, the vector of true parameters
is given by
$\mbox{\boldmath$p$}^{*}:=(3,0.8,0.5,0.6,0.4,0.5,0.02,0.06)\in{\mathop{\mathbb{R}}}^{8}$.
As in the previous case study, these parameters are first used to run a
simulation with sampling, initial, and final time equal to $0.01$, $0$, and
$0.1$, respectively. As shown in Fig. 10(a), the fire propagates faster in the
region where fuel A is present, it goes faster going uphill and slows down
going downhill.

The obtained fire front shapes are saved with a sample time $\Delta t=0.01$
and used as measurement fronts $\Gamma^{\rm meas}(t)$ to perform estimation of
the parameter vector $p$ using the procedure described in Section 3. The
initial values of the parameter vector $p$ are randomly chosen between
$\mbox{\boldmath$p$}_{\rm min}=(2,0.1,0.1,0.1,0.1,0.1,0.12,$ $0.12)$ and
$\mbox{\boldmath$p$}_{\rm max}=(4,3,3,3,3,1,0.001,0.001)$. The trends of the
components of $p$ and of the cost function $J$ during the optimization
procedure are shown in Fig. 11.

Likewise in the valley case study, after a transient behavior, the cost and
the components of the vector $p$ converge to stationary values, which
represent the optimal estimates. The final estimate of the parameter vector
$p$ obtained by the optimization algorithm is
$\hat{\mbox{\boldmath$p$}}=(2.38,0.73,0.45,0.87,0.36,$ $0.48,0.014,0.11)$,
while the optimal cost $J(\hat{\mbox{\boldmath$p$}})$ is equal to 17.
Optimization was performed in $5165$ s.

The fire front with $\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}$ propagates
as showcased in Fig. 10(b). Also for this case study, there is a good
correspondence between the simulation with real parameters (reported in Fig.
10(a)) and the one run with the identified parameters (see Fig. 10(b)), thus
confirming the effectiveness of the estimation procedure. In fact, we have
obtained an estimation error $e=0.20$ and a ratio between the value of the
optimal cost and the overall number of cells used to discretize the domain
equal to $r=0.017$.

###  4.4 Parameter Estimation: Troy Case Study

(a) (b)

Figure 13: Troy case study: 3D elevation map of the region where Troy fire
propagated (a) and initial condition at 14:05 (b). The initial fire front is
in red, the wind direction is in blue, the elevation level curves are in
black, fuel A is in light-yellow, and fuel B is in light-green.

(a) (b)

Figure 14: Troy case study: trends of the cost function (a) and of the
parameters (b) during the estimation procedure between 14:05 and 15:01.

(a) (b)

Figure 15: Comparison between simulation and real fronts of the Troy fire with
$\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}_{1}$. Fire evolution is
compared at 14:05, 14:25, 14:39, and 15:01 in (a) and at 15:01, 15:21, 15:42,
and 15:55 in (b). In both cases, the initial fronts at 14:05 or 15:01 are in
red, simulated fronts are in black, and measured fronts are in dotted green.

(a) (b)

Figure 16: Troy case study: trends of the cost function (a) and of the
parameters (b) during the estimation procedure between 15:01 and 15:55. Figure
17: Comparison between simulation and real images of the Troy fire with
$\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}_{2}$. Fire evolution is
compared at 15:01, 15:21, 15:42, and 15:55. The initial front at 15:01 is in
red, simulated fronts are in black, and measured fronts are in dotted green.

After testing the model on four known scenarios (see Section 4.1) and
verifying the effectiveness of the proposed parameter estimation approach in
two simulated case studies (see Sections 4.2 and 4.3), the proposed method is
applied to a real case study. Toward this end, we consider the 2002 Troy fire.
This fire started in Laguna Mountains in San Diego County, California, on June
19, 2002 and burnt about $1188$ acres of land [43]. The evolution of this fire
was recorded by Pacific Southwest Research Stations through measurements of
the emitted thermal-infrared light from an aerial view, between 13:30 and
15:55 [44]. From the observation of the aerial images (sketched in Fig. 12),
it is possible to notice a change of behavior occurring around 14:00. In fact,
the wildfire front spreads mainly towards East and South-East between 13:30
and 13:46, while it deviates towards North-East after 14:05.

We start considering the fire fronts at 14:05, 14:25, 14:39, and 15:01 for a
first round of test. Thus, the fire images captured by satellites are used to
perform parameter identification, setting the front at 14:05 as initial
condition. A squared area of $2879\times 2879$ m2 is taken as spatial domain,
i.e., we choose $\Omega=[0,2879]\times[0,2879]$. Such a region is discretized
with a regular grid made up of $101\times 101$ cells. Local map and elevation
data are obtained from Google Earth software and shown in Fig. 13(a).
Vegetation data are collected from the US Forest Service LANDFIRE website [45]
and then grouped together in two main fuel categories, named A and B as in the
previous two case studies. Fuel A includes short grass, timber grass,
chaparral, and brush, while fuel B is composed of short needle litter, long
needle or hardwood litter, timber litter, and understory. In order to apply
the method proposed in this paper, it is necessary to know wind speed and
direction during the whole fire evolution. Based on the observations of [36]
and on how the Troy fire evolved during the considered period, we fix wind
velocity to a constant value equal to 30 mph, while the wind direction is
supposed to be constant and equal to $63^{\circ}$ (being $0^{\circ}$ the North
direction and $90^{\circ}$ the East direction). Fig. 13(b) shows the fire
front at $t_{0}=\text{14:05}$, which, as said, is the initial condition of
this test, together with elevation level curves, fuel distribution, and wind
direction.

As previously pointed out, we first consider the fire fronts between time
instants $t_{0}=\text{14:05}$ and $t_{f}=\text{15:01}$ to perform parameter
estimation. The objective function $J$ is evaluated at 14:25, 14:39 and 15:01,
according to (8). The initial value of the parameter vector $p$ is randomly
chosen between $\mbox{\boldmath$p$}_{\rm
min}=(2,0.001,0.001,0.001,0.001,0.0001,0.0001,0.0001)$, and
$\mbox{\boldmath$p$}_{\rm max}$ = (4, 0.1, 0.1, 0.1, 0.1, 0.1, 0.12, 0.12).
The trends of the components of $p$ and of the cost function $J$ during the
optimization procedure are shown in Fig. 14. Optimization was performed in
$13213$ s.

Likewise in the valley and hill case studies, after a transient behavior, the
cost and the components of the vector $p$ converge to stationary values, which
represent the optimal estimates. The final estimate of the parameter vector
$p$ obtained by the optimization algorithm is $\hat{\mbox{\boldmath$p$}}_{1}=$
(3.91, 0.0028, 0.0635, 0.0463, 0.0313, 0.0056, 0.0980, 0.1030), while the
optimal cost is equal to 882. A simulation of fire propagation with
$\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}_{1}$ from $t_{0}=\text{14:05}$
to $t_{f}=\text{15:01}$ provides the fire front snapshots depicted in Fig.
15(a). The ratio $r$, computed as defined in (10), is equal to $0.0865$, while
it is not possible to compute the error $e$ defined in (9) since the value of
the true parameters $\mbox{\boldmath$p$}^{*}$ is not available.

To better evaluate the accuracy of parameter estimation, we introduce three
additional performance indexes taken from the reference literature (see [46]),
i.e., the Sørensen similarity index (denoted by SSI), the Jaccard similarity
coefficient (JSC, for short), and the Kappa statistics (named KS). They are
defined as follows:

$\text{SSI}(t,\mbox{\boldmath$p$}):=\dfrac{2|\Gamma^{\rm
meas}(t)\cap\Gamma(t,\mbox{\boldmath$p$})|}{|\Gamma^{\rm
meas}(t)|+|\Gamma(t,\mbox{\boldmath$p$})|},$ (11)
$\text{JSC}(t,\mbox{\boldmath$p$}):=\dfrac{|\Gamma^{\rm
meas}(t)\cap\Gamma(t,\mbox{\boldmath$p$})|}{|\Gamma^{\rm
meas}(t)\cup\Gamma(t,\mbox{\boldmath$p$})|},$ (12)
$\text{KS}(t,\mbox{\boldmath$p$}):=\dfrac{P_{a}(t,\mbox{\boldmath$p$})-P_{e}(t,\mbox{\boldmath$p$})}{1-P_{e}(t,\mbox{\boldmath$p$})},$
(13)

where $|\Gamma(t,\mbox{\boldmath$p$})|$ denotes the area of the region within
the front $\Gamma(t,\mbox{\boldmath$p$})$, as well as

$P_{a}(t,\mbox{\boldmath$p$}):=\frac{|\Gamma^{\rm
meas}(t)\cap\Gamma(t,\mbox{\boldmath$p$})|}{|\Omega|}+\frac{|\Omega\setminus(\Gamma^{\rm
meas}(t)\cup\Gamma(t,\mbox{\boldmath$p$}))|}{|\Omega|},$
$P_{e}(t,\mbox{\boldmath$p$}):=\frac{|\Gamma^{\rm
meas}(t)||\Gamma(t,\mbox{\boldmath$p$})|}{|\Omega|^{2}}+\frac{|\Omega\setminus\Gamma^{\rm
meas}(t)||\Omega\setminus\Gamma(t,\mbox{\boldmath$p$})|}{|\Omega|^{2}}.$

Table 1 contains the values of the performance indexes (11), (12), and (13)
computed on the basis of real images and results of simulations at 14:25,
14:39 and 15:01, with $\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}_{1}$.
Referring again to [46], we point out that the obtained values are large, thus
proving the effectiveness of the proposed approach.

Table 1: Values of performance indexes (11), (12), and (13) for the Troy case study obtained between 14:05 and 15:55 with $\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}_{1}$. | $t\\!=\\!\text{14:25}$ | $t\\!=\\!\text{14:39}$ | $t\\!=\\!\text{15:01}$ | $t\\!=\\!\text{15:21}$ | $t\\!=\\!\text{15:42}$ | $t\\!=\\!\text{15:55}$
---|---|---|---|---|---|---
$\text{SSI}(t,\mbox{\boldmath$p$})$ | 0.9416 | 0.9082 | 0.8439 | 0.9221 | 0.8190 | 0.7858
$\text{JSC}(t,\mbox{\boldmath$p$})$ | 0.9613 | 0.9362 | 0.8806 | 0.9396 | 0.8442 | 0.8078
$\text{KS}(t,\mbox{\boldmath$p$})$ | 0.9699 | 0.9519 | 0.9153 | 0.9595 | 0.9005 | 0.8801
Table 2: Values of performance indexes (11), (12), and (13) for the Troy case study obtained between 15:01 and 15:55, setting $\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}_{2}$. | $t\\!=\\!\text{15:21}$ | $t\\!=\\!\text{15:42}$ | $t\\!=\\!\text{15:55}$
---|---|---|---
$\text{SSI}(t,\mbox{\boldmath$p$})$ | 0.9424 | 0.8957 | 0.8987
$\text{JSC}(t,\mbox{\boldmath$p$})$ | 0.9564 | 0.9165 | 0.9178
$\text{KS}(t,\mbox{\boldmath$p$})$ | 0.9703 | 0.9450 | 0.9466

After estimating the parameters $\hat{\mbox{\boldmath$p$}}_{1}$, we perform a
further simulation run with the identified parameters
$\hat{\mbox{\boldmath$p$}}_{1}$ in order to evaluate the capability to
forecast the evolution of the fire front. In more detail, we let
$\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}_{1}$, as well as we set the
initial fire front shape equal to the real front at 15:01 and fixed the final
time at 15:55. The simulation results are then compared with the real fire
front shapes at 15:21, 15:42 and 15:55. The results are reported in Fig.
15(b). It is possible to observe an over-estimation of the fire propagation
towards North-East, which brings a small part of the fire front outside of the
spatial domain, and a good approximation of the fire evolution in the other
directions. The performance indexes (11), (12), and (13) of this simulation
run are shown again in Table 1. They are still large, even if we may observe a
reduction of their values, especially at 15:55.

The last result suggests that the accuracy of the forecasts decreases if we
consider a constant set of parameters over a large time span, possibly due to
changes in the wildfire regime. Thus, after postulating that a periodical re-
run of the parameter estimation procedure may increase the effectiveness of
estimation of the fire fronts, we performed a second round of tests by
estimating parameters using measurements of the fire front between
$t_{0}=\text{15:01}$ and $t_{f}=\text{15:55}$. The objective function $J$ is
evaluated at 15:21, 15:42 and 15:55, according to (8). The trends of the
components of $p$ and of the cost function $J$ during the optimization
procedure are shown in Fig. 16. Optimization was performed in $7655$ s. The
final estimate of the parameter vector $p$ obtained by the optimization
algorithm is $\hat{\mbox{\boldmath$p$}}_{2}=$ (3.97, 0.0135, 0.0321, 0.0223,
0.0054, 0.0994, 0.1168, 0.0840), while the optimal cost is equal to 958. The
value for the ratio $r$ in (10) is equal to $0.0939$. The verify the
effectiveness of the estimation, we let
$\mbox{\boldmath$p$}=\hat{\mbox{\boldmath$p$}}_{2}$ and simulate the wildfire
evolution from 15:01 to 15:55. Fig. 17 shows the obtained fire front shapes.
If we compare these results with those depicted in Fig. 15(b), we notice that
the fire fronts obtained with $\hat{\mbox{\boldmath$p$}}_{2}$ are
significantly more similar to the real fronts, with no over-estimation of the
fire evolution towards North-East. The larger values of the performance
indexes (11), (12), and (13) in Table 2 confirm this result from a
quantitative point of view. Therefore, we conclude that a periodical re-run of
the parameter estimation can improve the simulation accuracy, taking into
account possible changes of the fire regime.

##  5 Conclusions

In this paper, we have addressed the problem of modeling fire propagation
using level set methods to describe the front evolution over space and time
and an empirical model for the rate of spread. The model is composed of
several parameters depending on vegetation and characteristics of terrain,
which have been estimated by means of a procedure based on the minimization of
a least-squares cost accounting for the difference between measured and
predicted fronts. The proposed combination of modeling and identification has
been tested with both simulated data and real measurements of propagating
fronts. In particular, we have considered data of the 2002 Troy fire, and we
have obtained optimal estimates of the parameters of the propagating model
that have turned out to be effective in predicting the future evolution of the
fire front.

The successful results we have obtained suggest several future directions of
investigation. The first one is the use of other optimization techniques and
the adaptation over time of model parameters, which may be well-suited to
predicting fire propagation for the purpose of decision support to allocate
the firefighting resources and contain the fire. Moreover, we will extend our
approach to more complex models of the rate of spread than the one proposed in
[14, 36]. Another important subject of future effort will be the validation of
the considered fire propagation models with estimated parameters. Thus,
according to [47], we will test the identified models without variations of
some parameters against several real case studies in order to better detect
case-dependent and case-independent parameters as well as quantify the quality
of the overall estimation.

## References

  * [1] H. P. Hanson, M. M. Bradley, J. E. Bossert, R. R. Linn, L. W. Younker, The potential and promise of physics-based wildfire simulation, Environmental Science & Policy 3 (4) (2000) 161–172.
  * [2] R. Rothermel, A mathematical model for predicting fire spread in wildland fuels, Vol. INT-115, Intermountain Forest & Range Experiment Station, 1972.
  * [3] R. Burgan, R. Rothermel, BEHAVE: Fire behavior prediction and fuel modeling system - FUEL subsystem, Tech. rep., National Wildfire Coordinating Group (1984).
  * [4] T. L. Clark, M. A. Jenkins, J. Coen, D. Packham, A coupled atmosphere-fire model: Convective feedback on fire-line dynamics, Journal of Applied Meteorology 35 (6) (1996) 875–901.
  * [5] M. Mutlu, S. Popescu, K. Zhao, Sensitivity analysis of fire behavior modeling with LIDAR-derived surface fuel maps, Forest Ecology and Management 256 (2008) 289–294.
  * [6] I. Karafyllidis, A. Thanailakis, A model for predicting forest fire spreading using cellular automata, Ecological Modelling 99 (1) (1997) 87–97.
  * [7] I. A. French, D. H. Anderson, E. A. Catchpole, Graphical simulation of bushfire spread, Mathematical and Computer Modelling 13 (12) (1990) 67–71.
  * [8] T. Ghisu, B. Arca, G. Pellizzaro, P. Duce, An improved cellular automata for wildfire spread, Procedia Computer Science 51 (1) (2015) 2287–2296.
  * [9] L. Hernández Encinas, S. Hoya White, A. Martín del Rey, G. Rodríguez Sánchez, Modelling forest fire spread using hexagonal cellular automata, Applied Mathematical Modelling 31 (6) (2007) 1213–1227.
  * [10] A. Alexandridis, D. Vakalis, C. Siettos, G. Bafas, A cellular automata model for forest fire spread prediction: The case of the wildfire that swept through Spetses Island in 1990, Applied Mathematics and Computation 204 (2008) 191–201.
  * [11] M. Tonini, M. D’Andrea, G. Biondi, S. D. Esposti, A. Trucchia, P. Fiorucci, A machine learning-based approach for wildfire susceptibility mapping. The case study of the Liguria Region in Italy, Geosciences 10 (2020) 1–18.
  * [12] D. Muñoz-Esparza, B. Kosović, P. A. Jiménez, J. L. Coen, An accurate fire-spread algorithm in the weather research and forecasting model using the level-set method, Journal of Advances in Modeling Earth Systems 10 (4) (2018) 908–926.
  * [13] A. Mentrelli, G. Pagnini, Modelling and simulation of wildland fire in the framework of the level set method, Ricerche di Matematica 65 (2016) 523–533.
  * [14] V. Mallet, D. Keyes, F. Fendell, Modeling wildland fire propagation with level set methods, Computers & Mathematics with Applications 57 (7) (2009) 1089–1101.
  * [15] A. Bova, W. Mell, C. Hoffman, A comparison of level set and marker methods for the simulation of wildland fire front propagation, International Journal of Wildland Fire 25 (2015) 229–241.
  * [16] I. Kaur, A. Mentrelli, F. Bosseur, J.-B. Filippi, G. Pagnini, Turbulence and fire-spotting effects into wild-land fire simulators, Commun. Nonlinear Sci. Numer. Simul. 39 (2016) 300–320.
  * [17] J. Sethian, Level Set Methods and Fast Marching Methods, 2nd Edition, Vol. 3, Cambridge University Press, Cambridge, 1999.
  * [18] J. Sethian, P. Smereka, Level set methods for fluid interfaces, Annual Rev. Fluid Mechanics 35 (1) (2003) 341–372.
  * [19] S. Osher, R. Fedkiw, Level set methods and dynamic implicit surfaces, Vol. 153 of Applied Mathematical Sciences, Springer-Verlag, New York, 2003.
  * [20] R. Kimmel, Numerical Geometry of Images, Springer-Verlag, New York, 2004.
  * [21] M. Bernauer, R. Herzog, Optimal control of the classical two-phase Stefan problem in level set formulation, SIAM J. Sci. Comput. 33 (1) (2011) 342–363.
  * [22] I. Yang, C. Tomlin, Regularization-based identification for level set equations, in: Proc. IEEE Conf. on Decision and Control, 2013, pp. 1058–1064.
  * [23] A. Alessandri, P. Bagnerini, M. Gaggero, Optimal control of propagating fronts by using level set methods and neural approximations, IEEE Trans. Neural Networks and Learning Systems 30 (3) (2019) 902–912.
  * [24] A. Alessandri, P. Bagnerini, R. Cianci, M. Gaggero, Optimal propagating fronts using Hamilton-Jacobi equations, Mathematics 7 (11) (2019) 1–10, art.no. 1122\.
  * [25] A. Alessandri, P. Bagnerini, M. Gaggero, A. Rossi, State and observer-based feedback control of normal flow equations, Automatica 117 (2020) 1–7, art. no. 108980.
  * [26] D. P. Bertsekas, Nonlinear programming, 3rd Edition, Athena Scientific, 2016.
  * [27] C. Lautenberger, Wildland fire modeling with an Eulerian level set method and automated calibration, Fire Safety Journal 62 (2013) 289–298.
  * [28] M. Ambroz, K. Mikula, M. Frastia, M. Marcis, Parameter estimation for the forest fire propagation model, Tatra Mt. Math. Publ. 72 (2018) 1–22.
  * [29] A. Ervilha, J. Pereira, J. Pereira, On the parametric uncertainty quantification of the Rothermel’s rate of spread model, Applied Mathematical Modelling 41 (2017) 37–53.
  * [30] C. Zhang, M. Rochoux, W. Tang, M. Gollnera, J.-B. Filippi, A. Trouve, Evaluation of a data-driven wildland fire spread forecast model with spatially-distributed parameter estimation in simulations of the FireFlux I field-scale experiment, Fire Safety Journal 91 (2017) 758–767.
  * [31] M. Rochoux, R. S, D. Lucor, B. Cuenot, A. Trouve, Towards predictive data-driven simulations of wildfire spread - Part I: Reduced-cost ensemble Kalman filter based on a polynomial chaos surrogate model for parameter estimation, Nat. Hazards Earth Syst. Sci. 14 (2014) 2951–2973.
  * [32] J. Mandel, L. Bennethum, J. Beezley, J. Coen, C. Douglas, M. Kim, A. Vodacek, A wildland fire model with data assimilation, Mathematics and Computers in Simulation 79 (2008) 584–606.
  * [33] C. Zhang, A. Collin, P. Moireau, A. Trouve, M. Rochoux, State-parameter estimation approach for data-driven wildland fire spread modeling: Application to the 2012 RxCADRE S5 field-scale experiment, Fire Safety Journal 105 (2019) 286–299.
  * [34] T. Artes, A. Cencerrado, A. Cortes, T. Margalef, D. Rodriguez-Aseretto, T. Petroliagkis, J. San-Miguel-Ayanz, Towards a dynamic data driven wildfire behavior prediction system at European level, Procedia Computer Science 29 (2014) 1216–1226.
  * [35] M. Denham, K. Wendt, G. Bianchini, A. Cortes, T. Margalef, Dynamic data-driven genetic algorithm for forest fire spread prediction, Journal of Computational Science 3 (2012) 398–404.
  * [36] S.-E. Lo, A Fire Simulation Model for Heterogeneous Environments Using the Level Set Method, Ph.D. thesis, Claremont Graduate University (2012).  doi:10.5642/cguetd/72.
  * [37] C. Zhai, S. Zhang, Z. Cao, X. Wang, Learning-based prediction of wildfire spread with real-time rate of spread measurement, Combustion and Flame 215 (2020) 333–341.
  * [38] J. Hilton, C. Miller, A. Sullivan, C. Rucinski, Effects of spatial and temporal variation in environmental conditions on simulation of wildfire spread, Environmental Modelling & Software 67 (2015) 118–127.
  * [39] J. Hilton, C. Miller, J. Sharples, A. Sullivan, Curvature effects in the dynamic propagation of wildfires, International Journal of Wildland Fire 25 (2016) 1238–1251.
  * [40] L. Arroyo, C. Pascual, J. Manzanera, Fire models and methods to map fuel types: The role of remote sensing, Forest Ecology and Management 256 (2008) 1239–1252.
  * [41] I. M. Mitchell, The flexible, extensible and efficient toolbox of level set methods, Journal of Scientific Computing 35 (2-3) (2008) 300–329.
  * [42] P. Huabsomboon, Implicit level set methods for a fire spread model, in: Proc. Asian Simulation and Modeling Conf., 2007, pp. 127–133.
  * [43] E. H. Kim, K. Yang, M. N. Tran, Thermal-image-based wildfire spread simulation using a linearized model of an advection–diffusion–reaction equation, Simulation 88 (9) (2012) 1093–1115.
  * [44] Pacific Southwest Research Station, The USDA Website Fire Science (2007).
URL http://www.fireimaging.com/fires/2002/troy/index.html

  * [45] U.S. Department of Agriculture Forest Service and U.S. Department of the Interior, LANDFIRE, landscape fire and resource management planning tools (2019).
URL https://www.landfire.gov

  * [46] J.-B. Filippi, V. Mallet, B. Nader, Evaluation of forest fire models on a large observation database, Nat. Hazards Earth Syst. Sci. 14 (2014) 3077–3091.
  * [47] D. Morvan, Validation of Wildfire Spread Models, Springer International Publishing, Cham, 2019, pp. 1–6.