ADR 0008: Own integrators behind one trait, scipy as oracle

• Status: Accepted (2026-11-10)
• Date: 2026-11-10
• Deciders: Vincent Mader
• Related: Finding F5; M5, M6; E4, E9, E10, E13

Context

The thesis used scipy Radau through coag_py.react_0d, one call per macro step, with default tolerances. Studying sampled kernels needs control over stages, Newton iterations, Jacobian reuse, and step-size control (noise must not drive step rejection). Rust ODE crates exist, but none is known to offer Radau IIA with the hooks we need; this is re-checked in M6 before writing code.

Decision

Integrators are implemented in coag-integrate behind trait Integrator { fn step(&mut self, rhs: &dyn Rhs, t, y, dt) -> StepResult } with a separate step controller. Methods: explicit Euler, Dormand-Prince 5(4), backward Euler with Newton, linearized implicit Euler (Brauer 2008), one Rosenbrock-W method, Radau IIA order 5 following Hairer and Wanner (transformed stage system, simplified Newton, embedded error estimate). The Jacobian is requested from the Rhs in a representation chosen by config (ADR 0014). scipy Radau is used only as a test oracle (ADR 0013).

Alternatives considered

• scipy via Python bridge in production: no wasm, little control, Python overhead in the loop.
• Wrap Fortran RADAU5: fast and proven, but foreign build chain, no wasm, still limited control.
• An existing Rust crate: adopt if one passes the M6 parity tests and exposes Jacobian and stage hooks; otherwise use it as an extra oracle.

Consequences

• Significant implementation and testing effort in M6 (order tests, scipy parity, step statistics).
• Full control for the research questions RQ4 and RQ6.
• Note for write-ups: Radau IIA (3 stages) has order 5, not 4.