diff --git a/docs/development/hdf5-conventions.md b/docs/development/hdf5-conventions.md index 0fc0b4e04..c6f568931 100644 --- a/docs/development/hdf5-conventions.md +++ b/docs/development/hdf5-conventions.md @@ -67,6 +67,7 @@ Mechanism: writers stay table-driven — each writer keeps a `path => (; long_na - Complex quantities are stored as the native HDF5.jl compound type (readable by h5py as a compound dtype) — **never split into `*_real`/`*_imag` dataset pairs**. Sole sanctioned exception: `Input/RawInputs/ForcingTerms/amplitude_{real,imag}`, which mirrors the external forcing ingest-file format and keeps pre-existing snapshots replayable. - `NaN` is the not-computed sentinel in numeric datasets (e.g. auto-derived settings, rootless growth-rate entries). +- A **zero-extent array** is the not-computed sentinel for whole datasets that a given run never produces (e.g. `SingularSurfaces/ca_left`/`ca_right` on kinetic or galerkin-matched runs, the free-boundary energies when `vac_flag=false`, the on-demand derivative stores). Never write unpopulated (`undef`) memory. - Ragged (variable-length) data uses the flat-plus-`offsets` companion pattern (`offsets[k+1] - offsets[k]` = length of row `k`) rather than HDF5 VLEN types, e.g. `KineticForces//EnergyIntegrals/` and `Tearing/Diagnostics/*`. ## Back-compatibility policy diff --git a/docs/development/regression-harness.md b/docs/development/regression-harness.md index f06b6e356..e9246c962 100644 --- a/docs/development/regression-harness.md +++ b/docs/development/regression-harness.md @@ -98,14 +98,12 @@ regress --cases solovev_n1 --ref-range develop~10..develop - `--force` — re-run even if cached - `--verbose` — print GPEC subprocess output - `--no-instantiate` — skip `Pkg.instantiate()` (faster if deps are already resolved) -<<<<<<< HEAD - -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. -======= - `--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 +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. + ## Making source code the only variable `Manifest.toml` is untracked, so a worktree checked out at an old commit used to resolve whatever @@ -140,11 +138,9 @@ than leaving you to infer it from the numbers. Results cached before environment fingerprinting existed carry no environment and are therefore re-run once — those are exactly the entries whose provenance cannot be established. -Thread counts are recorded but **not** forced: the harness does not silently change how your runs -execute. If the two refs in a comparison ran under different thread counts, the report flags it. +If the two refs in a comparison ran under different thread counts, the report flags it. ## Exit status - `0` — every run completed (and, with `--fail-on-change`, nothing changed) - `1` — a run failed, or a quantity changed under `--fail-on-change` ->>>>>>> performance/regression-harness-worktree-reuse diff --git a/docs/src/stability.md b/docs/src/stability.md index 529304a90..d54459fbb 100644 --- a/docs/src/stability.md +++ b/docs/src/stability.md @@ -363,8 +363,10 @@ end ## Notes -- The standard path does not populate `delta_prime`; use `PerturbedEquilibrium.SingularCoupling` - for Δ' on the standard path (it reads `ca_l`/`ca_r` directly). +- The standard path does not populate `delta_prime`; the canonical Δ' is the STRIDE BVP + `SingularSurfaces/Delta_prime_matrix` from the parallel FM path. `ca_l`/`ca_r` are filled + only by ideal surface crossings (kinetic and galerkin-matched runs emit zero-extent + `ca_left`/`ca_right` sentinels). - The Riccati and parallel FM paths compute Δ' inline at each crossing, using the direct diagonal formula (no GR permutation). The result in `delta_prime_col[ipert_res, i]` equals `delta_prime[i]` to machine precision. diff --git a/regression-harness/cases/solovev_n1.toml b/regression-harness/cases/solovev_n1.toml index 1c625f22a..d7f4224e5 100644 --- a/regression-harness/cases/solovev_n1.toml +++ b/regression-harness/cases/solovev_n1.toml @@ -182,6 +182,16 @@ label = "pressure profile (checksum)" noise_threshold = 0 order = 71 +# Bitwise reproducibility guard for the ideal asymptotic coefficients (issue: kinetic +# runs used to dump uninitialized memory here; ideal values must stay byte-stable). +[quantities.ca_left] +h5path = "SingularSurfaces/ca_left" +type = "complex_matrix" +extract = "checksum" +label = "ca_left (checksum)" +noise_threshold = 0 +order = 72 + # Runtime (special: not from H5) [quantities.runtime] h5path = "" diff --git a/src/Analysis/ForceFreeStates.jl b/src/Analysis/ForceFreeStates.jl index 99dd60e2b..23d4db74b 100644 --- a/src/Analysis/ForceFreeStates.jl +++ b/src/Analysis/ForceFreeStates.jl @@ -327,6 +327,8 @@ function plot_delta_prime(h5path; save_path=nothing) end msing == 0 && return plot(; title="No singular surfaces found", legend=false) + # ca_left/ca_right are zero-extent sentinels on kinetic/galerkin-matched runs (never computed there). + isempty(ca_l) && return plot(; title="No asymptotic coefficients — ca_left/ca_right not computed for this run", legend=false) numpert_total = size(ca_l, 1) chi1 = 2π * psio diff --git a/src/ForceFreeStates/EulerLagrange.jl b/src/ForceFreeStates/EulerLagrange.jl index 88c1d9e4a..459b052db 100644 --- a/src/ForceFreeStates/EulerLagrange.jl +++ b/src/ForceFreeStates/EulerLagrange.jl @@ -641,6 +641,7 @@ function cross_ideal_singular_surf!( end # Get asymptotic coefficients after crossing rational surface odet.ca_r[:, :, :, ising] .= sing_get_ca(odet.u, ua, intr) + odet.ca_populated = true # Δ' is NOT computed for the standard path. The physical Δ' requires the solution # columns to be in the Riccati gauge (U₂=I), maintained only by Riccati renormalization. diff --git a/src/ForceFreeStates/ForceFreeStatesStructs.jl b/src/ForceFreeStates/ForceFreeStatesStructs.jl index 79735b895..589874073 100644 --- a/src/ForceFreeStates/ForceFreeStatesStructs.jl +++ b/src/ForceFreeStates/ForceFreeStatesStructs.jl @@ -532,6 +532,10 @@ and a small set of temporary matrices and factors used to compute singular-layer - `ca_l::Array{ComplexF64,4}` - Asymptotic coefficients just to the left of each singular surface with shape `(numpert_total, numpert_total, 2, msing)`. + - `ca_populated::Bool` - True once an ideal singular-surface crossing has filled `ca_l`/`ca_r`; kinetic and + galerkin-matched runs never populate them and leave this false, and the HDF5 writer then emits zero-extent + `ca_left`/`ca_right` datasets instead of unpopulated arrays. + - `edge_scan::EdgeScanState` - Edge dW scan state and results. Initialized as a disabled sentinel (N_edge=0) and replaced by `findmax_dW_edge!` when a scan runs. - `psifac::Float64` - Current normalized flux coordinate for the integrator. @@ -591,8 +595,9 @@ and a small set of temporary matrices and factors used to compute singular-layer u_store_el_basis::Bool = true du_store_populated::Bool = false crit_store::Vector{Float64} = Vector{Float64}(undef, numsteps_init) - ca_r::Array{ComplexF64,4} = Array{ComplexF64}(undef, numpert_total, numpert_total, 2, msing) - ca_l::Array{ComplexF64,4} = Array{ComplexF64}(undef, numpert_total, numpert_total, 2, msing) + ca_r::Array{ComplexF64,4} = zeros(ComplexF64, numpert_total, numpert_total, 2, msing) + ca_l::Array{ComplexF64,4} = zeros(ComplexF64, numpert_total, numpert_total, 2, msing) + ca_populated::Bool = false # Edge dW scan state and results (disabled sentinel when psiedge >= psilim, i.e. no edge scan) edge_scan::EdgeScanState = EdgeScanState(numpert_total, 0) diff --git a/src/ForceFreeStates/Riccati.jl b/src/ForceFreeStates/Riccati.jl index b1e29f0a7..274a71509 100644 --- a/src/ForceFreeStates/Riccati.jl +++ b/src/ForceFreeStates/Riccati.jl @@ -1265,6 +1265,7 @@ function _capture_right_crossing_data!(odet::OdeState, singp::SingType, sing_asy end end odet.ca_r[:, :, :, ising] .= sing_get_ca(odet.u, ua, intr) + odet.ca_populated = true end # STUB: per-surface ca-based Δ' (not physically valid; see SingType.delta_prime docstring). diff --git a/src/GeneralizedPerturbedEquilibrium.jl b/src/GeneralizedPerturbedEquilibrium.jl index 96f9f4b1a..ce78c1dde 100755 --- a/src/GeneralizedPerturbedEquilibrium.jl +++ b/src/GeneralizedPerturbedEquilibrium.jl @@ -1205,8 +1205,10 @@ function write_outputs_to_HDF5( out_h5["SingularSurfaces/rational_psi"] = [sing.psifac for sing in result.surfaces] out_h5["SingularSurfaces/rational_q"] = [sing.q for sing in result.surfaces] out_h5["SingularSurfaces/dqdpsi"] = [sing.q1 for sing in result.surfaces] - out_h5["SingularSurfaces/ca_left"] = diag !== nothing ? diag.ca_l : ComplexF64[] - out_h5["SingularSurfaces/ca_right"] = diag !== nothing ? diag.ca_r : ComplexF64[] + # Kinetic and galerkin-matched runs never populate ca_l/ca_r (only ideal surface + # crossings do); emit zero-extent sentinels instead of unpopulated arrays. + out_h5["SingularSurfaces/ca_left"] = (diag !== nothing && diag.ca_populated) ? diag.ca_l : zeros(ComplexF64, 0, 0, 0, 0) + out_h5["SingularSurfaces/ca_right"] = (diag !== nothing && diag.ca_populated) ? diag.ca_r : zeros(ComplexF64, 0, 0, 0, 0) if msing > 0 # Mode numbers at each surface (jagged — pad with 0 to max_modes width) diff --git a/src/HDF5Schema.jl b/src/HDF5Schema.jl index b8b0333c7..274f54f14 100644 --- a/src/HDF5Schema.jl +++ b/src/HDF5Schema.jl @@ -181,9 +181,11 @@ const MAIN_H5_ANNOTATIONS = [ "SingularSurfaces/D_I" => (; long_name="Mercier D_I evaluated at each rational surface", dims=("surface",), attach=(1 => "SingularSurfaces/rational_psi", 1 => "SingularSurfaces/rational_q")), "SingularSurfaces/ca_left" => - (; long_name="asymptotic large/small-solution coefficient matrices just left of each surface", dims=("mode", "solution", "large_small", "surface")), + (; long_name="asymptotic large/small-solution coefficient matrices just left of each surface (zero-extent when not computed — ideal crossings only)", + dims=("mode", "solution", "large_small", "surface")), "SingularSurfaces/ca_right" => - (; long_name="asymptotic large/small-solution coefficient matrices just right of each surface", dims=("mode", "solution", "large_small", "surface")), + (; long_name="asymptotic large/small-solution coefficient matrices just right of each surface (zero-extent when not computed — ideal crossings only)", + dims=("mode", "solution", "large_small", "surface")), "SingularSurfaces/E" => (; long_name="Glasser-Greene-Johnson coefficient E per surface", dims=("surface",), attach=(1 => "SingularSurfaces/rational_psi", 1 => "SingularSurfaces/rational_q")), "SingularSurfaces/F" => diff --git a/test/runtests_fullruns.jl b/test/runtests_fullruns.jl index ce6470385..a76515bf5 100644 --- a/test/runtests_fullruns.jl +++ b/test/runtests_fullruns.jl @@ -28,6 +28,10 @@ using HDF5 et = read(h5["ForceFreeStates/FreeBoundaryStability/eigenmode_energies"]) @test isfinite(real(et[1])) @test real(et[1]) > 0 + # Kinetic runs never populate the asymptotic ca coefficients; the writer must + # emit deterministic zero-extent sentinels, not uninitialized memory. + @test isempty(read(h5["SingularSurfaces/ca_left"])) + @test isempty(read(h5["SingularSurfaces/ca_right"])) end rm(joinpath(ex3, "gpec.h5"); force=true) true diff --git a/test/runtests_h5_schema.jl b/test/runtests_h5_schema.jl index 6b838ab53..c743d6a7d 100644 --- a/test/runtests_h5_schema.jl +++ b/test/runtests_h5_schema.jl @@ -95,6 +95,13 @@ end fwd = "ForceFreeStates/Solutions/ForwardIntegration" @test HDF5.API.h5ds_is_scale(h5["$fwd/psi"]) @test HDF5.API.h5ds_is_attached(h5["$fwd/q"], h5["$fwd/psi"], 0) + + # Ideal run with rational surfaces: the asymptotic ca coefficients are + # populated and finite (the kinetic/galerkin not-computed case emits + # zero-extent sentinels instead — asserted in runtests_fullruns.jl). + ca_l = read(h5["SingularSurfaces/ca_left"]) + @test !isempty(ca_l) + @test all(isfinite, ca_l) end end end