Input Files and Executable

Every GKX input is one runtime TOML deck, loaded by gkx.load_runtime_from_toml into a RuntimeConfig. The [grid], [time], [geometry] and [init] tables merge into GridConfig, TimeConfig, GeometryConfig and InitializationConfig; species, physics toggles, collisions, normalization, operator terms, output, quasilinear and parallel policy have their own tables (listed under Runtime sections). The same deck drives the executable and the Python API, for Cyclone, ETG, KBM or any other case, through one solver path.

Unified Runtime Schema

Resolution and time-step defaults

Velocity resolution. When neither the call nor the deck names Nl (Laguerre) and Nm (Hermite), the fallback is (Nl, Nm) = (12, 24) for linear runs and (4, 8) for nonlinear runs. Set both explicitly for any result you report.
ky layout. [grid] ky_layout = "half" (default) stores the 1 + Ny // 2 non-negative rows; "full" stores the two-sided ky axis of length Ny (gkx.core_ky_layout). Results do not depend on the choice: they agree to round-off, and "half" is faster and smaller, in the nonlinear step and in the heat-flux window gradient.
Time stepping. A deck that sets [time] dt runs rk2 at that fixed step unless it names another method. A deck without dt gets method = "rk4" with fixed_dt = false: the CFL controller chooses the step with cfl = 0.9. An explicit method or fixed_dt in the deck always wins.

Runtime precision

The default is float32 (complex64 states). JAX_ENABLE_X64=true, set before JAX is imported, selects double precision. The runtime’s own initial state is then built in complex128, so the Krylov basis, the eigenvector and the certification gate are all float64. This is what the shipped decks mean when their header says to run them under JAX_ENABLE_X64 for parity reproduction.

What x64 buys is the gate, not necessarily the eigenvalue. certifiable_residual_tolerance floors the eigenpair residual gate at 1e3 * eps(dtype), so a float32 run certifies against 1.19e-4 while a float64 run certifies against 1e-9. Use float64 when the number you are quoting has to carry a float64-scale residual.

A float64 run on a deck too coarse to support a certifiable eigenpair raises instead of returning an uncertified number. If that happens, raise the velocity resolution (Nl, Nm).

GKX_X64 is a test-harness variable, not a runtime switch: it only selects the precision banner in tests/conftest.py. Setting it without JAX_ENABLE_X64 does not change any run’s precision.

run_runtime_linear(cfg, initial_state=G128, ...) accepts an explicit seed and preserves its dtype. It requires shape (kinetic_species, Nl, Nm, grid.ky.size, grid.kx.size, grid.z.size) on the selected-ky runtime grid. Record the seed dtype and hash; casting a rounded f32 seed does not recover lost input precision. There is no TOML dtype option, and prepared linear solve(initial_state=...) raises NotImplementedError.

Minimal runtime TOML example

schema_version = 1

[[species]]
name = "ion"
charge = 1.0
mass = 1.0
density = 1.0
temperature = 1.0
tprim = 2.49
fprim = 0.8
kinetic = true

[physics]
linear = true
nonlinear = false
electrostatic = true
electromagnetic = false
adiabatic_electrons = true
tau_e = 1.0

[normalization]
contract = "cyclone"
diagnostic_norm = "rho_star"

[terms]
streaming = 1.0
mirror = 1.0
curvature = 1.0
gradb = 1.0
diamagnetic = 1.0
collisions = 1.0
hypercollisions = 1.0
end_damping = 1.0
apar = 0.0
bpar = 0.0
nonlinear = 0.0

[run]
ky = 0.3
Nl = 16
Nm = 48
solver = "auto"
fit_signal = "auto"

Nl = 16, Nm = 48 is the velocity resolution of the shipped Cyclone deck examples/01_linear_tokamak/case_full.toml.

Schema compatibility

Maintained runtime decks declare schema_version = 1 at the document root. Versionless GKX 1.8.2 decks remain readable as legacy schema 0 and are not rewritten in place. Resolved decks written by the equilibrium shorthand are schema 1. A non-integer or unsupported future version fails before numerical work and points to this migration contract; upgrade GKX before reading a newer schema. The removed Diffrax keys (use_diffrax, diffrax_solver, diffrax_adaptive, diffrax_rtol, diffrax_atol, diffrax_max_steps) in [time] are refused with the native replacement named in the error.

Quasilinear diagnostics

