A Simple Model for Wildland Fire Vortex–Sink Interactions

: A model is developed to explore ﬁre–atmosphere interactions due to the convective sink and vorticity sources in a highly simpliﬁed and idealized form, in order to examine their effect on spread and the stability of various ﬁre front geometries. The model is constructed in a cellular automata framework, is linear, and represents a background ﬂow, convective sink, and vortices induced by the ﬁre plume at every burning cell. We use standard techniques to solve the resulting Poisson equations with careful attention to the boundary conditions. A modiﬁed Bresenham algorithm is developed to represent convection. The three basic ﬂow types—large-scale background ﬂow, sink ﬂow, and vortex circulation—interact in a complex fashion as the geometry of the ﬁre evolves. Fire-generated vortex–sink interactions produce a range of ﬁre behavior, including unsteady spread rate, lateral spreading, and dynamic ﬁngering. In this simpliﬁed framework, pulsation is found associated with evolving ﬁre-line width, a ﬁre-front acceleration in junction ﬁres, and the breakup of longer initial ﬁre lines into multiple head ﬁres. Fuel is very simply represented by a single burn time parameter. The model fuel is uniform yet patchiness occurs due to a dynamic interaction of diffusive and convective effects. The interplay of ﬁre-induced wind and the geometry of the ﬁre front depends also on the fuel burn time.


