API

Types

BatchResult

BatchResult{T<:AbstractFloat}

Result from a single batch job.

BatchSimulationJob

BatchSimulationJob{T<:AbstractFloat}

A single simulation job for batch processing.

BuildingIgnitionResult

BuildingIgnitionResult{T<:AbstractFloat}

Result of a building ignition calculation.

CPUFireState

FireState{T<:AbstractFloat, A<:AbstractMatrix{T}}

Complete state of a fire simulation. Parameterized on element type T and array type A to support both CPU (Matrix{T}) and GPU arrays.

FireState{T}(ncols, nrows, cellsize; xllcorner=zero(T), yllcorner=zero(T), padding=2, band_thickness=5)

Create a new fire simulation state with specified precision.

Arguments

  • ncols: Number of grid columns
  • nrows: Number of grid rows
  • cellsize: Cell size in feet
  • xllcorner: X coordinate of lower-left corner (default 0.0)
  • yllcorner: Y coordinate of lower-left corner (default 0.0)
  • padding: Boundary padding for stencil operations (default 2)
  • band_thickness: Narrow band half-width (default 5)

CanopyGrid

CanopyGrid{T<:AbstractFloat}

Grid of canopy properties for each cell.

CanopyGrid{T}(ncols, nrows)

Create an empty canopy grid (no canopy).

CanopyGrid{T}(ncols, nrows, cbd, cbh, cc, ch)

Create a uniform canopy grid.

CanopyProperties

CanopyProperties{T<:AbstractFloat}

Canopy fuel properties for a cell.

ConstantWeather

ConstantWeather{T<:AbstractFloat}

Constant (spatially and temporally uniform) weather conditions.

ContainmentLine

ContainmentLine{T<:AbstractFloat}

A fire containment line (fireline, retardant line, etc.).

ContainmentLine{T}(cells; effectiveness=0.9, width=6.0, construction_time=0.0, resource_id=0)

Create a containment line from a vector of cells.

CrownFireResult

CrownFireResult{T<:AbstractFloat}

Results from crown fire calculation.

DiurnalConfig

DiurnalConfig{T<:AbstractFloat}

Configuration for diurnal (day/night) spread rate adjustment.

During the burn period (daytime), the adjustment factor is 1.0 (full spread). Outside the burn period (nighttime), spread rate is multiplied by overnight_factor.

Fields

  • forecast_start_hour: Simulation start hour in local time (0-24)
  • sunrise_hour: Local sunrise hour (default 6.0)
  • sunset_hour: Local sunset hour (default 20.0)
  • burn_period_center_frac: Fraction of daylight for burn period center (default 0.667)
  • burn_period_length: Duration of active burn period in hours (default 10.0)
  • overnight_factor: Spread rate multiplier outside burn period (default 0.1)

EllipticalSpread

EllipticalSpread{T<:AbstractFloat}

Elliptical fire spread parameters following Anderson (1982).

EnsembleConfig

EnsembleConfig{T<:AbstractFloat}

Configuration for ensemble simulation runs.

EnsembleMember

EnsembleMember{T<:AbstractFloat}

Results from a single ensemble member simulation.

EnsembleResult

EnsembleResult{T<:AbstractFloat}

Aggregated results from an ensemble of simulations.

FireState

FireState{T<:AbstractFloat, A<:AbstractMatrix{T}}

Complete state of a fire simulation. Parameterized on element type T and array type A to support both CPU (Matrix{T}) and GPU arrays.

FireState{T}(ncols, nrows, cellsize; xllcorner=zero(T), yllcorner=zero(T), padding=2, band_thickness=5)

Create a new fire simulation state with specified precision.

Arguments

  • ncols: Number of grid columns
  • nrows: Number of grid rows
  • cellsize: Cell size in feet
  • xllcorner: X coordinate of lower-left corner (default 0.0)
  • yllcorner: Y coordinate of lower-left corner (default 0.0)
  • padding: Boundary padding for stencil operations (default 2)
  • band_thickness: Narrow band half-width (default 5)

FuelModel

FuelModel{T<:AbstractFloat}

Complete fuel model with all computed Rothermel coefficients. Matches ELMFIRE’s FUELMODELTABLE_TYPE structure.

The 6 fuel size classes are:

  1. 1-hour dead
  2. 10-hour dead
  3. 100-hour dead
  4. Dynamic dead (transferred herbaceous)
  5. Live herbaceous
  6. Live woody

FuelModelArray

FuelModelArray{T<:AbstractFloat}

Dense struct-of-arrays representation of fuel model coefficients for GPU-friendly lookup. Replaces the Dict-based FuelModelTable with array indexing.

Fields are stored as [fuel_index, moisture_class_index] matrices where:

  • fuel_index maps from fuel ID via fuel_id_to_index
  • moisture_class_index = live_moisture_class - 29 (so class 30 → index 1, class 120 → index 91)

NTuple fields are stored as 3D arrays [fuel_index, moisture_class_index, component].

FuelModelArray(table::FuelModelTable{T}) -> FuelModelArray{T}

Convert a FuelModelTable to a dense FuelModelArray for GPU-friendly access.

Examples

table = create_standard_fuel_table(Float64)
arr = FuelModelArray(table)

FuelModelTable

FuelModelTable{T<:AbstractFloat}

2D table of fuel models indexed by (fuelid, livemoisture_class). Live moisture class ranges from 30 to 120 (representing 30% to 120% moisture).

GeoMetadata

GeoMetadata{T<:AbstractFloat}

Geospatial metadata for a raster dataset.

GeoRaster

GeoRaster{T<:AbstractFloat}

A raster dataset with geospatial metadata.

HamadaParameters

HamadaParameters{T<:AbstractFloat}

Parameters for the Hamada urban fire spread model.

The Hamada model describes fire spread between buildings based on separation distance, wind conditions, and building characteristics.

LandscapeData

LandscapeData{T<:AbstractFloat}

Complete landscape data for fire simulation including fuel, topography, and canopy.

NarrowBand

NarrowBand

Tracks the active cells near the fire front for efficient computation.

The narrow band contains cells within band_thickness cells of the fire front (where φ changes sign). Uses a BitMatrix for O(1) membership testing and a packed Vector for cache-friendly iteration.

ParallelConfig

ParallelConfig

Configuration for parallel execution.

PerturbationConfig

PerturbationConfig{T<:AbstractFloat}

Configuration for stochastic parameter perturbations.

RawFuelModel

RawFuelModel{T<:AbstractFloat}

Raw fuel model parameters as read from CSV file before coefficient calculation.

SimulationConfig

SimulationConfig{T<:AbstractFloat}

Configuration for full simulation with crown fire, spotting, and weather interpolation.

SpotFire

SpotFire{T<:AbstractFloat}

A spot fire ignition location.

SpotFireTracker

SpotFireTracker{T<:AbstractFloat}

Tracks pending and active spot fires during simulation.

SpottingParameters

SpottingParameters{T<:AbstractFloat}

Parameters controlling ember generation and transport.

SpreadRateDampeningConfig

SpreadRateDampeningConfig{T<:AbstractFloat}

Configuration for spread rate dampening.

Fields

  • mode::SpreadRateDampeningMode: Dampening strategy
  • max_spread_rate::T: Maximum spread rate for ABSOLUTE_CAP mode (ft/min, default 250.0)
  • linear_threshold::T: Wind speed threshold for LINEAR_DAMPENING mode (ft/min, default 400.0)

Examples

# No dampening (default)
cfg = SpreadRateDampeningConfig{Float64}()

# Cap spread rate at midflame wind speed
cfg = SpreadRateDampeningConfig{Float64}(WIND_SPEED_CAP)

# Hard cap at 200 ft/min
cfg = SpreadRateDampeningConfig{Float64}(ABSOLUTE_CAP, max_spread_rate=200.0)

