linstabpy 1.0.1

Install route and housekeeping release.

Install

pip install linstabpy
linstabpy install-stack      # builds complex PETSc/SLEPc + MUMPS with pip

install-stack needs a C compiler, gfortran, make and an MPI with compiler wrappers (sudo apt install build-essential gfortran libopenmpi-dev openmpi-bin). The conda route needs no compiler: conda env create -f https://raw.githubusercontent.com/openfluids/linstabpy/main/environment.yml.

Changes

• New: linstabpy install-stack builds the CPU solver stack into the current environment. Tested end to end on Linux: PETSc 3.24.2, SLEPc, MUMPS, full test suite 1791 passed, 0 failed.
• New: environment.yml installs linstabpy from PyPI, so the conda route needs no checkout. PETSc and SLEPc are pinned to 3.24.
• New: the source archive carries install.sh, environment.yml and .python-version.
• New: a solve without PETSc/SLEPc fails with a message that names the install routes.
• Changed: the minimum Python is 3.11. Pillow 13, which fixes the open security advisories, does not support 3.10.
• Fixed: the provenance record no longer carries the host name that built Open MPI, and drops the NUL padding Open MPI 5 adds.

Full list: CHANGELOG.md.

linstabpy 1.0.0

First stable release of linstabpy: stability and resolvent analysis of compressible viscous flows on curvilinear grids (SBP finite differences, PETSc/SLEPc, MUMPS).

Install

pip install linstabpy
linstabpy --version

The PyPI package has the CLI, config validation and the NumPy tools. Solving needs complex-scalar PETSc/SLEPc with MUMPS: clone the repository and run bash install.sh (see the README).

Changes that affect existing results

• Resolvent sign. The resolvent is now H(omega) = [(gamma - i*omega) I - L0]^-1, the response to a forcing exp(-i*omega*t). For W != 0 and beta != 0 the gain map is mirrored in omega, and forcing and response modes change in every case. Set solver.resolvent_sign_convention: matlab_parity to reproduce old results.
• MPI. All ranks now take the same walk-or-refactor and GMRES-deadline decisions. A multi-rank resolvent_map could deadlock before.
• Compatibility shims removed. CASE1/CASE2 and LINSTAB2D_* environment aliases, deprecated exports, solve_transpose (use solve_adjoint), the default CLI config and the default Re.
• Python 3.14 is the default of install.sh. PETSc and SLEPc stay at 3.24.2.

Full list, with the effect on old results: CHANGELOG.md.

Documentation: https://linstabpy.vercel.app/

v1.0.0rc11

[1.0.0rc11] - 2026-09-17

Fixed

• The dissipation advisory warned on every run without artificial
  dissipation, whatever the mesh. On a rectilinear mesh the bare operator
  cannot hold an odd-even (sawtooth) mode, so the warning fired against the
  correct setting for both released cases. The verdict now comes from
  mesh_advisories(), which measures the mesh that the run uses and reports
  the number in every case: the released Savarino and Song meshes measure
  4.8e-15 and 2.1e-12 against a gate of 1e-8 and get a note, and a bent mesh
  without dissipation still gets the warning. A config that does not state
  its dissipation choice at all is reported by config_advisories() as
  DISS-UNDECLARED.
• The base-flow check compared the largest value in the field with 1, so
  every compressible flow with a shock or a hot wall looked mis-scaled. It
  reported density 1.39 and 1.51 and temperature 1.81 and 1.83 for the two
  released base flows, whose freestream is exactly 1. The check now looks
  for the freestream itself, the part of the field where density, velocity
  and temperature are all 1, and records the share of the field it covers
  (1.9 % and 3.7 %). A field scaled with the speed of sound has no such
  region and is still reported, now with the values of the largest uniform
  region so the actual convention is visible.
• The MUMPS ICNTL entries went into the PETSc option database at startup on
  every run. A run with no MUMPS stage left them unused and PETSc ended with
  "There are options you set that were not used". They now go in only when a
  stage can use MUMPS, and the two cuDSS example configs name the card for
  every stage instead of naming MUMPS for stages they never run.
• examples/savarino/savarino_mapwalk.py wrote the real host name into the
  machine column of every line record, which put a private host name in
  three shipped run records. It follows the package policy now: the label is
  unknown unless LINSTABPY_HOSTNAME is set.
• examples/savarino/verify_norm_convention.py ended with a traceback when
  the Zenodo base flow and modes were not in the tree, which is the state of
  a fresh clone because their licence does not allow redistribution. It
  prints the two fetch commands and stops, like run_compare_points.py.

Changed

• The test that pins the committed beta = 40 farm record checks its gains,
  which a rerun reproduces, instead of the sum of its wall times, which is a
  property of the machine on the day.

Known issues

• The STRUMPACK LU returns a wrong solve about one run in three on this
  build. Measured at the shift-invert probe of the Blasius temporal
  problem, sixteen fresh processes on byte-identical input (c9, one A6000,
  PETSc 3.24.2): eleven runs gave a relative residual of 6.3e-13 to 8.7e-12
  and five gave 1.9e-5, 2.9e-5, 3.0e-5, 6.7e-2 and 1.1e-1. Nothing in the
  input varies: the probe right-hand side is seeded and the mesh is fixed.
  The accuracy probe rejects every bad solve, which is what it exists for,
  and the eigenvalue and resolvent paths refine with FGMRES and check the
  residual afterwards, so a bad factorization cannot reach a published
  number. The cause is not established. Prefer MUMPS or cuDSS for
  production work on this build until it is.

1.0.0rc10

Changed

