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 provide no vjp_p or second-order calls at M2, so they support only the adjoints of spec/operator-interface.md §8.3 (with respect to the initial state), which each implementation’s own tests check;
• 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 step sequences agree, and to rtol_steps_differ when they differ, in count or in any accepted step’s $t$ or $h$. Only states and step sequences are compared, not tape paths; 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.

4.4 Cases with data input (D87)

M2’s runtimes also fetch data and run leaves that read it (spec/operator-interface.md §8.4). Three problems with no M1 counterpart gate this, each with an exact solution, so their checks are against a reference rather than toy_equivalence. Together they exercise everything an esm data source can declare that reseact.esm uses, and the rest of esm-spec §8: data_decay the two-record bracket across a file boundary (this section), data_monthly calendar cadence, a custom substitution and mirrors (§4.5), data_points a points source with codes, a record filter and a discovered extent (§4.6). Their data files are small, checked in beside each model with the generator that writes them and the reference result files; CI regenerates them and compares decoded values (not bytes, which the HDF5 library stamps) and references.

data_decay (tests/conformance/toys/data_decay.esm). Model Decay, domain.temporal.reference_time 2016-01-01T00:00:00Z, index sets x (3) and rec (2):

• state u[x], $u(t_0) = 1$; one operator decay = {"components": "Decay"}, pointwise;
• Decay.F_k[rec, x], a data-fed parameter (default 1.0e-4) from source Rate: kind grid, url_template data/data_decay/rate_{date:%Y%m%d}.nc, temporal {"start": "2016-01-01T00:00:00Z", "frequency": "PT1H", "file_period": "P1D", "records_per_file": 24, "time_variable": "time"}, metadata.esio_format "netcdf", file_variable K, and update.from.records_per_sample 2 (spec/plan.md §4.1);
• Decay.s[x], a data-fed parameter (default 1.0) from source Scale, a static source with no temporal (const), file data/data_decay/scale.nc, file_variable S;
• the observed weight Decay.w $= \mathrm{frac}(t/3600)$, and $u_i’ = -s_i,\big((1 - w),F_k[1, i] + w,F_k[2, i]\big),u_i$.

Data. Two daily files, rate_20160101.nc and rate_20160102.nc, each with K(time, x) for the 24 hourly records of its day (CF time in hours since 2016-01-01 00:00:00), and scale.nc with S(x). All are netCDF-4, a few kilobytes, checked in under tests/conformance/toys/data/data_decay/ with their generator make_data.py, which also writes the reference result file. The values are $K_{h,i} = 10^{-4},(1 + i/2 + (h \bmod 5)/4)$ for hour $h$ counted from the origin and $i = 0, 1, 2$, and $S_i = 1 + i/4$. CI regenerates the files and compares their decoded values (not their bytes, which the HDF5 library stamps) and the reference.

Exact solution. The blend is linear in $t$ within each hour, so with $\kappa_i(t) = s_i,k_i(t)$, $u_i(t) = \exp!\big(-\int_{t_0}^{t} \kappa_i\big)$, the integral a sum of trapezoids $1800,s_i,(K_{h,i} + K_{h+1,i})$ over whole hours plus a partial one. The generator evaluates it in binary64 at each output time.

Cases. The interval is $[79200, 93600]$ (22:00 on 1 January to 02:00 on 2 January), so the run crosses a file boundary at midnight, and forcing epochs begin at 79200 ($t_0$), 82800, 86400, 90000 and 93600 ($t_f$, a refresh time read only by the stages at exactly $t_f$; spec/operator-interface.md §8.4.3), and the last sample needs the 03:00 record of 2 January, which rate_20160102.nc holds; output.saveat is 3600, and step.hit is ["forcing_refresh", "output"] (spec/plan.md §7.3):

The runtimes’ own tests add the check of spec/operator-interface.md §8.4.2 item 3 on this model: each buffer at $t_0$ equals EarthSciAST’s esm_problem forcing buffer for the same providers, bit for bit. They also run a variant whose initial condition reads data, ic(u[i]) = Decay.s[i], and check that the run starts from $u_i(t_0) = S_i$, the data value, not the default 1.0 (spec/operator-interface.md §8.4.5), and ends at $S_i$ times the exact solution above. An adaptive Tsit5 run also exercises the selection of spec/operator-interface.md §8.4.1: its stages with $c = 1$ land on each hourly boundary, where $w = 0$ and only the new epoch’s bracket gives the right rate. A Rosenbrock case waits for time-derivative documents (U17), since decay depends on $t$ through $w$.