Linear runtime runs can compute electrostatic quasilinear transport weights directly from the final linear state or Krylov eigenvector:

[quasilinear]
enabled = true
mode = "weights"                    # weights | saturated
saturation_rule = "none"            # none | mixing_length | lapillonne_2011
amplitude_normalization = "phi_rms" # phi_rms | phi_midplane | field_energy
kperp_average = "phi_weighted"
csat = 1.0
gamma_floor = 0.0
include_stable_modes = false
channels = ["es"]

The validated output level is linear weights plus optional uncalibrated electrostatic saturation rules. Electromagnetic quasilinear channels are rejected. The runtime writes *.quasilinear.summary.json and *.quasilinear_species.csv when [output].path or --out is set. Serial scan-runtime-linear runs also write *.quasilinear_spectrum.csv; combined-ky batched scans refuse quasilinear output. In scan spectra, the ky column is the requested scan coordinate and mode_ky records the signed selected grid-mode coordinate; they can differ for linked-boundary grids, so publication plots should use ky while audits can inspect mode_ky.

Equivalent executable flags are available for single-point runtime runs:

gkx run-runtime-linear \
  --config examples/01_linear_tokamak/case_full.toml \
  --quasilinear \
  --ql-mode saturated \
  --ql-saturation-rule mixing_length \
  --ql-csat 1.0 \
  --out tools_out/cyclone_quasilinear

Time-integration controls

Keys in [time] besides t_max, dt and method:

compressed_real_fft (default true): compressed real-FFT nonlinear bracket; false uses a full complex FFT for the nonlinear term.
sample_stride (record every N steps) and diagnostics_stride (compute streaming diagnostics every N steps). diagnostics = false (or --no-diagnostics) disables diagnostics entirely.
fixed_dt = false with cfl (default 0.9) and optional cfl_fac / dt_min / dt_max: CFL-controlled steps. When cfl_fac is omitted, rk3/sspx3 use 1.73, rk4 uses 2.82, and other methods 1.0. Adaptive runs record dt_t (per-sample step history) and dt_mean (average effective step). The nonlinear adaptive dt estimate combines the linear frequency cap with the instantaneous nonlinear cap. Leaving dt_max unset keeps dt_max = dt; raise it explicitly only when you trade strict comparison matching for throughput.
method = "rk4-pi" (nonlinear, fixed_dt = false): RK4 whose step is chosen from an embedded error estimate and capped at twice the CFL step; see Solvers. cfl_fac defaults to 2.82 as for rk4.
laguerre_nonlinear_mode: "grid" (reference quadrature, default) or "spectral" (spectral Jl without the quadrature transform) for the Laguerre handling in nonlinear brackets.
nonlinear_dealias (default true): 2/3 dealias masking. With masking on, nonlinear mode selection is dealias-aware: a requested ky removed by the mask is replaced by the nearest retained ky, and the executable prints the effective ky_sel/kx_sel. false keeps all configured modes for reference/debug runs.
state_sharding: "auto" or "ky" shards the packed state over multiple JAX devices along ky; "kx", "z", "l", "m" and "species" shard along that axis. Sharding-aware integration paths honor it; unsupported solver paths and one-device runs fall back to single-device execution. Treat it as a correctness-gated option unless the run also has a matched scaling artifact.

Diagonal nonlinear hypercollision splitting is enabled with collision_split = true. The conserving collision operator remains in the RHS because its field-particle correction is not diagonal. The collision_scheme key selects the hypercollision update: implicit (backward-Euler), exp (exact diagonal exponential), and sts/rkc aliases (treated as stabilized explicit/exponential updates for diagonal operators).

Geometry controls

