Rothermel Fire Spread Model

Mathematical foundation and implementation

The Rothermel surface fire spread model (Rothermel 1972) is the foundation of fire behavior prediction in Elmfire.jl. This page describes the mathematical equations and how they are implemented.

using Elmfire
using Plots

Overview

The Rothermel model predicts the rate of spread of a surface fire based on:

  • Fuel properties: Loading, particle size, moisture content, heat content
  • Environmental factors: Wind speed and terrain slope
  • Physical processes: Heat transfer, combustion chemistry, and energy balance

The model outputs several key fire behavior metrics:

Output Symbol Units Description
Rate of Spread \(R\) ft/min Fire front velocity
Reaction Intensity \(I_R\) BTU/ft²/min Heat release rate per unit area
Heat per Unit Area \(H_A\) BTU/ft² Total heat released
Fireline Intensity \(I\) BTU/ft/s Byram’s fireline intensity

The Rothermel Equation

The fundamental equation for rate of spread is:

\[ R = \frac{I_R \xi (1 + \phi_w + \phi_s)}{\rho_b \varepsilon Q_{ig}} \]

Where:

Symbol Description
\(I_R\) Reaction intensity (BTU/ft²/min)
\(\xi\) Propagating flux ratio
\(\phi_w\) Wind factor
\(\phi_s\) Slope factor
\(\rho_b\) Fuel bed bulk density (lb/ft³)
\(\varepsilon\) Effective heating number
\(Q_{ig}\) Heat of pre-ignition (BTU/lb)

Reaction Intensity

Reaction intensity represents the rate of heat release:

