Document and review spherical benchmark output/restart paths (write_timestep vs write_checkpoint)

[gthyagi]

Summary

This issue documents the current UW3 mesh/field output methods used in the spherical benchmark workflow and requests feedback on a restart-oriented follow-up design.

The immediate motivation is large spherical benchmark postprocessing, where:

• write_timestep() output works with read_timestep(...)
• but that reload path is expensive because it reconstructs fields from saved arrays and coordinates using a KDTree-based path
• exact restart-style reload needs PETSc restart metadata that is not written by write_timestep()

This issue is for documentation and design review only. No implementation is proposed until feedback is received.

@lmoresi could you please review the comparison below and comment on the proposed direction? I plan to wait for feedback before making implementation changes.

Current Method 1: mesh.write_timestep(...)

What it writes

Typical output:

• one mesh HDF5 file
• one HDF5 file per mesh variable
• optional XDMF file
• ParaView-compatible groups such as /vertex_fields

Typical file set:

• output.mesh.00000.h5
• output.mesh.Velocity.00000.h5
• output.mesh.Pressure.00000.h5
• output.mesh.00000.xdmf

Typical variable-file structure:

• /fields/<name>
• /fields/coordinates
• /vertex_fields/...

How reload works

read_timestep(...) does not perform an exact PETSc restart.

Instead it:

1. opens the variable HDF5 file with h5py
2. reads the saved field values from /fields/<name>
3. reads the saved coordinates from /fields/coordinates
4. builds a KDTree on the saved coordinates
5. reconstructs values on the current mesh variable by coordinate-based matching / interpolation

So this method works because the file contains field values and coordinates. It does not require PETSc section/order metadata.

Advantages

• Produces the standard benchmark output files already used by the repo
• Generates XDMF for ParaView / visualisation
• Works with read_timestep(...)
• Flexible when reloading onto a different mesh or ordering because the reload path is coordinate-based

Disadvantages

• Reload is expensive because read_timestep(...) uses h5py plus KDTree reconstruction
• At large MPI size this can become very memory-heavy
• Does not write the PETSc subDM section/order metadata needed for exact FE vector restart
• Direct PETSc vector load from the /fields/... datasets is not reliable for exact restart

Current Method 2: mesh.write_checkpoint(...)

What it writes

Typical output:

• one mesh HDF5 file
• one restart/checkpoint HDF5 file

Typical file set:

• restart.mesh.0.h5
• restart.checkpoint.00000.h5

The restart file includes PETSc restart metadata under paths like:

• /topologies/uw_mesh/dms/uw_mesh/order
• /topologies/uw_mesh/dms/uw_mesh/section/...
• /topologies/uw_mesh/dms/Velocity/order
• /topologies/uw_mesh/dms/Velocity/section/atlasDof
• /topologies/uw_mesh/dms/Velocity/section/atlasOff
• /topologies/uw_mesh/dms/Velocity/vecs/Velocity/Velocity
• /topologies/uw_mesh/dms/Pressure/order
• /topologies/uw_mesh/dms/Pressure/section/atlasDof
• /topologies/uw_mesh/dms/Pressure/section/atlasOff
• /topologies/uw_mesh/dms/Pressure/vecs/Pressure/Pressure

How reload should work

This is a restart-style path rather than a coordinate-based remap path.

Conceptually it should:

1. rebuild the distributed mesh / subDM layout
2. restore PETSc section/order metadata
3. load the saved FE vector into the distributed PETSc global vector
4. scatter to local FE storage

This is the metadata missing from write_timestep() for exact restart.

Advantages

• Writes the PETSc DM / subDM restart metadata needed for exact restart
• Stores vectors together with section/order information
• Better aligned with exact restart / exact FE field reconstruction
• Avoids relying on coordinate/KDTree reconstruction as the restart mechanism

Disadvantages

• Does not create XDMF output
• Current design places all checkpoint variables into one checkpoint file
• At very large scale, one shared checkpoint file may become an I/O hotspot
• The current restore path is less straightforward than the read_timestep(...) path

Practical Comparison

In the spherical benchmark workflow:

• write_timestep() is useful for visualisation and coordinate-based field reconstruction
• write_checkpoint() is the appropriate format for exact restart-style reload

So the two methods are related, but they serve different purposes:

• write_timestep() = visualisation / remap-oriented output
• write_checkpoint() = restart metadata / exact reload-oriented output

Request For Feedback