[geometry] drift_scale selects the benchmark-compatible drift convention (1.0, the default) or the alternate doubled-drift convention (2.0). The physical meaning of the runtime terms and geometry coefficients is detailed in Theory and Operators And Terms; the TOML layer here documents how those implemented models are selected and parameterized. For slab benchmarks, set model = "slab". Optional slab-specific keys are z0 (sets gradpar = 1/z0 when positive, matching the reference slab domain normalization) and zero_shat = true (forces the zero-shear slab metric gds2 = 1, gds21 = 0, gds22 = 1). It also accepts model = "imported-netcdf" with geometry_file = "external_geometry.nc" to run from imported sampled field-line geometry instead of the analytic s-alpha model. The imported file can be a tracked benchmark *.out.nc or a root-level *.eik.nc geometry file produced by the VMEC workflow. When that imported geometry is used with a linked boundary, GKX follows the file’s own theta range, jtwist/x0 geometry factor, and kxfac metadata instead of forcing the analytic s-alpha grid defaults. For direct VMEC workflows, the runtime also accepts model = "vmec". In that mode GKX calls the VMEC geometry helper to generate a matching *.eik.nc file on demand, then immediately reuses the same imported geometry path as the VMEC examples. Set vmec_file plus the flux-tube keys torflux, npol and optionally alpha. geometry_file can be used as an explicit output path for the generated *.eik.nc file, and geometry_helper_repo can point to a non-default helper checkout if needed. The VMEC path uses the booz_xform_jax backend, which is installed with GKX; GKX_BOOZ_XFORM_JAX_PATH or BOOZ_XFORM_JAX_PATH point to a source checkout instead. This is the recommended imported-geometry route for new stellarator cases. The shipped VMEC TOMLs point to locally generated files under examples/vmec; generate those files with vmex input.<case> before running the examples. vmec_file supports $ENV_VAR expansion, and relative paths are resolved against the TOML directory first. Command-line overrides are resolved from the shell working directory. Use geometry_helper_repo and geometry_helper_python when an imported geometry helper checkout or interpreter must be selected explicitly. When geometry_file is set for model = "vmec", GKX regenerates that target instead of reusing a stale file from an older VMEC conversion. For VMEC fix aspect runs, GKX follows the helper default contract and does not inject x0 from the runtime Lx. That keeps the generated *.eik.nc file aligned with the imported W7-X/HSX geometry output.

For Miller tokamak workflows, the runtime also accepts model = "miller". In that mode GKX uses its in-package Miller backend to generate a matching root-level *.eiknc.nc file, then re-enters the same imported geometry path used for VMEC eik.nc files. Set the Miller inputs directly in [geometry]: rhoc, q, s_hat, R0, optional R_geo, shift, akappa, akappri, tri, tripri, and betaprim. The generated Miller surface has inverse aspect ratio rhoc / R0. Set that ratio explicitly when reproducing a benchmark; the separate epsilon field belongs to analytic geometry metadata and does not replace the surface radius. geometry_file can be used as an explicit output path for the generated Miller *.eiknc.nc file. Existing output is reused by default; pass --force to gkx geometry miller to regenerate it. The geometry_helper_python and geometry_helper_repo settings apply to VMEC helper discovery, not to this in-package Miller backend.

Executable path overrides

The gkx executable accepts path overrides:

gkx run --config case.toml --vmec-file examples/vmec/wout_circular_tokamak.nc  # built by generate_wouts.sh
gkx run --config case.toml --geometry-file external_geometry.eik.nc
gkx run-runtime-nonlinear --config case.toml --init-file ~/restart.nc

These override paths expand ~ and environment variables and are resolved against the shell’s current working directory. TOML paths still resolve against the config file directory. --vmec-file replaces [geometry].vmec_file for VMEC-backed configs. --geometry-file replaces [geometry].geometry_file only; it does not change [geometry].model. For imported-geometry configs, that file is the imported EIK/NetCDF input. For model = "vmec", it remains the generated EIK target/cache path.

Solver and fit-signal keys

The [run] and [scan] sections accept solver and fit_signal keys:

solver = "auto" (default): choose time vs Krylov and fall back if needed
solver = "time": always use time integration
solver = "explicit_time": force the explicit single-mode time integrator used by controlled explicit-time comparisons
solver = "krylov": always use the eigensolver. Its default route is KrylovConfig(method="auto"): the Hermite-banded direct shift-invert (banded) when its block arrays fit in memory, else adaptive. Either returns an eigenpair only when its residual passes the certification gate (see Runtime precision). Other KrylovConfig.method values are adaptive, power, propagator, arnoldi, shift_invert (matrix-free GMRES inner solves) and sparse_shift_invert (assembles the sparse operator and factors it once on the host).
fit_signal = "auto" (default): pick phi vs density based on fit quality
fit_signal = "phi": use the electrostatic potential time trace
fit_signal = "density": use the density moment time trace

For large ky scans, scan-runtime-linear --batch-ky integrates all ky values in a single time integration pass (time integrator only) and then extracts the growth rates from the per-ky traces.

