Fix invalid symbolic atoms in DDt snapshot copies

[gthyagi]

Review Follow-Up (10 September)

Addressed lmoresi's four requests in 21a84519: cross-process disk regression, warnings for skipped state, documented by-reference limitations, and .array access. The numerical fix below applies to in-memory replay. Fresh-process disk replay still differs because symbolic history is skipped. The new disk regression records this as a strict expected failure, not a passing restart test; setup and exact-restoration failures remain ordinary failures. General symbolic reconstruction is outside this PR. The changelog also describes the broader folded-zero-slot hazard.

Problem

On current development (87091138), saving/restoring a legacy Diffusion
model with symbolic flux history and then solving raises AttributeError:
UWexpression has no is_finite because its sym property is invalid.
The reproduction fails in 4.91 s on the clean upstream build.

Generic copy.deepcopy reconstructs SymPy-derived UWexpression atoms without
their required wrapped state. These atoms are live references to the model's
expression/parameter machinery, not standalone values to reconstruct.

Fix

Give DDtSymbolicState a scoped deepcopy implementation: copy the mutable
history list, timestep metadata and matrix containers, but preserve symbolic
atoms through matrix copies. No change to the integrator, global SymPy
copying behavior, or disk format. Derived from mantle commit 11348b3f and
strengthened with an atom-identity/container-isolation regression.

Validation

Build/check	Result	Time
Clean development: solve after restore	Fails with invalid UWexpression	4.91 s
Clean PR: complete in-memory snapshot file	25 passed	9.96 s
Clean PR: new atom-copy and solve-after-restore tests, 8 MPI ranks	2 passed per rank	approximately 14.6 s

The clean PR was built in an isolated install using the existing Pixi
dependencies; its imported ddt.py was verified against this branch. No
mantle SUPG/geoid changes are included. Tests are small Level 1 regressions.

Original Limitation (Resolved by Follow-Up)

Update: commit d551b6fd resolves the continuation mismatch below by initializing symbolic flux history before compiling the initial residual. Clean-branch validation: 33 serial tests (11.38 s), 10 tests per rank on eight MPI ranks (17.31 s), maximum tight-tolerance replay difference 5.55e-16. The text below records the original investigation; see the follow-up comment for the root cause and changed legacy-solver behavior.

Original Finding

This fixes the exception, not full legacy Diffusion restart fidelity.
An additional experiment with Diffusion(order=2, theta=1), P1 unit-square
mesh (cellSize=0.3), initial T=x, diffusivity 0.05 and dt=0.01 advanced
three steps, captured a snapshot, advanced for a reference, restored, and
advanced again. With this fix, the two temperatures differ by about
0.00289 at maximum. Forced setup did not remove that discrepancy. Its
root cause is not established and is outside this scoped atom-copy fix.
The passing solver regression deliberately asserts a finite solve after
restore, not bit-identical continuation.

@lmoresi please review whether preserving the live symbolic atoms at this
state boundary is the appropriate contract. The separate continuation
discrepancy above remains a follow-up, not a claimed fix in this PR.

[gthyagi]

@lmoresi the previously disclosed numerical continuation mismatch is now fixed in d551b6f (developed and validated on the mantle branch first).

Root cause: checkpoint restoration was exact, including temperature, nodal histories, timestep history and startup coefficients. Legacy Diffusion compiled its first residual before populating symbolic flux-history slots. Zero placeholder terms disappeared from that kernel; restore rebuilt it with populated slots, changing the operator. In the order-2 constant-diffusivity probe, the old cached flux retained only the 5/12 term instead of all AM terms summing to one.

Fix: initialize uninitialized Symbolic flux history after assigning constitutive flux and before kernel compilation (five lines). No checkpoint format changes. This also corrects affected uninterrupted legacy Diffusion trajectories; composed AdvDiffusion, SUPG and SLCN are untouched.

Validation	Before	After
Original cellSize=0.3 probe, original tolerance	max replay difference 0.00289056	5.07e-10
New eight operator/replay regressions	6 failed, 2 passed	all pass
Clean PR serial snapshot + replay suite	-	33 passed, 11.38 s
Clean PR eight-rank focused suite	-	10 passed per rank, 17.31 s
Tight-tolerance maximum replay difference	fails assertions	<=5.55e-16