• A line of frequencies on one card is faster. For a serial cuDSS
  factorization, a refactorize now takes the CSR arrays of the new operator
  once, checks the nonzero pattern against the first operator, writes the
  values into the same PETSc matrix, and gives them to cuDSS. The adjoint
  factor gets K^H from a permutation of those values that is computed once
  per factorization, not from a sparse transpose at every frequency. On the
  450,000-unknown Savarino operator, one A6000, the median line point went
  from 32.1 s to 25.8 s over the eight frequencies after the first of the
  beta = 40 line, with gains identical to ten digits. The path is the default; set
  LINSTABPY_FAST_REFACTOR=0 to use the general path. A different nonzero
  pattern raises ValueError.
• The run summary reports the peak memory of the largest rank next to the
  sum over ranks.
• run_all_examples.sh prints SOFT_FAIL with the case names when the
  Martini parity comparison has soft failures. Before, it printed PASS.

Fixed

• examples/savarino/savarino_mapwalk.py compared each line gain with
  output/savarino_map2_200pt.csv, a farm computed before the filter-order
  fix. The farm_gain and rel_err columns of a line record then showed
  differences up to 4.4e-3 that are not in the line result. The default
  reference is now savarino_map3_200pt.csv, the committed map of record,
  in the script, in savarino_mapwalk.sbatch.template, which passed the old
  file explicitly, and in the test and docstrings that named it.

Added

• wiki/ architecture pages and .github/workflows/pages.yml. MkDocs
  builds the site. Pages are hand-written and hand-maintained.

Fixed

• The solvers/router/resolvent_map module docstring said every gain map
  needs SBP sawtooth dissipation. That claim is wrong on rectilinear
  meshes. The odd-even mode exists only on a mesh that is curvilinear in
  both directions, where the metric cross-terms mix grid-scale content.
  The docstring now states the mesh condition, points to
  driver.preflight.metric_cross_coupling() as the test, and scopes the
  tuned epsilon ~ 0.10 to the curved-wall mesh it came from. This aligns
  the last remaining text with the measured story already in
  config/defaults.py and driver/preflight.py.
• The Martini parity manifest still pointed Example 3 at a combined
  example_3_streaky_bl.npz that nothing writes since the example was split
  into the four paper campaigns, so every full parity run reported a hard
  failure on a missing artifact. The compare now walks the per-campaign
  NPZs, each against its own MATLAB bundle, with the same residual-first
  gates as before, and aligns the component-sliced mode subsets the
  repository bundles store (idx rows, saved-mode columns) before pairing.
  With the wiring repaired, the full parity run reports four passes, two
  documented structural soft failures, and no hard failures.

linstabpy 1.0.0rc9

Release-hardening pass before 1.0.0.

Changed

• Consolidated duplicated private helpers into shared modules: config lookup
  (cfg_get, as_bool) in utils/_solver_cfg, peak_rss_bytes() in
  utils/cost_metrics, the spatial Dirichlet helpers in
  boundary_conditions/dirichlet, the backend normalizer in
  linalg/backends, and one _matvec(operator, x, adjoint=...) in the
  solver router in place of a forward/adjoint pair. Net -142 lines with no
  behavior change; every module now imports the one canonical copy instead
  of carrying a private variant that could drift.

Fixed

• Four call sites in solvers/router/stable_sweep.py still called the
  removed _adjoint_matvec after the mat-vec merge, so PETSc-backend sweeps
  raised NameError. These paths only execute where petsc4py is installed,
  which is why a PETSc-less checkout did not catch it.
• Environment-absence tests (missing nvidia-smi, missing PETSc stack,
  missing repository .venv, empty PETSc options database) now simulate the
  absent component instead of asserting properties of the host machine. The
  full suite passes unmodified on CPU-only machines and on GPU nodes.
• The Martini export tools now locate the upstream MATLAB code through
  LINSTAB2D_ROOT only, and the one parity test that reads a raw MATLAB
  fixture bundle takes its directory from LINSTABPY_MATLAB_FIXTURES_NPZ;
  the machine-specific fallback paths are gone, and the public-provenance
  scanner now rejects those path forms.

Verification

The tagged commit passes the full test suite with zero failures on both
supported stacks: 1574 tests pass on a CPU-only environment (MKL BLAS,
MUMPS) and 1584 on a GPU environment (CUDA-enabled PETSc with STRUMPACK
plus the cuDSS backend, on an NVIDIA RTX A6000). Every remaining skip
names either an optional backend absent from that environment or a
reference file that is not distributed with the repository, and between
the two environments every backend-gated test runs at least once.

The six examples (the five Martini flows plus the Example 0 derivative
protocol) were recomputed from scratch for this release — every cached
NPZ deleted first — through ./run_all_examples.sh, and all six pass.
Against the shipped MATLAB reference bundles, the Blasius and jet
resolvent gains agree to 1.7e-11 and 3.3e-12 relative, the cavity spectra
agree to 6e-6 with eigenvalue residuals below 2e-8, and the Example 3
temporal and spatial eigenpairs satisfy their discrete problems to
residuals of 5e-9 or better. (The Example 3 parity check required a
manifest repair that landed just after the tag, in 05fc7ac and
7167080; the artifacts it verifies are the ones recomputed here.)

linstabpy v1.0.0rc8

Added --- cuDSS GPU sparse-direct backend

The largest performance change since the package was created. A new GPU
backend factors the linearized operator with NVIDIA cuDSS (through
nvmath-python) and wraps it in the same right-preconditioned FGMRES that
already owns the accuracy of the iterative path. On the 450,000-unknown
SWBLI resolvent operator (2.0--2.3e7 nonzeros, complex double) one
factorization drops from 41--55 s (STRUMPACK on an NVIDIA RTX A6000) to
5.7--9.0 s (cuDSS, double precision) on the same card. A complete
frequency-line resolvent point --- two refactorizations plus ten forward and
ten adjoint solves converged to 1e-10 --- costs 23 s, against 60--82 s for
the four-rank MUMPS CPU line mode and 130--180 s for the STRUMPACK GPU
path. The GPU changes from the slower arm into the fastest solver on the
node, in full double precision, with no accuracy loss.