For quasilinear spectra, use scan-runtime-linear --workers N instead of --batch-ky. This runs independent per-ky solves, computes the quasilinear state extraction for each mode, preserves serial spectrum ordering, and records the worker identity contract in the scan summary JSON. --batch-ky with [quasilinear] enabled = true raises ValueError; quasilinear scans require per-ky evaluation.

Executable usage

gkx examples/01_linear_tokamak/case_full.toml    # same as: gkx run --config ...
gkx scan-runtime-linear --config benchmarks/cases/etg_linear_scan.toml --out tools_out/runtime_etg_scan
gkx plot tools_out/runtime_etg_scan.scan.csv
gkx run-runtime-linear --config examples/01_linear_tokamak/case_full.toml --out tools_out/cyclone_runtime
gkx scan-runtime-linear --config benchmarks/cases/etg_linear_scan.toml --batch-ky
gkx run-runtime-nonlinear --config examples/03_nonlinear_tokamak/case_full.toml --sample-stride 5 --out tools_out/nonlinear_cyclone_diag.csv

examples/01_linear_tokamak/case_full.toml is a production-length run (32168 rk4 steps at Nl = 16, Nm = 48); see Quickstart for a ten-second demo. gkx CASE.toml dispatches to gkx run --config CASE.toml, which picks the linear or nonlinear command from [physics].

For run-runtime-nonlinear with fixed_dt = false, omit --steps unless you want a capped step count: the runtime then integrates in chunks until it reaches t_max (or saturates, below). With fixed_dt = true and no --steps, the step count is round(t_max / dt).

Diagnosed nonlinear runs stop at saturation by default ([time] run_to = "saturation"): the runtime integrates in chunks and, after each chunk, tests the accumulated heat-flux trace with spin-up excluded. It stops once the autocorrelation-corrected relative SEM of the windowed mean is at or below saturation_rel_sem, the retained window has at least 256 samples and spans at least 20 integrated autocorrelation times, any larger saturation_min_window is met, and the first and second halves agree within twice their combined SEM. The same drift-consistency guard is applied to the electrostatic field energy Wphi and gyrokinetic free energy Wg. These operational finite-sample safeguards provide neither a universal sequential- confidence guarantee nor qualification for irregular adaptive-time sampling. t_max remains the hard cap, so a run that never saturates integrates the full horizon. The chosen window and its mean +/- SEM are printed in the run summary and recorded under saturation in the summary JSON. Pass --no-until-saturated (or set run_to = "t_max") for a fixed horizon.

For single-point runtime commands, artifact output can be requested either with --out or directly in the runtime TOML:

[output]
path = "tools_out/runtime_case"

--out takes precedence over [output].path when both are provided. A plain prefix writes a JSON summary plus CSV sidecars; a nonlinear target ending in .out.nc (or another .nc suffix) writes the NetCDF bundle *.out.nc / *.big.nc / *.restart.nc instead. See Saved Outputs and Restart Files for every file, CSV column and NetCDF variable.

Python driver

# point CONFIG at the top of the script at benchmarks/cases/etg_linear_scan.toml
python examples/12_restart_and_analysis/run.py

The canonical KBM TOML is a controlled comparison input rather than a promoted standalone solve: its experimental shift-invert branch raises when the physical eigenpair-residual gate is not met. Reproduce the reviewed result with scripts/benchmarks/linear_benchmark.py kbm (retired; Retired generators) (removed from the tree; git show df2953026:scripts/benchmarks/linear_benchmark.py).

Runtime sections

load_runtime_from_toml reads these tables. Unknown keys in [grid], [time], [geometry] and [init] are ignored (the removed Diffrax [time] keys are refused). An unknown key in [[species]], [physics], [collisions], [normalization], [terms], [expert], [output], [quasilinear] or [parallel] raises TypeError. [run], [scan] and [fit] are read by the commands that use them.

schema_version (document root; see Schema compatibility)
[grid] (GridConfig), [time] (TimeConfig), [geometry] (GeometryConfig), [init] (InitializationConfig)
[[species]] (kinetic species definitions)
[physics] (electrostatic/electromagnetic, adiabatic/kinetic, linear/nonlinear)
[collisions] (collision and hypercollision controls)
[normalization] (contract key + optional overrides)
[terms] (term toggles used by modular RHS assembly)
[expert] (advanced fixed-mode controls for specialized workflows)
[output] (artifact path and restart controls)
[quasilinear] (quasilinear transport diagnostics)
[parallel] (parallelization policy; defaults to serial)
[damping_reference] (provenance written by migrate_end_damping_reference; see Explicit reference migration)
[run] / [scan] / [fit] (driver controls: ky, Nl, Nm, mu_vperp_max (nodal magnetic-moment grid, see Numerics), solver, method, dt, steps, the scan ky list, and the growth-rate window options auto_window, tmin, tmax, window_fraction, min_points, start_fraction, growth_weight, require_positive, min_amp_fraction, window_method, mode_method, fit_signal)

