Toy problems and the M1 conformance cases

Status: M1 draft, 2026-10-23; §4 (the .esm versions and the M2 end-to-end cases) added 2026-10-25. Normative for both implementations.

M1’s gate (architecture §14) is conformance tiers 2 and 3 on toy problems. The problems are written by hand in each language behind the operator interface (operator-interface.md ), not compiled from .esm documents (D53): at M1 Rust had no public reentrant right-hand side (U4, since released in EarthSciAST 0.5.0). .esm versions of the same problems join the suite at M2 (§4), planned and run end to end, and must give the same results.

Each problem below fixes its equations, operators, layout, parameters, initial state and interval exactly. Both implementations MUST implement every listed operator call (operator-interface.md §3) with hand-derived derivatives, not automatic differentiation, so the derivatives themselves are under test. The Julia implementation’s ForwardDiff Jacobian of rhs is a test oracle for them.

Indices are 0-based. “Cell” coordinates are $x_c = (c + \tfrac12)/\text{ncells}$. The names given to state components ($y$, $\mathrm{O_3}$, $A$, …) are the variable names a run config uses for them, for example in adjoint.functional; a problem with cells has one index set, cell. At M1 every toy leaf runs with dt_scope: global (tape.md covers no other scope until M3), whatever D17’s default for its locality. Toy operators carry no locality, so at M1 an imex partition (spec/run-config.md §4.5) names operators only, and a locality name in it is refused. A problem’s operators exist under the names given here, so a case’s run config names them in its schedule and has no operators key.

1. Problems

1.1 forced_decay

A linear, non-autonomous scalar problem with an exact solution. It tests convergence order for every method, the stage-time arguments $t_n + c_i h$, and Rosenbrock’s $f_t$ term.

• State: $u \in \mathbb{R}^1$, dense. Parameters $p = (k, \omega) = (2, 3)$. Interval $[0, 2]$, $u(0) = 1$.
• One operator, decay: $u’ = -k u + \sin(\omega t)$. depends_on_t = true, linear_in_u = true.
• Exact solution, with $A(t) = \dfrac{k \sin \omega t - \omega \cos \omega t}{k^2 + \omega^2}$: $u(t) = \big(u(0) - A(0)\big) e^{-kt} + A(t)$.

1.2 van_der_pol

Nonlinear, so the second-order terms of Rosenbrock adjoints are nonzero.

• State: $u = (y, z)$, dense. Parameter $p = (\mu)$. $u(0) = (2, 0)$.
• One operator, oscillator: $y’ = z$, $\ z’ = \mu (1 - y^2) z - y$. depends_on_t = false, linear_in_u = false.
• Variants: nonstiff, $\mu = 1$ on $[0, 6.3]$; stiff, $\mu = 1000$ on $[0, 1]$.

1.3 robertson

Classic stiff kinetics with a conservation law.

• State: $u = (y_1, y_2, y_3)$, dense. Parameters $p = (k_1, k_2, k_3) = (0.04,\ 3\times10^{7},\ 10^{4})$. Interval $[0, 40]$, $u(0) = (1, 0, 0)$.
• One operator, kinetics: $y_1’ = -k_1 y_1 + k_3 y_2 y_3$, $\ y_2’ = k_1 y_1 - k_3 y_2 y_3 - k_2 y_2^2$, $\ y_3’ = k_2 y_2^2$.
• Invariant: $y_1 + y_2 + y_3 = 1$.

1.4 photostationary_cells

A batch of independent cells with time-varying photolysis: the shape of reseact’s chemistry. It tests block_diagonal layout and batched block LU, and gives each cell different stiffness over time.