The next design step I am considering is to preserve the restart advantages of write_checkpoint() while avoiding a single large checkpoint file and keeping visualisation support available.

I am not starting implementation yet. I would like feedback first on whether that direction is reasonable, especially for large spherical runs on high MPI counts.

[gthyagi]

proposed new method

For large spherical benchmark runs, a cleaner long-term design would be to separate:

• mesh topology
• shared restart metadata
• per-variable restart data
• visualisation payload

A practical target layout would be:

• restart.mesh.0.h5
  • mesh topology / labels / mesh metadata
• restart.layout.00000.h5
  • shared mesh DM restart metadata
• restart.Velocity.checkpoint.00000.h5
  • Velocity section/order metadata
  • Velocity PETSc vector
  • optional /vertex_fields/... for visualisation
• restart.Pressure.checkpoint.00000.h5
  • Pressure section/order metadata
  • Pressure PETSc vector
  • optional /vertex_fields/... for visualisation

In this design:

• exact restart uses PETSc vectors plus section/order metadata
• visualisation uses /vertex_fields/...
• /fields/... would not be required

Why this layout is better

1. Avoids one giant checkpoint file

   The current write_checkpoint() design places all restart variables into a single checkpoint file. At very large scale, such as 1/128 spherical runs on 1000+ ranks, that can create a large shared-file I/O hotspot during both write and reload.

2. Avoids duplicating shared mesh restart metadata

   The mesh-level PETSc restart metadata should only be written once.

3. Keeps exact restart separate from visualisation

   Restart needs PETSc vector + section/order metadata, while visualisation needs projected array data that ParaView/XDMF can consume.

4. Better fit for discontinuous variables

   For discontinuous fields, coordinate-based reconstruction is even less reliable, because values may be element-local and not uniquely determined by spatial location alone.

5. Keeps visualisation lightweight and explicit

   If visualisation is wanted, the per-variable files can carry /vertex_fields/... without also carrying /fields/....

Main advantage

The proposed split layout gives:

• exact PETSc-based restart
• lighter and more modular per-variable files
• no duplication of shared mesh restart metadata
• optional visualisation support

without relying on the expensive KDTree-based read_timestep(...) path.

@lmoresi could you please share feedback on this proposal? I am not planning to implement it until feedback is received.

[gthyagi]

proposed new method update

One concrete shape for the proposed design is that write_checkpoint() would be modified, or a similar new function would be introduced, so that it:

• writes the mesh
• writes one HDF5 file per mesh variable
• writes the layout / metadata required for reconstruction and postprocessing

This keeps the restart metadata path explicit while avoiding one large shared checkpoint file for all variables.

This is still for discussion only. I am not planning to implement it until feedback is received.

[gthyagi]

PR #146 implements the checkpoint/reload path discussed in this issue:

• adds MeshVariable.read_checkpoint(...) for PETSc DMPlex section/vector reload
• updates Mesh.write_checkpoint(...) to support outputPath
• writes PETSc DMPlex HDF5 storage version 3.0.0 checkpoint files
• defaults checkpoint output to one HDF5 file per mesh variable
• preserves the combined-file option with separate_variable_files=False
• documents the distinction between write_timestep/read_timestep and write_checkpoint/read_checkpoint
• includes benchmark evidence for the large spherical postprocessing use case

The PR keeps write_timestep() as the visualisation / coordinate-remap path, and adds write_checkpoint() / read_checkpoint() as the restart/postprocessing path to avoid KDTree-based reload for large MPI runs.

See: #146

[lmoresi]

Thanks for writing this up clearly, Thyagi — the comparison matches what's in the code. A bit of history that explains the current shape, and then the direction we'd like the implementation to go.

Why the two paths exist. The original write_checkpoint was the PETSc section/order path. We rewrote it into write_timestep with KDTree reload because the section roundtrip was unreliable whenever partitioning changed, and sometimes even on plain restart when node numbering wasn't consistent. RBF reload is exact when points coincide and survives intentionally when they don't. That robustness is the thing we cannot give up. read_timestep also mirrors write_timestep 1:1 — symmetry is a deliberate design constraint.

Please confirm the PETSc situation first. Before we commit to maintaining both paths in parallel forever, can you verify what current PETSc does with a DMPlex section roundtrip when (a) the MPI count differs between save and load, (b) repartitioning has occurred, and (c) the same job restarts but produces inconsistent node numbering? If PETSc has fixed all three since we last checked, the calculus around how prominent the KDTree path needs to be may shift. A small standalone reproducer would be useful to settle this.