Notable runtime-only keys:

[time] damp_ends_rate: the field-line end-absorber rate $\nu_{end}$ in inverse simulation-time units, finite and nonnegative. Every route (time, Krylov, eigen, shift-invert, nonlinear) applies the same rate; no timestep enters. Omit it for the default, 30 thermal transits of the fastest kinetic species across the end layer. 0 turns the absorber off. See Parallel boundary and end damping.
[collisions] damp_ends_widthfrac: layer width as a fraction of the chain (default 0.125).
[physics] reduced_model: physics-family selector for runtime inputs. The maintained runtime supports full gyrokinetics via "gyrokinetic" and its full-GK aliases. Non-promoted reduced-model values fail closed with NotImplementedError instead of silently routing through the wrong equations.
[collisions] damp_ends_amp (removed): a positive value without damp_ends_rate on a non-periodic deck is refused, because GKX < 2.6 read it as damp_ends_amp/dt on the time route and as a rate elsewhere. damp_ends_amp = 0 still means no damping. damp_ends_scale_by_dt = true is refused too.
[collisions] hypercollisions_const / hypercollisions_kz: defaults are 0.0 / 1.0 (kz-proportional hypercollisions enabled by default, constant hypercollisions off). The parallel multiplier is spectral: $|\partial_z|G=\mathcal{F}^{-1}(|k_z|\mathcal{F}G)$, on each periodic line or linked chain. It is not pointwise multiplication of physical-z samples by the wavenumber array. GKX releases before 2.1.0 applied a pointwise multiplier on periodic lines: rerun periodic transport results from those versions that used nonzero kz damping and Hermite orders above two. Linked runs are unaffected.
[collisions] p_hyper_m: when omitted, the runtime path uses the resolution-aware default min(20, Nm/2) instead of a fixed exponent across Hermite resolutions.
[collisions] nu_hermite / nu_laguerre: coefficients entering the Lenard-Bernstein collision eigenvalue used by the modular collision kernel.
[collisions] nu_hyper / p_hyper: isotropic hypercollision amplitude and exponent.
[collisions] nu_hyper_l / nu_hyper_m / nu_hyper_lm and p_hyper_l / p_hyper_m / p_hyper_lm: Laguerre-only, Hermite-only, and mixed hypercollision channels. nu_hyper_l and nu_hyper_lm act only in the constant-coefficient branch – the $|k_z|$ kernel carries no Laguerre index – so a deck that declares either while hypercollisions_const = 0.0 (the default) is refused rather than run without the sink it asked for. Case validation also refuses either declaration when [physics] hypercollisions = false or [terms] hypercollisions = 0.0 disables the whole operator.
[collisions] nu_hyper_m_const: the constant branch’s Hermite rate. Omitted (the default) it reuses nu_hyper_m, which damps Hermite in both branches when both are on. Set it to 0.0 alongside hypercollisions_const = 1.0 to declare a pure Laguerre sink while the $|k_z|$ Hermite hypercollisions stay exactly as the reference code sets them. See Declared velocity-space regularization.
[collisions] D_hyper / p_hyper_kperp: perpendicular hyperdiffusion amplitude and exponent.
[normalization] flux_scale: multiplicative factor applied to heat/particle flux diagnostics (tracked comparison default 1.0).
[normalization] wphi_scale: multiplicative factor applied to Wphi diagnostics (default 1.0; no shipped configuration overrides it).
[normalization] rho_star: rescales the operator wavenumbers before building the drift/Bessel terms.
[normalization] omega_d_scale / omega_star_scale: multiplicative normalization factors for magnetic-drift and diamagnetic-drive terms.
[init] init_single with gaussian_init = false and init_single = false: initialize a random perturbation across the configured startup-loop bounds in (ky,kx).
[init] init_single with gaussian_init = true and init_single = false: initialize a Gaussian envelope across the same startup-loop bounds.
[init] init_single = true: initialize only the selected (ky,kx) mode. When combined with init_field = "phi" and gaussian_init = true, the selected electrostatic-potential mode is initialized with a Gaussian profile along the flux tube; this is the contract used by the W7-X zonal-flow response benchmark.
[init] init_field = "all": the runtime/TOML path uses Hermite/Laguerre moment-normalized scaling for this initializer, with reduced amplitudes for tpar (1/sqrt(2)) and qpar (1/sqrt(6)).
[init] init_field = "phi": initialize a requested electrostatic potential profile by inverting the same quasineutrality solve used during time advance. This is the appropriate contract for zonal-response literature tests that prescribe phi(t=0) rather than a density-moment perturbation; the masked ky=0, kx=0 gauge mode remains unavailable.
[init] init_electrons_only: if true in multispecies runs, initialize only electron species. If false (default), initialize all kinetic species.
[init] random_seed: RNG seed used for reference-compatible random initial conditions. The runtime follows a Linux glibc rand() sequence with positive-kx-major loop order and exact startup loop bounds, so random multi-mode perturbations are host-platform independent and reproduce the same seeded pattern on macOS and Linux.
[init] init_file: load a saved complex state from either the full-ky GKX layout or a packed positive-ky interchange layout.
[init] init_file_scale / init_file_mode: scale a loaded restart state and either replace the analytic seed (default) or add it to the fresh perturbation. This is the general runtime equivalent of restart-scaling and restart-with-perturb workflows.
[expert] fixed_mode with iky_fixed / ikx_fixed: keep one Fourier mode exactly frozen during nonlinear evolution; this is the eqfix-style contract used by the secondary benchmark.
[expert] source / phi_ext: runtime-only benchmark hooks for external electrostatic forcing. source = "phiext_full" with a small phi_ext injects the source into the solved phi field before the RHS is assembled, so it changes the actual evolution rather than only the saved diagnostics. The Merlo Case-III Rosenbluth-Hinton/GAM example uses the literature protocol instead: source = "default" plus an initial ion density perturbation.
[time] run_to: nonlinear stop policy. "saturation" (the default) ends a diagnosed nonlinear run as soon as the post-spin-up heat-flux window has converged, with t_max kept as the hard cap; "t_max" always integrates the full horizon. Linear runs and nonlinear runs with diagnostics = false ignore the key, because the stop test needs the streamed heat-flux trace. The executable overrides are --until-saturated and --no-until-saturated. Set "t_max" explicitly on any deck whose heat flux is not the observable it exists to produce – zonal-flow and secondary-instability responses, linear-response diagnostics – because a horizon chosen by a quantity that is zero by construction is not a horizon. Which shipped decks do that, and the measurement clearing the rest to run under the default, are held by the run_to audit in tests/release/test_release_gates.py.
[time] saturation_rel_sem: target relative standard error of the mean on the windowed heat flux (default 0.05, i.e. 5%). The SEM is corrected for autocorrelation using a first-negative-lag truncation estimate (function sokal_autocorrelation_time). It is not a calibrated sequential confidence guarantee. Fixed-horizon AR(1) coverage and a drifting-prefix rejection control are tested in tests/unit/quasilinear/test_quasilinear.py.
[time] saturation_min_window: minimum averaging-window span in time units before a run may stop. The effective requirement is the larger of this value and twenty integrated autocorrelation times; omitting it retains the correlation-time floor. Every accepted window also retains at least 256 samples. These are finite-sample safeguards calibrated on regularly sampled synthetic traces, not a guarantee for irregular adaptive-time sampling or arbitrary nonstationarity.
[time] nstep_restart: when writing a nonlinear NetCDF bundle, checkpoint every nstep_restart steps instead of waiting for the end of the run. This is useful for long adaptive runs and batch jobs.
[time] method = "sspx3": use the SSPx3 explicit scheme directly. This is the relevant explicit method for secondary and collisional-ETG benchmark families. Plain rk3 is the three-stage Heun-style scheme. The names euler, k10, rk3_heun and rk3_classic were retired and raise an error naming the replacement (rk2, rk4, rk3, rk3). Linear integration also accepts rk2, rk4, and the diagonal IMEX schemes imex and imex2 (ARS(2,2,2), second order).
[terms]: each key is a pure multiplicative operator weight: streaming, mirror, curvature, gradb, diamagnetic, collisions, hypercollisions, hyperdiffusion, end_damping, apar, bpar, nonlinear.