• Added linalg/cudss_factor.py: a cuDSS solver wrapper with the four
  operations the resolvent map needs: plan (analysis, once per sparsity
  pattern), factorize, refactorize (new values on the kept pattern,
  through reset_operands(a=...) --- nvmath copies a host operand to the
  device when it is set, not on each call), and per-right-hand-side solve
  (reset_operands(b=...) for the same reason). The cuDSS threading layer
  (libcudss_mtlayer_gomp) is discovered inside the installed NVIDIA wheel
  and enabled with the job's host threads; it halves the analysis time.
• Added factor_backend: cudss to the solver configuration: FGMRES in
  complex double, right-preconditioned by the cuDSS LU, tolerance owned by
  the iteration (ksp_rtol, default 1e-10). The direct MUMPS and STRUMPACK
  paths are unchanged.
• Adjoint solves build a second cuDSS factorization of the conjugate
  transpose. cuDSS 0.8 has no transposed or adjoint solve
  (CUDSS_CONFIG_SOLVE_MODE accepts only its default), so the dual factor
  is the design, not a heal: both factors share one analysis cost profile
  and refactorize together along a frequency line.
• Added factor_precision: single (opt-in, cudss only): the factorization
  runs in complex64 --- 0.7--1.2 s per factor, half the device memory ---
  while FGMRES stays in complex double and restores the solution to 1e-12.
  The option force-enables the two settings measured to be necessary:
  MC64-style matching (MAX_DIAG_PRODUCT) and double-precision row/column
  equilibration before the cast. Without matching the single-precision
  factor fails near zero frequency; with both, it converges at every tested
  point of the gain map.
• The adjoint self-test for this backend draws its probe vectors from the
  forcing-masked subspace instead of using dense random vectors. Random
  probes excite the boundary-penalty rows (row norms up to 5.5e12), where a
  single-precision factor cannot precondition; they would fail a
  configuration that is healthy on every right-hand side the resolvent map
  actually produces.
• Added the optional dependency group cudss (pip install "linstabpy[cudss]"): nvmath-python[cu12] >= 1.0, < 2, which brings the
  cuDSS 0.8 wheels. Requires a CUDA 12 runtime and an NVIDIA GPU supported
  by cuDSS. The core install is unchanged.

Performance --- measured (dual-socket Xeon Silver 4210R, MKL, one RTX A6000 48 GB)

All numbers from the 450,000-unknown SWBLI operator of the validation
study, on an otherwise idle node, with the solve tolerance at 1e-10.

• cuDSS double: analysis 3.3--3.4 s (threading layer on; 6.7--8.1 s
  without), factorization 5.7--9.0 s, 4.1--5.8e8 nonzeros in L+U, one solve
  26--30 ms, two resident factors (K and its adjoint) 15--20 GB device
  memory.
• Refactorization on the kept analysis costs a full numeric factorization
  (9.0--9.2 s); what a frequency line saves is the analysis. Timings are
  reproducible across independent runs to about 1 percent.
• Line-mode point (two refactorizations + 20 solves): 23 s double, 12--14 s
  single at a benign point --- the production sweep below shows single is
  slower overall. Reference points on the same operator and node: MUMPS CPU line
  mode 60--82 s (4 ranks), MUMPS factorization alone 24--33 s (10 threads),
  STRUMPACK GPU point 130--180 s.
• FGMRES iteration counts: 1--2 (double), 5--20 (single with matching and
  equilibration, at points where it converges; the 47-point sweep below found
  21 points where it exceeds the 400-iteration cap).
• Production campaign (two A6000s, one job per card, two concurrent): the
  47-point beta=40/45 subset runs at 57.1 s/pt as independent jobs (CPU record
  118.5 s/pt) and 34.0 s/pt in line mode (CPU 76.6 s/pt); the full 200-point
  map runs at 57.1 s/pt, 103.9 min wall against 83.6 min for five concurrent
  CPU jobs on twenty cores --- wall-clock is bound by device count, per-point
  cost halves. Run records committed as
  examples/savarino/savarino_cudss_*.csv.
• Single precision, measured end-to-end on the same 47 points: converged
  points average 207 s/pt against 58 s/pt double (the weaker factor needs the
  extra iterations), and 21 of 47 points --- the resonant band
  2 <= f <= 7 at beta=45, the response peak at beta=40, and the
  zero-frequency corner --- exceed 400 FGMRES iterations. Evaluated and
  rejected for production; the opt-in stays for the masked-RHS niche and as
  the memory route (halved factor storage).

Accuracy --- measured

• Solutions through the cuDSS-preconditioned FGMRES match the MUMPS double
  reference to 1.2e-12 relative; the discrete adjoint identity holds to
  2.5e-14--5.9e-12 with the dual factor.
• The full production campaign reproduces the committed map-of-record gains
  to every printed digit: all 200 map points and both 47-point arms (jobs and
  line mode) show zero difference at the CSV's ten-digit precision.
• The single-precision option reaches the same 1e-12 through the double
  iteration. Its scope is the production right-hand-side family (forcing-
  and response-masked vectors); on dense unmasked vectors the iteration can
  stall, which is why the option stays opt-in and the self-test uses masked
  probes.
• For a unit-norm right-hand side the attainable true-residual floor is
  about 1e-5, set by the boundary-penalty rows; MUMPS shows the same floor
  on the same operator. Solution error stays at 1e-12.

Fixed --- adjoint documentation

