Leftovers from the #1353 review: a commuting pair on one instant is refused, and a dropped fire is not asked #1425

Two things the review of #1353 (#1277) found and left. Neither returns a wrong number. Measured on main 15ba102.

1. A pair that commutes is refused across two batches

#1277's refusal of two fires on one instant whose times move differently does not ask whether the two commute, as #49's refusal within one batch does not. Within one batch that is right: the one shift is applied to every row, and the answer is wrong. Across two batches each fire gets its own shift, so a pair that commutes was computed right before #1353, and is now refused.

import bngsim, numpy as np
text = ("species A, B; A = 3; B = 3; q = 0; k = 1\nq' = k; A' = 1; B' = 1\n"
        "E0: at (time >= 2.5): A = A + 1\nE1: at (q >= 2.5): B = B/2\n")
grid = [0.0, 1.0, 2.0, 3.0]
for k in (1 - 1e-6, 1.0, 1 + 1e-6):
    m = bngsim.Model.from_antimony_string(text); m.set_param("k", k)
    r = bngsim.Simulator(m, method="ode").run(sample_times=grid, rtol=1e-10, atol=1e-12)
    print(k, np.asarray(r.species)[-1, :2])
bngsim.Simulator(bngsim.Model.from_antimony_string(text), method="ode",
                 sensitivity_params=["k"]).run(sample_times=grid, rtol=1e-10, atol=1e-12)

Either order gives A(3) = 7 and B(3) = 4.5 − t1/2 with t1 = 2.5/k, so dA/dk = 0 and dB/dk = 1.25. The sensitivity run is refused: "'E0' does not move and 'E1' does. Their order flips across the parameter, so the result jumps there and has no derivative (issues #49, #1277)". Nothing jumps here.

Disjoint writes are not enough to let a pair through. With B' = A in place of B' = 1, the two maps still commute, but E0's bump of A changes the rate B grows at before E1 halves it, and dB(3)/dt1 is −(A+1)/2 in one order and −A/2 in the other: a kink. The test the pair needs is the one #767 asks of an event beside a clock switch: with H1 E1's Jacobian and D(y) = H1·f(y) − f(E1(y)), the two orders agree where D(E0(x)) = H0·D(x), read at both states. Until then the refusal is safe, and only costs runs that were right.

2. A fire dropped from event_fires is not asked

event_fires keeps the last kFiresKept = 64 fires (src/cvode_simulator.cpp). The #1277 check (const double slack = 2.0 * event_root_reach(t_evt)) asks only the fires still kept. The #945 checks beside it refuse where a dropped fire could still be in reach (event_fires.size() >= kFiresKept && t_evt - event_fires.front().t <= reach); this one does not.

To miss, 65 batches have to land within 2·event_root_reach(t) (512 ulp of max(|t|, 1)) of one time. CVODE puts consecutive root batches at least a root tolerance, 100 ulp of |t| + |h|, apart, so at |t| ≥ 1 at most about five fit. Below |t| = 1 the slack stays at 512 ulp of 1 while the root tolerance shrinks with t and h, so a chain of 65 or more events re-arming each other one root tolerance apart could fit. No repro was made; found by reading. The fix is the #945 guard's shape: refuse where the oldest kept fire is within the slack and fires have been dropped.