Small Level 1 tests cover cold/partial/established history, orders 1-3, theta=0.5/1, changing timesteps, exact restored fields, and compiled-vs-live symbolic flux. Clean PR was independently rebuilt using existing dependencies; installed solvers.py and ddt.py match this branch. No A1 jobs or unrelated mantle commits are included. This is restart/operator consistency validation, not a claim of temporal accuracy for every legacy integration option.

Please review the initialization ordering as well as the original symbolic atom-copy contract.

[lmoresi]

Reviewed. We think this should go in — it fixes a crash and a real operator
inconsistency, and the root-cause note on the folded placeholder slots is
exactly right. Three findings, one of them a gap in the validation rather than
in the fix.

The validation is in-memory only, and the disk path drops the same state

psi_star holds sympy matrices, which disk_snapshot._serialise_field does not
handle. It falls through to the unserialisable branch, so a disk snapshot of a
legacy Diffusion carries:

python_state/Symbolic_20/psi_star__skipped = unserialisable list (len=2, first-type=MutableDenseMatrix)

Nothing warns at save time, and on restore the field is simply absent, so the
live object keeps whatever it happens to hold. Measured on a genuine
cross-process restart, same problem shape as the one in the PR (order 2,
theta 1, cellSize 0.3, diffusivity 0.05, three warm steps, then continue over
dt = 0.01 and 0.015):

	result
fields restored from disk	exact, max difference 0.0
restored DFDt metadata	history_initialised=True, n_solves=2, dt_history=[0.01, 0.01]
continuation vs reference	max difference 5.08e-04

So it looks intact and is not. The magnitude is small enough to read as a solver
tolerance and large enough to matter.

Caveat on that number: we measured it on development without this branch
applied. d551b6fd is solver-side, so it may reduce or remove this too. The ask
is not that we trust our number — it is that the PR add a cross-process disk
test alongside the in-memory ones, because that path is currently unexercised
for this flavour.

Cheap mitigation worth taking either way: warn at save time when a state field
is written as __skipped. The information is already in the file; nothing
surfaces it.

The folded-zero-slot failure is a class, not a one-off

The root cause here — a placeholder slot holding zero disappears from the
compiled kernel, so a later rebuild with populated slots is a different
operator — is the same mechanism behind the "rampable constant in exponent
position does not ramp" report in our planning notes. There, the atom is ramped
from zero, x**0 folds to 1 before collection, the atom never reaches
constants[], and it can never take effect afterwards. We reproduced both
halves: built at zero it is absent from constants_manifest and never moves;
built at 0.25 it manifests and ramps.

Worth naming in the changelog entry as the general hazard: an expression
compiled while a slot holds zero has that slot compiled out of existence.
Anywhere we initialise a symbolic or rampable slot to zero before compilation is
the same bug waiting.

On the contract question

Preserving atom identity is the right pragmatic call and the only thing that
works, but we would record it as a documented limitation rather than a contract.
Symbolic is the only DDt flavour whose state carries live objects; every other
one stores names and lets the data travel by the mesh-variable path. Holding the
atoms by reference means a snapshot does not isolate the caller from later
changes to them — fine for backstepping, not a general guarantee, and it does
not survive the disk boundary at all.

The principled alternative is to rebuild the symbolic form from the solver
rather than transport it, which is effectively what the other flavours do. That
is a bigger change and not this PR's job.

Minor

test_symbolic_flux_history_remains_valid_after_solver_restore uses the
deprecated access pattern:

with mesh.access(temperature):
    temperature.data[:, 0] = temperature.coords[:, 0]

The sibling test_1058_diffusion_restart.py uses .array correctly. Worth
making them agree.

Summary of asks

1. Add a cross-process disk restart test for the symbolic flavour.
2. Warn at save time when a state field is dropped as __skipped.
3. Note the by-reference limitation of the symbolic history in the docstring.
4. Use .array in the new in-memory test.

Underworld development team with AI support from Claude Code

[gthyagi]

@lmoresi addressed the four requests in 21a8451, developed first on the mantle branch and then applied alone to this PR.

1. Added a bounded fresh-process disk replay test (Level 2), with singleton or eight-rank workers selected by UW_DIFFUSION_TEST_RANKS. Save/reference and restore/continue are different interpreters.
2. Unsupported snapshot fields now emit RuntimeWarning with the full field path and skipped type; two small Level 1 warning tests cover unsupported scalar objects and matrix lists.
3. DDtSymbolicState now documents that copied containers still share live atoms: useful for backstepping, not mutation isolation or disk portability.
4. Updated the in-memory test to .array. The changelog also names the broader folded-zero symbolic slot hazard.

