From 387929b01b9a8ad39707533228fb943067b3241a Mon Sep 17 00:00:00 2001 From: Andrei Buzera Date: Mon, 26 Oct 2026 21:20:50 +0200 Subject: [PATCH] Bound learner covariance to stop RLS windup With forgetting 0.999 and no heat excitation (summer, long OFF periods), the heat-coefficient variance grew by 1/0.999 per sample without limit: about 8.7e12 after 104 simulated OFF days. The first ON samples could then swing coefficients to their bounds, and accumulated asymmetry approached the absolute 1e-8 restore tolerance, past which a restart silently reset the whole model. Symmetrize P and bound its diagonal by the prior after every update (row/column scaling keeps it positive definite). Restore uses a relative round-off tolerance and applies the same bound, so wound-up v0.2.1 models load instead of being discarded. Payload schema is unchanged. Co-Authored-By: Claude Opus 5.5 --- custom_components/heating_mpc/physics.py | 34 +++++++++++++++++++++--- tests/test_physics.py | 27 ++++++++++++++++++- 2 files changed, 57 insertions(+), 4 deletions(-) diff --git a/custom_components/heating_mpc/physics.py b/custom_components/heating_mpc/physics.py index ab5c620..e455819 100644 --- a/custom_components/heating_mpc/physics.py +++ b/custom_components/heating_mpc/physics.py @@ -12,18 +12,40 @@ from typing import Any STEP_MINUTES = 5 HEAT_LAG = 1.0 - exp(-STEP_MINUTES / 20.0) +# Prior covariance; also the upper bound that stops RLS windup (see _bound_covariance). +PRIOR_P = ((0.0005, 0.0, 0.0), (0.0, 0.8, 0.0), (0.0, 0.0, 0.08)) def bound(value: float, lo: float, hi: float) -> float: return max(lo, min(hi, value)) +def _bound_covariance(p: list[list[float]]) -> None: + """Keep P symmetric and no less certain than the prior, in place. + + With forgetting, unexcited directions (e.g. heat input during a summer with the + relay off) grow by 1/forgetting each step without limit. Scaling row and column i + together is a congruence transform, so positive definiteness is preserved. + """ + for i in range(3): + for j in range(i): + p[i][j] = p[j][i] = (p[i][j] + p[j][i]) / 2 + for i in range(3): + if p[i][i] > PRIOR_P[i][i]: + scale = (PRIOR_P[i][i] / p[i][i]) ** .5 + for j in range(3): + p[i][j] *= scale + for j in range(3): + p[j][i] *= scale + p[i][i] = PRIOR_P[i][i] + + @dataclass class Learner: """Recursive least squares for dT=a*(Tout-T)+b*effective_heat+c.""" theta: list[float] = field(default_factory=lambda: [0.002, 0.045, 0.0]) - p: list[list[float]] = field(default_factory=lambda: [[0.0005, 0., 0.], [0., 0.8, 0.], [0., 0., 0.08]]) + p: list[list[float]] = field(default_factory=lambda: [list(row) for row in PRIOR_P]) count: int = 0 error_sum: float = 0.0 on_count: int = 0 @@ -57,6 +79,7 @@ class Learner: for i in range(3): for j in range(3): self.p[i][j] = (self.p[i][j] - k[i] * px[j]) / forgetting + _bound_covariance(self.p) # Physical constraints; avoids unstable temperature predictions. self.theta = [ bound(self.theta[0], 0.0, 0.04), @@ -94,12 +117,17 @@ class Learner: return obj if not (0 <= theta[0] <= .04 and 0 <= theta[1] <= .25 and -.06 <= theta[2] <= .06): return obj + # Tolerate float round-off relative to magnitude (older versions could + # accumulate it during windup), then apply the current covariance bound. + for i in range(3): + for j in range(i): + if abs(p[i][j] - p[j][i]) > 1e-6 * max(1.0, abs(p[i][j]), abs(p[j][i])): + return obj + _bound_covariance(p) # A covariance must be symmetric positive definite (Cholesky check). lower = [[0.0] * 3 for _ in range(3)] for i in range(3): for j in range(i + 1): - if abs(p[i][j] - p[j][i]) > 1e-8: - return obj value = p[i][j] - sum(lower[i][k] * lower[j][k] for k in range(j)) if i == j: if value <= 0 or not isfinite(value): diff --git a/tests/test_physics.py b/tests/test_physics.py index 77555be..bebd8c3 100644 --- a/tests/test_physics.py +++ b/tests/test_physics.py @@ -90,7 +90,8 @@ class ModelTests(unittest.TestCase): def test_corrupt_model_restore_is_atomic(self): for changes in ({'theta': [-1, .1, 0]}, {'p': [[-1, 0, 0], [0, 1, 0], [0, 0, 1]]}, {'count': float('inf')}, {'count': '5'}, {'error_sum': float('nan')}, - {'off_count': 1}, {'theta': [float('nan'), 0, 0]}, {'p': None}): + {'off_count': 1}, {'theta': [float('nan'), 0, 0]}, {'p': None}, + {'p': [[1, .5, 0], [0, 1, 0], [0, 0, 1]]}): with self.subTest(changes=changes): data = m.Learner().as_dict() data.update(changes) @@ -109,6 +110,30 @@ class ModelTests(unittest.TestCase): with self.subTest(key=key), self.assertRaises(ValueError): m.optimize(**args) + def test_covariance_stays_bounded_without_heat_excitation(self): + learn = m.Learner() + t = 22.0 + for _ in range(30000): # about 104 days with the relay off + value = t + 0.003 * (15 - t) + self.assertTrue(learn.learn(t, value, 15, 0.0, command=False)) + t = value + for i in range(3): + self.assertLessEqual(learn.p[i][i], m.PRIOR_P[i][i]) + for j in range(3): + self.assertEqual(learn.p[i][j], learn.p[j][i]) + restored = m.Learner.from_dict(learn.as_dict()) + self.assertEqual(restored.count, 30000) + self.assertEqual(restored.p, learn.p) + + def test_restore_tolerates_round_off_and_bounds_wound_up_covariance(self): + # Shape of a v0.2.1 model after a summer: huge heat variance, tiny asymmetry. + data = m.Learner(count=40000, off_count=40000).as_dict() + data['p'] = [[0.03, 1.0, 1e-3], [1.0 + 1e-7, 8.7e12, 0.0], [1e-3, 0.0, 0.0012]] + restored = m.Learner.from_dict(data) + self.assertEqual(restored.count, 40000) + self.assertEqual(restored.p[1][1], m.PRIOR_P[1][1]) + self.assertEqual(restored.p[0][0], m.PRIOR_P[0][0]) + self.assertEqual(restored.p[0][1], restored.p[1][0]) if __name__ == '__main__': unittest.main()