Introduction
Sophisticated computational fluid dynamics and combustion models represent a myriad of effects and feedbacks inherent in the complex problem of wildland fire spread. Chemical and thermodynamic combustion processes, fuel structure and composition, terrain, and atmospheric boundary layer dynamics and thermodynamics are all ingredients to such models, depending on the situation being investigated. Canopy flow is recognized as a turbulent mixing layer with large eddies controlling vertical fluxes [1]; small fires within a canopy are subject to this control to a large degree, while more intense fires develop their own boundary layer flow. In simpler conditions, for instance, a grass fuel layer with no canopy, a different set of boundary layer turbulent structures have been used to describe interactions between boundary layer flow and fire spread [2]. A common thread is the convergence and divergence effects of buoyant convection and large eddies on nearsurface flow. Near-field turbulence due to buoyant plume effects develops from Kelvin-Helmholtz shear instabilities on the plume itself and is tightly linked to the surrounding boundary layer structure [3], giving rise to horizontal eddy flow at the surface. Many of these dynamically evolving turbulent structures are now being observed in wildland fires [4]. Here, our focus is on highly simplified models, with idealized representations of combustion rate and two key components of local, near-field flow-the convective inflow and vertical vorticity-in order to gain some insight into the larger-scale structure of wildland fire spread in the presence of local convectively driven inflow and vorticity sources.
Two important coupled fire-atmosphere effects occur even in one dimension. Babak et al. [5] studied 1D solutions of the advection-diffusion equation for heat in a wildland fire context, with an energy conserving approach. Fire-wind interaction is reflected, in this 1D case, by the heat advection term and a wind-dependent rate of combustion. They derived analytic traveling wave solutions to the spread rate in terms of the basic parameters. A characteristic form of spread rate versus wind speed was found for backing-fire spread rate and for head-fire spread rate. In this framework, direct coupling between the wind and the fire occurs in the combustion term, with stronger wind producing greater combustion rates. Although no new effects were found with this coupling, rates of spread were quite sensitive to the coupling, especially backing fire spread rates. An insightful 1D model was developed and applied by Morandini et al. [6] that coupled the vertical wind generated by the fire plume to the horizontal wind in the fuel bed. Note that the background wind used in this and similar models is scaled down to a much smaller value representing the heat transfer in the fuel bed. By representing the vertical flux of heat and mass in buoyant hot gases, they produced a corresponding mass convergence in the horizontal inflow. This vertical flow substantially modified the horizontal wind. Bebieva et al. [7] incorporated a similar divergent inflow, which they called the plume model component of their 1D model. They showed that the plume exerted a key control on the forward spread rate by setting the effective wind speed at the head of the fire.
Moving from 1D to 2D and 3D produces greater complexity in terms of fire-atmosphere interactions than the simple sink or plume effect. In 2D and 3D, the geometry and intensity of the fire will shape the wind field and produce fundamental effects on fire spread [8][9][10][11]. These pioneering numerical and experimental studies showed that the development of conical fire front shape is due to the location of convective updraft and convergence of flow at the fire line, and can produce single or multiple shapes depending on the length of the initial fire line. Moreover, at smaller scales, and for intense fires, a dynamic enhancement of fire front fingering can occur. This effect arises from the coupling of horizontal and vertical vorticity near the front. Both effects derive from the link between the convective columns above the fire and the fire front, and Clark et al. [8,12] and Miller et al. [13] infer the strength of the coupling as a function of a Froude number based on background wind, buoyancy, and fire-line width.
Coupled fire-atmosphere effects are embedded in the larger atmospheric boundary layer (ABL), itself a complex turbulent system, with various instabilities operating at many scales. The coupled effect of fire-induced flow and turbulence in the atmospheric boundary layer was investigated by Sun et al. [2]. Two fundamental ABL processes, buoyant or convective instability, and roll instabilities were studied. Interaction within a convectively unstable ABL was shown to strongly increase the fire front rate of spread (RoS), as did fire ignition line length. Convectively induced turbulent downdrafts behind the fire line are also responsible for increased RoS [2]. In a turbulent flow, the timing of these effects is not predictable and the random nature of interactions needs to be represented in models.
Another example of local fire-atmosphere interaction involves the development of small-scale streamwise vortices, and the along-front modulation of the fire by these vortices. This can lead to increases in the rate of spread through convective preheating by cross-front air flow [14]. Their results also suggest a competition between the head and lateral flank of the fire for upstream air, due to flow divergence over the burned area behind the fire.
Recent work by Thomas et al. [15] shows the small-scale structure of fine-scale firewind interactions on evolving fire fronts using the Weather Research and Forecasting (WRF) model and a simple fire combustion model. Their model also explicitly calculates the 3D flow field, with updraft and associated surface flow due to the fire heat source (they term this the pyro-convective process), and show that the resulting convergence exerts a strong control on fire line behavior. The role of both convectively generated inflow and the vorticity induced by fire-line structure are investigated, and the flow is shown to develop vortex pairs that enhance the forward advancement of fire-line perturbations.
While some of the larger-scale consequences of fire-atmosphere interactions are well known, observing the detailed flow interactions at the fire front scale is a challenging problem. Clements et al. [16] provide direct evidence for the influence of fire-induced winds on the movement of the fire front, depending on the intensity of the fire. The evolution of fire front geometry was suggested by Clements et al. [16] to lead to a positive feedback between the fire and background winds associated with greater inflow and spread rate.
Capturing some of these effects in models is another type of challenge. Simplified, reduced-physics, numerical models offer the advantage of isolating physical processes and rapid computations. The coupling between surface fire and wind on larger scales involves interactions between fire-generated convection and frontal motion (e.g., Hilton et al. [17]), background and induced local vorticity, and the tilting of the surrounding horizontal vorticity field. Sharples and Hilton [18] isolated and described such an effect in their model of flow over a ridge, with the associated vorticity-induced lateral fire spread. In this mechanism, the vorticity of the background flow is lifted by the fire plume, tilted, and modifies the horizontal flow. Our model, developed below, builds on these idealized studies and generalizes some of the same effects. Despite the contrived appearance of idealized models, they are an important and helpful step in the development of better physical understanding of fire behavior.
Our goal is to extract a minimal set of essential coupled dynamics in an idealized linear 2D model. Specifically, we will investigate the role of plume-generated sinks and vorticity on the evolution and geometry of the fire front. The stability of the flank fire, development of fingering, and dependence on basic parameters will be explored. Complex fire behavior is seen to arise owing to interactions of fire geometry and simple representations of fluid dynamic elements, leading to intricate wind flow patterns near the fire front and convective control on fire spread.