4.5 data_monthly: calendar cadence, a custom substitution and mirrors

tests/conformance/toys/data_monthly.esm mounts a component file inventory.esm (model Inv, with source Inventory and the data-fed parameter) through {"ref": "data_monthly/inventory.esm", "url_substitutions": {"sector": "point"}}, so the custom substitution is bound by the model at load (esm-spec §4.7, §8.2). Model Emis, domain.temporal.reference_time 2016-01-01T00:00:00Z, index set x (3):

• state c[x], $c(t_0) = 0$; one operator emis, pointwise; $c_i’ = E_i$;
• Emis.E[x], mapped by param_to_var from the data-fed Inv.E (which carries the default 1.0e-6, since the surviving parameter is the source’s declaration, EarthSciAST#634), from source Inventory: kind grid, url_template absent/inv_{sector}_{date:%Y%m}.nc (no such directory, so the read must fail over), mirrors ["inv_{sector}_{date:%Y%m}.nc"], temporal {"start": "2016-01-01T00:00:00Z", "frequency": "P1M", "file_period": "P1M", "records_per_file": 1, "time_variable": "time"}, metadata.esio_format "netcdf"; file_variable E_day, unit_conversion the expression $1/86400$ (per day to per second), a select taking cells 2–4 of the file’s 5 along its x dimension, one record per sample.
• Data. inv_point_201601.nc, inv_point_201602.nc, inv_point_201603.nc, each E_day(time, x) with one record, $E_{m,i} = 10^{-6},(m + i)$ for month $m$ = 1, 2, 3 and file cell $i$ = 1…5.
• Exact solution. $E$ is constant within each calendar month, so $c_i(t)$ is the sum over months of $E_{m,i+1}/86400$ times the seconds of the month inside $[t_0, t]$.
• Cases. Interval from 2016-01-20 to 2016-03-10 (time given as timestamps, D74), so the run crosses 1 February and 1 March 2016, a 31-day and a 29-day month, and a cadence of fixed length lands on neither; output.saveat 1 day; step.hit ["forcing_refresh", "output"]. m2.data_monthly.tsit5: one adaptive Tsit5 leaf; reference at rtol 1e-10 (the solution is piecewise linear and every month boundary is a tstop), cross_language with rtol_steps_differ.

4.6 data_points: a points source with codes, a record filter and an extent

tests/conformance/toys/data_points.esm. Model Src, index sets x (3) and stn sized by the metaparameter N_STN (default 4):

• state c[x], $c(0) = 0$; one operator sources, pointwise; $c_i’ = \sum_j \text{rate}_j,[\text{cell}_j = i]$, a faq sum over stn of ifelse on the code (it reads no state across cells, so the operator is pointwise);
• source Stations: kind points, no temporal (const), url_template data_points/stations.csv, metadata.esio_format "csv", record_filter {"require_finite": ["rate"]}, extent {"metaparameter": "N_STN"}; data-fed Src.rate[stn] (file_variable rate, default 1.0e-3) and Src.cell[stn] (file_variable cell, codes mapping the text codes "a", "b", "c" to 1, 2, 3, default 1.0).
• Data. stations.csv, six rows, two with a rate of NaN (not blank: EarthSciIO’s CSV readers never read a blank cell as NaN, so a blank would make the column text), so four survive the filter and the discovered extent equals the planned N_STN.
• Exact solution. $c_i(t) = t \sum_{j:,\text{cell}_j = i} \text{rate}_j$.
• Cases. Interval $[0, 3600]$, output.saveat 600. m2.data_points.tsit5: one adaptive Tsit5 leaf; reference at rtol 1e-12, cross_language. The runtimes’ tests add a variant run with N_STN bound to 5 at load, which must be refused before the first step for the extent mismatch (spec/operator-interface.md §8.4.2 item 4).

The Rust CSV reader this case needs arrived with EarthSciIO’s document providers (spec/plan.md §10 item 19).