• Corrected the explanation of why the GPU path needs a second
  factorization. PETSc's STRUMPACK interface registers no transposed solve,
  and MatSolveTranspose then silently substitutes the forward solve (the
  factor's symmetry flag is unknown, which the fallback treats as true).
  The measured "transpose defect" therefore comes from the PETSc interface,
  not from STRUMPACK's arithmetic. Earlier comments and the paper draft
  attributed it to STRUMPACK's transpose solve; code comments and
  documentation now state the fallback. The dual-factor construction was
  and remains the correct route on both GPU backends.

Fixed --- reproduction documentation

• A reproduction audit of the shipped tree closed four seams between the new
  backend and its documentation: the "Current Status" backend table in
  docs/COMPUTE_BACKENDS.md still said cuDSS was not wired (contradicting
  the new section in the same file); the install line named only the cudss
  extra when the iterative path also needs the petsc extra, and did not say
  that a standard CPU PETSc/SLEPc build is enough (the GPU work happens
  inside cuDSS); the Savarino example README did not mention the two shipped
  cuDSS configs (savarino_point_cudss.yaml, savarino_point_cudss32.yaml)
  or that the campaign generators drive them via --config; and both
  campaign generators' --config help text named only the STRUMPACK GPU
  route. All four are corrected.

Changed --- Savarino base flow fetched, no longer tracked

• The published Savarino base flow (reference/baseflow/*.dat, 10.9 MB) is no
  longer tracked: its Zenodo record (10.5281/zenodo.18854407) is open access
  but states no redistribution license, so the repository ships no file from
  it. fetch_zenodo.py baseflow downloads it (md5 checked against the Zenodo
  files API) and savarino_to_npy.py re-verifies the md5 before converting.
  The shipped reference/savarino_response_*.npz (repackaged from the same
  record) remain tracked pending the same license decision.

Added --- every published number now has its run record and recipe in the tree

A second reproduction audit, this time claim-by-claim against the manuscript,
brought in everything a reader needs to check or rerun the measurements:

• The STRUMPACK GPU campaign records (savarino_strumpack_farm_fp64.csv,
  savarino_strumpack_lines_fp64.csv) and the dedicated CPU 47-point
  references (savarino_cpu_farm_47pt.csv, savarino_cpu_lines_47pt.csv) —
  the previously unshipped halves of the backend comparison.
• The backend benchmarks (examples/savarino/benchmarks/): the five
  measurement CSVs behind the per-backend factorization times and the
  single-precision verdict, together with the harness scripts that produced
  them and the operator-dump tool they consume.
• savarino_campaign_walltimes.csv: wall-clock of every campaign arm from
  Slurm accounting, first job start to last job end, with job-id ranges. This
  corrected one published number: the cuDSS 200-point map wall-clock is
  103.9 min, not 100.2 (the earlier figure was measured from the submission
  loop's completion, while the first jobs were already running).
• savarino_mpi_scaling.csv: MPI strong scaling at the Mack point measured
  on the same node as every other row (190.5 / 154.6 / 148.2 s at 4/8/16
  ranks, gain identical to ten significant digits). This replaces an earlier
  scaling table whose hardware could not be traced to any machine or record.
• compare_backend_modes.py: prints the Chu-weighted correlation between the
  leading response modes of two runs of the same point, so the "GPU ...

linstabpy v1.0.0rc7

A documentation release, cut on top of rc6 the same day. It rewrites the
comments and docstrings across the whole tree and removes two internal
working documents whose job was done. No solver code changed.

What changed

Comments and docstrings, tree-wide. 155 Python files were edited:
50 in the package (src/linstabpy/), 72 in the test tree, 34 in the
example scripts, minus release metadata. Long multi-clause comment
paragraphs became short sentences, one idea per sentence, so each
statement reads in one pass. The edit is a reflow, not growth: about
2 900 lines removed and 2 900 added across 160 files.

Typical before/after, from src/linstabpy/linalg/factorize.py:

before — "The default path is KSPPREONLY + PCLU + MUMPS because the
release tests and reference examples are calibrated against exact
direct solves. The iterative path is available for large memory-bound
cases, but it is an explicit solver-policy choice rather than an
automatic fallback."

after — "The default path is KSPPREONLY + PCLU + MUMPS. The release
tests and reference examples use exact direct solves. The iterative
path is available for large memory-bound cases. It is an explicit
solver-policy choice rather than an automatic fallback."

ruff format over the test tree. The formatter had drifted on
tests/; this tag brings every test file back to the shared style.

Two working documents removed. MARTINI_PARITY_SUMMARY.md (210
lines) and examples/martini/AUDIT_PLAN.md were the running notes of
the Martini examples 0–5 parity audit. That audit closed in rc6, and
everything durable it produced — the reference data, the parity checks,
the per-case plotters, the shipped figures — is in the tree from rc6.
The notes that drove the audit are process material, not release
material, so they are gone.

Verification

The claim "no solver code changed" was checked mechanically, not
assumed: every changed .py file was parsed at both tags, docstrings
were stripped, and the syntax trees were compared (comments never reach
the syntax tree, so this compares exactly the executable code).

• 154 of 155 files: identical executable code at rc6 and rc7.
• The one exception is tests/regenerate_mumps_parity_reference.py,
  where a user-facing message string changed: the guard that refuses to
  regenerate the MUMPS baseline under ambient LINSTABPY_* variables
  now says "parity reference" instead of "parity fixture". The test
  that greps for that message was aligned in rc8.

Commits in this tag

• 46f1768 — remove the Martini audit plan and parity summary documents
• 793955f — comment pass over src/ and tests/
• cf7ae2c — ruff format over the test tree
• dd3671b — comment pass over the example scripts
• 0e0d036 — release commit: version to 1.0.0rc7, changelog entry

linstabpy v1.0.0rc6

The supplement examples ship complete: Martini 0-5, Song 32 deg, and
Savarino, with their parity audits, slim reference fixtures, and the
paper-style figures. Two operator corrections came out of that audit.
On a deformed y-symmetric mesh the metric terms now come from the full
eta operators, so the family second derivatives converge at the axis.
The Floquet shift defaults to the continuum chain rule; the LinStab2D
discrete form stays behind matlab_parity.

Added

• Each Savarino map-walk row records its commit, its Slurm job, and its machine.
• resolvent_analysis returns the factor, SVD, and residual stage clocks.
• The Savarino line campaign times a recycled sweep against a factor at each frequency.
• Song (2023) example complete (examples/song/, θ = 32° only): self-contained
  supplement pipeline from Song FV export through LinStabPy validation figures.
  • prepare_baseflow.py — Plot3D → nodal .npy (601→599 trim, interpolation, nondim).
  • plot_baseflow.py — mesh + mean-flow panels (examples.figure_style LaTeX fonts).
  • configs/run_fd4_n30.yaml — FD4, n = 30, β = 19.2, no operator dissipation (ε = 0).
  • song.py — EVP driver; default config → configs/run_fd4_n30.yaml.
  • report_matched20.py, plot_spectrum_matched20.py, plot_eigenmode_lead.py —
    Hungarian matched-20 spectrum and lead eigenmode vs Song reference.
  • song_plot_utils.py, reference_paths.py; shipped reference under reference/.
  • Shipped validation figures under figures/ (baseflow mesh/fields, matched-20
    spectrum, lead eigenmode; FD4 n = 30, ε = 0).
  • Local-only outputs: inputs/baseflow/, output/fd4_n30/ (gitignored).
  • Validated (FD4 n = 30, ε = 0): 20/20 Song modes within |Δω| ≤ 0.01; lead |Δω| ≈ 0.002.
• Martini (2024) examples 0–5 — presentation + parity closure (examples/martini/):
  audit CLOSED 2026-08-18 (AUDIT_PLAN.md; steps A→E, user sign-off all six cases).
  • Per-case plotters: example_{0..5}_plots.py plus shared
    plot_compare_utils.py, split_martini_fixtures.py,
    run_all_martini_examples_and_plots.py — thesis-style mosaics (LinStab2D |
    linstabpy | peak-normalized % error), spectra, and baseflow panels.
  • Slim MATLAB reference fixtures under examples/martini/LinStab2DResults/
    (compact u+v / u+v+w, complex128); bundled linstab2D_scripts/ and
    regenerate_martini_fixtures.sh for MATLAB → NPZ regen.
  • Bundled baseflow .mat files under examples/martini/baseflows/
    (BL_BF.mat, Step_BF.mat, M09jet.mat, Paraboloid_BF.mat).
  • Shipped presentation PNGs next to each driver (baseflow, spectrum, mode
    mosaics for Ex. 1–5; validation triptych for Ex. 0).
• Savarino (2026) example — reorganized supplement reproduction (examples/savarino/,
  θ = 30.8° SWBLI resolvent validation, Figs. S4–S6):
  • Layout: configs/ (YAMLs + grids), reference/ (author baseflow + packed
    response modes), output/ (run records, figures, local intermediates).
  • Configs (configs/): map_base.yaml, mack_f2_b45.yaml,
    liftup_f1em4_b163.yaml, point.sbatch.template, savarino_map317_points.csv
    (317-pt β–f grid for dense-map campaigns).
  • Shipped reference: Tecplot base flow under reference/baseflow/; packed
    reference/savarino_response_{mack,liftup}.npz (primitive ρ′, u′, v′, w′, T′).
  • Shipped results (output/): savarino_map2_200pt.csv (200-pt gain-map run
    record), savarino_map317pt.csv (317-pt minimal dense subset from PhD
    map_minimal_keep_author200), validation figures under output/figures/
    (mode comparisons + 200-/317-pt gain maps).
  • collect_campaign.py --from-npz — aggregate NPZ + .metrics.json trees into
    the run-record CSV layout (for PhD dense-map campaigns without Slurm logs).
  • gen_campaign.py --grid {author200,dense317} — Slurm/local campaign generator
    for 200- or 317-point grids.
  • Plotters aligned with PhD Savarino2026 style:
    • plot_gain_map.py — tensor/PCHIP interpolation, σ₁² colorbar with ×10⁴
      on top, Mack / lift-up peak markers (plot_resolvent_map_savarino_style.py
      parity); optional --show-points / --baseline-summary for diagnostics.
    • compare_modes.py — Song eigenmode layout, peak-normalised ±1 colourbar,
      ymax 0.06 (Mack) / 0.1 (lift-up) (plot_savarino_modes_linstabpy_template.py
      parity).

Changed

• Recycled GMRES restarts. The restart is a Krylov length. It is never the system size.
• certified_rsvd logs every check when one check fails. One seed skips the multi-seed check.
• Floquet shift default is now continuum (LINSTABPY_FLOQUET_SHIFT,
  numerics.floquet_shift): DZ += iα I, D2RZ += iα DR, D2Z += 2iα DZ − α² I,
  the chain rule of (∂z + iα). The LinStab2D algebra (D2RZ += iα I and D2Z
  built from the already-shifted DZ, i.e. −3α² I) stays available as
  matlab_parity; the Martini Ex. 3 YAMLs and examples/sawtooth pin it.
  Plane Poiseuille at Re = 10⁴, α = 1 (Orszag 1971): continuum converges to
  ω = 0.23752649 + 0.00373967i (|Δω| 3.8e-4 → 1.8e-6 for ny 128 → 512);
  matlab_parity stalls at |Δω| ≈ 3.2e-4, a growth rate 8.5 % low, at every
  resolution. Any unpinned Floquet run changes with this default.
• Replaced the legacy Song figure/pack/matrix scripts with the pipeline above;
  see examples/song/README.md for run order.
• song.py creates the output NPZ parent directory before save and reads
  numerics.operator_dissipation from YAML (disabled, ε = 0 for this test).
• song_plot_utils.load_song_leading_mode_npz transposes the shipped reference
  eigenmode NPZ to LinStabPy (ny, nx) layout for plotting.
• .gitignore re-ignores regenerated example outputs (examples/**/figures/,
  examples/**/output/, Song inputs/baseflow/, logs) after the example PNG
  whitelist; Song figures/ and reference NPZ are shipped; raw Plot3D exports
  under inputs/raw/ stay local.
• Song plot scripts require LaTeX via figure_style.check_plot_environment;
  see examples/song/README.md for system packages (texlive-latex-extra, …).
• examples.figure_style: check_plot_environment(require_latex=…) with install
  hint, optional dpi / pad_inches on save(); tests in tests/test_figure_style.py.
• Grouped Song (2023) reference assets under examples/song/reference/:
  baseflow/song_32deg.npz, eigenmode/song_32deg_leading_mode.npz,
  spectrum/eigenvalues.dat, spectrum/omega_i.csv. Added reference_paths.py
  and updated configs, figure scripts, and tests.
• Martini fixture layout and git policy:
  • Shipped NPZ: examples/martini/LinStab2DResults/*.npz only (slim MATLAB
    reference). Root-level example_*.npz from linstabpy drivers are
    gitignored (regenerated locally; several exceed the 60 MiB tracked-file
    guard).
  • Raw MATLAB .mat exports stay outside Git in martini_examples/ (gitignored
    scratch; recreated by regenerate_martini_fixtures.sh).
  • Run logs under examples/martini/_run_logs/ are gitignored.
  • Renamed examples/martini/data/ → baseflows/.
• Martini plot/compare layer (repo mosaics ≤ 0.49% peak-norm vs campaign;
  Ex. 4 cavity audit closed — residual tail is MATLAB eigs convergence, not
  operator/BC):
  • plot_compare_utils.py: used_ind 0-based gridding, matlab_plot_grid_shape,
    compact_resolvent_modes_uvw, peak-normalized mosaic triptychs.
  • example_4_plots.py / split_martini_fixtures.py: slim MATLAB u+v+w compare
    with joint complex scaling; Example 4 YAML pins Pr: 0.7 and
    zero_dmudT: true for LinStab2D parity.
  • compare_martini_examples.py and tests/test_martini_parity_compare.py
    updated for LinStab2DResults/ fixture paths and compact mode layout.
• Moved library validation CLIs from examples/martini/ to tests/tools/
  (assembly_parity.py, assembly_metrics_diff.py, lnse_equivalence_harness.py,
  matlab_source_parity.py, symbolic_compressible_nse.py,
  coefficients_reference_base.py, NSE_DERIVATION.md). See
  tests/tools/README.md for rationale and owner review checklist.
• Removed unused maintainer-only scripts from examples/martini/:
  check_parity_matrices.py, compare_operator_construction.py,
  perf_memory_operator_parity.py, operator_matrix.json, branch_matrix.json,
  martini_c64_max_rel.txt, martini_examples.sbatch.
• Savarino example layout and git policy:
  • Moved point YAMLs, campaign template, and 317-pt grid into configs/; run
    records into output/ (savarino_map2_200pt.csv was at example root).
  • reference_paths.py: GAIN_MAP_CSV, GAIN_MAP_317_CSV, figure paths under
    output/figures/; compare outputs under output/compare/.
  • .gitignore ships only curated Savarino outputs: both map CSVs and
    output/figures/*.png. Gitignored: output/baseflow/*.npy (from
    savarino_to_npy.py), output/compare/**/*.npz (from run_compare_points.py),
    reference/zenodo/. Whitelist repeated after root output/ rule.
  • run_compare_points.py, gen_campaign.py, collect_campaign.py, mode_oracle.py
    updated for configs/ paths.
  • See examples/savarino/README.md for minimal workflow and maintainer scripts
    (verify_norm_convention.py, mode_oracle.py, fetch_zenodo.py,
    pack_reference_modes.py).

Fixed

• Mesh.deform on a y-symmetric mesh with a non-separable map: the dJ/Q
  metric terms of the symm/asymm family operators are now differentiated with
  the full Dy/D2y (the coordinate map has no η parity), as LinStab2D
  DeformMesh.m does. Before, the family stencils were applied to the
  coordinates themselves and D2y_asymm/D2y_symm did not converge at the axis
  rows (relative error growing like 1/h under refinement: 7.7 → 61 for
  n = 16 → 128 on a polynomial map; 3.5 → 29 on the Z² map). After: 2.2e-4 → 3.9e-7
  (order 3, asymm) and 2.1e-3 → 3.1e-5 (order 2, symm). Only `y_symmetry...

linstabpy v1.0.0rc4

[1.0.0rc4] - 2026-08-13

The Example 2 jet resolvent matched the MATLAB LinStab2D gains. The axis
coefficient used the true radius. The spatial filter sat inside the
forcing and response masks.

A forcing region composed with the boundary-condition mask.
input_region.space: conservative_momentum mapped conservative momentum
forcing onto the primitive residual.

Case files became YAML. One loader read a packed field. The tree dropped
JSONC, import shims, scripts/, and the unused Dockerfile.

Fixed

• compute_cylindrical_coefficients no longer floored R at the smallest
  positive spacing. A grid that contains r = 0 raised ValueError.
• ResolventOperator applied the spatial filter inside the input and output
  masks, matching MATLAB resolvent.m.
• The Example 2 leading gain moved from 411.72 to the MATLAB 413.72 reference.
• A forcing or response region composed with the boundary-condition mask.
• mask_bc_forcing stayed true by default. A region no longer forced
  constraint rows.
• A region that selected only constrained rows raised.
• mask_bc_forcing: false still forced those rows on purpose.
• The Floquet shift also updated D2RZ. The default mode matched LinStab2D.
• The spatial EVP finite-differenced L in kx, not in the Floquet exponent.
• Example 3 applied spatial Dirichlet after Lc + omega Lw.

Added

• input_region.space: conservative_momentum mapped conservative momentum
  forcing onto the primitive residual with a mass-weighted L2 norm.
• Momentum forcing induced a consistent temperature forcing
  T' = -kappa(kappa-1) Ma^2 (ubar · u').
• load_config loaded YAML case files.
• load_case_fields read a packed npy directory or one npz.
• examples/song and examples/savarino shipped the reproduction chains.
• examples/sawtooth added the Lui 2022 field on a packed stride-2 300 x 340 mesh.
• Four paper Example 3 YAML files reused example_3_streaky_bl.py.
• examples/martini/matlab_baseflow.py held the Martini MATLAB field constructors.
• build_runtime_payload grouped the runtime NPZ keys before save_results.
• install.sh --diagnose printed the host preflight and installed nothing.
• tests/tracked_file_size.py and tests/public_provenance.py replaced scripts/.

Changed

• Example figures probed TeX by writing a PNG. They used usetex when that worked.
• A failed usetex save retried with Computer Modern mathtext.
• A missing TeX stack no longer crashed the figure.
• Case files under examples/ became YAML. The loader rejected a leftover .jsonc suffix.
• Shipped YAML files now comment each key. validation_cases.yaml is YAML.
• Comments and the tutorial now say YAML, not JSONC.
• Sponge presets dropped the case1_ prefix.
• Moved the Martini, Song, and Savarino cases to examples/martini,
  examples/song, and examples/savarino.
• Packed the Song reference fields as float64 npz. The 32 deg field uses
  velocity / c_inf and temperature / T_inf.
• Packed Martini mode arrays as complex64. The largest file is
  example_3_Streaky_BL_results_temporal_modes.npz.
• Raised the tracked-file size limit from 50 MiB to 60 MiB for that file.

Removed

• load_jsonc, load_hamada_case, results_runtime, and the solver import shims.
• The run/ campaign kit and the broken Docker PETSc CI job.
• The mumps_icntl_7 boolean alias. Use mumps_ordering: metis.
• The unused Dockerfile that did not configure PETSc.
• The scripts/ folder and diagnose.sh.

linstabpy v1.0.0rc3

[1.0.0rc3] - 2026-08-10

Changed

• BREAKING (default): numerics.operator_dissipation.apply_to_resolvent now
  defaults to true instead of false, matching apply_to_stability. Previously,
  enabling dissipation gave a damped eigenvalue problem and a bare resolvent
  operator — a combination that runs to completion and returns grid-scale sawtooth
  (spurious high-wavenumber streaks, gains orders of magnitude too large) with nothing
  in the log to indicate it. Configs that deliberately want an undamped resolvent must
  now set apply_to_resolvent: false explicitly. No config shipped with the project
  changes behaviour: every one that runs a resolvent either disables dissipation or sets
  the switch explicitly. Both defaults now live in one place,
  config.defaults.DISSIPATION_APPLY_DEFAULTS, so the driver and the preflight checks
  cannot disagree about them.
• solver.factor_backend is now honoured as a fallback for all solver stages
  (linsolve, resolvent, and the shift-invert eigen path). It was previously read by
  nothing: every resolver receives a sub-config, so a key at solver level silently did
  nothing — one symptom being a run intended for GPU that executed entirely on CPU. It
  sits at the bottom of the precedence chain (per-stage key > environment >
  solver.factor_backend > built-in default), so no configuration that resolved before
  resolves differently now. Route-specific defaults (for example strumpack_fgmres
  requiring strumpack) are unaffected.

Added

• Preflight and post-run advisories (driver.preflight), reported in the run log. Each
  advisory names a concrete remedy rather than only a symptom. Preflight covers:
  dissipation enabled but not reaching a stage that is about to run; no dissipation at
  all; dissipation enabled with epsilon = 0; a spanwise wavenumber inside the band
  where grid-scale sawtooth dominates while undamped; a base flow that is not normalised
  to the canonical freestream (rho = u = T = 1); a shadowed factor_backend; and a
  strided (reduced-resolution) mesh. Post-run checks flag a leading gain far above the
  physical scale for these cases and a postcheck residual above tolerance.

• A run banner and closing summary in the log, reporting effective values —
  including the per-stage dissipation epsilon actually applied — so a run's real
  configuration and the trustworthiness of its result are visible without loading the
  output archive.

• Translation of opaque factorization failures into actionable messages. A MUMPS
  workspace exhaustion (INFOG(1) = -9) previously surfaced only as
  SLEPc EPS.setUp failed; it now names the cause and the remedy (ICNTL(14)).

• A band-marking backend advisory. resolvent_map runs a mandatory band-marking
  eigenvalue problem, and that eigenproblem reads solver.eigen — not
  solver.resolvent. A configuration that moved only the resolvent onto a scalable
  backend therefore left the eigenproblem on the default direct solver, where at full
  resolution it does not fail but simply never finishes; the run looks like a slow
  resolvent sweep for hours. The advisory names the three ways out (solver.factor_backend,
  an explicit solver.eigen.factor_backend, or analysis.skip_band_marking).

• A source-wide lint (tests/test_mpi_collective_lint.py) rejecting MPI collectives
  inside rank-0-only regions, in both the if rank == 0: and if rank != 0: return
  forms. Such a call deadlocks every multi-rank run while remaining invisible to the
  serial test suite, which runs at size 1. The lint ships with positive controls, since
  a lint never shown to fire cannot be distinguished from one that cannot fire.

Fixed

• The single-point resolvent path now honours input_region and output_region, and no
  longer forces the response mask to equal the forcing mask. driver.run._stage_run_resolvent
  built its own boundary-condition-only diagonal and then set C_mat = B_mat.copy(), so two
  things went wrong at once and neither said anything: a config naming a forcing region ran
  with no region at all, because only solvers.router.resolvent_map ever read those keys; and
  even a correct region was copied onto the response, which turns an asymmetric problem
  (restricted forcing, free response) into a different one. Both paths now call the same
  _forcing_response_mats helper, so one config key means one thing everywhere, and C
  defaults to the boundary-condition mask rather than to B. Runs that set no region are
  unchanged: B and C are still the boundary-condition mask, asserted byte-for-byte by a
  regression test. scripts/validation/song_2023/song_2023.py builds its masks separately and
  cannot honour a region, so it now raises instead of ignoring the key.
• MPI is now initialised at MPI_THREAD_MULTIPLE on every import path, restoring the
  STRUMPACK factorization backend. Whichever library calls MPI_Init first fixes the
  thread level for the process, and it can only be requested at init — never raised
  afterwards. mpi4py requests MPI_THREAD_MULTIPLE, PETSc only
  MPI_THREAD_FUNNELED. utils.petsc_init.init_petsc_from_argv already imported
  mpi4py first, but seven modules imported petsc4py/slepc4py at module scope
  and bypassed it, so whichever of them was imported first decided the level. When
  PETSc won, MatGetFactor_aij_strumpack() failed with PETSc error 57 and SLATE requires MPI_THREAD_MULTIPLE, taking the whole STRUMPACK/GPU route out of service —
  8 of 14 failures in the GPU test venv had this single cause, and none of them
  reported anything about MPI. The guard is a plain import mpi4py.MPI, chosen over a
  first-party helper module because isort orders third-party imports above first-party
  ones: a linstabpy.* guard could not stay above petsc4py, and ruff --fix would
  have silently reordered the fix away. mpi4py sorts above petsc4py on its own.
  tests/test_mpi_thread_level.py enforces the ordering as a source-wide AST lint plus
  a runtime thread-level assertion that skips loudly without petsc4py.
• The freestream-normalisation check and the MUMPS error translation were defined but
  never called: both shipped with full unit-test coverage and zero call sites, so
  neither ran in production. They are now invoked from the driver and the eigensolver
  respectively, and a test asserts that every advisory producer is reachable from the
  runtime — unit tests establish that a guard works, not that anything calls it.
• ASM/BJacobi sub-preconditioner settings on the inexact shift-invert path are now
  actually applied. They were configured before eps.setUp(), where PETSc cannot yet
  materialise sub-solvers (ierr=58), and the error was swallowed — so sub_pc_type
  and pc_factor_levels were discarded. PETSc's default sub-PC is also ilu, which
  masked the type being dropped, but the fill level was not masked: a config asking for
  ILU(2) ran with ILU(0). The settings are now applied immediately after eps.setUp(),
  failures raise instead of downgrading, and the realized type and fill level are logged
  rather than assumed. This mirrors what linalg.factorize already did for the linsolve
  path. Only st_solver: "inexact" reached the defect — no shipped configuration selects
  it — so no previously reported result changes.
• MATLAB-parity tests are opt-in. This repository ships no MATLAB reference data,
  so the seven tests that compare against it called pytest.skip when it was absent —
  which is every machine but the two the data lives on. A run reporting
  "1361 passed, 13 skipped" read as near-total coverage without disclosing that seven
  of those skips could never be satisfied by anything in the repository. They now carry
  a matlab_parity marker, deselected by default, and are requested explicitly with
  LINSTABPY_MATLAB_PARITY=1 pytest -m matlab_parity. When requested, a missing
  reference is a failure, not a skip: a comparison you asked for that quietly does
  not run is indistinguishable from one that passed. A default run now reports
  6 skips, every one an absent optional backend — a true statement about the machine.
• Subprocess timeouts in the test suite are now load-aware everywhere, via a single
  shared helper. One test file already scaled its timeout with load average; its two
  siblings used hard limits of 120 s and 60 s, leaving as little as 1.4x headroom over
  the test's own runtime. On a shared or busy machine those expired and surfaced as
  test failures, indistinguishable from a numerical disagreement — a contributor
  running the suite with a browser open would see phantom errors. These timeouts guard
  against a hang; they are not performance assertions, so the bases are deliberately
  generous. The helper documents a limitation found in practice: load average counts
  runnable processes and does not detect memory-bandwidth starvation, so a co-resident
  memory-bound job slowed a test run ~3x while the computed multiplier was 1.17x.
• production_case_template.jsonc, the file the documentation tells a new user to copy
  for a production run, shipped with no operator_dissipation block at all. A first
  production run made from it was sawtooth-contaminated with nothing in the log to say
  so. It now enables dissipation and documents that epsilon must be re-tuned at the
  target resolution rather than inherited.
• The Streaky-BL MATLAB parity test guarded only one of the two references it needs, so
  restoring the parity bundle while the base-flow data directory was still absent turned
  a clean skip into a hard error that resembled a numerical disagreement. The guard now
  checks every reference and announces which one is missing.
• PETSc is now imported lazily by linalg.backends, linalg.factorize,
  linalg.shifted, and solvers.router.stable_sweep, matching the existing
  utils.petsc_core pattern. Importing these modul...