Docs - DOCS - Rewrite the user documentation around the pipeline, its physics and its outputs

.claude/agents/regression-guardian.md

@@ -35,7 +35,7 @@ Useful invocations:
 - `--list-cases` — enumerate available cases. **Always run this first** to get the current set;
   cases are added over time, so never assume the list. (At time of writing: `diiid_n1`,
   `diiid_n1_riccati`, `diiid_slayer_n1`, `efit_fixedbdy_separatrix`, `gal_resistive_diiid`,
-  `gal_resistive_pe`, `ggj_ray_q500i`, `ggj_reference`, `solovev_n1`, `solovev_multi_n`,
+  `ggj_ray_q500i`, `ggj_reference`, `solovev_n1`, `solovev_multi_n`,
   `solovev_kinetic_nuzero`, `tj_n1`, `tj_kinetic_calculated`, `tj_kinetic_multiion`,
   `tj_kinetic_ntv`.)
 - `--cases <case[,case]> --refs <refA,refB>` — compare two refs (branches/commits), e.g.

README.md

@@ -1,12 +1,13 @@
 # GPEC — Generalized Perturbed Equilibrium Code
 
-GPEC is a Julia implementation of the [GPEC suite](https://github.com/PrincetonUniversity/GPEC) for magnetohydrodynamic (MHD) analysis of fusion plasmas. The Fortran suite (DCON, VACUUM, STRIDE, GPEC, PENTRC and SLAYER) has long been used to model the response of tokamak plasmas to small non-axisymmetric magnetic perturbations, from error-field correction and resonant magnetic perturbation optimization to kink and tearing stability. This implementation consolidates them into one code base of cooperating physics modules that share data in memory behind a unified, modern interface, which removes the hand-offs between separate codes and makes parameter scans, uncertainty quantification and integrated modeling practical.
+This Julia implementation of GPEC continues the legacy of the Fortran [GPEC suite](https://github.com/PrincetonUniversity/GPEC) for magnetohydrodynamic (MHD) analysis of fusion plasmas. The Fortran suite (DCON, VACUUM, STRIDE, GPEC, PENTRC and SLAYER) has long been used to model the response of tokamak plasmas to small non-axisymmetric magnetic perturbations, from error-field correction and resonant magnetic perturbation optimization to kink and tearing stability. This implementation consolidates them into one code base of cooperating physics modules that share data in memory behind a unified, modern interface, which removes the hand-offs between separate codes and makes parameter scans, uncertainty quantification and integrated modeling practical.
 
-It performs a three-stage analysis pipeline:
+One input deck runs the stages it names, in order, into one output file:
 
-1. **Equilibrium reconstruction** from experimental data (EFIT, CHEASE) or analytical models — solves the Grad-Shafranov equation, computes flux surfaces, safety factor q(ψ), and builds field quantity splines.
-2. **Ideal MHD stability analysis** via Newcomb's criterion (DCON-style ODE integration) — identifies singular surfaces where q = m/n, computes force-free eigenfunctions ξ(ψ,θ), evaluates the tearing stability parameter Δ' at each rational surface, and checks fixed/free boundary stability.
-3. **Perturbed equilibrium** — quantifies the self-consistent plasma response to external magnetic perturbations (RMPs, error fields, correction coils). Key outputs include resonant flux, island half-widths, and the Chirikov overlap parameter — the principal diagnostics for ELM suppression, error field correction, and disruption avoidance.
+1. **Equilibrium** — reads an EFIT or CHEASE reconstruction, an IMAS data dictionary or an analytic model, maps it to straight-field-line flux coordinates, and finds the rational surfaces q = m/n.
+2. **Stability** — ideal MHD stability by Newcomb's criterion (the DCON integration, or the STRIDE Riccati and RDCON Galerkin formulations): the force-free states ξ(ψ,θ), the free-boundary kink energies with the VACUUM response, and the tearing stability matrix Δ′.
+3. **Perturbed equilibrium** — the plasma response to external non-axisymmetric fields (RMPs, error fields, correction coils): resonant fields, island half-widths and the Chirikov overlap parameter.
+4. **Torque, error fields and tearing** — the neoclassical toroidal viscosity torque (PENTRC), coil-tolerance and locking-risk analysis, and tearing growth rates and critical resonant fields from inner-layer models (SLAYER, GGJ).
 
 Full documentation: [https://openfusiontoolkit.github.io/GPEC/dev/](https://openfusiontoolkit.github.io/GPEC/dev/)
 

docs/DOC_STANDARD.md

@@ -13,32 +13,26 @@ render it to the website. It is a contributor-facing reference, read here in the
 repository — which is why `developer_notes.md` links to it by GitHub URL rather
 than as a site-internal page.
 
-The `Ballooning and Mercier Local Stability` and `Inner Layer Module` pages are
+The `Stability` and `Local Stability` pages are
 the reference exemplars; new or reworked module pages should follow their shape.
 
 ## Module-page template
 
 A module reference page (`docs/src/<module>.md`) has these sections, in order.
 Not every section is large for every module, but the skeleton is the same.
 
-1. **Overview** — what the module computes and where it sits in the
-   equilibrium → stability → perturbed-equilibrium pipeline.
-2. **Governing equations** — the physics the module discretizes, written as
-   numbered display equations. Each equation or block cites the paper it comes
-   from (the PDFs in `docs/resources/`), e.g. *(GWP2016 Eq. 11)*. Citing
-   equation numbers is required, not optional — it is what makes the
-   implementation auditable against theory.
-3. **Numerical method** — the discretization and algorithm: what is expanded in
-   what basis, how the linear/eigenvalue/BVP problem is assembled and solved,
-   and the *design choices* that matter (why this contour, this ordering, this
-   precision). This is the part that distinguishes a methods page from an API
-   dump.
-4. **Validation & benchmarks** — evidence the method is correct: a table of
-   pinned values cross-checked against an independent code or analytic result,
-   and figures (see below) showing convergence, invariance, or agreement.
-5. **Practical usage** — the configuration keys and the public entry points a
-   user calls, with a short worked example.
-6. **API Reference** — an `@autodocs` block for the module so every exported
+1. **What it computes** — the physics of the stage in at most two paragraphs,
+   and where it sits in the pipeline (a link back to the Workflow page).
+2. **Governing equations** — the equations the module discretizes, as display
+   equations, each citing its paper by "Author (year)" with a link to the
+   citations page; equation numbers where the paper is in `docs/resources/`.
+3. **Running it** — the deck section(s), as a snippet copied from a shipped
+   example that is named, and the matching scripting-API call.
+4. **Outputs** — a table of the `gpec.h5` datasets the stage writes, by full
+   path, with meaning and units.
+5. **Numerical notes** — the discretization and algorithm choices that matter
+   to a user, kept short, with the validation figures (see below).
+6. **API reference** — an `@autodocs` block for the module so every exported
    docstring is rendered and `checkdocs=:exports` stays satisfied.
 
 ## Figures

docs/development/regression-harness.md

@@ -205,3 +205,42 @@ If the two refs in a comparison ran under different thread counts, the report fl
 
 - `0` — every run completed (and, with `--fail-on-change`, nothing changed)
 - `1` — a run failed, or a quantity changed under `--fail-on-change`
+
+## Pinned Δ′ and grid convergence
+
+The ideal Δ′ diagonal pinned in `test/runtests_parallel_integration.jl` and tracked by the harness is a
+single-grid snapshot that catches unintended changes, not a converged value. A `psi_accuracy` scan of the
+DIII-D-like SLAYER deck on the automatic grid (`implied` is the knot count `implied_knot_count` asks for
+after the refined solve):
+
+| psi_accuracy | knots | implied | dpm[1,1] | dpm[2,2] | dpm[3,3] | gamma 2/1 (Hz) |
+|---|---|---|---|---|---|---|
+| 2.0e-3   | 172 | –    | 6.850 | -5.783 | -16.308 | 165.3 |
+| 1.0e-3   | 205 | –    | 6.387 | -5.873 | -16.867 | 154.3 |
+| 5.0e-4   | 252 | –    | 8.003 | -6.105 | -16.413 | 192.9 |
+| 2.5e-4   | 307 | 522  | 8.567 | -6.235 | -16.529 | 206.2 |
+| 1.25e-4  | 372 | 666  | 8.146 | -6.285 | -16.649 | 196.2 |
+| 6.25e-5  | 457 | 916  | 8.911 | -6.477 | -16.510 | 214.4 |
+| 3.125e-5 | 559 | 1226 | 9.514 | -6.260 | -16.485 | 228.7 |
+
+The same deck on a fixed `core_edge_packed` grid:
+
+| mpsi | dpm[1,1] | dpm[2,2] | dpm[3,3] | gamma 2/1 (Hz) |
+|---|---|---|---|---|
+| 128  | 5.208 | -8.666 | -19.157 | 126.2 |
+| 256  | 7.097 | -5.187 | -15.729 | 171.2 |
+| 512  | 8.489 | -6.217 | -16.517 | 204.4 |
+| 1024 | 8.932 | -6.298 | -16.529 | 214.9 |
+| 2048 | 8.914 | -6.385 | -16.468 | 214.5 |
+
+- Convergence is per surface: the q = 4 `dpm[3,3]` settles near -16.5 on both grid families, while the q = 2
+  `dpm[1,1]` converges to about 8.91 only on the fixed grid and keeps moving on the automatic one.
+- The SLAYER growth rate is linear in Δ′ here (`gamma/dpm[1,1]` ≈ 24.1 within 0.5%), so a converged Δ′ gives a
+  converged growth rate.
+- The automatic grid stops after its second pass, and the knot count it implies grows faster than the grid in
+  use as `psi_accuracy` tightens, so tightening moves it further from self-consistency.
+
+A plateau criterion (scan `psi_accuracy` and the truncation controls, pin the stationary value per surface)
+would be a real convergence statement; it needs either the fixed grid or knot refinement iterated to a fixed
+point, and is not implemented. Until then, read a pinned diagonal as a sign and order-of-magnitude signal and
+re-pin it when the equilibrium grid changes.

docs/make.jl

@@ -28,15 +28,15 @@ makedocs(;
         "Workflow" => "workflow.md",
         "Scripting API" => "api.md",
         "Conventions Reference" => "conventions.md",
-        "API Reference" => [
-            "Vacuum" => "vacuum.md",
+        "Modules" => [
             "Equilibrium" => "equilibrium.md",
-            "Stability Analysis" => "stability.md",
+            "Vacuum" => "vacuum.md",
+            "Stability" => "stability.md",
             "Galerkin Solver" => "galerkin.md",
-            "Ballooning Local Stability" => "ballooning.md",
-            "KineticForces" => "kinetic_forces.md",
+            "Local Stability" => "ballooning.md",
             "Forcing Terms" => "forcing_terms.md",
             "Perturbed Equilibrium" => "perturbed_equilibrium.md",
+            "Kinetic Forces" => "kinetic_forces.md",
             "Error Fields" => "error_fields.md",
             "Tearing" => "inner_layer.md",
             "Analysis" => "analysis.md",

docs/src/analysis.md

@@ -1,14 +1,18 @@
-# Analysis Module
+# Analysis
 
-The Analysis module provides post-processing and visualization utilities for GPEC simulation outputs.
+`Analysis` holds the plotting and post-processing functions for GPEC results, grouped by the stage they read; each works on the in-memory result or on `gpec.h5`.
 
 ## Submodules
 
-- `ForceFreeStates`: Plotting functions for DCON-style ideal MHD stability results
-- `Equilibrium`: Plotting functions for equilibrium objects
-- `CoilForcing`: Plotting functions for coil geometry and normal field spectra
-- `PerturbedEquilibrium`: Plotting functions for perturbed equilibrium results
-- `PerturbedEquilibriumModes`: Data helpers to convert modal GPEC output to (ψ, θ) and (ψ, θ, φ) grids
+- `ForceFreeStates`: the eigenmode energies, the `crit` scan and the displacement profiles of a stability solve
+- `Equilibrium`: flux surfaces, profiles and the radial grid of an equilibrium
+- `CoilForcing`: coil geometry and the normal-field spectra it drives
+- `PerturbedEquilibrium`: the plasma response, the resonant fields, island widths and Chirikov parameter against ψ_N or q
+- `PerturbedEquilibriumModes`: helpers that convert the modal output to (ψ, θ) and (ψ, θ, φ) grids for real-space plots
+- `ErrorFields`: coil sensitivities, tolerance distributions, locking risk and correction limits (see [Error Fields](error_fields.md))
+
+`examples/DIIID-like_ideal_example/analyze_example.jl`, `examples/Solovev_ideal_example/analyze_example.jl` and `examples/DIIID-like_error_field_example/analyze_example.jl` are worked examples of each; the s–α maps of local stability use `salpha_reference`, `compute_ballooning_stability!` and `scan_delta_prime_map`.
+Figure conventions for spectrum plots are in `docs/development/plotting.md` of the repository, and `write_imas` on the [Utilities](utilities.md) page exports results to an IMAS data dictionary.
 
 ```@docs
 GeneralizedPerturbedEquilibrium.Analysis