• State: block_diagonal, ncells = 16, n = 3; species order $(\mathrm{NO_2}, \mathrm{NO}, \mathrm{O_3})$, species-major (operator-interface.md §2). Parameters $p = (j_0, k, \omega) = (10^{-2},\ 4.4\times10^{-4},\ 2\pi/3600)$. Interval $[0, 3600]$. Initial state in every cell: $(10, 5, 40)$.
• One operator, photochemistry. With $\varphi_c = 2\pi c/\text{ncells}$ and $j_c(t) = \tfrac12 j_0 \big(1 + \sin(\omega t + \varphi_c)\big)$, and $r_c = j_c,\mathrm{NO_2} - k,\mathrm{NO},\mathrm{O_3}$ in cell $c$: $\mathrm{NO_2}’ = -r_c$, $\ \mathrm{NO}’ = r_c$, $\ \mathrm{O_3}’ = r_c$. depends_on_t = true, linear_in_u = false.
• Invariants per cell: $\mathrm{NO_2} + \mathrm{NO}$ and $\mathrm{NO_2} + \mathrm{O_3}$.

1.5 advection_reaction

Two operators that do not commute, for splitting order, subcycle, and adjoints through sequence and strang.

• Periodic 1-D grid, ncells = 64, $\Delta x = 1/64$. Species $(A, B)$, species-major. Parameters $p = (v, k_1, k_2) = (1, 50, 25)$. Interval $[0, 0.5]$. $A(0) = 1 + \tfrac12 \sin(2\pi x_c)$, $B(0) = 0$.
• Operator transport (dense structure): first-order upwind advection of $A$ only, $A_c’ = -v,(A_c - A_{c-1})/\Delta x$ with $A_{-1} \equiv A_{63}$. $B$ does not move. linear_in_u = true.
• Operator reaction (block_diagonal, n = 2): $A’ = -k_1 A + k_2 B$, $\ B’ = k_1 A - k_2 B$. linear_in_u = true.
• Invariant: $\sum_c (A_c + B_c)$.
• Because $B$ does not advect, the operators do not commute, so Lie splitting has a first-order and Strang a second-order splitting error.

1.6 source_decay_map

A leaf and a map in sequence, with an exact discrete solution. It tests map, bind_dt and the map adjoint.

• State: $u \in \mathbb{R}^1$, dense. Parameters $p = (a, r) = (1, 0.5)$. Interval $[0, 4]$, $u(0) = 0$.
• Operator source: $u’ = a$.
• Map operator decay: $u \leftarrow u, e^{-r \Delta}$, with $\Delta$ the length of the map node’s interval (spec/run-config.md §4.4).
• With schedule sequence[leaf source (ForwardEuler, dt "step"), map decay] and macro step $\Delta$, the exact discrete solution is $u_{n+1} = (u_n + a\Delta), e^{-r\Delta}$.

2. Cases

../tests/conformance/m1/cases.json lists the cases. It is generated by make_cases.py beside it; CI checks the two agree. A case is a run config (spec/run-config.md) whose model is toy:<problem> or toy:<problem>/<variant>, plus what the run config does not carry: the interval (tspan, which the harness passes as the caller’s interval, spec/run-config.md §1.1), output times, the functional for gradient cases, and the checks. The file validates against schemas/cases.schema.json .

model: "toy:..." is legal only in conformance cases. Implementations accept it only from the test harness.

2.1 Checks

A check fails loudly with the measured quantity and the threshold.

3. Result files

For cross_language and for references, an implementation writes one result file per case, validated by schemas/result.schema.json :

The values in this example are illustrative:

{
  "case": "robertson.Rodas4.adaptive",
  "producer": "julia",
  "t": [0, 10, 40],
  "u": [[1, 0, 0], [0.84136, 1.6e-05, 0.15862], [0.71583, 9.2e-06, 0.28416]],
  "stats": {"accepted": 63, "rejected": 2, "rhs_calls": 410},
  "functional": 0.28416,
  "gradient": {"u0": [0.1, 0.2, 0.3], "p": [-1.2, 3.4e-08, -5.6e-05]}
}