Runtime parallelization controls

The [parallel] section is parsed by RuntimeParallelConfig and is serial by default:

[parallel]
strategy = "serial"
axis = "ky"
batch_size = 2
num_devices = 2
strict_identity = true
profile = false
backend = "auto"

Accepted strategies are "serial", "batch", "combined_ky", "device_batch", "pmap", "pjit", "shard_map", "state", and "velocity". "batch-ky" is an alias for "combined_ky" (the combined-k_y time-integration scan path), and "none"/"off" are aliases for "serial". Quasilinear scan artifacts require per-k_y evaluation. auto = true resolves to strategy = "shard_map" with axis = "species_hermite".

For runtime k_y scans, strategy = "batch" with axis = "ky" selects the production independent-worker path. num_devices is interpreted as the requested worker count when the executable call does not explicitly pass workers; backend = "thread" or "process" selects the executor, and "auto" keeps the executable default. The runtime output records whether the worker policy came from the TOML file or from explicit executable arguments, together with the requested/effective worker counts and the ordering-preservation identity contract.

A k_y scan reads strategy at two levels, and each strategy belongs to exactly one of them. "batch", "combined_ky", "device_batch", "pmap" and "pjit" are scan-level: they say how the k_y points are dispatched, the scan honours them itself, and each worker then runs an ordinary single-k_y solve with strategy = "serial". Their reported growth rates are bit-identical to a serial scan. "velocity", "state" and "shard_map" are solver-level: they say how one solve is decomposed, so they reach the per-k_y solve unchanged and are honoured or refused there on the solver’s own terms. Only "velocity" has a linear RHS route today, so a k_y scan requesting "state" or "shard_map" raises rather than quietly running unsharded.