\[ I_R = \Gamma' w_n h \eta_M \eta_s \]

Where: - \(\Gamma'\) = optimum reaction velocity (1/min) - \(w_n\) = net fuel loading (lb/ft²) - \(h\) = heat of combustion (BTU/lb) - \(\eta_M\) = moisture damping coefficient - \(\eta_s\) = mineral damping coefficient

Moisture Damping

The moisture damping coefficient reduces fire spread as fuel moisture increases:

\[ \eta_M = 1 - 2.59r + 5.11r^2 - 3.52r^3 \]

Where \(r = M/M_{ex}\) is the ratio of fuel moisture to moisture of extinction.

# Visualize moisture damping
r_values = 0:0.01:1
eta_m = [Elmfire.moisture_damping(r) for r in r_values]

plot(r_values, eta_m,
    xlabel = "Moisture Ratio (M/Mₑₓ)",
    ylabel = "Moisture Damping (ηₘ)",
    title = "Rothermel Moisture Damping Coefficient",
    linewidth = 2,
    legend = false,
    ylims = (0, 1.1)
)

Wind and Slope Factors

Wind accelerates fire spread through convective heat transfer:

\[ \phi_w = C \left(\frac{\beta}{\beta_{op}}\right)^{-E} U^B \]

Slope accelerates uphill spread through radiant heat transfer:

\[ \phi_s = 5.275 \beta^{-0.3} \tan^2\theta \]

Where: - \(\beta\) = packing ratio (dimensionless, fuel volume / bed volume) - \(\beta_{op}\) = optimal packing ratio (dimensionless) - \(U\) = midflame wind speed (ft/min) - \(\tan^2\theta\) = squared tangent of slope angle (dimensionless) - \(B\), \(C\), \(E\) = empirical coefficients (dimensionless, functions of SAV ratio)

Implementation in Elmfire.jl

SpreadResult Type

The surface_spread_rate function returns a SpreadResult struct:

struct SpreadResult{T<:AbstractFloat}
    velocity::T       # Spread rate (ft/min)
    vs0::T            # Base spread rate without wind/slope (ft/min)
    ir::T             # Reaction intensity (kW/m²)
    hpua::T           # Heat per unit area (kJ/m²)
    flin::T           # Fireline intensity (kW/m)
    phiw::T           # Wind factor
    phis::T           # Slope factor
end

Computing Spread Rate

# Create fuel table and get fuel model
fuel_table = create_standard_fuel_table(Float64)
fm = get_fuel_model(fuel_table, 1, 60)  # Short grass at 60% live moisture

# Define conditions
M1 = 0.06    # 1-hour moisture (6%)
M10 = 0.08   # 10-hour moisture (8%)
M100 = 0.10  # 100-hour moisture (10%)
MLH = 0.60   # Live herbaceous moisture (60%)
MLW = 0.90   # Live woody moisture (90%)
wind = 440.0 # Midflame wind speed (ft/min ≈ 5 mph)
slope = 0.0  # tan²(slope angle)

# Calculate spread rate
result = Elmfire.surface_spread_rate(fm, M1, M10, M100, MLH, MLW, wind, slope)

println("Fire Spread Results:")
println("  Rate of spread: $(round(result.velocity, digits=1)) ft/min")
println("  Base rate (no wind/slope): $(round(result.vs0, digits=1)) ft/min")
println("  Reaction intensity: $(round(result.ir, digits=1)) kW/m²")
println("  Heat per unit area: $(round(result.hpua, digits=0)) kJ/m²")
println("  Fireline intensity: $(round(result.flin, digits=1)) kW/m")
println("  Wind factor (φw): $(round(result.phiw, digits=2))")
println("  Slope factor (φs): $(round(result.phis, digits=2))")
Fire Spread Results:
  Rate of spread: 103.3 ft/min
  Base rate (no wind/slope): 4.6 ft/min
  Reaction intensity: 156.4 kW/m²
  Heat per unit area: 1030.0 kJ/m²
  Fireline intensity: 540.6 kW/m
  Wind factor (φw): 21.43
  Slope factor (φs): 0.0

Effect of Wind Speed

# Spread rate vs wind speed for different fuel models
fuel_ids = [1, 3, 4, 8]  # Short grass, tall grass, chaparral, timber litter
fuel_names = ["Short grass", "Tall grass", "Chaparral", "Timber litter"]
wind_speeds = 0.0:50.0:1000.0  # ft/min

p = plot(
    xlabel = "Midflame Wind Speed (ft/min)",
    ylabel = "Rate of Spread (ft/min)",
    title = "Wind Speed Effect on Fire Spread",
    legend = :topleft
)

for (i, fuel_id) in enumerate(fuel_ids)
    fm = get_fuel_model(fuel_table, fuel_id, 60)
    rates = [Elmfire.surface_spread_rate(fm, 0.06, 0.08, 0.10, 0.60, 0.90, w, 0.0).velocity
             for w in wind_speeds]
    plot!(p, wind_speeds, rates, label=fuel_names[i], linewidth=2)
end
p

Effect of Fuel Moisture

# Spread rate vs 1-hour fuel moisture
fm = get_fuel_model(fuel_table, 1, 60)  # Short grass
moistures = 0.02:0.01:0.15  # 2% to 15%
wind = 440.0

rates = [Elmfire.surface_spread_rate(fm, m, m+0.02, m+0.04, 0.60, 0.90, wind, 0.0).velocity
         for m in moistures]

plot(moistures .* 100, rates,
    xlabel = "1-Hour Fuel Moisture (%)",
    ylabel = "Rate of Spread (ft/min)",
    title = "Fuel Moisture Effect on Fire Spread (Short Grass)",
    linewidth = 2,
    legend = false
)

Effect of Slope

# Spread rate vs slope angle
fm = get_fuel_model(fuel_table, 1, 60)
slopes_deg = 0.0:5.0:45.0
slopes_tan2 = tan.(deg2rad.(slopes_deg)).^2

rates_no_wind = [Elmfire.surface_spread_rate(fm, 0.06, 0.08, 0.10, 0.60, 0.90, 0.0, s).velocity
                 for s in slopes_tan2]
rates_with_wind = [Elmfire.surface_spread_rate(fm, 0.06, 0.08, 0.10, 0.60, 0.90, 440.0, s).velocity
                   for s in slopes_tan2]

plot(slopes_deg, [rates_no_wind rates_with_wind],
    xlabel = "Slope Angle (degrees)",
    ylabel = "Rate of Spread (ft/min)",
    title = "Slope Effect on Fire Spread",
    label = ["No wind" "5 mph wind"],
    linewidth = 2
)

Elliptical Fire Spread

Real fires don’t spread uniformly in all directions. Under the influence of wind or slope, fires spread in an elliptical pattern (Anderson 1982).

Length-to-Breadth Ratio

Direction of Maximum Spread

Rothermel gives the spread rate along a single direction. ELMFIRE resolves that direction by treating the wind and slope factors as vectors, so a fire under a crosswind on a slope heads somewhere between the two:

\[ \begin{aligned} \phi_x &= \alpha\left[\phi_w \sin(\theta_w - \pi) + \phi_s \sin(\theta_a - \pi)\right] \\ \phi_y &= \alpha\left[\phi_w \cos(\theta_w - \pi) + \phi_s \cos(\theta_a - \pi)\right] \\ |V_{DMS}| &= V_{s0}\left(\alpha + |\phi|\right) \end{aligned} \]

\(\theta_w\) is the direction the wind blows from, \(\theta_a\) the topographic aspect, and \(\alpha\) the acceleration factor. Wind and slope that oppose each other partially cancel; a scalar \(1 + \phi_s + \phi_w\) would let them add.

# Wind from the west, east-facing slope: the two factors oppose
phimag, dms_x, dms_y = Elmfire.spread_direction(2.0, 3.0, 90.0, 270.0)
println("opposed:  |phi| = ", round(phimag, digits = 3))

# West-facing slope: they reinforce
phimag, dms_x, dms_y = Elmfire.spread_direction(2.0, 3.0, 270.0, 270.0)
println("aligned:  |phi| = ", round(phimag, digits = 3),
        "  heading = (", round(dms_x, digits = 2), ", ", round(dms_y, digits = 2), ")")
opposed:  |phi| = 1.0
aligned:  |phi| = 5.0  heading = (1.0, 0.0)

Because \(|\phi|\) carries both effects, inverting Rothermel’s Eq. 47 recovers an effective mid-flame wind speed that includes slope. That is what drives the ellipse shape below.

\[ U_{mf,e} = \left(\frac{|\phi|}{C(\beta/\beta_{op})^{-E}}\right)^{1/B} \]

Velocities from this construction are parallel to the local slope. slope_projection_factors maps them onto the horizontal map grid:

\[ \frac{U_x}{U_{x,\parallel}} = 1 - |\sin\theta_a|(1 - \cos\gamma), \qquad \frac{U_y}{U_{y,\parallel}} = 1 - |\cos\theta_a|(1 - \cos\gamma) \]

The ellipse shape is characterized by the length-to-breadth ratio:

\[ \frac{L}{B} = 0.936 e^{0.2566U} + 0.461 e^{-0.1548U} - 0.397 \]

Where \(U\) is the effective mid-flame wind speed, scaled from ft/min by ELMFIRE’s WSMFEFF_LOW_MULT (5.07955e-3). elliptical_spread takes that wind speed in ft/min; multiply mi/h by 88 to convert.

The effective wind speed is not the raw wind: it is recovered from the combined wind and slope factor by inverting Rothermel’s Eq. 47 (see effective_wind_speed), so terrain elongates the ellipse too.

# Length-to-breadth ratio vs wind speed
wind_mph = 0:0.5:20
lb_ratios = [Elmfire.elliptical_spread(100.0, w * 88).length_to_breadth for w in wind_mph]

plot(wind_mph, lb_ratios,
    xlabel = "Effective Mid-Flame Wind Speed (mi/h)",
    ylabel = "Length-to-Breadth Ratio",
    title = "Fire Shape vs Wind Speed",
    linewidth = 2,
    legend = false,
    ylims = (0, 10)
)

EllipticalSpread Type

struct EllipticalSpread{T<:AbstractFloat}
    head::T           # Head fire spread rate (ft/min)
    back::T           # Backing fire spread rate (ft/min)
    flank::T          # Flanking fire spread rate (ft/min)
    eccentricity::T   # Ellipse eccentricity
    length_to_breadth::T  # Length to breadth ratio
end

Computing Elliptical Spread

# Head fire velocity from Rothermel
head_velocity = 50.0  # ft/min
wsmfeff = 10.0 * 88   # 10 mi/h effective mid-flame wind, in ft/min

es = Elmfire.elliptical_spread(head_velocity, wsmfeff)

println("Elliptical Fire Spread:")
println("  Head fire rate: $(round(es.head, digits=1)) ft/min")
println("  Backing fire rate: $(round(es.back, digits=1)) ft/min")
println("  Flanking fire rate: $(round(es.flank, digits=1)) ft/min")
println("  Eccentricity: $(round(es.eccentricity, digits=3))")
println("  Length/Breadth ratio: $(round(es.length_to_breadth, digits=2))")
Elliptical Fire Spread:
  Head fire rate: 50.0 ft/min
  Backing fire rate: 1.7 ft/min
  Flanking fire rate: 9.3 ft/min
  Eccentricity: 0.933
  Length/Breadth ratio: 2.78

Visualizing Elliptical Spread

# Plot fire ellipse shape for different wind speeds
wind_speeds = [0.0, 5.0, 10.0, 15.0]
p = plot(aspect_ratio=:equal, xlabel="Distance (ft)", ylabel="Distance (ft)",
         title="Fire Spread Ellipses (1 hour)", legend=:topright)

for wind in wind_speeds
    es = Elmfire.elliptical_spread(50.0, wind * 88)

    # Generate ellipse points (60 minutes of spread)
    t = 60.0  # minutes
    theta = range(0, 2π, length=100)

    # Semi-axes from spread rates
    a = (es.head + es.back) / 2 * t  # semi-major axis
    c = (es.head - es.back) / 2 * t  # focus distance from center

    if es.eccentricity > 0.001
        b = a * sqrt(1 - es.eccentricity^2)  # semi-minor axis
    else
        b = a
    end

    # Ellipse centered at origin, then shifted so the ignition point sits on the
    # REAR focus: the head runs a + c = head * t ahead, the back -a + c = -back * t
    # behind. Shifting the other way mirrors the fire and puts the long axis upwind.
    x = a .* cos.(theta) .+ c
    y = b .* sin.(theta)

    plot!(p, x, y, label="$(Int(wind)) mph", linewidth=2)
end
p

Derivation of Constants

The Rothermel model constants were derived through a combination of theoretical analysis, laboratory experiments, and field observations. This section documents the origin of each constant.

Physical Constants

These constants describe fundamental properties of wildland fuels (Rothermel 1972, Table II; Albini 1976):

Constant Value Units Description Source
\(\rho_p\) 32 lb/ft³ Particle density (ovendry wood) Rothermel (1972)
\(S_T\) 0.055 lb/lb (mass fraction) Total mineral content (mass fraction) Rothermel (1972)
\(S_e\) 0.01 lb/lb (mass fraction) Effective (silica-free) mineral content Rothermel (1972)
\(h\) 8000 BTU/lb Low heat of combustion (typical) Rothermel (1972)

The particle density of 32 lb/ft³ is representative of ovendry wood and was measured by Byram et al. (1952). The mineral content values are averages for forest fuels from laboratory analysis.

Fixed Surface-Area-to-Volume Ratios

Fuel Class SAV (1/ft) Source
10-hour dead 109 Rothermel (1972), Table I
100-hour dead 30 Rothermel (1972), Table I

These values are fixed because the timelag fuel classes (10-hr, 100-hr) are defined by particle diameter, which directly determines SAV for cylindrical fuel elements.

Empirically-Derived Coefficients

The following coefficients were derived by fitting experimental data from laboratory fire spread tests (Rothermel 1972, pp. 16-26):

Wind Coefficients (A, B, C, E)

Coefficient Formula Derivation
\(A\) \(133.0 / \sigma^{0.7913}\) Reaction velocity exponent (Eq. 38)
\(B\) \(0.02526 \times \sigma^{0.54}\) Wind velocity exponent (Eq. 49)
\(C\) \(7.47 \times e^{-0.133 \times \sigma^{0.55}}\) Wind coefficient (Eq. 48)
\(E\) \(0.715 \times e^{-0.000359 \times \sigma}\) Wind/packing exponent (Eq. 50)

These were determined by regression analysis of 31 laboratory fires in ponderosa pine needle fuel beds with SAV ratios from 1,140 to 2,000 ft⁻¹ and wind speeds from 0 to 9 mi/h (Rothermel 1972, Table IV).

Optimal Packing Ratio

\[ \beta_{op} = 3.348 \times \sigma^{-0.8189} \quad \text{(dimensionless)} \]

Where \(\sigma\) is the characteristic surface-area-to-volume ratio (ft⁻¹). The coefficient 3.348 has units of ft⁻⁰·⁸¹⁸⁹.

Derived from the relationship between packing ratio and reaction velocity that maximizes heat release (Rothermel 1972, Eq. 37). The optimal packing provides the best balance between fuel availability and oxygen access.

Propagating Flux Ratio

\[ \xi = \frac{\exp[(0.792 + 0.681\sqrt{\sigma})(0.1 + \beta)]}{192 + 0.2595\sigma} \]

This ratio represents the fraction of reaction intensity that heats adjacent fuel particles. Derived from no-wind fire spread experiments (Rothermel 1972, Eq. 42).

Reaction Velocity

\[ \Gamma' = \Gamma'_{max} \left(\frac{\beta}{\beta_{op}}\right)^A \exp\left[A\left(1 - \frac{\beta}{\beta_{op}}\right)\right] \]

where:

\[ \Gamma'_{max} = \frac{\sigma^{1.5}}{495 + 0.0594\sigma^{1.5}} \]

The maximum reaction velocity equation (Rothermel 1972, Eq. 36) was fitted to data from Fons (1946) and Anderson (1969) for cribs and fuel beds.

Residence Time

\[ \tau_r = \frac{384}{\sigma} \quad \text{(min)} \]

Where \(\sigma\) is in ft⁻¹ and the coefficient 384 has units of ft⁻¹·min.

The residence time represents how long the flaming front persists at a given point. The coefficient 384 was derived from laboratory observations of flame duration vs. fuel particle size (Rothermel 1972, Eq. 77; Anderson 1969).

Damping Coefficients

Moisture Damping

\[ \eta_M = 1 - 2.59r + 5.11r^2 - 3.52r^3 \quad \text{(dimensionless)} \]

where \(r = M/M_{ex}\) is the ratio of fuel moisture content \(M\) (lb water/lb dry fuel) to moisture of extinction \(M_{ex}\) (lb water/lb dry fuel). Both are typically expressed as fractions (e.g., 0.06 = 6% moisture).

This polynomial was fitted to laboratory data showing the reduction in reaction intensity as fuel moisture increases (Rothermel 1972, Eq. 29). The curve passes through (0, 1) at zero moisture and (1, 0) at extinction moisture.

Mineral Damping

\[ \eta_s = 0.174 \times S_e^{-0.19} \quad \text{(dimensionless)} \]

where \(S_e\) is the effective (silica-free) mineral content (lb mineral/lb fuel).

Accounts for the heat-absorbing effect of fuel minerals. Derived from laboratory combustion experiments (Rothermel 1972, Eq. 30). With the standard effective mineral content of 0.01, this gives \(\eta_s \approx 0.42\).

Heat Transfer Constants

Heat of Pre-ignition

\[ Q_{ig} = 250 + 1116M \quad \text{(BTU/lb)} \]

Where \(M\) is fuel moisture content (lb water/lb dry fuel, e.g., 0.06 for 6%).

The heat required to raise fuel to ignition temperature. The constant 250 BTU/lb represents the sensible heat to raise dry fuel to ignition (~600°F), while the coefficient 1116 BTU/lb accounts for the latent heat of vaporizing fuel moisture (Rothermel 1972, Eq. 12).

Effective Heating Number

\[ \varepsilon = \exp\left(\frac{-138}{\sigma}\right) \quad \text{(dimensionless)} \]

Where \(\sigma\) is surface-area-to-volume ratio (ft⁻¹) and the coefficient 138 has units of ft⁻¹.

The fraction of a fuel particle that is heated to ignition temperature before the flame front passes. Smaller particles (higher SAV) heat more completely (Rothermel 1972, Eq. 14). The coefficient 138 ft⁻¹ was calibrated from experimental data.

Slope Factor

\[ \phi_s = 5.275 \times \beta^{-0.3} \times \tan^2\theta \quad \text{(dimensionless)} \]

Where \(\beta\) is packing ratio (dimensionless) and \(\tan^2\theta\) is passed directly to the function (dimensionless).

The slope factor accounts for radiant heat transfer to upslope fuels. The coefficients 5.275 and -0.3 are dimensionless and were derived from laboratory experiments on inclined fuel beds (Rothermel 1972, Eq. 52). The dependence on \(\tan^2\theta\) reflects the geometric increase in view factor with slope angle.

Elliptical Spread (Anderson 1982)

The length-to-breadth ratio equation:

\[ \frac{L}{B} = 0.936 e^{0.2566U} + 0.461 e^{-0.1548U} - 0.397 \]

was derived by Anderson (1982) from aerial observations of 28 wildfires in Australia and the United States. The data showed that fire shape depends primarily on effective wind speed \(U\) (mi/h), with fires becoming more elongated at higher wind speeds.

Viewing Computed Coefficients

fm = get_fuel_model(fuel_table, 1, 60)

println("Rothermel Coefficients (Fuel Model 1 - Short Grass):")
println("  A = $(round(fm.A, digits=4))")
println("  B = $(round(fm.B, digits=4))")
println("  C = $(round(fm.C, digits=4))")
println("  E = $(round(fm.E, digits=6))")
println("")
println("Derived Properties:")
println("  Overall SAV (σ): $(round(fm.SIG_overall, digits=0)) 1/ft")
println("  Bulk density (ρb): $(round(fm.rhob, digits=4)) lb/ft³")
println("  Packing ratio (β): $(round(fm.beta, digits=5))")
println("  Optimal packing (βop): $(round(fm.betaop, digits=5))")
println("  Reaction velocity (Γ'): $(round(fm.gammaprime, digits=3)) 1/min")
println("  Propagating flux (ξ): $(round(fm.xi, digits=4))")
println("  Residence time (τr): $(round(fm.tr, digits=3)) min")
println("  Mineral damping (ηs): $(round(fm.etas, digits=4))")
Rothermel Coefficients (Fuel Model 1 - Short Grass):
  A = 0.2087
  B = 2.0712
  C = 0.0001
  E = 0.203524

Derived Properties:
  Overall SAV (σ): 3500.0 1/ft
  Bulk density (ρb): 0.034 lb/ft³
  Packing ratio (β): 0.00106
  Optimal packing (βop): 0.00419
  Reaction velocity (Γ'): 14.201 1/min
  Propagating flux (ξ): 0.0578
  Residence time (τr): 0.11 min
  Mineral damping (ηs): 0.4174

References

Primary Sources

  1. Rothermel, R.C. (1972). A mathematical model for predicting fire spread in wildland fuels. USDA Forest Service Research Paper INT-115. Ogden, UT. 40 p.
    • The foundational paper defining the spread rate equation and deriving all major coefficients.
  2. Anderson, H.E. (1982). Aids to determining fuel models for estimating fire behavior. USDA Forest Service General Technical Report INT-122. Ogden, UT. 22 p.
    • Defines the 13 standard fuel models and elliptical fire spread equations.
  3. Albini, F.A. (1976). Estimating wildfire behavior and effects. USDA Forest Service General Technical Report INT-30. Ogden, UT. 92 p.
    • Provides physical constants and fuel property values used in fire modeling.

Historical Sources (cited by Rothermel)

  1. Byram, G.M., Sauer, F.M., Fons, W.L., & Arnold, R.K. (1952). Thermal properties of forest fuels. USDA Forest Service, Division of Fire Research. Berkeley, CA.
    • Source of particle density (32 lb/ft³) measurements.
  2. Fons, W.L. (1946). Analysis of fire spread in light forest fuels. Journal of Agricultural Research 72(3): 93-121.
    • Early fire spread experiments used to calibrate reaction velocity.
  3. Anderson, H.E. (1969). Heat transfer and fire spread. USDA Forest Service Research Paper INT-69. Ogden, UT. 20 p.
    • Laboratory data on heat transfer and residence time.

Modern Interpretations

  1. Andrews, P.L. (2018). The Rothermel surface fire spread model and associated developments: A comprehensive explanation. USDA Forest Service General Technical Report RMRS-GTR-371. Fort Collins, CO. 121 p.
    • Comprehensive modern explanation of the Rothermel model with derivations and examples.
  2. Scott, J.H. & Burgan, R.E. (2005). Standard fire behavior fuel models: A comprehensive set for use with Rothermel’s surface fire spread model. USDA Forest Service General Technical Report RMRS-GTR-153. Fort Collins, CO. 72 p.
    • Expanded set of 40 fuel models (Scott/Burgan models) for diverse fuel types.

See Also