Files
MeanField/experiments/plateau_diagnostic_README.md

11 KiB
Raw Permalink Blame History

Plateau diagnostics

This opt-in experiment runs the current production context, operator, normalization, preconditioner and Newton implementation. It adds observations and uncommitted frozen trials. It does not select a replacement numerical algorithm.

Build from the repository root:

cmake --build cmake-build-release-homebrew --target plateau_diagnostic -j 2

Run a calibrated spherical family with G=M=R_reference=1:

./cmake-build-release-homebrew/plateau_diagnostic \
  --index 1.5 --angular-momentum 0 --no-central --max-newton 8 \
  --output /absolute/path/to/new-results

The output directory must not already exist. The mesh defaults to the repository's sandbox.smesh; --mesh selects another mesh. The driver is restricted to one MPI rank. An MPI runtime must be available even for a single-rank run.

For each index, the driver obtains K and central density from the Lane–Emden solution so the spherical reference mass and radius are one. Index one uses its analytic constants. Other indices use the production radial integrator, so this calibration is seed-consistent, not an independent exact reference. Angular momentum sets a physical target, not a fraction of critical rotation. --no-central omits the central phase for controlled comparisons; it changes the model specification. The central phase remains enabled by default to reproduce earlier runs. Pass --no-central for the n=1 and n=1.5 campaign models; retain the optional phase for the n=3 homology-control experiment.

Default solver controls are restart 40, 200 linear iterations, linear relative tolerance 0.03, Newton relative tolerance 1e-8 and absolute tolerance zero. Override with --restart, --max-linear, --linear-tol, --rtol, --atol, and --max-newton. Exit 0 means the requested Newton criterion was met, 3 means it was not met (or a seed-only inspection was requested), and 2 means setup/diagnostics failed. A nonconverged result still saves the retained accepted state and its physical measurements. newton.csv also records linear_initial_residual, linear_rhs_norm, and linear_initial_over_rhs to measure the incoming correction warm start. The initial residual is ||b - J p_initial||; the ratio is zero for an exactly solved zero RHS and infinity for a nonzero initial residual with zero RHS. Missing solve reports produce NaN.

Frozen investigations

Append --probes to inspect the retained state after a solve. This compares cold FGMRES windows 40 and 160 at the same iteration budget and independently recomputes both true residuals. These are direct mfem::FGMRESSolver calls with empty Krylov history, not calls through the configured production linear backend; they do not test a modified backend or its recycling behavior. It also exports normalization consistency checks, the scalar Schur matrix and its coordinate labels, inverse-diagonal quality, and full-Jacobian actions on preconditioner-completed scalar directions.

To inspect a previously saved physical root vector without repeating Newton:

./cmake-build-release-homebrew/plateau_diagnostic \
  --index 1.5 --angular-momentum 0 --no-central --max-newton 0 \
  --probe-restart 160 \
  --probe-state /absolute/path/to/previous/accepted_state.txt \
  --output /absolute/path/to/new-frozen-results

Use the same mesh, compiled field orders, index, K/central-density calibration, and model flags. The snapshot is prepared as an uncommitted trial inside the existing diagnostic hook; the context's accepted state is unchanged. This is a frozen inspection facility, not a Newton continuation/restart facility. The saved vector is dimension-checked, but does not contain a mesh/model fingerprint.

Set --probe-restart 40 or --probe-restart 160 to select the correction used for positivity, finite-difference, and nonlinear-trial probes. The default is 40 for reproducibility; both cold solves still run, and the choice does not change the production Newton --restart setting. metadata.txt records probe_restart and frozen_probe_backend=direct_mfem_fgmres_cold; the selected correction is saved as probe_normalized_correction.txt. The probes compute the exact sampled enthalpy step bounds for both production closure and pressure-force quadrature. Every trial is restored and none is committed. Finite differences are one-sided; interpret their epsilon trend and absolute error, particularly for tiny Jacobian actions. global_correction_norm in the block table is the whole-vector norm.

Physical measurements

The observer integrates over stellar material only, at independent quadrature orders 14 and 18. It reconstructs physical positions and Piola-transformed gravity flux. For zero surface pressure it measures

  • W = 1/2 integral rho Phi, W_force = -integral rho x dot gradPhi;
  • rotational energy T and both P(h) and P(rho) pressure integrals;
  • energy and force virials abs(2T + W + 3 integral P) / abs(W_energy);
  • mass, angular momentum, center of mass, Bernoulli and EOS mismatch, hydrostatic gradient mismatch, sampled density/enthalpy signs, and geometry quality.

Pressure is evaluated with the actual polytropic index, including noninteger n. Signed finite DG density is retained in mass/gravity/rotation integrals. Any negative density sample makes the complete P(rho) pressure and virial undefined (NaN); P(h) remains available for positive enthalpy. The observer reports the negative-density mass and volume fractions instead of silently clipping them. NaN in P(rho) alone is not a quantitative measure of the global physical error. Quadrature samples are not a global positivity certificate.