The velocity-space RHS route is diagnostic: strategy = "velocity", axis = "hermite", and backend = "streaming_only", backend = "streaming_electrostatic", or backend = "electrostatic_linear_slices". backend = "auto" selects the electrostatic-slices route when the active terms satisfy that gate; otherwise the runtime raises instead of silently falling back to an unvalidated path. These backends are accepted only by gkx.linear_rhs_parallel_cached. The first two require all non-streaming linear terms disabled. The electrostatic-slices backend allows only streaming, mirror, curvature, grad-B, and diamagnetic-drive weights; collisions, linked boundaries, electromagnetic terms, and nonlinear terms are refused. Velocity RHS routes are limited to single-species periodic 5D electrostatic states.

For full runtime TOML files this velocity-space route is exposed only through the fixed-step linear executable path with fit_signal = "phi". Density-assisted automatic fitting is serial only, so requesting velocity parallelization there raises an error.

For independent scan, sensitivity, and UQ workloads, use gkx.batch_map for JAX-array maps and gkx.independent_map for file-backed Python tasks such as calibration rows or leave-one-out UQ holdouts. Require a serial identity artifact before using timing results in a publication claim. The helpers preserve ordering, so diagnostics such as growth rates, frequencies, quasilinear weights, and covariance summaries can be carried through one parallel map.

Runtime output and restart controls

The [output] section controls runtime artifact layout and restart behavior:

path: artifact target. Use a plain prefix such as tools_out/runtime_case for JSON/CSV sidecars, or *.out.nc for a nonlinear NetCDF restart bundle.
restart: force loading from restart_from_file or from the derived sibling *.restart.nc next to path. Raise an error if the restart file is missing.
restart_if_exists: opportunistically resume from an existing restart file without requiring one to be present.
save_for_restart: write the *.restart.nc checkpoint when a nonlinear NetCDF bundle is requested.
restart_to_file / restart_from_file: explicit checkpoint paths when the default sibling naming is not desired.
restart_with_perturb: combine the loaded restart state with a fresh analytic seed instead of fully replacing it. Internally this maps onto init_file_mode = "add".
restart_scale: multiplicative scale applied to the loaded restart state.
append_on_restart: append continued diagnostic history to the existing *.out.nc file instead of replacing it.
warm_start: carry a converged state from one point of a repeated workload to the next instead of restarting every point cold. It is opt-in (default false) and has no effect on a single run, which has no predecessor. See Warm-started scans below.
resolved_diagnostics: materialize mode-resolved spectra in nonlinear diagnostic history. It defaults to true for publication and restart artifacts. Set it to false when only scalar time traces such as heat and particle fluxes, free energy, field energy, growth rate, and frequency are required. The compact path preserves those scalar channels while avoiding the larger spectral-history arrays.
nsave: checkpoint cadence fallback, in steps, for nonlinear NetCDF bundles when time.nstep_restart is not set.
plots: draw the run’s figures beside its output once it finishes. It defaults to true; set it to false (or pass --no-plots) on batch surfaces that only want the data. See Saved Outputs and Restart Files for the figure set each run kind writes.