# Linear dampening above 300 ft/min midflame wind
cfg = SpreadRateDampeningConfig{Float64}(LINEAR_DAMPENING, linear_threshold=300.0)

SpreadRateDampeningMode

SpreadRateDampeningMode

Configurable strategy for limiting fire spread rates at extreme wind speeds.

The Rothermel model’s power-law wind factor (phi_w = C * wsmf^B) overpredicts spread rates at high wind speeds. These modes provide physically motivated dampening.

Values

  • NO_DAMPENING: Status quo — only Rothermel’s 0.9*IR wind limit applies
  • WIND_SPEED_CAP: Andrews et al. (2013) — cap ROS at midflame wind speed (ft/min)
  • ABSOLUTE_CAP: Hard cap at a user-specified maximum spread rate (ft/min)
  • LINEAR_DAMPENING: Transition phi_w from power-law to linear above a threshold

SpreadResult

SpreadResult{T<:AbstractFloat}

Results from the Rothermel surface fire spread rate calculation.

SuppressionResource

SuppressionResource{T<:AbstractFloat}

A fire suppression resource (crew, engine, etc.).

SuppressionState

SuppressionState{T<:AbstractFloat}

Complete state of fire suppression activities.

SuppressionState{T}(ncols::Int, nrows::Int)

Create an empty suppression state.

ThreadLocalState

ThreadLocalState{T<:AbstractFloat, R<:AbstractRNG}

Thread-local simulation state to avoid data races.

WUIBuilding

WUIBuilding{T<:AbstractFloat}

A building in the WUI that can be affected by wildfire.

WUIGrid

WUIGrid{T<:AbstractFloat}

Grid of buildings for WUI fire simulation.

WUIGrid{T}(buildings::Vector{WUIBuilding{T}}, ncols::Int, nrows::Int)

Create a WUI grid from a vector of buildings.

WUIGrid{T}(ncols::Int, nrows::Int)

Create an empty WUI grid.

WeatherGrid

WeatherGrid{T<:AbstractFloat}

A grid of weather values at a single time.

WeatherGrid{T}(ncols, nrows, cellsize; xllcorner=zero(T), yllcorner=zero(T))

Create an empty weather grid.

WeatherGrid{T}(weather::ConstantWeather{T}, ncols, nrows, cellsize)

Create a uniform weather grid from constant weather conditions.

WeatherInterpolator

WeatherInterpolator{T<:AbstractFloat}

Handles bilinear spatial and linear temporal interpolation of weather data to simulation grid and time.

WeatherInterpolator(
    weather_series::WeatherTimeSeries{T},
    sim_ncols::Int, sim_nrows::Int,
    sim_cellsize::T,
    sim_xllcorner::T = zero(T),
    sim_yllcorner::T = zero(T)
)

Create a weather interpolator for the given simulation grid.

WeatherTimeSeries

WeatherTimeSeries{T<:AbstractFloat}

A time series of weather grids.

WeatherTimeSeries{T}(grids, times)

Create a weather time series from a vector of grids and times.

WeatherTimeSeries{T}(weather::ConstantWeather{T}, ncols, nrows, cellsize, duration)

Create a constant weather time series.

Functions

acceleration_factor

acceleration_factor(t::T, tau::T) -> T

Compute the fire acceleration factor at time t with time constant tau (minutes).

The acceleration factor ramps from 0 at ignition to 1 at steady state following 1 - exp(-t/τ). This models the observed behavior that fires initially spread slower than the Rothermel steady-state prediction.

When tau ≤ 0.5 (minutes), acceleration is disabled and returns 1.0.

active_mask

active_mask(nb::NarrowBand, nx::Int, ny::Int) -> BitMatrix

Create a dense BitMatrix indicating which cells are in the active narrow band. Useful for GPU kernels that need a mask instead of an index list.

Examples

mask = active_mask(state.narrow_band, size(state.phi)...)