2D Model
In this section, we develop a highly idealized 2D model to investigate fire-atmosphere interaction. The key physical atmospheric processes of flow divergence and vorticity are represented in a simplified but physically meaningful context. The main ingredients are a background uniform wind, a potential flow, and a vortical flow. The potential flow is due to wind sink on complex fire front geometry, and the vortical flow is due to lifting of horizontal background vortex lines by the buoyant plume and local horizontal shear induced by the fire's buoyancy source, analogous to shear flow around a bluff object. A simple combustion term with plausible burn times is used, and turbulent diffusion is incorporated representing atmospheric eddies colliding with the fuel bed.
Potential flow has been used to infer the flow field and convergence of fire lines and plumes (see [9,19] and references therein). Our model incorporates the potential flow and solenoidal flow similar to Sharples and Hilton [18]. However, unlike their model, which used a level-set method for a specific ridge geometry, we solve for the flow generated by each burning elemental cell in a cellular automata (CA) model of combustion. Thus, we apply the notion of a tilting vortex to small-scale flow around individual fire elements, which sum to produce the total effect on the fire front.

Solution Method
We discretize the domain with a two-dimensional rasterized grid. Each cell is assigned one of three states: (1) the cell is actively burning; (2) the cell is not burning but has remaining fuel; (3) the cell has no remaining fuel. Depending on the state of the system, we define a two-dimensional velocity field, v(x), which is the sum of three components. The first component is a background velocity that we always take to be the uniform velocity U BG (0, 1). The other two velocity components are fire-induced, and following Sharples and Hilton [18], we decompose it as the sum of an irrotational and a solenoidal velocity field As was performed by Hilton et al. [17], we impose mass conservation by requiring that where w is the vertical velocity of the buoyant plume. Letting ∂w ∂z be a positive constant at combusting cells and zero otherwise results in a velocity field ∇ψ that Hilton et al. called the pyrogenic potential [17].
The solenoidal velocity field introduces vorticity, and η solves where ω S represents vorticity sources. Fire-induced vorticity is an active area of investigation, and we focus on a simplified model for vorticity that is only present in the two-dimensional plane. At the flanks of the fire, we impose a vorticity source that is clockwise at the left flanks and counterclockwise at the right flanks. These vortical forces are placed along the entire flanks, including between broken fire lines ( Figure 1). We also consider vortical flows that are counterclockwise on the left flanks and clockwise on the right flanks. This vorticity occurs when the shear direction above the fire is reversed by strong fire regions, and the result is what Clark et al. referred to as dynamic fingering [8].
Computing ψ and η requires solving the two-dimensional Poisson equation with a known forcing function. We use the standard second-order central difference formula to discretize the Laplacian, and we impose a constant Neumann boundary condition. The resulting sparse linear system is solved with Matlab's backslash operator, which applies an LU solver. For the pyrogenic Neumann boundary condition, a constant value is chosen so that This means that the normal flow into Ω is uniform, which is a valid assumption as long as the fire is sufficiently far from ∂Ω. Since the vorticity sources integrate to 0 by construction, the boundary condition for Equation (3) is Neumann and homogeneous. Once both ψ and η are calculated, we compute the resulting velocity field (1) using standard centered difference formulas. After defining the velocity field v(x), new cells must be ignited in the direction of the wind velocity. We make use of Bresenham's line algorithm [20]. Bresenham's line algorithm is a line drawing algorithm that determines cells of a rasterized geometry that are between a start and end point. We slightly modify the algorithm since we only have a start point (a burning cell) and a direction (the wind velocity). Since fire does not spread through burnt or combusting cells, we further modify Bresenham's line algorithm so that it stops igniting cells along the line once it reaches a cell that is either combusting or has no remaining fuel.

Without Vorticity Term
With Vorticity Term

