ARAGOG entropy solver stalls when crossing the first liquidus phase boundary

[timlichtenberg]

Symptom

test_zalmoxis_aragog_calliope_two_timesteps fails within its 2-timestep / 1e3 yr budget because the interior cannot cross the first liquidus phase boundary. In three full validation runs, the CVODE ΔΦ_global cap rootfn fires immediately at t=8.607e2 yr (Phi_global(start)=0.9878), after 54 rootfn evaluations, with the accepted step covering 0.00% of the 45 yr interval requested.

Mechanism

The CVODE rootfn _PhiCapRootFunction.evaluate() (aragog/src/aragog/solver/entropy_solver.py:547-577) returns g[0] as the minimum of up to three margins: the melt-fraction cap, the temperature cap, and the entropy cap. CVODE stops the solve as soon as any one margin reaches zero.

Near a phase boundary a cell's melt fraction reaches phi=1 and its derivative saturates, so the melt-fraction cap goes blind (comment at entropy_solver.py:524-530). When that happens, the temperature or entropy cap becomes the active margin instead, and the melt-fraction cap has no further effect on when the solve stops.

Diagnostic evidence

Three real validation runs grew only the melt-fraction cap, from its nominal value up to 8x (0.1 to 0.8), while the temperature and entropy caps stayed fixed. The stall time, evaluation count, and starting Phi_global were identical across all three runs (t=8.607e2 yr, 54 evals, Phi_global(start)=0.9878). Under a min() combination, an 8x change in one term with zero effect on the result means that term is not the active margin: the temperature or entropy cap is binding, not the melt-fraction cap.

Where a fix would act

Step acceptance in AragogRunner._solve_with_retry (src/proteus/interior_energetics/aragog.py:2268) accepts any terminal-event stop with dt_actual > 0, with no lower bound. A cell that advances by a near-zero fraction of the requested interval is accepted as a completed step rather than retried with a different strategy. A fix needs to act on this acceptance path, or grow the temperature and entropy caps jointly with the melt-fraction cap, not the melt-fraction cap alone.

Scope

This is pre-existing and predates #789; #789 does not introduce it and does not fix it. test_zalmoxis_aragog_calliope_two_timesteps is marked xfail(strict=False) pending this fix.

[timlichtenberg]

Diagnosis of the interior phase-boundary stall

Recording what the stall is, and what can and cannot be its cause, from direct measurement on a reproduced run. The geometry below is measured on the current solver code from the pristine as-stalled snapshot.

Symptom

A coupled run aborts at the first solidus crossing with InteriorStalledError (interior status 21). Over the last 20 outer steps the interior advances a fraction of a percent of the requested time. The abort reproduces on the production configuration and on a truncated one at the same state (t = 2.745e5 yr, global melt fraction about 0.594), so it is not a timestep-truncation artifact. It reproduces on aragog main, so it is not release-specific.

Binding constraint

