Make decayEffect a continuous function of sample timestamp (#33)

* Make decayEffect a continuous function of sample timestamp

Reformulates decayEffect using a closed-form quadratic in time-since-sample
rather than accumulating step-by-step from the floored simulation boundary.
This makes the effect value at any future absolute timestamp independent of
which delta-sized simulation bucket the sample's startDate falls into.

For samples aligned to delta boundaries the two formulations are
mathematically identical. For unaligned samples (the common case with real
CGM streams) the new formulation removes a small discontinuity that the
old code exhibited at bucket boundaries.

Adds LoopMathTests covering continuity across a delta boundary. Existing
fixture-calibrated tests are re-pinned to the new values; per-prediction
drift is on the order of 0.1 mg/dL.

Ports the LoopMath change from LoopKit/LoopKit#556 by Moti Nisenson-Ken to
the LoopAlgorithm package, where decayEffect now lives.

* Space to kick off tests

---------

Co-authored-by: Pete Schwamb <pete@pete-mbp-eth.maplect.net>

Author: ps2

Sources/LoopAlgorithm/LoopMath.swift

@@ -192,19 +192,25 @@ extension GlucoseValue {
         let glucoseUnit = LoopUnit.milligramsPerDeciliter
         let velocityUnit = GlucoseEffectVelocity.perSecondUnit
 
-        // The starting rate, which we will decay to 0 over the specified duration
-        let intercept = rate.doubleValue(for: velocityUnit) // mg/dL/s
-        let decayStartDate = startDate.addingTimeInterval(delta)
-        let slope = -intercept / (duration - delta)  // mg/dL/s/s
+        let firstChange = rate.doubleValue(for: velocityUnit) * delta // mg/dL/s * s = mg/dL
+        let secondChange = firstChange * (1 - delta / (duration - delta))
+
+        // Solve for f(t) = a*t^2 + b*t + c, where t is relative to self.startDate.
+        // f(0) = c
+        // f(delta) - c = firstChange = a*delta^2 + b*delta
+        // f(2*delta) - c = firstChange + secondChange = 4*a*delta^2 + 2*b*delta
+        // --> firstChange - secondChange = 2*a*delta^2
+        let c = quantity.doubleValue(for: glucoseUnit)
+        let a = (secondChange - firstChange) / (2 * delta * delta) // mg/dL/s^2
+        let b = (firstChange + secondChange - 4 * a * delta * delta) / (2 * delta) // mg/dL/s
 
         var values = [GlucoseEffect(startDate: startDate, quantity: quantity)]
-        var date = decayStartDate
-        var lastValue = quantity.doubleValue(for: glucoseUnit)
+        var date = startDate.addingTimeInterval(delta)
 
         repeat {
-            let value = lastValue + (intercept + slope * date.timeIntervalSince(decayStartDate)) * delta
+            let time = min(duration, date.timeIntervalSince(self.startDate))
+            let value = a * time * time + b * time + c
             values.append(GlucoseEffect(startDate: date, quantity: LoopQuantity(unit: glucoseUnit, doubleValue: value)))
-            lastValue = value
             date = date.addingTimeInterval(delta)
         } while date < endDate
 

Tests/LoopAlgorithmTests/Fixtures/carbs_with_isf_change_recommendation.json

@@ -1,5 +1,5 @@
 {
   "manual" : {
-    "amount" : 10.546890782709953
+    "amount" : 10.52269112701204
   }
-}
+}