Repo - BUGFIX! - Make reruns self-contained and flag vanished harness outputs

.gitignore

@@ -44,3 +44,6 @@ profiling/
 
 # Per-developer package preferences (e.g. the precompile_workload switch)
 LocalPreferences.toml
+
+# Python bytecode from the CI convention scripts
+__pycache__/

docs/development/hdf5-conventions.md

@@ -29,7 +29,7 @@ These rules govern `gpec.h5` (and any future GPEC-produced HDF5 output); harness
 
 ## Inputs live only under `Input/`
 
-`Input/gpec_toml_raw` stores the deck that reproduces the run: each result's save adds the sections stating the problem, algorithm and controls it was computed with (`Utilities.record_inputs!`), so a deck run and a script run of the same workflow record the same deck. `Input/RawInputs/` stores the raw equilibrium/forcing/coil data and tolerance file — together they make `gpec.h5` a self-contained rerun snapshot (`Rerun.jl` reconstructs every control struct from them; the writer and rerun reader carry cross-reference comments marking the mirrored path pair). **Never echo TOML flags or control-struct values into any other group** — every group outside `Input/` is derived output. (The former `kinetic/` and `slayer/settings/` echoes were removed under this rule.)
+`Input/gpec_toml_raw` stores the deck that reproduces the run: each result's save adds the sections stating the problem, algorithm and controls it was computed with (`Utilities.record_inputs!`), so a deck run and a script run of the same workflow record the same deck. `Input/RawInputs/` stores the raw equilibrium/forcing/coil data, the kinetic profiles, the wall file and the tolerance file — together they make `gpec.h5` a self-contained rerun snapshot (`Rerun.jl` reconstructs every control struct from them; the writer and rerun reader carry cross-reference comments marking the mirrored path pair). **Never echo TOML flags or control-struct values into any other group** — every group outside `Input/` is derived output. (The former `kinetic/` and `slayer/settings/` echoes were removed under this rule.)
 
 ## Schema
 
@@ -38,7 +38,7 @@ Top level (11 groups):
 | Group | Contents |
 |---|---|
 | `Info/` | Run metadata: `git_version`, mode-number ranges (`mpert`, `mlow`, …, `mn_index`), `psilim`, `qlim`, `Runtimes/` (per-stage wall-clock seconds) |
-| `Input/` | Rerun snapshot: `gpec_toml_raw`, `RawInputs/{Equilibrium, ForcingTerms, Coils/<name>}` |
+| `Input/` | Rerun snapshot: `gpec_toml_raw`, `RawInputs/{Equilibrium, ForcingTerms, Coils/<name>, Kinetic, Wall, ErrorFields}` |
 | `Equilibrium/` | Scalars (`beta_N`, `q_axis`, `q_95`, `I_p`, …) plus `Profiles/` (1-D on `psi`: 2piF, mu0p, dVdpsi, q) and `Geometry/` (2-D on `psi`×`theta`: rcoords, offset, nu, jac) |
 | `ForceFreeStates/` | `Solutions/ForwardIntegration/` (u-solutions), `Solutions/GalerkinIntegration/` (closed ξ profiles in the shared layout, the gal surface list, `gal_basis_output`-gated `Basis/`), `write_el_matrices`-gated `EulerLagrangeMatrices/{Ideal,Kinetic}`, `FreeBoundaryStability/`, `EdgeScan/` |
 | `LocalStability/` | Mercier `D_I`, resistive interchange `D_R`, `ballooning_Delta_prime` on `psi`; the ballooning α boundary on `ballooning_psi` |

docs/development/regression-harness.md

@@ -72,6 +72,11 @@ singular q values            [4 elements]             [4 elements]             0
 Summary: 23 unchanged, 3 missing/N/A
 ```
 
+Each row's status is `OK` (within the case's noise threshold), `** CHANGED **` (moved past it),
+`** CHANGED (missing) **` (on the earlier ref, gone from the later one), `** CHANGED (type) **` (stored
+with a different shape or type), `NEW` (only on the later ref), `N/A` (on neither ref) or `FAILED`.
+Every `CHANGED` status counts toward `--fail-on-change`; `NEW` and `N/A` do not.
+
 **Compare your uncommitted working tree against develop:**
 ```bash
 regress --cases solovev_n1 --refs develop,local