Model Configuration
The configuration of the CA model is designed to represent scales associated with real wildland fires. While, in principle, the model can represent any size domain over any period of time, the nature of the representation and parameterization of fire spread will differ across the scales chosen. Here, the cell size L = ∆x = ∆y = 1 m, and the time evolution occurs over minutes-hours. The domain is rectangular, flat, and hundreds of meters on a side. Other choices are possible but this relatively small scale allows for a reasonable separation between elements mostly on fire and those that are not, as well as rapid computations (seconds on a laptop using Matlab). If, for instance, the width of the fire were much less than the grid size, a model would simply need to represent the bulk spread rate across the cell rather than the details of spread from cell to cell. If, on the other hand, the cell size were much smaller than the width of the fire, then the representation of the combustion and turbulent fluxes would be more complicated, with nondiffusive effects.
We use a simple burnout model where cells extinguish after a fixed burn time τ b . This is the only fuel parameter in the model. This time is not a function of the cell size but rather reflects, at least in an idealized form, the exponential decay time of a simple fuel evolution model dF/dt = −F/τ b . This parameter depends on the fuel and can be estimated from observations [21]; although the distributions of this quantity can be diverse, mean values of a few seconds to minutes are plausible in light to moderate fuels.
To model diffusion, we introduce a stochastic term that allows cells to ignite neighboring cells at the time step before the original cells extinguish. A cell with remaining fuel that shares an edge with a cell that is about to extinguish is ignited with probability p ig . Cells with remaining fuel that share a corner are ignited with the same probability, but scaled by √ 2 to account for the additional distance. The diffusion-based ignition occurs independently of the convection-based ignition. Model runs were made with no wind to define a reasonable range for p ig . A final choice of p ig = 0.4 was made based on a balance between the smoothness of the solutions and a desire to represent random fine-scale structure. Smaller choices provide a smooth and uniform spread, while larger choices lead to more of a random-walk character of the fire spread. As ρ ig approaches 1, the spread becomes square, reflecting the initial condition. The rate of spread under diffusion is related to the burn time (τ b ) since the model waits for the cell to burnout before igniting a neighboring cell. With the diffusion D = p ig ∆x 2 /∆t, and a diffusive spread rate of √ D/τ b , then in dimensional terms RoS ∼ p 1/2 ig ∆x/τ b . Indeed, the rate of spread is relatively insensitive to p ig in the range we explored.
We can use basic plume scaling to provide guidance for the representation of convection across the CA grid. Each cell, when on fire, produces a buoyancy of b = −g∆ρ/ρ, where g is gravity and ∆ρ/ρ is the density change (negative) due mainly to the high temperatures of the hot gases. The net buoyancy source of each cell is B c = bL 2 , and for L = 1 m, this is about 10-100 m 3 /s 2 . Fuel burn time together with the cell buoyancy source B c provides a scale for the buoyant plume velocity w b = (B c /τ b ) 1/3 . With a choice of 10 s for burn time, this scale is about w b = 1-2 m/s. A fundamental parameter is the ratio of mean background wind U BG to the plume velocity U BG /w b . Horizontal wind u s is also generated by the fire itself, due to entrainment, in proportion to the vertical velocity scale w b . The associated inflow toward the plume can be a large perturbation of the background flow (see e.g., [22]). This inflow, due to the convective updraft, is modeled here as the potential sink flow.
The convection step represents the advection of hot gases over cells due to wind. This effect is stronger for higher wind relative to the plume velocity. The case U BG w b implies that the flames lie more or less flat along the ground for a large distance leading to a rate of spread that can approach a large fraction of the wind velocity itself in light fuels. Here, in order to represent the convection effect, the fire is carried downwind a distance proportional to the total wind velocity as defined above-the sum of background, potential, and vorticity sources. Thus, a simple form for this distance is L f = c|u total |∆t. The convection parameter c represents the net effect of wind in the fuel bed. Since the time step ∆t also satisfies a consistency relation ∆t < L/w b , and c is order 1, the distance L f is typically not much greater than L over a single time step. Values of c between 0.25 and 1 were investigated, with c = 0.75 selected as the reference value. Note that, due to the sink and vorticity flows, the total wind velocity near the fire is generally quite different from the background flow (see Figure 1); the presence of sink and vorticity sources disrupts the simple relation that c has on the overall rate of spread when only background flow is present. No jumps are allowed, so the ignition stops if a burning cell or a burnt cell is encountered in the direction of the velocity vector. This constraint could be relaxed to allow for spotting.
Finally, the burn time τ b in some runs is increased or decreased according to the strength of the background wind U BG , to represent a linear wind-dependent combustion rate or an inversely dependent decay time scale, τ b (w b /U BG ). Thus, a higher background wind speed will reduce the burn time and increase the rate of spread.
Vorticity sources from both background shear (positive) and a fire-induced reversed wind shear (negative) are considered. An estimate of the strength of a background vorticity source relative to the convective sink term can be derived by separating the vorticity equation into a background and perturbation and approximating the vorticity equation source as the background. This represents the source due to the effect of the buoyant plume, lifting the background horizontal vorticity filaments to generate a vertical component [8,18], η ∼η, whereη is the background horizontal vorticity.
Sharples and Hilton [18] used this approach to derive an explicit idealized form for a background flow over a ridge and a vorticity source term parallel to the ridge. If the wind near the surface varies over the height scale H, the background shear U BG /H with or without a ridge may be taken as the scale for the vertical component of the vorticity source lofted by the plume. Then, the approximate vorticity balance produces a scale for the ratio of the vorticity to sink source terms where η is the vertical component of vorticity, is the aspect ratio H/L, H is the boundary layer shear height scale, and L is the horizontal plume width scale. The boundary layer scale H and plume entrainment scale are taken to be similar. In our runs, the scaled vorticity source varies from zero to about one. However, the strength of the vorticity source and sink are varied separately to investigate their distinct effects.