• producer is julia, rust, or reference (a high-accuracy OrdinaryDiffEq solve, for cases without an exact solution).
• t is $t_0$ followed by the case’s output times: the times of the case’s saveat list that lie in $(t_0, t_f]$, in increasing order and each once, or $t_f$ alone when the case has no saveat. A listed $t_0$ is not repeated, a listed time outside $[t_0, t_f]$ is dropped, and $t_f$ is not added to a list that omits it; so saveat [0, 1, 2] with tspan [0, 2] gives t = [0, 1, 2]. u holds the full state at each of them, in the layout of §1. Numbers are written in shortest round-trip form.
• functional and gradient are present for gradient cases only.
• stats.accepted and stats.rejected sum over all leaves and macro steps.

Result files under tests/conformance/m1/reference/ are checked in. They are produced by a Julia script with OrdinaryDiffEq (Rodas5P for stiff and Vern9 for non-stiff problems, reltol = abstol = 1e-13); for advection_reaction the reference integrates the unsplit sum of both operators, so order checks measure the splitting error. The script and the files arrive with the Julia M1 implementation. Julia and Rust results are produced by CI and compared there, not checked in.

4. The toy problems as .esm documents (M2)

M2’s runtimes execute plans end to end (D80). The gate replays the M1 cases on .esm versions of the problems: the planner plans each case, the runtime builds its operators from the plan (spec/operator-interface.md §8) and runs it, and the result must match the same implementation’s run of the M1 case on the hand-written problem. Julia and Rust both run them from M2; both build their right-hand sides through the public reentrant rhs_program (U4, EarthSciAST 0.5.0; D53).

4.1 Documents

The documents are in ../tests/conformance/toys/ , esm 1.3.0, each validated by EarthSciAST in conformance tier 0 (spec/plan.md §0.2). Each writes §1’s equations with §1’s parameter values, initial state and interval (solver.time_span), and evaluates to the same right-hand side up to the order of floating-point operations:

source_decay_map is planned only from M4b, when map operators arrive (U6); at M2 its document is validated and its planning is a refusal fixture (unsupported_map_operator).

4.2 Cases

../tests/conformance/m2/cases.json is generated from the M1 cases by make_cases.py beside it; CI checks the two agree. An M2 case is the M1 case, changed as follows:

• id is the M1 case’s id prefixed with m2., so the two suites’ result files (§3) never carry the same case;
• model is the problem’s document, relative to the cases file, and operators defines the problem’s operators by the selectors of §4.1, so the M1 schedule’s operator names resolve;
• every leaf’s control.dt_scope is global, as at M1, where the planner would resolve pointwise leaves to cell (D17);
• adjoint is removed and the only check is toy_equivalence (§4.3): planned operators do not yet provide vjp_p or second-order calls (spec/operator-interface.md §8.3);
• layout lists the variables of §4.1’s last column.

Left out: source_decay_map (maps), and every case with a Rosenbrock leaf on a time-dependent problem (forced_decay, photostationary_cells), since planned operators have no ft at M2 (spec/plan.md §10 item 4; U17 restores them). That leaves 41 of the 51 M1 cases.

4.3 The toy_equivalence check

A runtime writes an M2 case’s result file (§3) with u in the case’s layout: the listed variables’ cells, each variable column-major, in the listed order. That is §1’s layout. The check compares this result file with the same implementation’s result file for the M1 case named by against, by the rule of the cross_language check (§2.1) with the check’s rtol and rtol_steps_differ: states agree to rtol relative to each component’s largest magnitude when the accepted step counts agree, and to rtol_steps_differ when they differ. Only states and step counts are compared; an M1 gradient case’s functional and gradient have no counterpart in its M2 run. The two runs evaluate the same equations in a different order of floating-point operations (the hand-written right-hand side against EarthSciAST’s compiled one, and species order inside a block_lu block), so step sequences may differ in adaptive runs, as between languages.