[wshlavacek]

Closes #1179. This re-vendors NFsim on top of #1447's carry 0021, which put #1419's direct edit to the vendored shim into the patch series; without it the refresh would have reverted that edit.

What was wrong

NFsim read a table function indexed by time only at its x values and after each reaction (#892), and held that value in between. A linear table keeps moving between its x values, so wherever little fired the run lagged the ramp, with no error. From the issue's cold start (100 X(p~1) leaving at a rate that ramps from 0 at t = 0 to 10 at t = 1), nf kept all 100 through t = 1, where the exact mean is 1.75. Main refused such a model on method="nf" instead; that refusal is unreleased.

The fix: carry 0020, amended

Where a linear table indexed by time ramps, System::stepTo now integrates the hazard, as the issue proposed and as bngsim's SSA (#719) and RuleMonkey 3.11 do:

• Placing the firing. Between firings the state is fixed, so the propensity changes with the clock alone. The next firing comes where the propensity, integrated from the current time, reaches an Exp(1) draw. The reaction is picked from the propensities at that point.
• The integral. It uses error-controlled 5-point Gauss-Lobatto panels, ported from RuleMonkey's next_clocked_event_time, which shares bngsim's SSA scheme, with the same tolerances and the same refusal of a singular rate. Each node is a propensity refresh with the clock moved there. No panel is longer than 1/1000 of the ramp.
• Breakpoints. No panel crosses one. A panel that ends on a breakpoint reads an ulp short of it, so a step table that jumps there is read on the side being integrated.
• Draws across calls. A stepTo call that stops inside a ramp leaves the rest of its draw to the next call. A breakpoint drops it, as it already drops a waiting time; that is exact because the waiting time is memoryless.
• Outside every ramp. The propensities are constant between firings, and stepTo runs exactly as before. Same-seed trajectories of every model without a ramping linear table are unchanged.

The carry is landed on bngsim/vendor, patches 0020 and 0021 are re-exported, and the tree is re-vendored. vendor_nfsim.py --verify-clean --strict passes with 21 patches.

The refusal (refuse_nfsim_linear_time_tables and its scan) is removed, along with its unreleased changelog entry. The #892 entry, the user guide and the carry's docs now describe the integration.

Evidence

The issue's cold start, 200 seeds, against exact means:

t	exact	before	now
0.25	73.4	100	73.3 ± 0.3
0.5	29.6	100	29.9 ± 0.3
1.0	1.75	100	1.82 ± 0.09

New tests in test_nfsim_time_table_breakpoints.py compare against an ODE solve of one molecule's probability. That gives the exact mean, because every molecule is independent. They cover:

• the four starts the step-table test uses (cold, near-cold, frequent, and the fixture as written), on outputs at and between the x values
• a file table listed in each order
• a saturating rate law that reads the table through another function, so the quadrature is not exact and its error control has to work
• a step table whose jump falls inside the ramp
• an NfsimSession run in segments that end inside the ramp

With the ramp branch switched off, 11 of these 13 fail. The two that pass are the fixture as written, which fires often enough that the old held value was only about 3 SE off.

test_nfsim.py runs the two linear fixtures again, and checks that a file table reads its method: the linear one interpolates at t = 0.5 and the step one holds.

Cost

Each firing inside a ramp costs a few refreshes of every propensity instead of one. On the fixture with 10,000 molecules, a firing took 10.9 µs with the linear table against 10.2 µs with the step table, and the same within noise over a long tail after the ramp. That model has two rules; a model with many rules pays proportionally more per firing, but only while a linear table ramps.

Testing

The NFsim cluster, the network-free checks, the vendoring and the CI-bound tests pass locally (529). The pre-push suite passed. The Windows NFsim cluster grows to 329, so its floor moves to 321.