@@ -102,7 +107,7 @@ regress --cases solovev_n1 --ref-range develop~10..develop
 - `--no-instantiate` — skip `Pkg.instantiate()` (faster if deps are already resolved)
 - `--no-pin-manifest` — let each ref resolve its own package set (see below)
 - `--allow-env-mismatch` — reuse cached results produced in a different environment
-- `--fail-on-change` — exit non-zero when any tracked quantity changed
+- `--fail-on-change` — exit non-zero when any tracked quantity changed, disappeared or changed storage type
 
 GPEC subprocesses run with `-t auto` (all cores) so GPEC's threaded kernels are active; set `GPEC_REGRESS_THREADS=1` to force single-threaded runs. Tracked quantities are thread-count independent, and the count each run actually used is recorded in its environment fingerprint (shown in the report's `env:` lines). Thread count is deliberately not part of the cache key, so `Runtime (s)` rows cached from single-threaded runs are not comparable to threaded ones — re-baseline with `--force` if runtime tracking matters.
 

docs/src/workflow.md

@@ -184,7 +184,7 @@ All results are written to a single HDF5 file (default: `gpec.h5`). The top-leve
 | Group | Contents |
 |---|---|
 | `Info/` | Run metadata: git version, mode-number ranges, ψ limit, `Runtimes/` (per-stage wall-clock seconds) |
-| `Input/` | Self-contained rerun snapshot: the recorded deck, raw equilibrium/forcing/coil inputs |
+| `Input/` | Self-contained rerun snapshot: the recorded deck, raw equilibrium/forcing/coil inputs, kinetic profiles and wall file |
 | `Equilibrium/` | Equilibrium scalars (`beta_N`, `q_axis`, `q_95`, …), 1-D profiles (`Profiles/`), 2-D geometry (`Geometry/`) |
 | `ForceFreeStates/` | Stability solve: `Solutions/{ForwardIntegration,GalerkinIntegration}`, `EulerLagrangeMatrices/` (opt-in: `write_el_matrices = true`), `FreeBoundaryStability/`, `EdgeScan/` |
 | `LocalStability/` | Mercier D_I, resistive interchange D_R, ballooning Δ' profiles |

regression-harness/src/extractor.jl

@@ -155,21 +155,26 @@ function _json_element_abs(x)::Float64
 end
 
 """
-Compare two extracted values. Returns (abs_diff, rel_diff, status).
-For arrays, returns max element-wise absolute diff.
+Compare two extracted values, `q1` from the baseline ref and `q2` from the newer one. Returns
+(abs_diff, rel_diff, status). For arrays, returns max element-wise absolute diff. A quantity
+missing from both refs is `N/A`; one that disappears or changes storage type is a change, and
+one that only the newer ref has is `NEW`.
 """
 function compare_values(q1::NamedTuple, q2::NamedTuple)
     vtype = q1.value_type
 
-    # Handle missing values
-    if vtype == "missing" || q2.value_type == "missing"
+    if vtype == "missing" && q2.value_type == "missing"
         return (NaN, NaN, "N/A")
+    elseif q2.value_type == "missing"
+        return (NaN, NaN, "CHANGED (missing)")
+    elseif vtype == "missing"
+        return (NaN, NaN, "NEW")
     end
 
     # Storage type differs between refs (e.g. a quantity whose schema changed from a
     # vector to a scalar across the two commits); the values are not directly comparable.
     if vtype != q2.value_type
-        return (NaN, NaN, "N/A (type changed)")
+        return (NaN, NaN, "CHANGED (type)")
     end
 
     threshold = q1.noise_threshold

regression-harness/src/reporter.jl