Results
Initial results are obtained with head-fire simulations without any fire-induced vorticity. These produce well-known spread patterns but, in addition, the model provides a useful, simplified framework for flow conditions surrounding the fire. We then consider junction fires with several different opening angles. The distribution of first arrival time and associated rate-of-spread reveal detailed fire-atmosphere interaction due to the evolving fire line geometry and convective sink. Finally, we consider the effects of both positive and negative vorticity on a head fire.
The Two important parameters are the strengths of the fire-induced sinks and vorticity. For each example, we report the dimensionless quantity where U BG is the background speed. The scalars ψ and η satisfy (2) and (3), respectively. The maximum is taken over all spatial points for an initial condition corresponding to a spot fire.

Head Fire
We start by considering a head fire with four different background speeds and a fixed fire-induced sink velocity. The simulated wind speeds and dimensionless variables are reported in Table 1, and the burn time is τ b = 15 s. A fire line with width 1 m is initialized at y = 34 m and spans the middle 66 m of the computational domain. In Figure 2, we plot the resulting fire line at t = 48 s for three of the background speeds. The blue curves are streamlines induced by the background and fire-induced velocities. At low wind speeds, the streamlines cannot pass through the combustion region, and the result is that the fire propagates mostly under diffusion. In contrast, at larger wind velocities, the streamlines can travel through the combustion region, and the fire propagates mostly under convection. The streamlines provide information not just about the flow but also the flow geometry: they illustrate the concept of a flow "barrier" around the fire at the surface that divides the air entering the plume from the air that makes it around the fire. The region of entrainment can extend around the fire some distance downstream due to the convective sink. A concentration of streamlines can occur when the wind overcomes the sink and is able to focus the flow (hence smoke) on the downwind side in a narrow stream (e.g., 8 m/s case, Figure 2). Note that, although streamlines are also curves of the constant stream function in this quasi-steady limit, they do not represent a uniform contour interval because some streamlines end at a sink.  Table 1 and the burnout time is τ b = 15 s.
In Figure 3, we plot a heat map of the first arrival time. Then, in Figure 4a, we plot the first arrival time along the line x = 100 m, which is in the center of the ignition pattern. Given this first arrival time t(y), we can define a rate of spread as RoS(y) = d dt (y(t))) = 1 t (y) (7) and this is plotted in Figure 4b. Note that a small amount of smoothing is applied to both t(y) and RoS(y). The head-fire cases show rates of spread that fall within the experimental values determined from field experiments by Johnston et al. [23] and Paugam et al. [24]. Model rates of spread show, in addition, a plausible dependence on wind speed according to accepted empirical models and observations (see e.g., Cruz and Alexander [25]). Our model produces periods of faster and slower rates of spread, and two such instances are indicated by the black dots in Figure 4b. The configurations of the system at these locations are plotted in Figure 5 (slow rate of spread) and Figure 6 (fast rate of spread). In the slow case, the fire line near x = 100 m is sufficiently wide that the convective sink does not allow heat to be advected away from the fire, as is indicated by the streamlines, and this persists for 15 s. The resulting rate of spread is about 0.12 m/s. It is only once the fire line thins due to burnout that the convective sink weakens, and the head is finally able to propagate forward. In contrast, in the fast case, the fire line near x = 100 m is thinner, and the result is that the fire propagates two cells forward in 6 s, resulting in a larger rate of spread of 0.14 m/s. The rate of spread may also be described as a time series following the fire front, and time-series analysis becomes applicable. For the central line in a head-fire case, the timeseries look similar to the results versus distance (Figure 4b), and periodograms suggest multiple-frequency components. If instead the focus is on the fastest moving part of the front, this may not occur along a line, and different techniques are needed. The spread rate (7) may be written in a vector form to map locations of rapid spread. The more complex ignition and spread patterns produce highly variable rates of spread, and integral measures become necessary to describe the overall progress of the fire.
When the convective sink is allowed to vary with the wind, as it does in real fires due to a wind-dependent combustion reaction rate, an important additional feedback occurs. We investigated this feedback with a simple proportionality between the sink term and background wind. This can be done directly, or via the burn time τ b , with the sink proportional to the inverse burn time. We chose the former to isolate the effect; some runs were also made in the latter-but in this case, there is not a single time scale for combustion and the interpretation in terms of this scale is more difficult. The effects described for the fixed sink case are still present, but the "breakthrough" of the background wind across the fire front producing fast spread or "runs" is a more sensitive function of the background wind, since the sink inflow increases along with the background wind. Thus, the rate of spread at higher wind speed can be restrained by the intense inflow due to the sink term. In most cases, the breakthrough events were not regular, as might be anticipated in the fixed burn time cases, but intermittent, producing strongly nonstationary behavior.