add_raw_model! {#add_raw_model!}

add_raw_model!(table::FuelModelTable{T}, raw::RawFuelModel) where {T}

Add a raw fuel model to the table and compute coefficients for all live moisture classes (30-120).

add_resource! {#add_resource!}

add_resource!(state::SuppressionState{T}, resource::SuppressionResource{T})

Add a suppression resource to the state.

add_spot_fires! {#add_spot_fires!}

add_spot_fires!(tracker::SpotFireTracker{T}, fires::Vector{SpotFire{T}})

Add new spot fires to the pending queue.

aggregate_ensemble_statistics! {#aggregate_ensemble_statistics!}

aggregate_ensemble_statistics!(result::EnsembleResult{T})

Compute aggregate statistics from ensemble members.

apply_containment! {#apply_containment!}

apply_containment!(fire_state::FireState{T}, suppression_state::SuppressionState{T})

Apply containment effects to the fire state by reducing spread velocities at contained cells.

apply_spread_rate_dampening

apply_spread_rate_dampening(
    velocity::T, wsmf::T, vs0::T, phis::T, phiw::T,
    phiwterm::T, B::T,
    config::SpreadRateDampeningConfig{T}
) -> T

Apply spread rate dampening to a computed fire spread velocity.

Called after surface_spread_rate() and before elliptical_spread().

Arguments

  • velocity: Spread rate from Rothermel model (ft/min)
  • wsmf: Midflame wind speed (ft/min)
  • vs0: Base spread rate without wind/slope (ft/min)
  • phis: Slope factor
  • phiw: Wind factor (power-law)
  • phiwterm: Fuel model wind coefficient (C * (beta/betaop)^(-E))
  • B: Fuel model wind exponent (0.02526 * sigma^0.54)
  • config: Dampening configuration

assign_resource! {#assign_resource!}

assign_resource!(
    suppression::SuppressionState{T},
    resource_id::Int,
    targets::Vector{Tuple{Int,Int}}
)

Assign a resource to construct lines at target locations.

building_ignition_probability

building_ignition_probability(
    building::WUIBuilding{T},
    heat_flux::T,
    exposure_time::T
) -> T

Calculate the probability of building ignition given heat flux and exposure time.

Uses a dose-response model based on heat flux accumulation.

Arguments

  • building: Building properties
  • heat_flux: Incident heat flux (kW/m²)
  • exposure_time: Duration of exposure (minutes)

Returns

  • Ignition probability (0-1)

calculate_tanslp2

calculate_tanslp2(slope_degrees::T) -> T

Calculate tan²(slope) from slope in degrees.

check_convergence

check_convergence(result::EnsembleResult{T}; threshold::T = T(0.001)) -> Bool

Check if the ensemble has converged based on burn probability changes.

combined_fireline_intensity

combined_fireline_intensity(
    surface_result::SpreadResult{T},
    crown_result::CrownFireResult{T},
    fm::FuelModel{T}
) -> T

Calculate the combined fireline intensity (kW/m) for surface and crown fire.

combined_spread_rate

combined_spread_rate(
    surface_result::SpreadResult{T},
    crown_result::CrownFireResult{T}
) -> T

Calculate the combined spread rate considering both surface and crown fire.

For active crown fires, returns the crown fire spread rate. For passive crown fires, returns surface spread rate (crown fire doesn’t add to ROS). For surface fires only, returns the surface spread rate.

compute_burn_probability

compute_burn_probability(members::Vector{EnsembleMember{T}}, ncols::Int, nrows::Int) -> Matrix{T}

Compute burn probability from ensemble members.

compute_cfl_timestep

compute_cfl_timestep(
    ux::AbstractMatrix{T},
    uy::AbstractMatrix{T},
    active_cells::AbstractVector{CartesianIndex{2}},
    dx::T,
    dt_current::T;
    target_cfl::T = T(0.9),
    dt_max::T = T(10)
) -> T

Compute the adaptive timestep based on the CFL condition.

CFL = umax * dt / dx < target_cfl

Implementation follows elmfire_level_set.f90:2010-2027.

compute_fuel_model

compute_fuel_model(raw::RawFuelModel{T}, live_moisture_class::Int) -> FuelModel{T}

Compute all Rothermel coefficients from raw fuel model data.

For dynamic fuel models, the live herbaceous fuel is partitioned between dead and live classes based on the live moisture class (30-120).

Implementation follows elmfire_init.f90:786-880.

compute_mean_arrival_time

compute_mean_arrival_time(members::Vector{EnsembleMember{T}}, ncols::Int, nrows::Int) -> Tuple{Matrix{T}, Matrix{T}}

Compute mean and standard deviation of time of arrival from ensemble members. Returns (meantoa, stdtoa).

compute_normal

compute_normal(
    phi::AbstractMatrix{T},
    ix::Int, iy::Int,
    dx::T
) -> Tuple{T, T}

Compute the unit normal vector to the level set at cell (ix, iy). The normal points in the direction of increasing φ (outward from fire).

Uses central differences.

compute_num_embers

compute_num_embers(
    params::SpottingParameters{T},
    flin::T,
    crown_fire_type::Int,
    dt::T,
    cellsize::T;
    rng::AbstractRNG = Random.default_rng()
) -> Int

Calculate the number of embers to generate for a burning cell.

compute_radiative_heat_flux

compute_radiative_heat_flux(flin::T, distance::T, flame_length::T) -> T

Compute radiative heat flux (kW/m²) at a given distance from the fire.

Based on the solid flame model with view factor calculations.

Arguments

  • flin: Fireline intensity (kW/m)
  • distance: Distance from fire front (m)
  • flame_length: Flame length (m)

Returns

  • Heat flux in kW/m²

compute_slope_aspect

compute_slope_aspect(dem::GeoRaster{T}) -> Tuple{Matrix{T}, Matrix{T}}

Compute slope (degrees) and aspect (degrees) from a DEM using a 3x3 moving window.

Aspect is in degrees clockwise from north (0-360), with flat areas set to 0.

compute_slope_aspect(elevation::Matrix{T}, cellsize::T) -> Tuple{Matrix{T}, Matrix{T}}

Compute slope and aspect from an elevation matrix.

compute_view_factor

compute_view_factor(flame_height::T, distance::T) -> T

Compute the view factor between a vertical rectangular flame and a point.

Uses the approximation for a vertical radiating surface.

Arguments

  • flame_height: Height of the flame (m)
  • distance: Horizontal distance to the target point (m)

construct_containment_line! {#construct_containment_line!}

construct_containment_line!(
    state::SuppressionState{T},
    resource::SuppressionResource{T},
    start_ix::Int, start_iy::Int,
    target_ix::Int, target_iy::Int,
    dt::T,
    cellsize::T,
    t::T
) -> Tuple{Vector{Tuple{Int,Int}}, T}

Construct a containment line from start to target (or as far as possible in dt).

Arguments

  • state: Suppression state (modified in place)
  • resource: Resource constructing the line
  • start_ix, start_iy: Starting grid coordinates
  • target_ix, target_iy: Target grid coordinates
  • dt: Time step (minutes)
  • cellsize: Grid cell size (ft)
  • t: Current simulation time (minutes)

Returns

  • Tuple of (cellsconstructed, lengthconstructed)

create_building_grid

create_building_grid(
    ncols::Int, nrows::Int,
    building_spacing::Int,
    building_footprint::Int,
    origin_x::Int, origin_y::Int;
    construction_type::Symbol = :wood
) -> Vector{WUIBuilding{T}}

Create a regular grid of buildings for WUI simulation.

Arguments

  • ncols, nrows: Grid dimensions
  • building_spacing: Spacing between buildings (cells)
  • building_footprint: Building footprint size (cells)
  • origin_x, origin_y: Starting position for building grid
  • construction_type: Default construction type for buildings

create_constant_interpolator

create_constant_interpolator(
    weather::ConstantWeather{T},
    sim_ncols::Int, sim_nrows::Int,
    sim_cellsize::T
) -> WeatherInterpolator{T}

Create a weather interpolator for constant weather conditions.

create_grid_mapping

create_grid_mapping(
    weather_grid::WeatherGrid{T},
    sim_ncols::Int, sim_nrows::Int,
    sim_cellsize::T,
    sim_xllcorner::T,
    sim_yllcorner::T
) -> Tuple{Vector{T}, Vector{T}}

Create mapping from simulation grid to weather grid fractional coordinates.

Returns (xcolweather, yrowweather) vectors of fractional weather grid positions for bilinear interpolation. Values are clamped to [1, ncols] and [1, nrows].

create_standard_fuel_table

create_standard_fuel_table(::Type{T}) -> FuelModelTable{T}

Create a fuel model table with the standard FBFM 1-13 models plus non-burnable.

create_thread_local_states

create_thread_local_states(template::FireState{T}, n::Int) -> Vector{ThreadLocalState{T, MersenneTwister}}

Create n thread-local states from a template.

critical_fireline_intensity

critical_fireline_intensity(cbh::T, foliar_moisture::T) -> T

Calculate the critical fireline intensity (kW/m) needed to initiate crown fire.

Uses Van Wagner (1977) criterion: I_crit = (0.01 * CBH * (460 + 26 * FMC))^1.5

Arguments

  • cbh: Canopy base height (m)
  • foliar_moisture: Foliar moisture content (%, e.g., 100 for 100%)

Returns

  • Critical fireline intensity (kW/m)

crown_spread_rate

crown_spread_rate(
    canopy::CanopyProperties{T},
    flin_surface::T,
    ws20::T,
    m1::T,
    vs0::T;
    crown_fire_adj::T = one(T),
    spread_rate_limit::T = T(1000),
    critical_canopy_cover::T = T(0.4)
) -> CrownFireResult{T}

Calculate crown fire spread rate using the Cruz (2005) model.

Arguments

  • canopy: Canopy fuel properties
  • flin_surface: Surface fireline intensity (kW/m)
  • ws20: 20-ft wind speed (mph)
  • m1: 1-hour dead fuel moisture (fraction)
  • vs0: Base surface spread rate (ft/min)
  • crown_fire_adj: Crown fire adjustment factor (default 1.0)
  • spread_rate_limit: Maximum crown fire spread rate (ft/min)
  • critical_canopy_cover: Minimum canopy cover for active crown fire

Returns

  • CrownFireResult with crown fire type, spread rate, and related values

Implementation follows elmfire_spread_rate.f90:136-211.

diurnal_adjustment

diurnal_adjustment(config::DiurnalConfig{T}, t::T) -> T

Compute diurnal spread rate adjustment factor at simulation time t (minutes).

Returns 1.0 during the burn period (daytime) and config.overnight_factor outside it.

elliptical_spread

elliptical_spread(velocity::T, effective_windspeed::T) -> EllipticalSpread{T}

Calculate elliptical fire spread dimensions from the head fire rate and effective windspeed.

Uses Anderson (1982) length-to-breadth ratio: LB = 0.936 * exp(0.2566 * U) + 0.461 * exp(-0.1548 * U) - 0.397

Where U is effective windspeed in mi/h.

find_time_indices

find_time_indices(wts::WeatherTimeSeries{T}, t::T) -> Tuple{Int, Int, T}

Find the indices and interpolation weight for time t.

Returns (ilo, ihi, f) where the interpolated value is: value = (1-f) * grids[i*lo] + f * grids[i*hi]

generate_spot_fires

generate_spot_fires(
    ix0::Int, iy0::Int,
    flin::T,
    ws20::T,
    wd20::T,
    crown_fire_type::Int,
    params::SpottingParameters{T},
    cellsize::T,
    ncols::Int, nrows::Int,
    time_now::T,
    burned::BitMatrix;
    weather_interp = nothing,
    use_sardoy::Bool = false,
    rng::AbstractRNG = Random.default_rng()
) -> Vector{SpotFire{T}}

Generate potential spot fire locations from a burning cell.

Arguments

  • ix0, iy0: Source cell grid indices
  • flin: Fireline intensity (kW/m)
  • ws20: 20-ft wind speed (mph)
  • wd20: Wind direction (degrees, FROM convention)
  • crown_fire_type: 0 = surface, 1 = passive crown, 2 = active crown
  • params: Spotting parameters
  • cellsize: Cell size (ft)
  • ncols, nrows: Grid dimensions
  • time_now: Current simulation time (minutes)
  • burned: Matrix of burned cells
  • weather_interp: Optional WeatherInterpolator for wind-field ember advection
  • use_sardoy: Use Sardoy model for distribution parameters
  • rng: Random number generator

Returns

  • Vector of potential spot fire locations (may be empty)

geo_to_grid

geo_to_grid(x::T, y::T, metadata::GeoMetadata{T}) -> Tuple{Int, Int}

Convert geographic coordinates to grid indices.

get_active_cells

get_active_cells(nb::NarrowBand) -> AbstractVector{CartesianIndex{2}}

Get a view of active cell indices for iteration. Zero-allocation.

get_burned_area

get_burned_area(state::FireState{T}) -> T

Get the total burned area in square feet.

get_burned_area_acres

get_burned_area_acres(state::FireState{T}) -> T

Get the total burned area in acres.

get_canopy_properties

get_canopy_properties(grid::CanopyGrid{T}, ix::Int, iy::Int) -> CanopyProperties{T}

Get canopy properties for a specific cell.

get_exceedance_probability

get_exceedance_probability(result::EnsembleResult{T}, threshold_acres::T) -> T

Calculate the probability that burned area exceeds a threshold.

get_fire_perimeter

get_fire_perimeter(state::FireState) -> Vector{Tuple{Int, Int}}

Get the grid coordinates of cells on the fire perimeter (burned cells adjacent to unburned).

get_fuel_model

get_fuel_model(table::FuelModelTable, fuel_id::Int, live_moisture::AbstractFloat) -> FuelModel

Get the fuel model for the given fuel ID and live moisture content. Live moisture is clamped to [30, 120] and rounded to nearest integer.

get_fuel_model(table::FuelModelTable, fuel_id::Int, live_moisture_class::Int) -> FuelModel

Get the fuel model for the given fuel ID and live moisture class (30-120).

get_percentile_fire

get_percentile_fire(result::EnsembleResult{T}, percentile::T) -> Union{Nothing, EnsembleMember{T}}

Get the ensemble member closest to a given percentile of burned area.

get_ready_ignitions! {#get_ready_ignitions!}

get_ready_ignitions!(tracker::SpotFireTracker{T}, time_now::T) -> Vector{Tuple{Int,Int}}

Get spot fires that are ready to ignite and remove them from pending queue.

Returns grid indices of cells to ignite.

get_suppression_statistics

get_suppression_statistics(suppression::SuppressionState{T}) -> NamedTuple

Get summary statistics for suppression activities.

get_weather_at

get_weather_at(interp::WeatherInterpolator{T}, ix::Int, iy::Int, t::T) -> NamedTuple

Get bilinearly interpolated weather values at simulation grid cell (ix, iy) and time t.

Returns named tuple with fields: ws, wd, m1, m10, m100, mlh, mlw

get_wui_statistics

get_wui_statistics(wui_grid::WUIGrid{T}) -> NamedTuple

Get summary statistics for the WUI simulation.

grid_to_geo

grid_to_geo(ix::Int, iy::Int, metadata::GeoMetadata{T}) -> Tuple{T, T}

Convert grid indices to geographic coordinates (cell center).

grid_to_padded

grid_to_padded(state::FireState, ix::Int, iy::Int) -> Tuple{Int, Int}

Convert grid coordinates to padded array coordinates.

half_superbee

half_superbee(r::T) -> T

Superbee flux limiter function (half-superbee variant). Returns ψ(r)/2 where ψ is the Superbee limiter.

Used for second-order upwind schemes to limit gradients and prevent oscillations.

Implementation matches elmfire_level_set.f90:1533-1541.

hamada_spread_probability

hamada_spread_probability(
    source::WUIBuilding{T},
    target::WUIBuilding{T},
    params::HamadaParameters{T},
    wind_speed::T,
    wind_direction::T,
    cellsize::T
) -> T

Calculate the probability of fire spreading between buildings using the Hamada model.

Arguments

  • source: Building that is on fire
  • target: Potential target building
  • params: Hamada model parameters
  • wind_speed: Wind speed (mph)
  • wind_direction: Wind direction (degrees, FROM)
  • cellsize: Grid cell size (ft)

ignite! {#ignite!}

ignite!(state::FireState{T}, ix::Int, iy::Int, t::T)

Ignite a cell at grid coordinates (ix, iy) at time t. Sets up the level set as a signed distance function near the ignition point.

ignite_circle! {#ignite_circle!}

ignite_circle!(state::FireState{T}, center_x::Int, center_y::Int, radius_cells::T, t::T)

Ignite all cells within a circle of given radius (in grid cells) centered at (centerx, centery).

ignite_point! {#ignite_point!}

ignite_point!(state::FireState{T}, x::T, y::T, t::T)

Ignite a cell at world coordinates (x, y) at time t.

initialize_circular_fire! {#initialize_circular_fire!}

initialize_circular_fire!(
    phi::AbstractMatrix{T},
    center_x::Int, center_y::Int,
    radius_cells::T,
    dx::T
)

Initialize the level set with a circular fire of given radius centered at (centerx, centery).

initialize_phi! {#initialize_phi!}

initialize_phi!(phi::AbstractMatrix{T}, ignition_points::Vector{Tuple{Int,Int}}, dx::T)

Initialize the level set field φ with signed distance from ignition points.

φ < 0 inside the fire (burned) φ > 0 outside the fire (unburned)

interpolate_wind_direction

interpolate_wind_direction(wd1::T, wd2::T, f::T) -> T

Interpolate wind direction, handling the 0°/360° wrap-around.

isnonburnable

isnonburnable(fm::FuelModel) -> Bool

Check if a fuel model is non-burnable (fuel model 256 or similar).

level_set_step! {#level_set_step!}

level_set_step!(
    phi::AbstractMatrix{T},
    phi_old::AbstractMatrix{T},
    ux::AbstractMatrix{T},
    uy::AbstractMatrix{T},
    active_cells::AbstractVector{CartesianIndex{2}},
    dt::T,
    dx::T
)

Perform one complete level set propagation step using RK2.

  1. Copy phi to phi_old
  2. Perform RK2 stage 1
  3. Perform RK2 stage 2

limit_gradients

limit_gradients(
    phi::AbstractMatrix{T},
    ux::T, uy::T,
    ix::Int, iy::Int,
    rcellsize::T
) -> Tuple{T, T}

Calculate flux-limited gradients ∂φ/∂x and ∂φ/∂y at cell (ix, iy).

Uses the Superbee flux limiter with upwind differencing based on the velocity direction (ux, uy).

Implementation follows elmfire_level_set.f90:2040-2118.

Arguments

  • phi: Level set field (with boundary padding)
  • ux, uy: Velocity components at this cell
  • ix, iy: Cell indices (1-based, assuming 2-cell padding)
  • rcellsize: Reciprocal of cell size (1/dx)

Returns

  • (dphidx, dphidy): Flux-limited gradients

load_fuel_models

load_fuel_models(::Type{T}, filepath::AbstractString) -> FuelModelTable{T}

Load fuel models from a CSV file and compute all coefficients with specified precision.

The file should be in ELMFIRE fuel_models.csv format.

load_fuel_models(::Type{T}) -> FuelModelTable{T}

Load fuel models from the default ELMFIRE fuel_models.csv bundled with the package.

lognormal_params

lognormal_params(mean::T, normalized_variance::T) -> Tuple{T, T}

Convert mean and normalized variance to lognormal μ and σ parameters.

The normalized variance is defined as variance / mean².

moisture_damping

moisture_damping(r::T) -> T

Calculate the moisture damping coefficient ηm. Formula: ηm = 1 - 2.59r + 5.11r² - 3.52r³

Where r = M/M_ex (moisture ratio to extinction moisture).

Implementation follows Rothermel (1972) equation.

padded_to_grid

padded_to_grid(state::FireState, px::Int, py::Int) -> Tuple{Int, Int}

Convert padded array coordinates to grid coordinates.

parallel_map

parallel_map(f, items::AbstractVector; n_threads::Int = 0)

Apply function f to each item in parallel.

parallel_reduce

parallel_reduce(f, g, items::AbstractVector; init, n_threads::Int = 0)

Apply function f to each item and reduce with g.

Arguments

  • f: Map function to apply to each item
  • g: Reduce function to combine results (must be associative)
  • items: Items to process
  • init: Initial value for reduction

parse_fuel_model_line

parse_fuel_model_line(::Type{T}, line::AbstractString) -> RawFuelModel{T}

Parse a single line from the fuel_models.csv file with specified precision.

CSV format: fuelid, name, dynamic, W01hr, W010hr, W0100hr, W0herb, W0woody, SIG1hr, SIGherb, SIGwoody, delta, mexdead, hoc

Example line: 1,FBFM01,.FALSE.,0.034,0,0,0,0,3500,9999,9999,1,12,8000

perturb_ignition

perturb_ignition(ix::Int, iy::Int, radius::T, ncols::Int, nrows::Int, rng::AbstractRNG) -> Tuple{Int, Int}

Perturb ignition location within a given radius.

perturb_weather

perturb_weather(weather::ConstantWeather{T}, config::PerturbationConfig{T}, rng::AbstractRNG) -> ConstantWeather{T}

Apply stochastic perturbations to weather conditions.

plan_direct_attack

plan_direct_attack(fire_state::FireState{T}) -> Vector{Tuple{Int,Int}}

Plan a direct attack on the current fire perimeter.

plan_indirect_attack

plan_indirect_attack(
    fire_state::FireState{T},
    weather::ConstantWeather{T},
    buffer_distance::Int
) -> Vector{Tuple{Int,Int}}

Plan an indirect attack line ahead of the fire front.

Arguments

  • fire_state: Current fire state
  • weather: Weather conditions for wind direction
  • buffer_distance: Distance ahead of fire front (cells)

Returns

  • Vector of grid cells for the planned line

read_dem

read_dem(::Type{T}, path::String) -> GeoRaster{T}

Read a Digital Elevation Model raster.

read_fuel_raster

read_fuel_raster(::Type{T}, path::String) -> Tuple{Matrix{Int}, GeoMetadata{T}}

Read a fuel model raster (integer fuel IDs) and return the matrix with metadata.

read_geotiff

read_geotiff(::Type{T}, path::String) -> GeoRaster{T}

Read a GeoTIFF file and return a GeoRaster with the specified precision.

Arguments

  • T: Element type for the raster data (e.g., Float64, Float32)
  • path: Path to the GeoTIFF file

read_landscape

read_landscape(::Type{T}, fuel_path, dem_path;
               cbd_path=nothing, cbh_path=nothing, cc_path=nothing, ch_path=nothing
) -> LandscapeData{T}

Read a complete landscape dataset from GeoTIFF files.

Arguments

  • T: Element type for the data
  • fuel_path: Path to fuel model ID raster
  • dem_path: Path to DEM raster
  • cbd_path: Optional path to canopy bulk density raster (kg/m³)
  • cbh_path: Optional path to canopy base height raster (m)
  • cc_path: Optional path to canopy cover raster (fraction 0-1)
  • ch_path: Optional path to canopy height raster (m)

resample_to_match

resample_to_match(source::GeoRaster{T}, target_metadata::GeoMetadata{T}) -> GeoRaster{T}

Resample a raster to match the resolution and extent of another. Uses nearest-neighbor interpolation.

reset! {#reset!}

reset!(state::FireState{T})

Reset a FireState to initial conditions for reuse in ensemble simulations.

rk2_step! {#rk2_step!}

rk2_step!(
    phi::AbstractMatrix{T},
    phi_old::AbstractMatrix{T},
    ux::AbstractMatrix{T},
    uy::AbstractMatrix{T},
    active_cells::AbstractVector{CartesianIndex{2}},
    dt::T,
    dx::T,
    stage::Int
)

Perform one stage of the 2-stage Runge-Kutta time integration.

Stage 1: φ^(1) = φ^n - dt(ux∂φ/∂x + uy∂φ/∂y) Stage 2: φ^(n+1) = 0.5(φ^n + φ^(1) - dt(ux∂φ/∂x + uy*∂φ/∂y))

Implementation follows elmfire_level_set.f90:1957-1983.

run_batch_simulations

run_batch_simulations(
    jobs::Vector{BatchSimulationJob{T}},
    state_template::FireState{T},
    fuel_ids::AbstractMatrix{Int},
    fuel_table::FuelModelTable{T},
    slope::AbstractMatrix{T},
    aspect::AbstractMatrix{T};
    canopy::Union{Nothing, CanopyGrid{T}} = nothing,
    config::SimulationConfig{T} = SimulationConfig{T}(),
    show_progress::Bool = true
) -> Vector{BatchResult{T}}

Run a batch of independent simulations in parallel.

This is useful for running multiple scenarios with different ignition points or weather conditions.

run_ensemble! {#run_ensemble!}

run_ensemble!(
    config::EnsembleConfig{T},
    state_template::FireState{T},
    fuel_ids::AbstractMatrix{Int},
    fuel_table::FuelModelTable{T},
    weather::ConstantWeather{T},
    slope::AbstractMatrix{T},
    aspect::AbstractMatrix{T},
    ignition_ix::Int,
    ignition_iy::Int,
    t_start::T,
    t_stop::T;
    canopy::Union{Nothing, CanopyGrid{T}} = nothing,
    show_progress::Bool = true,
    callback::Union{Nothing, Function} = nothing
) -> EnsembleResult{T}

Run a Monte Carlo ensemble of fire simulations.

Arguments

  • config: Ensemble configuration
  • state_template: Template FireState to clone for each simulation
  • fuel_ids: Matrix of fuel model IDs
  • fuel_table: Fuel model lookup table
  • weather: Base weather conditions (will be perturbed)
  • slope: Slope in degrees
  • aspect: Aspect in degrees
  • ignition_ix, ignition_iy: Grid coordinates of ignition point
  • t_start, t_stop: Simulation time range (minutes)
  • canopy: Optional canopy grid for crown fire
  • show_progress: Show progress bar (default true)
  • callback: Optional callback(member_id, state) called after each member completes

run_ensemble_threaded! {#run_ensemble_threaded!}

run_ensemble_threaded!(
    config::EnsembleConfig{T},
    state_template::FireState{T},
    fuel_ids::AbstractMatrix{Int},
    fuel_table::FuelModelTable{T},
    weather::ConstantWeather{T},
    slope::AbstractMatrix{T},
    aspect::AbstractMatrix{T},
    ignition_ix::Int,
    ignition_iy::Int,
    t_start::T,
    t_stop::T;
    canopy::Union{Nothing, CanopyGrid{T}} = nothing,
    parallel_config::ParallelConfig = ParallelConfig(),
    show_progress::Bool = true,
    callback::Union{Nothing, Function} = nothing
) -> EnsembleResult{T}

Run a Monte Carlo ensemble using multiple threads.

Arguments

Same as run_ensemble! with additional:

  • parallel_config: Configuration for parallel execution

sample_lognormal

sample_lognormal(mu::T, sigma::T, rng::AbstractRNG) -> T

Sample from a lognormal distribution with given μ and σ.

sardoy_parameters

sardoy_parameters(ws::T, flin::T) -> NTuple{4,T}

Calculate spotting distance distribution parameters using the Sardoy (2008) model.

Arguments

  • ws: 20-ft wind speed (mph)
  • flin: Fireline intensity (kW/m)

Returns

  • (mu_dist, sigma_dist, mu_spanwise, sigma_spanwise): Lognormal parameters for downwind and crosswind distributions

simulate! {#simulate!}

simulate!(
    state::FireState{T},
    fuel_ids::AbstractMatrix{Int},
    fuel_table::FuelModelTable{T},
    weather::ConstantWeather{T},
    slope::AbstractMatrix{T},
    aspect::AbstractMatrix{T},
    t_start::T,
    t_stop::T;
    dt_initial::T = one(T),
    target_cfl::T = T(0.45),
    dt_max::T = T(10),
    spread_rate_adj::T = one(T),
    accel_time_constant::T = zero(T),
    diurnal::Union{Nothing, DiurnalConfig{T}} = nothing,
    callback::Union{Nothing, Function} = nothing
)

Run the fire simulation from tstart to tstop.

Arguments

  • state: Fire state (modified in place)
  • fuel_ids: Matrix of fuel model IDs for each cell
  • fuel_table: Fuel model lookup table
  • weather: Weather conditions
  • slope: Slope in degrees for each cell
  • aspect: Aspect direction in degrees for each cell
  • t_start: Start time (minutes)
  • t_stop: Stop time (minutes)
  • dt_initial: Initial timestep (minutes, default 1.0)
  • target_cfl: Target CFL number (default 0.45)
  • dt_max: Maximum timestep (minutes, default 10.0)
  • spread_rate_adj: Spread rate adjustment factor (default 1.0, multiplies base spread rate)
  • accel_time_constant: Fire acceleration time constant (minutes, default 0 = disabled). When > 0.5, spread rate is multiplied by 1 - exp(-t/τ), modeling fire buildup from ignition to steady-state.
  • diurnal: Optional DiurnalConfig for day/night spread rate adjustment (default nothing = disabled)
  • callback: Optional callback function(state, t, dt, iteration) called each timestep

simulate_full! {#simulate_full!}

simulate_full!(
    state::FireState{T},
    fuel_ids::AbstractMatrix{Int},
    fuel_table::FuelModelTable{T},
    weather_interp::WeatherInterpolator{T},
    slope::AbstractMatrix{T},
    aspect::AbstractMatrix{T},
    t_start::T,
    t_stop::T;
    canopy::Union{Nothing, CanopyGrid{T}} = nothing,
    config::SimulationConfig{T} = SimulationConfig{T}(),
    dt_initial::T = one(T),
    target_cfl::T = T(0.45),
    dt_max::T = T(10),
    spread_rate_adj::T = one(T),
    accel_time_constant::T = zero(T),
    diurnal::Union{Nothing, DiurnalConfig{T}} = nothing,
    callback::Union{Nothing, Function} = nothing,
    rng::AbstractRNG = Random.default_rng()
)

Run full fire simulation with crown fire, spotting, and weather interpolation.

Arguments

  • state: Fire state (modified in place)
  • fuel_ids: Matrix of fuel model IDs for each cell
  • fuel_table: Fuel model lookup table
  • weather_interp: Weather interpolator for spatially/temporally varying weather
  • slope: Slope in degrees for each cell
  • aspect: Aspect direction in degrees for each cell
  • t_start: Start time (minutes)
  • t_stop: Stop time (minutes)
  • canopy: Optional canopy properties grid (required if crown fire enabled)
  • config: Simulation configuration (crown fire, spotting settings)
  • dt_initial: Initial timestep (minutes, default 1.0)
  • target_cfl: Target CFL number (default 0.45)
  • dt_max: Maximum timestep (minutes, default 10.0)
  • spread_rate_adj: Spread rate adjustment factor (default 1.0, multiplies base spread rate)
  • accel_time_constant: Fire acceleration time constant (minutes, default 0 = disabled). When > 0.5, spread rate is multiplied by 1 - exp(-t/τ), modeling fire buildup from ignition to steady-state.
  • diurnal: Optional DiurnalConfig for day/night spread rate adjustment (default nothing = disabled)
  • callback: Optional callback function(state, t, dt, iteration) called each timestep
  • rng: Random number generator for stochastic processes

Returns

  • spot_tracker: SpotFireTracker with any remaining pending spot fires

simulate_full_uniform! {#simulate_full_uniform!}

simulate_full_uniform!(
    state::FireState{T},
    fuel_id::Int,
    fuel_table::FuelModelTable{T},
    weather::ConstantWeather{T},
    slope_deg::T,
    aspect_deg::T,
    t_start::T,
    t_stop::T;
    canopy_cbd::T = zero(T),
    canopy_cbh::T = zero(T),
    canopy_cc::T = zero(T),
    canopy_ch::T = zero(T),
    config::SimulationConfig{T} = SimulationConfig{T}(),
    kwargs...
)

Run full simulation with uniform conditions across the domain.

simulate_gpu! {#simulate_gpu!}

simulate_gpu!(
    state::FireState{T},
    fuel_ids::AbstractMatrix{Int},
    fuel_array::FuelModelArray{T},
    weather::ConstantWeather{T},
    slope::AbstractMatrix{T},
    aspect::AbstractMatrix{T},
    t_start::T,
    t_stop::T;
    kwargs...
)

GPU-accelerated fire simulation. Requires loading KernelAbstractions and Adapt:

using KernelAbstractions, Adapt

Uses KernelAbstractions.jl kernels for velocity calculation, CFL reduction, and RK2 level set integration. The narrow band is managed on the CPU, with an active mask uploaded to GPU each timestep.

See also: simulate!, simulate_gpu_uniform!

simulate_gpu_uniform! {#simulate_gpu_uniform!}

simulate_gpu_uniform!(
    state::FireState{T},
    fuel_id::Int,
    fuel_array::FuelModelArray{T},
    weather::ConstantWeather{T},
    slope_deg::T,
    aspect_deg::T,
    t_start::T,
    t_stop::T;
    kwargs...
)

GPU-accelerated simulation with uniform fuel, slope, and aspect. Requires loading KernelAbstractions and Adapt.

See also: simulate_gpu!

simulate_uniform! {#simulate_uniform!}

simulate_uniform!(
    state::FireState{T},
    fuel_id::Int,
    fuel_table::FuelModelTable{T},
    weather::ConstantWeather{T},
    slope_deg::T,
    aspect_deg::T,
    t_start::T,
    t_stop::T;
    kwargs...
)

Run simulation with uniform fuel, slope, and aspect across the domain.

simulate_with_suppression! {#simulate_with_suppression!}

simulate_with_suppression!(
    state::FireState{T},
    suppression::SuppressionState{T},
    fuel_ids::AbstractMatrix{Int},
    fuel_table::FuelModelTable{T},
    weather::ConstantWeather{T},
    slope::AbstractMatrix{T},
    aspect::AbstractMatrix{T},
    t_start::T,
    t_stop::T;
    dt_initial::T = one(T),
    target_cfl::T = T(0.9),
    dt_max::T = T(10),
    callback::Union{Nothing, Function} = nothing
)

Run fire simulation with active suppression.

surface_spread_rate

surface_spread_rate(
    fm::FuelModel{T},
    M1::T, M10::T, M100::T,
    MLH::T, MLW::T,
    wsmf::T,
    tanslp2::T;
    adj::T = one(T)
) -> SpreadResult{T}

Calculate surface fire spread rate using the Rothermel (1972) model.

Arguments

  • fm: Fuel model with pre-computed coefficients
  • M1: 1-hour dead fuel moisture (fraction, e.g., 0.06 for 6%)
  • M10: 10-hour dead fuel moisture (fraction)
  • M100: 100-hour dead fuel moisture (fraction)
  • MLH: Live herbaceous moisture (fraction)
  • MLW: Live woody moisture (fraction)
  • wsmf: Mid-flame wind speed (ft/min)
  • tanslp2: tan²(slope angle)
  • adj: Adjustment factor (default 1.0)

Returns

  • SpreadResult with velocity, reaction intensity, heat per unit area, and fireline intensity

Implementation follows elmfire_spread_rate.f90:13-132.

surface_spread_rate(
    fm::FuelModel{T},
    moisture::NTuple{5,T},
    wsmf::T,
    tanslp2::T;
    adj::T = one(T)
) -> SpreadResult{T}

Convenience method accepting moisture as tuple (M1, M10, M100, MLH, MLW).

surface_spread_rate_flat

surface_spread_rate_flat(
    fuel_arr::FuelModelArray{T}, fuel_index::Int, mc_index::Int,
    M1::T, M10::T, M100::T, MLH::T, MLW::T,
    wsmf::T, tanslp2::T, adj::T
) -> Tuple{T, T, T}

Compute surface fire spread rate from a FuelModelArray using dense array indexing. This variant is designed for GPU kernels where Dict-based FuelModel lookup is not available.

Returns (velocity, vs0, flin) as a tuple instead of a SpreadResult struct.

Arguments

  • fuel_arr: Dense fuel model array
  • fuel_index: Row index in fuel_arr (from fuel_arr.fuel_id_to_index[fuel_id])
  • mc_index: Moisture class index (live_moisture_class - 29)
  • M1, M10, M100, MLH, MLW: Fuel moisture fractions
  • wsmf: Mid-flame wind speed (ft/min)
  • tanslp2: tan²(slope angle)
  • adj: Spread rate adjustment factor

tag_band! {#tag_band!}

tag_band!(nb::NarrowBand, center::CartesianIndex{2}, nx::Int, ny::Int, padding::Int=2)

Add cells within band_thickness of center to the active set.

transport_ember

transport_ember(
    x0::T, y0::T,
    ws20::T,
    wd20::T,
    spotting_distance::T,
    cellsize::T,
    ncols::Int, nrows::Int;
    weather_interp::Union{Nothing, WeatherInterpolator{T}} = nothing,
    time_now::T = zero(T),
    rng::AbstractRNG = Random.default_rng()
) -> Tuple{Int, Int, T}

Transport a single ember from source location to landing location.

When weather_interp is provided, the ember is advected step-by-step through the spatially varying wind field. Otherwise, uses source-cell wind only.

Arguments

  • x0, y0: Source location (in grid units, 1-based)
  • ws20: 20-ft wind speed at source (mph)
  • wd20: Wind direction at source (degrees, FROM convention)
  • spotting_distance: Transport distance (m)
  • cellsize: Cell size (ft)
  • ncols, nrows: Grid dimensions
  • weather_interp: Optional weather interpolator for wind-field advection
  • time_now: Current simulation time (minutes), used with weather_interp
  • rng: Random number generator

Returns

  • (ix, iy, dist): Landing grid indices and actual distance traveled

untag_isolated! {#untag_isolated!}

untag_isolated!(
    nb::NarrowBand,
    phi::AbstractMatrix{T},
    burned::Union{AbstractMatrix{Bool}, BitMatrix},
    padding::Int = 0
)

Remove cells from the active set that are:

  1. Already burned (φ ≤ 0) and have no unburned neighbors
  2. Far from the fire front and won’t affect propagation

Note: The active set uses padded coordinates if padding > 0, but the burned matrix uses grid coordinates. This function handles the conversion.

Uses swap-and-compact on the active list for zero-allocation removal.

update_suppression_state! {#update_suppression_state!}

update_suppression_state!(
    suppression::SuppressionState{T},
    fire_state::FireState{T},
    dt::T,
    t::T
)

Update suppression state: move resources, construct lines, etc.

update_wui_state! {#update_wui_state!}

update_wui_state!(
    wui_grid::WUIGrid{T},
    fire_state::FireState{T},
    weather::ConstantWeather{T},
    t::T,
    dt::T,
    rng::AbstractRNG
) -> Vector{BuildingIgnitionResult{T}}

Update the WUI grid based on fire state and check for building ignitions.

Arguments

  • wui_grid: WUI grid (modified in place)
  • fire_state: Current fire simulation state
  • weather: Weather conditions
  • t: Current time (minutes)
  • dt: Time step (minutes)
  • rng: Random number generator

Returns

  • Vector of new building ignitions during this time step

velocity_at_angle

velocity_at_angle(es::EllipticalSpread{T}, theta::T) -> T

Calculate fire spread rate at angle theta from the wind direction. theta = 0 is head fire (downwind), theta = π is backing fire (upwind).

velocity_components

velocity_components(
    es::EllipticalSpread{T},
    wind_direction_rad::T,
    normal_x::T,
    normal_y::T
) -> Tuple{T, T}

Calculate velocity components (ux, uy) using the Richards (1990) ellipse equation.

Uses semi-major and semi-minor axes derived from head/back velocities and L/B ratio to compute the exact elliptical velocity at any angle relative to the wind direction.

Arguments

  • es: Elliptical spread parameters (head, back, eccentricity, lengthtobreadth)
  • wind_direction_rad: Wind direction in radians (meteorological convention)
  • normal_x, normal_y: Unit normal vector to the fire front

Returns

  • (ux, uy): Velocity components in grid coordinates (ft/min)

wind_adjustment_factor

wind_adjustment_factor(fuel_bed_depth::T) -> T

Calculate unsheltered wind adjustment factor (WAF) from 20-ft wind to mid-flame wind.

Uses: WAF = 1.83 / ln((20 + 0.36H) / (0.13H)) where H is the fuel bed depth in feet.

wind_adjustment_factor(fuel_bed_depth::T, canopy_cover::T, canopy_height::T) -> T

Calculate sheltered wind adjustment factor (WAF) under forest canopy.

Uses the Albini & Baughman (1979) sheltered formula when canopy cover and height are non-trivial, otherwise falls back to the unsheltered formula.

Arguments

  • fuel_bed_depth: Fuel bed depth (ft)
  • canopy_cover: Canopy cover fraction (0-1)
  • canopy_height: Canopy height (m)

write_fire_perimeter

write_fire_perimeter(path::String, state::FireState, metadata::GeoMetadata;
                     format::Symbol=:geojson)

Write the fire perimeter to a GeoJSON or other format.

Arguments

  • path: Output file path
  • state: Fire simulation state
  • metadata: Geospatial metadata for coordinate conversion
  • format: Output format (:geojson supported)

write_geotiff

write_geotiff(path::String, data::AbstractMatrix{T}, metadata::GeoMetadata{T})

Write a matrix to a GeoTIFF file with the given metadata.

write_geotiff(path::String, raster::GeoRaster{T})

Write a GeoRaster to a GeoTIFF file.

Constants

ABSOLUTE_CAP

SpreadRateDampeningMode

Configurable strategy for limiting fire spread rates at extreme wind speeds.

The Rothermel model’s power-law wind factor (phi_w = C * wsmf^B) overpredicts spread rates at high wind speeds. These modes provide physically motivated dampening.

Values

  • NO_DAMPENING: Status quo — only Rothermel’s 0.9*IR wind limit applies
  • WIND_SPEED_CAP: Andrews et al. (2013) — cap ROS at midflame wind speed (ft/min)
  • ABSOLUTE_CAP: Hard cap at a user-specified maximum spread rate (ft/min)
  • LINEAR_DAMPENING: Transition phi_w from power-law to linear above a threshold

BTUPFT2MIN_TO_KWPM2

Float64 <: AbstractFloat <: Real

64-bit floating point number type (IEEE 754 standard). Binary format is 1 sign, 11 exponent, 52 fraction bits. See bitstring, signbit, exponent, frexp, and significand to access various bits.

This is the default for floating point literals, 1.0 isa Float64, and for many operations such as 1/2, 2pi, log(2), range(0,90,length=4). Unlike integers, this default does not change with Sys.WORD_SIZE.

The exponent for scientific notation can be entered as e or E, thus 2e3 === 2.0E3 === 2.0 * 10^3. Doing so is strongly preferred over 10^n because integers overflow, thus 2.0 * 10^19 < 0 but 2e19 > 0.

See also Inf, NaN, floatmax, Float32, Complex.

Float64(x [, mode::RoundingMode])

Create a Float64 from x. If x is not exactly representable then mode determines how x is rounded.

Examples

julia> Float64(pi, RoundDown)
3.141592653589793

julia> Float64(pi, RoundUp)
3.1415926535897936

See RoundingMode for available rounding modes.

FT_TO_M

Float64 <: AbstractFloat <: Real

64-bit floating point number type (IEEE 754 standard). Binary format is 1 sign, 11 exponent, 52 fraction bits. See bitstring, signbit, exponent, frexp, and significand to access various bits.

This is the default for floating point literals, 1.0 isa Float64, and for many operations such as 1/2, 2pi, log(2), range(0,90,length=4). Unlike integers, this default does not change with Sys.WORD_SIZE.

The exponent for scientific notation can be entered as e or E, thus 2e3 === 2.0E3 === 2.0 * 10^3. Doing so is strongly preferred over 10^n because integers overflow, thus 2.0 * 10^19 < 0 but 2e19 > 0.

See also Inf, NaN, floatmax, Float32, Complex.

Float64(x [, mode::RoundingMode])

Create a Float64 from x. If x is not exactly representable then mode determines how x is rounded.

Examples

julia> Float64(pi, RoundDown)
3.141592653589793

julia> Float64(pi, RoundUp)
3.1415926535897936

See RoundingMode for available rounding modes.

LINEAR_DAMPENING

SpreadRateDampeningMode

Configurable strategy for limiting fire spread rates at extreme wind speeds.

The Rothermel model’s power-law wind factor (phi_w = C * wsmf^B) overpredicts spread rates at high wind speeds. These modes provide physically motivated dampening.

Values

  • NO_DAMPENING: Status quo — only Rothermel’s 0.9*IR wind limit applies
  • WIND_SPEED_CAP: Andrews et al. (2013) — cap ROS at midflame wind speed (ft/min)
  • ABSOLUTE_CAP: Hard cap at a user-specified maximum spread rate (ft/min)
  • LINEAR_DAMPENING: Transition phi_w from power-law to linear above a threshold

M_TO_FT

Float64 <: AbstractFloat <: Real

64-bit floating point number type (IEEE 754 standard). Binary format is 1 sign, 11 exponent, 52 fraction bits. See bitstring, signbit, exponent, frexp, and significand to access various bits.

This is the default for floating point literals, 1.0 isa Float64, and for many operations such as 1/2, 2pi, log(2), range(0,90,length=4). Unlike integers, this default does not change with Sys.WORD_SIZE.

The exponent for scientific notation can be entered as e or E, thus 2e3 === 2.0E3 === 2.0 * 10^3. Doing so is strongly preferred over 10^n because integers overflow, thus 2.0 * 10^19 < 0 but 2e19 > 0.

See also Inf, NaN, floatmax, Float32, Complex.

Float64(x [, mode::RoundingMode])

Create a Float64 from x. If x is not exactly representable then mode determines how x is rounded.

Examples

julia> Float64(pi, RoundDown)
3.141592653589793

julia> Float64(pi, RoundUp)
3.1415926535897936

See RoundingMode for available rounding modes.

NO_DAMPENING

SpreadRateDampeningMode

Configurable strategy for limiting fire spread rates at extreme wind speeds.

The Rothermel model’s power-law wind factor (phi_w = C * wsmf^B) overpredicts spread rates at high wind speeds. These modes provide physically motivated dampening.

Values

  • NO_DAMPENING: Status quo — only Rothermel’s 0.9*IR wind limit applies
  • WIND_SPEED_CAP: Andrews et al. (2013) — cap ROS at midflame wind speed (ft/min)
  • ABSOLUTE_CAP: Hard cap at a user-specified maximum spread rate (ft/min)
  • LINEAR_DAMPENING: Transition phi_w from power-law to linear above a threshold

PI

Float64 <: AbstractFloat <: Real

64-bit floating point number type (IEEE 754 standard). Binary format is 1 sign, 11 exponent, 52 fraction bits. See bitstring, signbit, exponent, frexp, and significand to access various bits.

This is the default for floating point literals, 1.0 isa Float64, and for many operations such as 1/2, 2pi, log(2), range(0,90,length=4). Unlike integers, this default does not change with Sys.WORD_SIZE.

The exponent for scientific notation can be entered as e or E, thus 2e3 === 2.0E3 === 2.0 * 10^3. Doing so is strongly preferred over 10^n because integers overflow, thus 2.0 * 10^19 < 0 but 2e19 > 0.

See also Inf, NaN, floatmax, Float32, Complex.

Float64(x [, mode::RoundingMode])

Create a Float64 from x. If x is not exactly representable then mode determines how x is rounded.

Examples

julia> Float64(pi, RoundDown)
3.141592653589793

julia> Float64(pi, RoundUp)
3.1415926535897936

See RoundingMode for available rounding modes.

PIO180

Float64 <: AbstractFloat <: Real

64-bit floating point number type (IEEE 754 standard). Binary format is 1 sign, 11 exponent, 52 fraction bits. See bitstring, signbit, exponent, frexp, and significand to access various bits.

This is the default for floating point literals, 1.0 isa Float64, and for many operations such as 1/2, 2pi, log(2), range(0,90,length=4). Unlike integers, this default does not change with Sys.WORD_SIZE.

The exponent for scientific notation can be entered as e or E, thus 2e3 === 2.0E3 === 2.0 * 10^3. Doing so is strongly preferred over 10^n because integers overflow, thus 2.0 * 10^19 < 0 but 2e19 > 0.

See also Inf, NaN, floatmax, Float32, Complex.

Float64(x [, mode::RoundingMode])

Create a Float64 from x. If x is not exactly representable then mode determines how x is rounded.

Examples

julia> Float64(pi, RoundDown)
3.141592653589793

julia> Float64(pi, RoundUp)
3.1415926535897936

See RoundingMode for available rounding modes.

WIND_SPEED_CAP

SpreadRateDampeningMode

Configurable strategy for limiting fire spread rates at extreme wind speeds.

The Rothermel model’s power-law wind factor (phi_w = C * wsmf^B) overpredicts spread rates at high wind speeds. These modes provide physically motivated dampening.

Values

  • NO_DAMPENING: Status quo — only Rothermel’s 0.9*IR wind limit applies
  • WIND_SPEED_CAP: Andrews et al. (2013) — cap ROS at midflame wind speed (ft/min)
  • ABSOLUTE_CAP: Hard cap at a user-specified maximum spread rate (ft/min)
  • LINEAR_DAMPENING: Transition phi_w from power-law to linear above a threshold