With the binding-cap diagnostic (FormingWorlds/aragog#35), every root-function fire at the stall (51 of 51) is the entropy step cap. The melt-fraction cap and the temperature cap never fire. Growing the melt-fraction cap does not help, because that cap is not the active constraint. In the released build the three caps share one root function that always prints the melt-fraction label, which can misattribute the stop (FormingWorlds/aragog#30).

What cannot be the cause (ruled out by measurement)

• Integrator stiffness collapse. CVODE net failures are 0 to 1 and internal steps 14 to 27 per solve, and every solve succeeds. The internal step is tiny only because the entropy cap arrests it, not because the integrator is failing.
• An unresolved mush band, at this stall state. At the stall (global melt fraction about 0.59) the mush band at the base is about 6766 J/kg/K, which is 3.36 times wider than one radial cell, and 75 of 79 nodes lie inside it. So the arrest here is not from an unresolved mush band. A sharper single-node front develops later in the freeze and is a separate regime, described in Solidus-crossing near-discontinuity produces single-node solidification-front cliffs and 1-2% energy non-conservation aragog#26.
• The CV_TOO_MUCH_WORK budget-exhaustion path. That is the caps-off behavior (see Partial-melt window exhausts the hardcoded CVODE step budget at production tolerances aragog#27). With the caps on, the stop is the entropy-cap arrest.

The candidate mechanism

The mantle freezes from the base upward, so the front sits at the core-mantle boundary. Two single-cell features there are the candidate drivers of the entropy-cap arrest:

• A solidus corner: the melt fraction rises 0.00, 0.22, 0.52 over the first three cells, a slope discontinuity at the solidus from the hard clip in phi(S), where the latent heat is released. The maximum single-cell entropy jump is about 2016 J/kg/K.
• A rheology cliff: the viscosity drops about 11 orders of magnitude (log10 from 22.0 to 11.1) over the same three cells. The mixing-length flux depends strongly on viscosity, so this is the sharpest flux-stiffness feature.

Either or both drive the front-cell entropy rate high enough to saturate the entropy cap over near-zero time. The cap arrests the step, the outer time-step target collapses, and the stall detector aborts.

Boundary of the claim: which of the two features is the primary driver is a structural inference pending a decomposition test (smooth the corner only, widen the rheology transition only, then both), which also gives the per-cell entropy rate inside the rejected step.

Figures

The interior profile at the stall, and a zoom on the base-of-mantle front. The mush is broadly resolved (75 of 79 nodes lie in it), but the base holds a solidus corner and an order-11 viscosity cliff, each within a single cell.

Relation to the aragog solidus-cliff issue

This stall is an early front state, at the base, on the same hard-clip solidus mapping. A sharper single-node terminal front develops later in the freeze and is a separate regime, described and resolution-tested in FormingWorlds/aragog#26. The two are the same front at different freeze times: the single-cell entropy jump is smaller here (about 2016 J/kg/K) than at the terminal cliff. A durable fix must cross both.

Direction

Regularize the front: smooth the solidus corner and, if needed, widen the rheology transition. The related energy non-conservation is covered in FormingWorlds/aragog#26. The CVODE-statistics logging that enabled these measurements is FormingWorlds/aragog#32.

Related issues

This stall, the solidus cliff it sits on, and the energy-budget diagnostic at the same boundary are one cluster of the interior freezing front, and are tracked together:

• the interior phase-boundary stall, this issue: ARAGOG entropy solver stalls when crossing the first liquidus phase boundary #840
• the solidus-crossing cliff and the 1 to 2 percent energy non-conservation: Solidus-crossing near-discontinuity produces single-node solidification-front cliffs and 1-2% energy non-conservation aragog#26
• the per-call over-count of the core-mantle-boundary flux transient in the energy integrals: Energy budget at the freezing front: over-counted boundary transients and the 1-2% non-conservation step aragog#31

[timlichtenberg]

Cap-based and front-smoothing approaches do not cross; the fix is mesh-level

Following the diagnosis above: neither a per-cell step cap nor a front-smoothing change makes the interior cross this front.

Decisive test: the repro config was run from t=0 with a window-aware entropy cap armed. The cap armed at the liquidus event but reported no cell inside the two-phase window on every arming, the temperature cap was the binding constraint rather than the entropy cap, and the run stalled at the same point as with the cap off (near the liquidus, global melt fraction about 0.99). With no in-window cell to act on, the cap is a no-op here.

The root cause is that the two-phase band is thinner than one radial cell at the front, so a cell crosses the boundary in one step and is never mid-transition. A per-cell cap has nothing to bound, and constitutive front-smoothing is refuted at every width. This is the sub-cell cliff of FormingWorlds/aragog#26 , appearing at the first liquidus as well as the terminal solidus.

The real fix resolves the sub-cell front in the interior solver (a sub-cell interface or front-tracking treatment, an entropy-space remap, or an enthalpy-method reformulation), tracked in FormingWorlds/aragog#26 and being scoped before implementation.

[timlichtenberg]

The arrest is the CVODE cap root function itself. A comparison against SPIDER, which integrates the same staggered entropy formulation, is informative: SPIDER installs no step cap and no CVODE root function (its optional rollback is disabled by default), and it crosses phase boundaries on CVODE's own local error control. The stall is therefore a consequence of the cap in the Aragog path, not a limitation of the shared physics or discretisation.

On a single-column reproducer of the crossing, removing the arrest lets the solver step through the boundary. Both accepting the capped state and continuing, and removing the cap so the integrator steps on its own error control, cross the front. The other Aragog-specific integrator settings are not required to cross on that reproducer.

A sub-cell face-coefficient average was tested as a candidate and does not resolve the arrest: the controlling quantity is the per-cell entropy step the cap acts on, not the face transport coefficient.

Two fixes are under evaluation in isolated branches: (1) accept the capped state and continue, and (2) align the integrator setup with SPIDER. Validation on a real case that exhibits the stall is still required before either is adopted.

[timlichtenberg]

Root-cause and reproduction update, with a warning about a diagnostic path that misled the investigation.

The arrest is the step caps catching a real face-level drain

With the aragog step caps on (the zalmoxis defaults phi 0.1, T 100 K, S 100 J/kg/K), a cap root fires at the crossing and the solve then advances a negligible fraction of its allotted time, which trips the stall detector: the arrest reported here. With the caps off, the same crossing completes as a bottom-cell entropy drain (the cliff of FormingWorlds/aragog#26): the bottom cell drops to the solidus in one step while its neighbours do not move, F_cmb collapses to about 0.16 W/m2, and the run continues. It is the same event with two outcomes, set by whether the caps are on.

Reproduced on main lineage

A fresh t = 0 run on aragog main (01b728d, release integrator settings, caps off, the 1.5 M_earth config) reaches the crossing at about 1.15e5 yr with Phi_global 0.60 and shows the drain there. This matches the production arrest cases (the campaign grid arrests at 2.745e5 yr, Phi 0.594, on the 26.07.04 release wheel): same melt fraction and phenomenology, a different planet and config.

Warning, so the path is not repeated

Diagnostic runs on the relaxed "caps-off" CVODE settings (order 5, phase-boundary max_step clamp removed, atol_nd floored at 1e-6) show an earlier cliff at 2 to 3 kyr with Phi about 1. That early event is a numerical artefact of those settings, not the physical crossing, and it does not reproduce on release settings. A discriminator of four runs across branch and integrator config showed the cliff timing is set by the CVODE settings, not the aragog branch. Any closure fix must be validated against the Phi about 0.6 event on release settings.

Full mechanism and the discriminator table are in FormingWorlds/aragog#26; the integrator-budget angle is in FormingWorlds/aragog#27.

[timlichtenberg]

Closing with the landing.

As established above, the arrest is the step caps catching the bottom-cell drain of FormingWorlds/aragog#26; with the caps off the same crossing completes as the drain. Both halves are now addressed on main (6fd76e2).

The drain is removed by FormingWorlds/aragog#38 (aragog 26.09.06): the viscosity in the gravitational separation velocity can be the bulk viscosity blended through the rheological transition, so separation locks up with the bulk rheology and the bottom cell no longer drains ahead of its neighbours. #848 exposes this as interior_energetics.aragog.separation_viscosity with "mixture" as the coupled default, #849 lowers the grain size default to 1e-3 m, #847 raises the aragog floor. Validation numbers are in FormingWorlds/aragog#26: no collapse at 80 or 160 levels with the caps off, zero cell-drain violations, resolution-independent.

The cap semantics are settled in #850: an explicit 0.0 is rejected at load with a message, an omitted key takes the schema default (which Zalmoxis promotes to a positive cap), -1.0 disables a cap, a positive value is used as given, and the effective caps are logged once per process so a run states whether it ran capped.

Boundary of what was checked: the validation runs were caps-off. The exact configuration of this report, caps on at the Zalmoxis defaults, was not re-run on the landed stack. On main the coupled test test_zalmoxis_aragog_calliope_two_timesteps passes in the nightly for all six cases; the failure in the report was observed on the accretion branch (#789), which needs the new aragog floor when it rebases. If the arrest reappears with the caps on, reopen with the cap log line that names the cap that fired.