Junction Ignition Patterns
We next consider a junction ignition pattern with three different opening angles: 30 • , 60 • , and 90 • . In Figure 7, we compute the first arrival for each of these ignition patterns. In all three examples, the background wind is U = 4 m/s, U sink = 3.63 × 10 −2 , U vort = 0, and the burn time is τ b = 15 s. By applying Equation (7), we use the first arrival time to define a rate of spread along the line that bisects the two flanks, which is also parallel to the background wind direction ( Figure 8). As expected, we see that smaller initial junction angles result in faster rates of spread.
These results are consistent with observations of junction fire merging and enhanced rates of spread from field experiments [26], suggesting that the idealized model plausibly represents junction fire behavior. Both the typical magnitude of spread and the spread acceleration in time are obtained, with a strong dependence on opening angle.  Figure 7. (b) The rate of spread, defined using Equation (7), along the vertical line bisecting the ignition pattern.

Lateral Spread in the Presence of Positive Vorticity
We examine the effect of vorticity on a head fire, starting with the effect of positive vorticity in a weak convective sink environment. For each background wind speed U BG , the dimensionless sink and vortical velocities are U sink = 3.63 × 10 −3 and U vort = 3.63 × 10 −1 , respectively. Therefore, the convective sink velocity is 100 times smaller than the vortical velocity. For all simulations, we use the burn time of τ b = 15 s. As discussed in Section 3.1, a head fire is ignited at y = 34 m, and it is initially 66 m-wide and 1 m-thick. In Figure 9, we plot the first arrival time for four different wind speeds. This first arrival time is calculated along the horizontal line going through the ignition line, and is plotted in Figure 10a. Even though the spread includes a random term that breaks symmetry, the rates of spread of the two flanks are comparable. Therefore, we only plot the lateral rate of spread of the right flank in Figure 10b.
We observe that introducing positive vorticity promotes lateral spread. In the case of a stronger convective sink, a significant lateral spread is not observed due to the restraining effect of inflow. A weak sink might be expected in patchy or wet fuel conditions, or as a fire burns down. The rate of spread increases as wind increases, with greater variability at higher wind speed.  Heat maps of the first arrival time are in Figure 9. (b) The rate of spread, defined using (7), along the horizontal line passing through the ignition pattern. The rate of spread is only reported for the right-moving flank, and the left-moving flank has a similar trend.