@@ -114,11 +114,11 @@ function report_two_ref_comparison(db::SQLite.DB, case_spec::CaseSpec,
 
     if info1 === nothing
         println("ERROR: No results for ref 1 ($(ref1.name))")
-        return (n_ok=0, n_changed=0, n_missing=0, n_failed=1)
+        return (n_ok=0, n_changed=0, n_missing=0, n_new=0, n_failed=1)
     end
     if info2 === nothing
         println("ERROR: No results for ref 2 ($(ref2.name))")
-        return (n_ok=0, n_changed=0, n_missing=0, n_failed=1)
+        return (n_ok=0, n_changed=0, n_missing=0, n_new=0, n_failed=1)
     end
 
     failed1 = !info1.success
@@ -135,6 +135,7 @@ function report_two_ref_comparison(db::SQLite.DB, case_spec::CaseSpec,
     n_ok = 0
     n_changed = 0
     n_missing = 0
+    n_new = 0
 
     # Helper: pick the value-cell text for one ref/quantity, accounting for
     # whole-ref failure (which should render as "FAILED" rather than "N/A").
@@ -179,14 +180,16 @@ function report_two_ref_comparison(db::SQLite.DB, case_spec::CaseSpec,
         q1 = q1_all[qname]
         q2 = q2_all[qname]
         abs_diff, rel_diff, status = compare_values(q1, q2)
-        st = status == "CHANGED" ? "** CHANGED **" : status
+        st = startswith(status, "CHANGED") ? "** $status **" : status
         push!(rows, [spec.label, format_value(q1), format_value(q2),
             format_diff(abs_diff, rel_diff, status, q1.value_type), st])
 
         if status == "OK"
             n_ok += 1
-        elseif contains(status, "CHANGED")
+        elseif startswith(status, "CHANGED")
             n_changed += 1
+        elseif status == "NEW"
+            n_new += 1
         else
             n_missing += 1
         end
@@ -231,10 +234,11 @@ function report_two_ref_comparison(db::SQLite.DB, case_spec::CaseSpec,
     parts = String[]
     n_changed > 0 && push!(parts, "$n_changed changed")
     n_ok > 0 && push!(parts, "$n_ok unchanged")
+    n_new > 0 && push!(parts, "$n_new new")
     n_missing > 0 && push!(parts, "$n_missing missing/N/A")
     println("Summary: ", join(parts, ", "))
     println()
-    return (n_ok=n_ok, n_changed=n_changed, n_missing=n_missing,
+    return (n_ok=n_ok, n_changed=n_changed, n_missing=n_missing, n_new=n_new,
         n_failed=count(identity, (failed1, failed2)))
 end
 
@@ -269,6 +273,7 @@ function report_multi_ref(db::SQLite.DB, case_spec::CaseSpec,
     n_ok = 0
     n_changed = 0
     n_missing = 0
+    n_new = 0
 
     for spec in case_spec.quantities
         qname = spec.name
@@ -316,13 +321,15 @@ function report_multi_ref(db::SQLite.DB, case_spec::CaseSpec,
             q_prev = haskey(all_qs[end-1], qname) ? all_qs[end-1][qname] : nothing
             if q_last !== nothing && q_prev !== nothing
                 abs_diff, rel_diff, status = compare_values(q_prev, q_last)
-                st = status == "CHANGED" ? "** CHANGED **" : status
+                st = startswith(status, "CHANGED") ? "** $status **" : status
                 push!(row, format_diff(abs_diff, rel_diff, status, q_prev.value_type))
                 push!(row, st)
                 if status == "OK"
                     n_ok += 1
-                elseif contains(status, "CHANGED")
+                elseif startswith(status, "CHANGED")
                     n_changed += 1
+                elseif status == "NEW"
+                    n_new += 1
                 else
                     n_missing += 1
                 end
@@ -384,11 +391,12 @@ function report_multi_ref(db::SQLite.DB, case_spec::CaseSpec,
         parts = String[]
         n_changed > 0 && push!(parts, "$n_changed changed")
         n_ok > 0 && push!(parts, "$n_ok unchanged")
+        n_new > 0 && push!(parts, "$n_new new")
         n_missing > 0 && push!(parts, "$n_missing missing/N/A")
         println("Summary (last vs prev): ", join(parts, ", "))
     end
     println()
-    return (n_ok=n_ok, n_changed=n_changed, n_missing=n_missing,
+    return (n_ok=n_ok, n_changed=n_changed, n_missing=n_missing, n_new=n_new,
         n_failed=count(identity, failed_mask))
 end
 
@@ -434,7 +442,7 @@ function show_quantity_history(db::SQLite.DB, case_name::String, qty_name::Strin
         else
             abs_diff, rel_diff, status = compare_values(prev_q, q)
             diff_str = format_diff(abs_diff, rel_diff, status, q.value_type)
-            st = status == "CHANGED" ? "** CHANGED **" : status
+            st = startswith(status, "CHANGED") ? "** $status **" : status
             push!(rows, [short, date_str, val_str, diff_str, st])
         end
         prev_q = q

regression-harness/test/runtests_compare.jl

@@ -0,0 +1,18 @@
+using HDF5, JSON, Test
+
+include(joinpath(@__DIR__, "..", "src", "types.jl"))
+include(joinpath(@__DIR__, "..", "src", "extractor.jl"))
+
+# A stored quantity as get_quantities returns it.
+stored(value_type; value_real=nothing) =
+    (label="q", value_real=value_real, value_int=nothing, value_text=nothing, value_type=value_type, noise_threshold=1e-10)
+
+@testset "A quantity that disappears or changes type is a change" begin
+    real1, real2 = stored("real"; value_real=1.0), stored("real"; value_real=2.0)
+    @test compare_values(real1, real1)[3] == "OK"
+    @test compare_values(real1, real2)[3] == "CHANGED"
+    @test compare_values(stored("missing"), stored("missing"))[3] == "N/A"
+    @test compare_values(real1, stored("missing"))[3] == "CHANGED (missing)"
+    @test compare_values(stored("missing"), real1)[3] == "NEW"
+    @test compare_values(real1, stored("integer"))[3] == "CHANGED (type)"
+end

src/Equilibrium/KineticProfiles.jl

@@ -120,25 +120,27 @@ Read the GPEC HDF5 kinetic schema. Only `psi` is required; other standard datase
 optional and surfaced as-is. Non-standard `n_*` datasets become named per-species densities.
 """
 function _read_kinetic_h5(path::AbstractString; group::AbstractString="/")
-    h5open(path, "r") do f
-        g = group == "/" ? f : f[group]
-        haskey(g, "psi") || error("kinetic HDF5 '$path' group '$group': missing required dataset 'psi'")
-        rd(name) = haskey(g, name) ? Float64.(vec(read(g[name]))) : nothing
-        prov = ""
-        ats = attributes(g)
-        haskey(ats, "provenance") && (prov = string(read(ats["provenance"])))
-        # Non-standard datasets named `n_*` (e.g. "n_D", "n_T") are named per-species density
-        # profiles for explicit multi-ion input; other names are reserved for future schema fields.
-        extras = Dict{String,Vector{Float64}}()
-        for k in keys(g)
-            (k == "psi" || any(p -> first(p) == k, _KINETIC_H5_PROFILES) || !startswith(k, "n_")) && continue
-            v = read(g[k])
-            v isa AbstractArray && (extras[k] = Float64.(vec(v)))
-        end
-        return KineticProfileData(; psi=Float64.(vec(read(g["psi"]))),
-            (Symbol(name) => rd(name) for (name, _) in _KINETIC_H5_PROFILES)...,
-            species_densities=(isempty(extras) ? nothing : extras), provenance=prov)
+    h5open(f -> _read_kinetic_group(group == "/" ? f : f[group], "kinetic HDF5 '$path' group '$group'"), path, "r")
+end
+
+# The kinetic profiles in the open HDF5 file or group `g`; `label` names it in errors.
+function _read_kinetic_group(g, label::AbstractString)
+    haskey(g, "psi") || error("$label: missing required dataset 'psi'")
+    rd(name) = haskey(g, name) ? Float64.(vec(read(g[name]))) : nothing
+    prov = ""
+    ats = attributes(g)
+    haskey(ats, "provenance") && (prov = string(read(ats["provenance"])))
+    # Non-standard datasets named `n_*` (e.g. "n_D", "n_T") are named per-species density
+    # profiles for explicit multi-ion input; other names are reserved for future schema fields.
+    extras = Dict{String,Vector{Float64}}()
+    for k in keys(g)
+        (k == "psi" || any(p -> first(p) == k, _KINETIC_H5_PROFILES) || !startswith(k, "n_")) && continue
+        v = read(g[k])
+        v isa AbstractArray && (extras[k] = Float64.(vec(v)))
     end
+    return KineticProfileData(; psi=Float64.(vec(read(g["psi"]))),
+        (Symbol(name) => rd(name) for (name, _) in _KINETIC_H5_PROFILES)...,
+        species_densities=(isempty(extras) ? nothing : extras), provenance=prov)
 end
 
 """
@@ -177,30 +179,38 @@ with `read_kinetic_file`. The group carries `schema_version` and `provenance` at
 function write_kinetic_h5(path::AbstractString, data::KineticProfileData;
     group::AbstractString="/", schema_version::AbstractString="1.0",
     provenance::AbstractString=data.provenance)
-    h5open(path, "w") do f
-        g = group == "/" ? f : create_group(f, group)
-        function put(name, v, units)
-            v === nothing && return
+    h5open(f -> write_kinetic_group!(group == "/" ? f : create_group(f, group), data; schema_version, provenance), path, "w")
+    return path
+end
+
+"""
+    write_kinetic_group!(g, data; schema_version="1.0", provenance=data.provenance) -> g
+
+Write the profiles of `data` into the open HDF5 file or group `g` in the GPEC kinetic schema, as
+[`write_kinetic_h5`](@ref) does for a new file; the composition and closure fields are not written.
+"""
+function write_kinetic_group!(g, data::KineticProfileData; schema_version::AbstractString="1.0", provenance::AbstractString=data.provenance)
+    function put(name, v, units)
+        v === nothing && return
+        g[name] = collect(Float64, v)
+        attributes(g[name])["units"] = units
+    end
+    put("psi", data.psi, "normalized poloidal flux")
+    for (name, units) in _KINETIC_H5_PROFILES
+        put(name, getfield(data, Symbol(name)), units)
+    end
+    # Named per-species densities (multi-ion input); dataset names must be `n_*` to round-trip.
+    if data.species_densities !== nothing
+        for (name, v) in data.species_densities
+            startswith(name, "n_") ||
+                error("species density dataset `$name` must be named `n_*` to be readable as a per-species profile")
             g[name] = collect(Float64, v)
-            attributes(g[name])["units"] = units
+            attributes(g[name])["units"] = "m^-3"
         end
-        put("psi", data.psi, "normalized poloidal flux")
-        for (name, units) in _KINETIC_H5_PROFILES
-            put(name, getfield(data, Symbol(name)), units)
-        end
-        # Named per-species densities (multi-ion input); dataset names must be `n_*` to round-trip.
-        if data.species_densities !== nothing
-            for (name, v) in data.species_densities
-                startswith(name, "n_") ||
-                    error("species density dataset `$name` must be named `n_*` to be readable as a per-species profile")
-                g[name] = collect(Float64, v)
-                attributes(g[name])["units"] = "m^-3"
-            end
-        end
-        attributes(g)["schema_version"] = schema_version
-        attributes(g)["provenance"] = provenance
     end
-    return path
+    attributes(g)["schema_version"] = schema_version
+    attributes(g)["provenance"] = provenance
+    return g
 end
 
 """

src/ForceFreeStates/InputSnapshot.jl

@@ -13,3 +13,16 @@ function save_equilibrium_ingest!(out_h5, equil)
     end
     return nothing
 end
+
+"""
+    save_input_files!(out_h5, equil, wall_settings)
+
+Record the external files a run read, so a rerun and `equilibrium_from_h5` replay them without the originals:
+the equilibrium's kinetic profiles under `Input/RawInputs/Kinetic/` (the kinetic HDF5 schema; the composition is
+restated by the deck) and a file wall's text as `Input/RawInputs/Wall/wall_file_raw`.
+"""
+function save_input_files!(out_h5, equil, wall_settings)
+    equil.kinetic === nothing || Equilibrium.write_kinetic_group!(create_group(out_h5, "Input/RawInputs/Kinetic"), equil.kinetic)
+    wall_settings isa Vacuum.FileWall && (out_h5["Input/RawInputs/Wall/wall_file_raw"] = read(wall_settings.file, String))
+    return nothing
+end

src/ForceFreeStates/Output.jl

@@ -52,6 +52,7 @@ function save_to_hdf5(
     end
 
     save_equilibrium_ingest!(out_h5, equil)
+    save_input_files!(out_h5, equil, result.wall_settings)
     Utilities.record_inputs!(out_h5, result)
 
     # Write derived run parameters