For direct restart control outside the [output] helper path, the generic [init] init_file / init_file_scale / init_file_mode keys remain the lower-level mechanism.

Warm-started scans

[output] restart reuses a state written to disk by an earlier run. [output] warm_start is its in-memory sibling for work that repeats inside one invocation: a converged eigenmode at one ky is an excellent initial guess at a neighbouring ky, so the next point can start from the answer instead of rediscovering it from the analytic seed. How much that is worth depends entirely on the solver, and Performance measures both cases before you decide.

gkx scan-runtime-linear therefore walks the requested ky list in monotone order – so that neighbours follow neighbours regardless of the order you asked for – and seeds each point from the previous point’s final state. On the eigensolver paths that state is the converged eigenvector and becomes the Krylov starting vector; on the time-integration paths it is the final distribution function. Results are written back into the order you requested, so the reported arrays are unchanged in shape, ordering, and meaning.

The carried state is rescaled by an exact power of two before reuse. Both the eigensolve and the growth-rate fit are invariant to the overall scale of the initial condition, so this loses nothing, and it stops a state that has grown by exp(gamma T) from overflowing. A state that comes back non-finite or zero is refused and the next point starts cold.

Turn it on with --warm-start, with [scan] warm_start = true, or with [output] warm_start = true; the flag wins over [scan], which wins over [output]. It is opt-in rather than the default because it is not free and not always identical – see Performance for the measurements. GKX also declines to warm start on its own when carrying state would break a contract you asked for:

--workers greater than one, or [parallel] strategy = "batch", because independent workers must stay independent;
solver = "explicit_time", which cannot return a final state;
--batch-ky and [parallel] strategy = "combined_ky", which already integrate every ky in one pass and so have nothing to carry.

Two things to know before switching it on. First, on a fixed horizon a warm-seeded time integration reports a different growth rate than a cold one, because it has removed a startup transient the cold run still contains; the two agree only once the horizon is long enough that both are converged, which is also why the GX-parity decks are left cold. Second, an eigensolver started inside one branch can converge to that branch. On a case whose branch structure you do not already know, run the scan once with --no-warm-start and compare point by point before trusting the warm numbers.

run_runtime_parameter_scan takes the same idea as an opt-in warm_start=True: the previous point’s state is reused only while the parameter moves by no more than warm_start_max_step (relative, default 0.25) from the point before it. That is weaker than continuation=True, which follows one branch and requires every point to return a state.

For the explicit equations attached to these controls, see Operators And Terms.

Explicit reference migration

For a resolved legacy deck with a known execution route:

from gkx.workflows.runtime.toml import migrate_end_damping_reference

case = migrate_end_damping_reference(
    "legacy.toml", route="fixed_linear", dt_step=0.02,
)  # actual reference step, not an assumed input default
case.to_toml("fixed_rate.toml")

Use route="timestep_free" without dt_step for a reference RHS that received no solver timestep. The helper does not run a simulation or overwrite the source. It records [damping_reference]: source path/SHA-256, route, amplitude, input timestep, actual reference step and old scaling flag. The record survives gkx.load / case.to_toml; no source file is needed to reload it. A changed rate with stale provenance is rejected at runtime: clear damping_reference when intentionally defining a new model. Migration is host-side bookkeeping; differentiate the resulting LinearParams.damp_ends_rate, not the file-reading helper.

Let $s=1$ for the old true scaling flag and $s=0$ otherwise:

\[\begin{split}\nu_{ref} = \frac{A}{\Delta t_{input}^{s}} \begin{cases} 1/\Delta t_{step}, & \text{fixed linear route},\\ 1, & \text{timestep-free RHS}. \end{cases}\end{split}\]

This preserves the isolated RHS and RK stage map at the reference step. Subsequent timestep refinement holds $\nu_{ref}$ fixed; it does not preserve the legacy per-step model at every timestep. Adaptive linear sources, conflicting run.dt overrides and already-converted decks are rejected. Resolve any command-line timestep override into the source deck first. This is provenance and algebraic equivalence, not benchmark recertification.