normalized_unbordered_residual removes the central phase's artificial source from the hydrostatic equation. A small bordered residual can otherwise hide a pointwise modification at the central vertex. Virial success and Newton stopping success are deliberately separate measurements.

Reproducible serial campaigns

run_plateau_campaign.py requires a committed source revision and an already configured CMake build directory. It rebuilds the diagnostic target, records the source commit, build settings, binary and mesh hashes, and copies the executable into the result directory before running cases serially. Every exact command, exit code, log, state and compact summary is retained. A model exit code 3 is an observed nonconvergence; campaign exit 0 means the observations completed, not that every model converged. Output directories must be new.

python3 experiments/run_plateau_campaign.py \
  --source . --build cmake-build-release-homebrew --mesh sandbox.smesh \
  --cases experiments/plateau_control_cases.json \
  --output /absolute/path/to/new-campaign

Use --cmake if CMake is not on PATH. Build products, raw states and full result folders should remain outside source commits; retain compact result tables and reproduction instructions with each hypothesis branch. The controls omit the central phase at n=1 and n=1.5. A separate n=3 control retains it to fix homology; the production model API does not insert it implicitly.

Check a completed endpoint with python3 experiments/check_plateau_physics.py /absolute/path/to/result-folder. This gates strict recorded solver convergence and all eight finite virial errors at 1e-5, while reporting local field errors separately. See check_plateau_physics_README.md for the exit codes and interpretation limits.

Pressure obstruction and shape curvature (September 17)

These switches belong only to this experiment executable; no production model or solver options are added.

--capture-iterations 2,7,11 copies the selected zero-based Newton iteration's physical/normalized base, residual, and exact first-trial candidate from observer spans into capture_N/. It reconstructs the normalized direction as (candidate - base) / first_trial_alpha. capture.txt records the step, observed rejection, reconstruction error, and a heuristic subtraction-roundoff scale. This is not a direct copy of the solver's internal correction. Prefer an onset capture whose first step is large for derivative measurements; retain the exact candidate for unambiguous replay even when a late direction loses digits. No operator or preconditioner is evaluated from the observer, and no accepted state is changed. Compare non-timing trajectory columns and final-state hashes against the control to verify that capture is observational.

A focused frozen investigation bypasses both cold Krylov solves:

./cmake-build-release-homebrew/plateau_diagnostic \
  --index 1 --angular-momentum 0.1 --no-central --max-newton 0 \
  --probe-state /absolute/path/to/capture_2/physical_base.txt \
  --probe-direction /absolute/path/to/capture_2/normalized_direction.txt \
  --replay-candidate /absolute/path/to/capture_2/normalized_candidate.txt \
  --pressure-probe --shape-curvature \
  --output /absolute/path/to/new-focused-results

The state input is physical; direction and candidate inputs are normalized in the same recorded normalization. Mesh, model, field layout, and scaling must match the capture. If a sibling normalized_base.txt exists, its agreement with the newly normalized physical base is checked before probing. This check is not a complete model fingerprint. The runner hashes each explicit input plus every text/CSV/JSON output, including nested capture files. Exact replay logs root rejection reason, stage and EOS code; the EOS code has meaning only for thermodynamic rejections.

frozen_input_checks.csv records the current F, Jp+F defect, captured normalization agreement, and final restoration error. The two methodology notes describe the pressure-bound crossings and common-step vector Taylor decomposition. Every candidate remains uncommitted. The flags may be used separately to limit work. These are diagnostic outcomes, not new equilibrium-convergence or virial claims.

--surface-trace exports the actual compiled stellar-surface carrier row indices, physical base enthalpy and direction coefficients. At each production-rule limiting point it separates the full enthalpy polynomial into contributions from surface coefficients and all remaining coefficients. The former is a polynomial contribution at an interior point, not the value of the physical surface trace there. This inspection adds no trial preparations, Jacobian evaluations or linear solves, and restores the saved correction. A surface coefficient sign or zero trace alone does not certify interior polynomial positivity.

plateau_supplemental_cases_2026-09-17.json replays one onset candidate strictly between the closure and pressure enthalpy bounds and compares surface traces at onset and late plateau. The candidate is generated from captured normalized base and direction at the arithmetic midpoint of those two measured bounds; its derivation and hash are retained with the supplemental input.

The surface trace output also records actual FE nodal h and dh values on each limiting element, obtained through full physical-field scatters and evaluations at GetNodes(). Shared elements are emitted once per quadrature family. This allows sign comparisons at arbitrary alpha without confusing coefficient signs with interior polynomial positivity. The union-bound column is an affine prediction at the boundary, not a separately prepared candidate.

The completed campaign and evidence limits are recorded in PRESSURE_AND_SHAPE_CAMPAIGN_2026-09-17.md. Portable offline analyses are in plateau_analysis/; pass the retained artifact root explicitly. The Q3 follow-up checks a finite interior grid against the exported polynomial and does not certify positivity everywhere.