The disk limitation is confirmed on this PR, including d551b6f. All restored fields and recorded history metadata agree exactly, but continuation differs by 0.00275634 after the first step and 0.00350097 over both steps, identically in serial and eight ranks. Symbolic matrices were skipped while history_initialised was restored true. Thus the in-memory fix does not close disk fidelity.

The disk regression is a strict expected failure, narrowly restricted to numerical continuation mismatch after all ranks' exact restoration checks. Setup errors, incorrect restored data, missing outputs or missing warnings fail normally. A future accurate replay becomes strict XPASS, requiring this limitation marker to be removed. I have not implemented general symbolic reconstruction or unsafe expression deserialization; that remains separate work as you suggested.

Independently rebuilt focused PR	Result	Time
Snapshot + disk-format + in-memory replay + new tests, serial	62 passed, 1 disk xfail	25.11 s
Warning tests and fresh-process disk test, 8-rank workers	2 passed, 1 disk xfail	23.89 s

The description now explicitly distinguishes fixed in-memory replay from uncorrected disk replay. No production benchmark jobs were needed.

[lmoresi]

Review from the 2026-09 style/clutter audit. Findings against the Style Charter, scoped to lines this PR adds.

All three are in the new worker tests/parallel/ptest_1058_diffusion_disk_restart.py.

1. The test compares through .npz sidecars.

np.savez(f"{phase}_{uw.mpi.rank}.npz", **values)

§11: avoid .npz or .csv logging when we have parallel output that is
properly flushed. It is also a standing ruling in this project that we compare
through the checkpoint, never through npz sidecars. That bites hardest here,
because the thing under test is the checkpoint round-trip: the sidecar is
doing the comparison the checkpoint is supposed to be trusted for, so a defect
in what the checkpoint writes can pass while the sidecars agree.

2. phase = sys.argv[1].

§11 names option parsing as the canonical case: use uw.Params / Param, with
-uw_name value overrides through PETSc options — not ad-hoc sys.argv.

3. Numbering. ptest_1058_* joins test_1058_dropped_meshvariable_aux_layout.py
and test_1058_solve_report.py. Three unrelated suites on one number (§8).

On the fix itself. The in-memory path looks right. The question this worker
is well placed to answer, and does not: does the disk path carry the same
state? A cross-process restart is the case where psi_star__skipped would show
up as a wrong answer while the metadata still looks intact. If the sidecars are
what is being compared, that would not be visible.

Underworld development team with AI support from Claude Code

[gthyagi]

@lmoresi implemented the style-review requests in 71c8a12 (developed first on mantle, dd34120).

• Removed npz output and comparisons entirely. The worker now writes model checkpoints at the saved/restored state and at both continuation steps. The parent reads wrapper metadata and follows external_file references into the flushed PETSc checkpoint bulk files. It checks exact initial vector values and layouts, exact saved Python-state metadata, and explicitly asserts psi_star__skipped contains MutableDenseMatrix. Continued checkpoint vectors are compared only after those checks.
• Replaced positional sys.argv parsing with uw.Params and -uw_phase write/resume.
• Renamed in-memory tests to test_1074_diffusion_restart.py, and disk tests/worker to test_1079_diffusion_disk_restart.py / ptest_1079_diffusion_disk_restart.py. Checked these numbers against mantle and fetched upstream development; the old files are removed.

The checkpoint files themselves confirm the disk state gap: restored fields/metadata match exactly, but symbolic history is skipped and continuation differs by 0.0035009713383 in both serial and eight-rank runs. The strict numerical xfail remains explicit; missing data, missing skip markers, bad restored values, or changed layouts are ordinary failures. No disk fidelity claim is added.

Validation	Result	Time
Mantle checkpoint-based disk tests, serial	2 passed, 1 known xfail	16.90 s
Mantle checkpoint-based disk tests, eight-rank workers	2 passed, 1 known xfail	26.98 s
Mantle complete in-memory snapshot/replay tests	33 passed	11.65 s
Focused PR disk + replay tests, serial	10 passed, 1 known xfail	22.79 s

Test-only changes: no new solver build needed. The PR installed checkpoint and solver source still matches the checkout byte-for-byte.