Fix shadow observer lifecycle, learning and prediction correctness

custom_components/heating_mpc/physics.py

@@ -26,25 +26,32 @@ class Learner:
     p: list[list[float]] = field(default_factory=lambda: [[0.0005, 0., 0.], [0., 0.8, 0.], [0., 0., 0.08]])
     count: int = 0
     error_sum: float = 0.0
+    on_count: int = 0
+    off_count: int = 0
 
     def predict_delta(self, temp: float, outside: float, heat: float) -> float:
         a, b, c = self.theta
         return a * (outside - temp) + b * heat + c
 
-    def learn(self, temp_before: float, temp_after: float, outside: float, effective_heat: float) -> None:
+    def learn(self, temp_before: float, temp_after: float, outside: float, effective_heat: float,
+              duration_seconds: float = 300, command: bool | None = None) -> bool:
         if not all(isfinite(x) for x in (temp_before, temp_after, outside, effective_heat)):
-            return
-        y = temp_after - temp_before
+            return False
+        if not 0 <= effective_heat <= 1:
+            return False
+        if not isfinite(duration_seconds) or not 270 <= duration_seconds <= 330:
+            return False
+        y = (temp_after - temp_before) * 300 / duration_seconds
         if abs(y) > 0.8:  # sensor glitch/door opened/restarted device
-            return
+            return False
         x = [outside - temp_before, effective_heat, 1.0]
         predicted = sum(t * v for t, v in zip(self.theta, x))
         residual = y - predicted
         px = [sum(self.p[i][j] * x[j] for j in range(3)) for i in range(3)]
         forgetting = 0.999
         denominator = forgetting + sum(x[i] * px[i] for i in range(3))
-        if denominator <= 0:
-            return
+        if not isfinite(denominator) or denominator <= 0:
+            return False
         k = [value / denominator for value in px]
         self.theta = [t + k[i] * residual for i, t in enumerate(self.theta)]
         for i in range(3):
@@ -58,30 +65,58 @@ class Learner:
         ]
         self.count += 1
         self.error_sum += abs(residual)
+        if command is not None:
+            if command:
+                self.on_count += 1
+            else:
+                self.off_count += 1
+        return True
 
     @property
     def mae(self) -> float | None:
         return self.error_sum / self.count if self.count else None
 
     def as_dict(self) -> dict[str, Any]:
-        return {"theta": self.theta, "p": self.p, "count": self.count, "error_sum": self.error_sum}
+        return {"theta": self.theta, "p": self.p, "count": self.count, "error_sum": self.error_sum,
+                "on_count": self.on_count, "off_count": self.off_count}
 
     @classmethod
     def from_dict(cls, source: dict[str, Any]) -> "Learner":
         obj = cls()
         try:
+            if not isinstance(source, dict):
+                return obj
             theta = [float(v) for v in source["theta"]]
             p = [[float(v) for v in r] for r in source["p"]]
             if len(theta) != 3 or len(p) != 3 or any(len(r) != 3 for r in p):
                 return obj
             if not all(isfinite(v) for v in theta + [v for row in p for v in row]):
                 return obj
-            obj.theta, obj.p = theta, p
-            obj.count = max(0, int(source.get("count", 0)))
-            obj.error_sum = max(0., float(source.get("error_sum", 0)))
-        except (KeyError, TypeError, ValueError):
-            return cls()
-        return obj
+            if not (0 <= theta[0] <= .04 and 0 <= theta[1] <= .25 and -.06 <= theta[2] <= .06):
+                return obj
+            # 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):
+                            return obj
+                        lower[i][j] = value ** .5
+                    else:
+                        lower[i][j] = value / lower[j][j]
+            counts = [source.get(key, 0) for key in ("count", "on_count", "off_count")]
+            if any(type(v) is not int or v < 0 for v in counts):
+                return obj
+            count, on_count, off_count = counts
+            error_sum = float(source.get("error_sum", 0))
+            if not isfinite(error_sum) or error_sum < 0 or on_count + off_count > count:
+                return obj
+            return cls(theta, p, count, error_sum, on_count, off_count)
+        except (KeyError, TypeError, ValueError, OverflowError):
+            return obj
 
 
 def simulate_step(temps: dict[str, float], outside: float, lag: float,
@@ -112,9 +147,18 @@ def optimize(
         raise ValueError("Invalid temperature")
     if any(not isfinite(t) for t in temps.values()):
         raise ValueError("Invalid room temperature")
-    if any(float(weights.get(room, 0)) <= 0 for room in temps):
+    if any(not isfinite(float(weights.get(room, 0))) or float(weights.get(room, 0)) <= 0 for room in temps):
         raise ValueError("All rooms need a positive weight")
+    if not isfinite(effective_heat) or not 0 <= effective_heat <= 1:
+        raise ValueError("Invalid heating lag")
+    if any(not isfinite(v) or v < 0 for v in (
+        seconds_since_change, minimum_on_minutes, minimum_off_minutes,
+        switching_cost, heating_cost, comfort_band,
+    )):
+        raise ValueError("Invalid optimizer parameter")
     total_weight = sum(float(weights[room]) for room in temps)
+    if not isfinite(total_weight):
+        raise ValueError("Invalid total weight")
     norm = {room: float(weights[room]) / total_weight for room in temps}
     blocks = horizon_minutes // 15
     nsteps_per_block = 3
@@ -156,7 +200,7 @@ def optimize(
             elapsed_s += 15 * 60
         if feasible and cost < best:
             best = cost
-            result = (bits, at_60 or dict(state), dict(state))
+            result = (bits, at_60 or {}, dict(state))
     if result is None:
         raise RuntimeError("No feasible schedule")
     plan, temps_60m, temps_end = result