Target. One per-variable file layout, with the PETSc DM section/order metadata written as additive HDF5 groups inside the existing `.mesh...h5` files. The current `/fields/...` and `/vertex_fields/...` groups stay exactly where they are (so XDMF / ParaView / the KDTree reader keep working unchanged). Section + ordering + the per-DM vector go into a sibling group, using the same `sectionView` / `globalVectorView` calls that `write_checkpoint` already uses today — just issued through one viewer per variable file rather than one combined viewer. That removes the single-file I/O hotspot you flagged and gives every variable file enough metadata for an exact PETSc restart.

One reader with internal dispatch. Keep `read_timestep` as the single public entry point. Internally it should try the fast PETSc path first, detect mismatch (different MPI count, partitioning, node numbering), and fall back to KDTree with a clear log message describing why. The user can refuse the fast path entirely with a flag — something like `read_timestep(..., fast=False)` or `force_kdtree=True` (your preference). Symmetry with `write_timestep` is preserved; the two-mode internal logic is owned by our code, not buried in PETSc behaviour.

Things to avoid:

• Don't replace `write_timestep`'s file layout — augment it. Old files written by the current code must keep working with the KDTree reader, and new files must keep working with ParaView / XDMF.
• Don't let the fast/slow choice be invisible. When the reader falls back from fast to KDTree it must say so in the log with the reason — silent partition mismatches were exactly what motivated the original rewrite, and the user needs to know which path actually ran.
• Don't expose two reader entry points — keep the single `read_timestep` surface and put the dispatch behind a flag the user can flip.
• Don't change the swarm story in this issue. The swarm-proxy variables written today are visualisation only. Proper swarm checkpointing is a separate, also-RBF-based path that we already have plans for and is out of scope here — please leave swarm output untouched in this PR.
• Don't ship without a roundtrip test for the fast path. There's currently zero coverage of `write_checkpoint` + `load_from_h5_plex_vector` — please land that test alongside the writer change so we have a regression net.

`write_checkpoint` itself can become a thin alias to `write_timestep(..., create_xdmf=False)` once the section metadata is being written by `write_timestep`. We can decide on the deprecation timeline separately.

If this direction looks right to you — and once we've seen what current PETSc does on the partitioning-change cases above — go ahead and put up a draft PR. We can sort out the flag naming and the fallback log wording in review.

Underworld development team with AI support from Claude Code

[lmoresi]

Update: the gap this issue identified has landed

The core concern here — 'exact restart-style reload needs PETSc restart metadata that is not written by write_timestep()' — is now resolved on development:

• Use PETSc-native XDMF visualization topology #218 added mesh.write_timestep(..., petsc_reload=True), which writes the PETSc DMPlex section/vector metadata (under /topologies/uw_mesh/dms/...) consumed by MeshVariable.read_checkpoint() for a PETSc-native reload — i.e. the restart-oriented direction this issue proposed. The three modes are now selectable:
  • create_xdmf=True, petsc_reload=False — ParaView/XDMF + KDTree-remap read_timestep
  • create_xdmf=False, petsc_reload=True — exact PETSc-native restart
  • create_xdmf=True, petsc_reload=True — both
• Follow-ups from #218: route snapshot subsystem off deprecated write_checkpoint + XDMF cleanups #252 deprecated write_checkpoint in favour of write_timestep(..., petsc_reload=True) (the snapshot subsystem migration is tracked there).
• The comparison/reload paths are documented in docs/developer/subsystems/checkpoint-output-and-reload-methods.md (already mentions petsc_reload, the output table, and read_timestep vs read_checkpoint).

So both the document and the design-direction halves of this issue are effectively satisfied: the restart metadata is written, selectable per-mode, and documented.

@lmoresi — unless you want a further design pass on the snapshot/restart end-state (that broader work is in the snapshot v1.1 thread), this can be closed as addressed by #218/#252. Happy to close it on your say-so, or leave it as the umbrella for the remaining snapshot design discussion.

[lmoresi]

Closing as addressed: write_timestep(petsc_reload=True) (#218) writes the PETSc restart metadata this issue identified as missing, write_checkpoint deprecated in its favour (#252), and the paths are documented in checkpoint-output-and-reload-methods.md. Broader snapshot/restart design continues in the snapshot v1.1 thread.