Formation of Heads by Negative Vorticity
As a final example, we consider a head fire that is sufficiently intense that the background vorticity field is set by the inflow or the dynamic fingering case. Therefore, the induced vortex flow rotates counterclockwise on left flanks and clockwise on right flanks.
In this situation, the vorticity promotes fingering, some of which develop into new head fires. Using our model, we are able to observe this behavior under various conditions. Somewhat heavier fuels and a strong sink so that fingers do not rapidly diverge from the main fire line destabilize the front and promote multiple head fires. In the case of a moderate wind speed U BG = 4 m/s, a strong negative vortical velocity U mathrmvort = 4.54 × 10 −1 , and a sink velocity U sink = 7.27 × 10 −2 , the head breaks into roughly four smaller head fires ( Figure 11). Here, we also increase the burn time to τ b = 25 s, and with the presence of the random diffusion term, this increases the probability of a finger forming. Finally, the length of the initial fire line is increased to 135 m so that additional heads have the potential to form. In Figure 11, we plot the first arrival time for the first 300 s. We see that at early times, near t = 100 s, the initially uniform fire line splits into about six separate head fires (green contour) that merge to form four major head fires at later times. Figure 11. The first arrival time, in seconds, of a head fire. The background wind is U BG = 4 m/s, the sink velocity is U sink = 7.27 × 10 −2 , the vortical velocity is U vort = 4.54 × 10 −1 , and the burn time is τ b = 25 s. The negative vorticity results in dynamic fingers that then turn into separate head fires. The green curve is the fire line at t = 100 s.

Discussion and Conclusions
We have developed a simplified model to explore fire-atmosphere interactions, and their effect on spread. The model couples the fire and atmosphere by parameterizing a fire-induced velocity field that includes both convective sinks and vortices. This velocity field is coupled to the fire spread by using a modification of Bresenham's line algorithm. Finally, turbulent diffusion is included with a stochastic term, and a constant burn time is used as a highly idealized fuel model.
Using the model, we investigate several important feedback effects. Foremost among these is the feedback between the fire-induced flow and the spread of the fire front into complex shapes, producing a new flow field and altering the spread. The spread is quantified as a local rate of spread that is defined in terms of the first arrival time. These rates of spread are not constant and display multiple time scales of variation. We also demonstrate that our model reproduces known qualitative behaviors, including dynamic fingering and the development of multiple heads from an initially uniform head fire. It is of interest that discontinuous fire fronts are realizable even in this simplified model. We also note that a dynamical burn patchiness arises in many of the runs. Since our fuel model is homogeneous, this is not due to patchiness in fuels (e.g., [27]) or moisture. Diffusion is necessary for this flow-generated patchiness, but its extent also depends on the sink and vorticity-induced spread effects.
A large simplification of our model is that we keep the background wind constant. This is not necessary, and both time and spatially variable background flow are possible-for instance, to represent large-scale gusts or topographic effects. The underlying assumption here is that the components of the wind generated by the fire evolve on a faster time scale than the background wind, the key parameter being the ratio of mean background wind U BG to the plume velocity w b . The orientation of the vortices would need to change with a variable background wind direction, since their generation mechanism is linked to background shear. Moreover, the model is linear; hence, there is no interaction between the different components of the wind in terms of their dynamics. In a nonlinear model, a major effect would be the advection of vorticity; hence, vortex filaments would become concentrated by the sink. On the other hand, a nonlinear model would have to halt the nonlinear steepening or vorticity concentration at some point with diffusion; the vorticity effect, as modeled here, can be viewed as an idealized representation of the larger-scale fire-induced vortex flow elements in that limit.
Fuel effects are not explored here but deserve greater emphasis. The simple fuel burnout time used here is not adequate to represent the vast range of interactions that a three-dimensional fuel structure can impose on the boundary. Porosity and wind speed reduction factors, wake effects, open and closed canopies, and moisture content all provide important factors in fire evolution.
In the future, we expect the model to be of use for research, educational, and training purposes to explore the basic physical processes operating in the surface fire environment. A growing approach to fire prediction involves integrating data with physical models. Without physical models, statistical or machine learning approaches tend to be highly specialized to particular environments or situations and require large amounts of data to learn patterns. Even simple physical models can provide crucial constraints and physically based parameters to improve fire prediction.