From 343845b58bad0ed414be93bd47daacccdbd4b55b Mon Sep 17 00:00:00 2001 From: Antonio Leblanc Date: Tue, 10 Nov 2026 22:19:35 -0300 Subject: [PATCH 1/3] farsite: ellipse at flank and back, as in the reference --- CHANGELOG.md | 3 ++ src/propagation/FarsitePropagationModel.cpp | 34 +++++++++++++++++++++ tests/python/test_farsite_reference.py | 11 +++++-- 3 files changed, 45 insertions(+), 3 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index be2751e..fadac33 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -39,6 +39,9 @@ history, so they summarise each release rather than list every change. The ### Fixed +- The Farsite model now matches the FARSITE reference (mbedward/farsite) at + the flanks and back of the fire, not only at the head. This changes Farsite + results. ([#NNN]) - The Farsite model reads the wind in m/s and the slope as a tangent, like the other Rothermel models, and slope only speeds up spread uphill. It used to read the wind as if it were in mph and double the slope, so downhill spread diff --git a/src/propagation/FarsitePropagationModel.cpp b/src/propagation/FarsitePropagationModel.cpp index 78b0cf6..247e429 100644 --- a/src/propagation/FarsitePropagationModel.cpp +++ b/src/propagation/FarsitePropagationModel.cpp @@ -29,6 +29,8 @@ namespace libforefire size_t idx_slope; size_t idx_wind; + size_t idx_windU; + size_t idx_windV; size_t idx_m_ones; size_t idx_m_liveh; @@ -88,6 +90,8 @@ namespace libforefire /* defining the properties needed for the model */ idx_slope = registerProperty("slope"); idx_wind = registerProperty("normalWind"); + idx_windU = registerProperty("windU"); + idx_windV = registerProperty("windV"); idx_m_ones = registerProperty("moist.ones"); idx_m_liveh = registerProperty("moist.liveh"); @@ -197,6 +201,18 @@ namespace libforefire slope = s / c; windspd /= c; } + // FARSITE spreads each point as an ellipse along the wind, so the rate + // is computed with the whole wind and then taken at the angle between + // the node's normal and the wind. + double windU = valueOf[idx_windU]; + double windV = valueOf[idx_windV]; + double windNorm = adjustementWind * std::hypot(windU, windV); + double cosWind = 1.0; + if (windNorm > 1e-9) + cosWind = std::max(-1.0, std::min(1.0, windspd / windNorm)); + else + windNorm = windspd; + windspd = windNorm; // moinstures variables double m_ones = valueOf[idx_m_ones]; double m_liveh = valueOf[idx_m_liveh]; @@ -492,6 +508,24 @@ namespace libforefire rate_of_spread = rateo * (1 + phiew); } + // Ellipse: FARSITE's Mechanix::ellipse (Anderson 1983 length-to-breadth, + // Alexander 1985 head-to-back), wind in mph. The rate normal to the + // front at angle theta from the wind is Richards (1990): + // c cos(theta) + sqrt(a^2 cos^2(theta) + b^2 sin^2(theta)). + double windMph = std::max(windspd, 0.0) / 0.44704; + double lb = 0.936 * std::exp(0.1147 * windMph) + + 0.461 * std::exp(-0.0692 * windMph) - 0.397; + lb = std::max(1.0, std::min(8.0, lb)); + double part = std::sqrt(lb * lb - 1.0); + double head = rate_of_spread; + double back = head * (lb - part) / (lb + part); + double ea = (head + back) / 2.0; + double ec = (head - back) / 2.0; + double eb = ea / lb; + double sin2 = 1.0 - cosWind * cosWind; + rate_of_spread = ec * cosWind + + std::sqrt(ea * ea * cosWind * cosWind + eb * eb * sin2); + rate_of_spread = rate_of_spread/60.0; diff --git a/tests/python/test_farsite_reference.py b/tests/python/test_farsite_reference.py index fef3d1d..6f354b7 100644 --- a/tests/python/test_farsite_reference.py +++ b/tests/python/test_farsite_reference.py @@ -23,8 +23,9 @@ is in how ForeFire derives the normal wind and slope it hands over. The test fails when the head rate of spread differs beyond REL_TOL in either -comparison. Flank, back and length-to-breadth are only reported: ForeFire has -no ellipse, so a node not facing the wind spreads at the no-wind rate. +comparison, or the flank or back one does on flat ground. On a slope they are +only reported: the model sees the slope along each node's normal, not the +whole slope vector FARSITE turns its ellipse with. Every case runs in its own interpreter, because a second ``ForeFire()`` in a process inherits the first one's parameters (see test_moisture_invariants.py). @@ -267,9 +268,13 @@ def main(argv=None): "head(model)": rel(h["ROS"], mo["head"]), "head(pipeline)": rel(h["ROS"], pi["head"]), } + if c["slope"] == 0: + errs["flank"] = rel(flank["ROS"], pi["flank"]) + errs["back"] = rel(back["ROS"], pi["back"]) # The spline normals do not point at the head: report that gap only. checked = [k for k in errs - if not (k == "head(pipeline)" and args.scheme == "splines")] + if not (k in ("head(pipeline)", "flank", "back") + and args.scheme == "splines")] bad = [k for k in checked if errs[k] > REL_TOL] if bad: failures.append((c, bad)) From 8f0da62515d2960df2815e9e8e6246f4ae41076c Mon Sep 17 00:00:00 2001 From: Antonio Leblanc Date: Tue, 10 Nov 2026 22:19:49 -0300 Subject: [PATCH 2/3] changelog: link #270 --- CHANGELOG.md | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index fadac33..a155c88 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -41,7 +41,7 @@ history, so they summarise each release rather than list every change. The - The Farsite model now matches the FARSITE reference (mbedward/farsite) at the flanks and back of the fire, not only at the head. This changes Farsite - results. ([#NNN]) + results. ([#270]) - The Farsite model reads the wind in m/s and the slope as a tangent, like the other Rothermel models, and slope only speeds up spread uphill. It used to read the wind as if it were in mph and double the slope, so downhill spread @@ -346,3 +346,4 @@ repository stays easy to return to. [#172]: https://github.com/forefireAPI/forefire/issues/172 [#265]: https://github.com/forefireAPI/forefire/issues/265 [#242]: https://github.com/forefireAPI/forefire/pull/242 +[#270]: https://github.com/forefireAPI/forefire/pull/270 From 0cda9320cce78ac49a15d9d97696399f1d163e57 Mon Sep 17 00:00:00 2001 From: Antonio Leblanc Date: Tue, 10 Nov 2026 22:28:38 -0300 Subject: [PATCH 3/3] farsite: note the slope approximation --- CHANGELOG.md | 4 ++-- src/propagation/FarsitePropagationModel.cpp | 6 ++++-- 2 files changed, 6 insertions(+), 4 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index a155c88..e4606f1 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -40,8 +40,8 @@ history, so they summarise each release rather than list every change. The ### Fixed - The Farsite model now matches the FARSITE reference (mbedward/farsite) at - the flanks and back of the fire, not only at the head. This changes Farsite - results. ([#270]) + the flanks and back of the fire, not only at the head. On slopes it is an + approximation. This changes Farsite results. ([#270]) - The Farsite model reads the wind in m/s and the slope as a tangent, like the other Rothermel models, and slope only speeds up spread uphill. It used to read the wind as if it were in mph and double the slope, so downhill spread diff --git a/src/propagation/FarsitePropagationModel.cpp b/src/propagation/FarsitePropagationModel.cpp index 247e429..9f476a6 100644 --- a/src/propagation/FarsitePropagationModel.cpp +++ b/src/propagation/FarsitePropagationModel.cpp @@ -187,7 +187,7 @@ namespace libforefire //----------------------------------------------------------------------------- // 2) Read actual values from valueOf //----------------------------------------------------------------------------- - double windspd = adjustementWind * valueOf[idx_wind]; // wind times cos(slope) here, m/s + double windspd = adjustementWind * valueOf[idx_wind]; // along the normal, times cos(slope), m/s double slope = valueOf[idx_slope]; // Except with splines, ForeFire puts the slope in the normal's z, so // models get sin(slope) and the wind times cos(slope). FARSITE wants the @@ -203,7 +203,9 @@ namespace libforefire } // FARSITE spreads each point as an ellipse along the wind, so the rate // is computed with the whole wind and then taken at the angle between - // the node's normal and the wind. + // the node's normal and the wind. On a slope this is an approximation: + // the node only sees the slope along its normal, so the ellipse follows + // the wind alone. double windU = valueOf[idx_windU]; double windV = valueOf[idx_windV]; double windNorm = adjustementWind * std::hypot(windU, windV);