ErrorFields Module

The ErrorFields module quantifies how sensitive a perturbed equilibrium's resonant drive is to the placement of each coil set. It is the foundation of an error-field tolerance assessment: once one plasma solve exists, the question "how much resonant field does a misaligned coil produce?" is linear in the coil's motion, so the module linearizes it once and everything downstream (tolerance Monte Carlo, locking risk, allowable-tolerance scans) becomes linear algebra on a small table.

What is computed

For every named coil set of a coil-forced run, the module evaluates the root-area-weighted control-surface spectrum $\tilde{b}$ of the set as built and its central-difference derivatives with respect to the six rigid-body degrees of freedom: Cartesian shifts $(\Delta x, \Delta y, \Delta z)$ in metres and rotations $(\theta_x, \theta_y, \theta_z)$ about the machine axes in degrees. The spectra are placed on the run's $(m, n)$ ordering and conformed with the control-surface operator, so they are exactly the vectors the resonant coupling matrix and its singular vectors act on (see Dominant resonant-coupling mode).

The overlap of a coil set with singular mode $k$ of the coupling matrix, normalized by the axis toroidal field, is the dimensionless

\[\delta = \frac{V_k^{\mathrm{H}}\,\tilde{b}}{B_{T0}},\]

and because the projection is linear, the derivatives of $\tilde{b}$ project to the derivatives of $\delta$. A rigid shift $(\Delta x, \Delta y)$ therefore moves the overlap by $S_x \Delta x + S_y \Delta y$ with complex $S_x = \partial\delta/\partial\Delta x$ and $S_y = \partial\delta/\partial\Delta y$, which is exact for any coil shape; for an axisymmetric coil $S_y = \pm i S_x$ and the response reduces to a single magnitude with a free phase, the model the OMFIT tolerance tool used.

The stored primitive is the linearization of the spectrum, not a scalar, so the resonant surfaces retained, the singular mode, and the normalization are all post-hoc choices.

Running it

Add an [ErrorFields] section to a deck whose [ForcingTerms] uses forcing_data_format = "coil" and whose [PerturbedEquilibrium] computes the singular coupling:

[ErrorFields]
fd_step_shift_m = 1e-3          # Central-difference step for the rigid shifts [m]
fd_step_tilt_deg = 0.1          # Central-difference step for the rigid tilts [degrees]
rotation_center = "conductor"   # Tilt pivot: each conductor's own centre ("conductor") or the whole set's ("set")
write_outputs_to_HDF5 = true    # Write ErrorFields/CoilSensitivities/ to the output file
verbose = false                 # Log per-coil-set progress and linearity diagnostics
exclude_coils = []              # Coil sets that are not error-field sources (the NTV efc_coils are excluded automatically)

The stage runs after the perturbed equilibrium and writes ErrorFields/CoilSensitivities/: the spectra field_as_designed, shift_sensitivity_per_m, tilt_sensitivity_per_deg, a finite-difference curvature diagnostic per tap, the current pattern the spectra were evaluated at, and a DominantMode/ summary projected onto the run's full-window dominant mode (delta_as_designed, its shift and tilt sensitivities, their direction-averaged in-plane magnitudes, and the in-plane shift and tilt that would cancel delta_as_designed).

examples/DIIID-like_error_field_example/ is a complete case: the DIII-D-like equilibrium with the C-coil as the nominal n = 1 source and the eighteen DIII-D F coils added as single-filament hoops at the centroids of their winding packs (from OpenFUSIONToolkit's TokaMaker DIIID_geom.json). An axisymmetric hoop drives no n = 1 field as built, so each F coil's error-field content is entirely its sensitivity to misalignment; analyze_example.jl ranks the coils by error field per millimetre of shift and per tenth of a degree of tilt, over all rational surfaces and over the edge only.

The tilt pivot matters for multi-filament winding packs: "conductor" rotates each filament about its own arc-length centre, the Fortran coil_read convention inherited by the OMFIT tolerance tool, while "set" rotates the pack rigidly about its common centre, which is what an engineering axis-line tolerance constrains. Single-conductor sets give identical results either way.

Tolerance input

Manufacturing tolerances are engineering data, reviewed and versioned independently of any run, so they live in their own TOML file named by tolerance_file in [ErrorFields] (a path relative to the run directory). The run validates the file, echoes its text into Input/RawInputs/ErrorFields/tolerance_toml_raw for replay, and leaves the device numbers where they belong: outside the repository. The sampling and Monte Carlo stages that consume these tolerances follow in later releases; this release fixes the format.

# Fallbacks for keys a coil block leaves out (every key optional)
[ErrorFields.defaults]
shift_sigma_mm = 0.0             # Gaussian uncertainty added to the sampled shift [mm]
tilt_sigma = 0.0                 # Gaussian uncertainty added to the sampled tilt, in tilt_units
tilt_units = "deg"               # "deg", or "m" for a rim displacement converted through the nominal radius
radial_shape = "hollow"          # Radial sampling density on the disk: flat, uniform_area, hollow, or ring
tolerance_model = "additive"     # "additive" (independent shift and tilt) or "cylinder" (correlated axis line)
cylinder_half_height_m = 0.0     # Cylinder half-height for the cylinder model [m]
current_factor = 1.0             # Current as built relative to the run's; scales the coil's overlap and sensitivities

# One block per coil set that moves on its own
[[ErrorFields.coil]]
name = "PF1U"                    # Coil set name; must match a [[ForcingTerms.coil_set]] of the run
shift_tol_mm = 0.5               # Radius of the in-plane displacement disk the coil centre may lie in [mm]
tilt_tol = 0.019                 # Tilt tolerance, in tilt_units
tolerance_model = "cylinder"     # Axis line confined to a cylinder; shift and tilt drawn together
cylinder_half_height_m = 1.5     # Cylinder half-height; sets the tilt reachable at the given radius [m]

# Coil sets that move together: one shared draw per sample
[[ErrorFields.coherent_group]]
name = "upper_pf_brace"          # Label of the coherently moving group
members = ["PF1U", "PF2U"]       # Coil sets sharing the draw
shift_tol_mm = 2.0               # Coherent shift amplitude shared by the group [mm]
tilt_tol = 0.0                   # Coherent tilt amplitude shared by the group, in tilt_units
phase_group = "braces"           # Groups naming the same phase_group share the random direction
rotation_center_z_m = 0.0        # Height of the pivot of the group's rigid rotation on the machine axis [m]
tilt_lever_arm_m = 3.2           # Lever arm converting a tilt_tol in metres to the group's rotation angle [m] (optional)

[ErrorFields.correctability]
uncorrectable_coils = ["EFCC_U"] # Coil sets the error-field correction cannot reduce
efc_factor = 2.0                 # Divisor applied to correctable contributions under error-field correction

[ErrorFields.other_field]
magnitude = 8.2e-6               # Overlap budget of sources not attributed to any coil
sigma = 0.0                      # Gaussian uncertainty on that budget
radial_shape = "ring"            # Radial sampling density of the budget magnitude

The shift tolerance is one number, the radius of the disk the coil centre may sit in; the direction is sampled. A tilt may be given in degrees or, as legacy tolerance tables do, as a rim displacement in metres, which tilt_tolerance_deg converts through the coil set's arc-length-weighted major radius exactly as the coil loader's tilt_in_meters does. A coherent group's tilt is a rigid rotation of all its members about a pivot on the machine axis at rotation_center_z_m, so a member at height z also shifts laterally by (z − z_pivot)·θ. Groups sharing a phase_group label draw the same random direction with independent amplitudes. A group's tilt given in metres is converted through tilt_lever_arm_m when the group has one (a reference structure's lever arm is its own, not any member coil's radius), and otherwise through the members' common major radius. A group is correctable only if none of its members is listed as uncorrectable. Coil sets of the run without a tolerance block contribute their as-designed overlap only. A coil's current_factor scales its as-designed overlap and its sensitivities together, since the field is linear in current: 0 switches a coil off, 0.5 runs it at half current, a negative value reverses it; this is the knob for a risk against the current in an uncorrectable coil.

Every key is checked against the schema, so a misspelled key is an error rather than a silent default. read_tolerance_toml parses a file into a ToleranceSet, validate_tolerances checks its names against a run's coil sets, and update copies a set, a coil, a group or the unattributed budget with named fields replaced, which is how a programmatic scan over a quantity the file fixes (a sigma, a budget, a current factor) is written:

ts = EF.read_tolerance_snapshot("gpec.h5")
half = EF.update(ts; coils=[c.name == "F6A" ? EF.update(c; current_factor=0.5) : c for c in ts.coils])
tight = EF.update(ts; other_field=EF.update(ts.other_field; magnitude=1e-5))

Sampling a tolerance

A tolerance is one number, the radius of the disk the coil centre may lie in; the direction is random and the radial density is a RadialDistribution: Flat (uniform in radius, the OMFIT tool's flat), UniformArea (uniform over the disk), Hollow (peaked toward the edge, the OMFIT default), Ring (always on the edge), or PowerLaw(p). sample_disk draws a point as Δx + iΔy, sample_uncertainty adds the Gaussian uncertainty on where the coil actually sits, and the two tolerance models combine them: sample_additive draws a shift (metres) and a tilt (degrees) independently, while sample_cylinder confines the coil's axis line to a cylinder of radius R and half-height z_top, drawing its two endpoints and deriving the correlated midplane shift and lean. Tilts are the rotation angles θx + iθy about the machine axes in the sense apply_transforms uses. Every sampler takes the random generator and optional fixed directions, so coherent groups sharing a direction pass one phase and draw their own radii.

rng = Random.Xoshiro(1)
Δ = EF.sample_disk(rng, 0.5e-3, EF.Hollow())          # a point in a 0.5 mm disk, m
Δ, θ = EF.sample_cylinder(rng, 0.5e-3, 1.5, EF.randpow(EF.Flat()))   # correlated shift [m] and tilt [deg]

Tolerance Monte Carlo

With a tolerance_file named, the run samples every coil set's misalignment within its tolerance and histograms the dominant-mode overlap |δ|: per sample and coil set the overlap moves by S·(Δ + u) + T·(θ + v) for the coil's own draw (Δ, θ) and Gaussian placement uncertainties (u, v); coherent groups add one shared draw per group, with the lateral shift a rigid rotation about the group pivot gives each member; the unattributed budget adds a random direction. Because the overlap is linear in the misalignments, a million samples take about a second and no field is recomputed. Two histograms are written to ErrorFields/MonteCarlo/: the intrinsic |δ| and the corrected one, in which every correctable term (coil sets and groups not listed as uncorrectable, and the unattributed budget) is divided by efc_factor. Batches are seeded individually, so results are bit-identical for any thread count, and their spread is the statistical error bar of anything derived from them.

A coil set the run energizes but that is not an error-field source — a correction array named in [ErrorFields.NTV] efc_coils, or anything listed in [ErrorFields] exclude_coils — is swept like the others, so its sensitivities and couplings are tabulated, but it is left out of the as-designed error field, the Monte Carlo, the worst-case bound and the scans, in the run and in every post-hoc entry point that rebuilds the table from the file. A tolerance that names it is an error. without_coils(table, names) is the same operation on a table in memory, and excluded_coil_names("gpec.h5") says which sets a run left out.

[ErrorFields]
tolerance_file = "tolerances.toml"      # Manufacturing-tolerance TOML, relative to the run directory

[ErrorFields.MonteCarlo]
nsample = 1000000               # Samples per batch
nbatch = 10                     # Independent batches; their spread is the statistical error bar
seed = 1                        # Base seed; batch b uses Xoshiro(hash((seed, b)))
nbins = 300                     # Histogram bins, linear on [0, delta_max]
delta_max = 0.0                 # Upper histogram edge; 0 = 1.5 × the worst-case alignment bound
tolerance_scale = 1.0           # Multiplies every shift and tilt tolerance (for tolerance scans)
coil_subset = []                # Coil sets whose tolerances are sampled; empty = all
scale_subset = []               # Coil sets and groups that tolerance_scale applies to; the rest hold their tolerance (empty = all)
scale_map = {}                  # Extra multiplier per coil set or group name, e.g. {F6A = 2.0}

coil_subset and scale_subset answer different questions. coil_subset samples only the named coils and gives every other coil zero tolerance, which isolates one coil's contribution. scale_subset applies tolerance_scale to the named coils and groups and leaves every other name at its own tolerance, which is the product question of a design review: how far can the tolerance of one class of coils be relaxed while the rest of the machine holds its own? Names in either that are not coil sets of the run or groups of the tolerance file are errors. Neither changes which random numbers are drawn, so runs with the same seed differ only through the tolerances they scale.

The run's histogram is the full-window, dominant-mode summary. Any other window or mode, a tolerance scale, or a coil subset is a post-hoc re-run of the same kernel:

mc = EF.run_monte_carlo("gpec.h5"; psi_low=0.5, tolerance_scale=2.0, coil_subset=["PF1U", "PF2U"])
mc.abs_delta_pdf, mc.abs_delta_bin_edges           # intrinsic |δ| density
mc.abs_delta_efc_pdf                     # corrected
mc.abs_delta_sampled_mean, mc.abs_delta_total_as_designed

Locking risk and allowable tolerance

An overlap distribution becomes a locking risk through the empirical ITPA penetration-threshold scalings (citations): the n = 1 fits of Logan et al., Plasma Phys. Control. Fusion 62, 084001 (2020) and of Bursch et al., Plasma Phys. Control. Fusion (2026), doi:10.1088/1361-6587/aea7d6, and the n = 2 fits of Logan et al., Nucl. Fusion 60, 086010 (2020): δ_thresh = 10^α_c · n_e^α_n · B_T^α_B · R_0^α_R · (β_N/l_i)^α_β · I_p^α_I, with n_e in 10¹⁹ m⁻³, B_T in T, R_0 in m and I_p in MA. Each fit is chosen by year, dataset and fit; only the 2026 fits carry the current term:

Fit (`year`, `dataset`, `fit`)	α_c	α_n	α_B	α_R	α_β	α_I
`2026`, `"O,L"`, `"OLS"` (Eq. 7)	−4.31 ± 0.09	0.77 ± 0.08	0.19 ± 0.09	1.88 ± 0.16	0.25 ± 0.08	−0.97 ± 0.08
`2026`, `"O,L"`, `"WLS"` (Eq. 8)	−4.26 ± 0.09	0.56 ± 0.08	0.30 ± 0.10	1.57 ± 0.15	0.13 ± 0.06	−1.01 ± 0.07

The 2026 fits use only ohmic and L-mode discharges of conventional tokamaks (C-Mod, DIII-D, EAST, JET, J-TEXT and KSTAR; no NSTX or COMPASS), so they suit conventional-aspect-ratio designs in linear ohmic confinement or L-mode. Sampling the fitted exponents within their standard errors turns the threshold into a distribution; its cumulative distribution is the probability that an overlap δ locks, and the locking probability of the assembled machine is 100 ∫ pdf(δ) P(lock|δ) dδ over the Monte Carlo bins, per batch. The operating point is an [ErrorFields.scenario] table: density must be given (it is not an equilibrium output); field, major radius, βN, li and the plasma current default from the equilibrium.

[ErrorFields.scenario]
n_e = 5.0                       # Electron density for the threshold scaling [1e19 m^-3]
# i_p = 1.2                     # Plasma current magnitude [MA]; defaults to the equilibrium's

[ErrorFields.Risk]
year = 2020                     # Publication year of the threshold fit: 2020 (n = 1 and n = 2) or 2026 (n = 1)
dataset = "O,L"                 # ITPA dataset of the fit: 2020 "O,L" or "O,L,H" (n = 1), "O,L", "O,L,-C", "O,L,N" (n = 2); 2026 "O,L"
fit = "WLS"                     # Fitting method: "OLS", "DSOLS", or "WLS" (2026: "OLS" or "WLS")
distribution = "normal"         # How the fit exponents are sampled: "normal", "flat", or "normal_truncated"
nsample_threshold = 1000000     # Threshold samples
seed = 1                        # Seed of the threshold sampling
scan_scales = [0.25, 0.5, 1.0, 2.0, 4.0]   # Tolerance multipliers of the allowable-tolerance scan (empty: no scan)

ErrorFields/Risk/ holds the threshold density and P(lock|δ) on the Monte Carlo grid, the locking probability of the intrinsic and corrected distributions (with per-batch values), the as-designed risk, and the sharp-threshold risk; ErrorFields/Risk/ToleranceScan/ the risk against tolerance scale. The scan is the stored quantity; the allowable tolerance for a target risk is a post-hoc inversion, and every window or fit choice is re-evaluated from the file:

scan = EF.ToleranceScan("gpec.h5")
EF.allowable_tolerance(scan, 1.0)                  # tolerance multiplier at 1 % locking risk
EF.allowable_tolerance(scan, 1.0; corrected=true)  # with error-field correction
risk = EF.locking_risk("gpec.h5"; n_e=5.0, psi_low=0.7, risk_ctrl=EF.RiskControl(; dataset="O,L,H"))
risk26 = EF.locking_risk("gpec.h5"; n_e=5.0, risk_ctrl=EF.RiskControl(; year=2026, fit="OLS"))
scan2 = EF.tolerance_scan("gpec.h5"; n_e=5.0, scales=[0.5, 1, 2, 4], coil_subset=["F6A", "F7A"])
scan3 = EF.tolerance_scan("gpec.h5"; n_e=5.0, scales=[0.5, 1, 2, 4], scale_subset=["F6A", "F7A"])  # scan the outboard pair, hold the rest

A threshold scaling is fitted for one toroidal mode number, so the in-run risk and every file-path entry point refuse a run that spans several (single_toroidal_mode); project onto one n's coupling first.

A number quoted at the target risk level should be shown free of sampling and binning bias. risk_convergence repeats the Monte Carlo over a list of sample counts and a list of bin counts, holding the other at its control value, and returns the locking probability with its batch spread for each: the spread must fall as 1/√nsample, and the value must not move with nbins by more than the spread. The spread alone cannot see a binning bias, which is why the bin sweep is there.

conv = EF.risk_convergence(table, ts, coil_sets, sc, scen; nsamples=[10^5, 10^6, 10^7], nbins_list=[100, 300, 3000])
conv.locking_probability_percent_by_nsample, conv.locking_probability_spread_percent_by_nsample
conv.locking_probability_percent_by_nbins

Plots and coil-array phasing

Analysis.ErrorFields plots everything above from gpec.h5, and every function takes a list of label => path pairs so coil-design revisions overplot on one axis:

AEF = GeneralizedPerturbedEquilibrium.Analysis.ErrorFields
AEF.plot_coil_sensitivities(["rev A" => "revA/gpec.h5", "rev B" => "revB/gpec.h5"]; quantity=:shift)
AEF.plot_tolerance_pdf("gpec.h5"; corrected=true)
AEF.plot_locking_risk("gpec.h5"; target_percent=1.0)     # marks the allowable scale
AEF.plot_threshold_scaling("gpec.h5")
AEF.plot_dominant_mode_spectrum("gpec.h5")
AEF.plot_error_field_summary("gpec.h5"; save_path="error_field_summary.png")

Four more views answer the questions the assessment's single numbers hide. plot_overlap_phasors lays each coil's complex overlap head to tail, so a coil that adds to the error field is told apart from one that cancels another's; plot_tolerance_budget stacks the terms of the worst-case bound abs_delta_worst_case (EF.worst_case_terms) per coil, group and unattributed budget, so the tolerance the budget is spent on is visible; plot_linearity_residuals shows the curvature the linear sensitivity model neglects, at the finite-difference step or rescaled to each coil's own tolerance; and quantity=:fraction on the sensitivity bars is the resonant share of each coil's own spectrum, which separates a coil that is small from one that drives the wrong harmonics.

AEF.plot_overlap_phasors("gpec.h5")                                  # one panel per source
AEF.plot_tolerance_budget("gpec.h5"; sort=:total, top=10)            # largest budget terms first
AEF.plot_linearity_residuals("gpec.h5"; at=:tolerance)               # curvature over the tolerance range
AEF.plot_coil_sensitivities("gpec.h5"; quantity=:fraction)           # resonant fraction, in percent
AEF.plot_coil_sensitivities("gpec.h5"; quantity=:as_designed, yscale=:log10)   # stems, never bars, on a log axis
AEF.plot_tolerance_pdf("gpec.h5"; show_batches=true)                 # batch spread and the clamped fraction
EF.worst_case_terms("gpec.h5").total                                 # = ErrorFields/MonteCarlo/abs_delta_worst_case

When several independently powered coil arrays share the job of correcting the error field, the relative phases of their current patterns decide how much dominant-mode field they can drive per ampere-turn. phasing_map evaluates |Σ_k δ_k e^{iφ_k}| per kilo-ampere-turn and the resonant fraction of the applied field on a grid of the N − 1 relative phases from the stored nominal spectra, a closed form with no optimizer, and plot_phasing_map draws it (a line for two arrays, a contour for three):

pmap = EF.phasing_map("gpec.h5", ["EFCC_L", "EFCC_M", "EFCC_U"]; psi_low=0.5)
EF.extreme_phasing(pmap)                        # best |δ| per kAt and the phases giving it
AEF.plot_phasing_map(pmap; quantity=:resonant_fraction_percent)

NTV limits of error-field correction

A correction coil cancels the dominant-mode overlap at C_c per kilo-ampere-turn, but the non-resonant remainder of its field drives a neoclassical toroidal viscosity (NTV) torque T·I² that a perfect correction does not remove. With a torque budget T_0 and the threshold taken to fall in proportion to the torque spent, the current that corrects an intrinsic overlap δ_EF solves δ_EF − C_c I = s δ_thresh (1 − T_residual I²/T_0). This is the model used for SPARC by Logan et al., Nucl. Fusion (2026), doi:10.1088/1741-4326/ae6086, and for ARC by Leuthold et al., J. Plasma Phys. 92, E49 (2026), doi:10.1017/S0022377826101421. It holds only while the plasma still rotates, so the current is capped at I_max = √(T_0/T_residual), and the largest correctable overlap is whichever comes first: the peak of the quadratic, or the overlap cancelled at I_max. max_correctable_overlap returns that limit (with_ntv) together with the two published zero-rotation limits, C_c √(T_0/T_residual) (SPARC, residual_only) and C_c √(T_0/T_full) (ARC, torque_only). [ErrorFields.NTV] names the correction arrays; the run evaluates each one's C_c, resonant fraction, and NTV torque per kAt² for its whole field and for its field with the dominant mode projected out — two plasma-response evaluations of the unit-current spectrum followed by the kinetic torque, which needs a [KineticForces] section — and writes ErrorFields/NTV/. The spectrum is evaluated on the run's [ForcingTerms] boundary grid and normalized by the magnitude of the array's ampere-turns, |nw| × max|I|, so the winding sense and current pattern of the deck stay the phase reference; the torques are stored with their sign, and the limits consume the budget with their magnitude (the sign depends on the rotation and on conventions, so a negative torque is never read as no torque). Torque budget, threshold and safety factor are analysis choices:

[ErrorFields.NTV]
efc_coils = ["d3d_c"]           # Coil set names of the correction arrays to evaluate; never error-field sources
method = "fgar"                 # KineticForces torque method (must be enabled in [KineticForces])

couplings = EF.read_efc_couplings("gpec.h5")
curve = EF.efc_current_curve(couplings[1]; delta_threshold=1.4e-4, torque_budget=4.0)
EF.max_correctable_overlap(couplings[1]; delta_threshold=1.4e-4, torque_budget=4.0)
AEF.plot_efc_ntv_limits("gpec.h5"; torque_budget=4.0)   # threshold from the run's Risk/ group
# The SPARC paper's uncertainty: ±50 % on the NTV torque, T_0 = 4 ± 2 N·m
AEF.plot_efc_ntv_limits("gpec.h5"; torque_budget=4.0, torque_rtol=0.5, budget_rtol=0.5)
AEF.plot_efc_ntv_limits("gpec.h5"; torque_budget=4.0, profiles=true)   # adds T(ψ) integrated from the axis

The torque budget is the plasma's intrinsic torque. The SPARC paper uses 4 N·m for SPARC and the ARC paper 5–20 N·m for ARC; pick a value for the machine being analysed.

Designing the correction

The NTV limits say how much correction current an array may carry; correction_requirement says how much it needs, and on which surfaces the correction falls short. For one array the dominant-mode answer is the ratio −δ_source/δ_array; for several arrays the dominant-mode condition alone does not fix the currents (one complex equation in several unknowns leaves a family of solutions in which any two arrays can cancel each other, of which the result names the member with the least total ampere-turns squared); the question worth a function is the joint one over every rational surface: the complex factors on all arrays that minimize the resonant field C·(b̃_source + Σ_k f_k b̃_k) in the least-squares sense, the field left on each surface after that and after each array's dominant-mode correction alone, and how much of the source's spectrum each array can see. Each current comes twice: as a complex factor on the array's current as given, and in kilo-ampere-turns through the array's ampere_turns_kat (|winding multiplier| × max|I| / 1000, the same normalization as the EFC couplings), which is the form to compare arrays whose patterns were given at different currents and to set beside the NTV-limited current. Any CoilOverlap can be a source or an array, so an as-built assembly (combine_overlaps) can be corrected with arrays that were never part of the run. Over the tolerance samples the needed current is a distribution, and uncorrectable_probability is the fraction of the Monte Carlo's overlap beyond what the NTV-limited array can correct at all: the explicit-correction counterpart of the corrected locking probability's efc_factor model. The Monte Carlo leaves the correction arrays out of its overlap, so its as-designed current and a correction_requirement whose source combines the error-field coils describe the same field.

req = EF.correction_requirement("gpec.h5", "F_coils_as_built", ["d3d_c", "iu", "il"])
req.current_factor_least_squares, req.least_squares_rank      # unique when the rank is the number of arrays
req.current_factor_dominant, req.current_factor_dominant_minimum_norm
req.current_least_squares_kat                                   # the same in kilo-ampere-turns, comparable across arrays
abs.(req.resonant_field_least_squares_t) ./ abs.(req.resonant_field_source_t)   # what is left, per surface
AEF.plot_correction_requirement(req)
mc = EF.MonteCarloResult("gpec.h5"); c = EF.read_efc_couplings("gpec.h5")[1]
EF.uncorrectable_probability(mc, c; delta_threshold=1.4e-4, torque_budget=4.0)
AEF.plot_needed_current(mc, array; coupling=c, delta_threshold=1.4e-4, torque_budget=4.0)

Analysis after the run

Window the coupling to any range of rational surfaces and project onto any singular mode without re-running anything, from memory or from the file:

using GeneralizedPerturbedEquilibrium
EF = GeneralizedPerturbedEquilibrium.ErrorFields

# From the file: the coupling is rebuilt, windowed to 0.5 ≤ ψ_N ≤ 1, and projected onto mode 1
table = EF.sensitivity_table("gpec.h5"; psi_low=0.5)
table.delta_as_designed            # complex overlap of each coil set
table.abs_delta_shift_per_mm       # direction-averaged |∂δ/∂Δ| per millimetre of shift, per coil set
table.abs_delta_rim_per_mm         # the same for tilt, as rim displacement at the coil's major radius
table.tilt_sensitivity_per_deg[1, :]               # ∂δ/∂θx per degree, per coil set

# In memory, from the run's returned state (no I/O)
rc  = PerturbedEquilibrium.ResonantCoupling(run.pe, run.ffs)
dom = PerturbedEquilibrium.dominant_coupling(rc; psi_low=0.5)
table = EF.sensitivity_table(run.coil_sensitivities, dom; mode=1)

How far the linear model holds

Every sweep above assumes the overlap is linear in the misalignments out to the tolerance edge. The run checks that only at the finite-difference step (the stored *_linearity_residual ratios, plot_linearity_residuals). linearity_check displaces each coil set with a tolerance block by ±scale × tolerance along each in-plane shift and tilt axis, rotates each coherent group as one body about its pivot, recomputes the overlap of the displaced geometry and reports the relative error of the linear prediction per row. Include the largest scale a tolerance scan will use, since that is the regime it relies on; the cost is one Biot–Savart pass per row.

lin = EF.linearity_check("gpec.h5"; scales=(1.0, 2.0, 8.0))
maximum(filter(!isnan, lin.relative_error))      # the worst row
AEF.plot_linearity_check(lin)                     # stems per coil, log axis
AEF.plot_linearity_residuals("gpec.h5"; at=:step, yscale=:log10)   # the stored step residuals, readable below 1 %

A coil set whose field at the run's toroidal mode is round-off (an axisymmetric hoop at n = 1) has no resonant fraction: coil_overlaps, combine_overlaps and plot_coil_sensitivities(quantity=:fraction) report NaN for it rather than the ratio of two noise vectors, which would read as a fraction between 1 % and 30 %; the floor is 1e-8 of the largest field norm among the coil sets judged together (resonant_fraction_percent).

Comparing coil revisions

Assessing a new coil design against an existing run costs one Biot-Savart pass per coil set and no plasma solve: the equilibrium, the control surface and the resonant coupling all come from the stored file. Gather them once with ResonantDriveContext and hand it whatever geometry you like.

ctx = EF.ResonantDriveContext("gpec.h5")                    # nzeta_coil=…, dat_dir=… override the deck
old = EF.coil_overlaps(ctx, ForcingTerms.load_coil_sets(old_cfg, 1))
new = EF.coil_overlaps(ctx, ForcingTerms.load_coil_sets(new_cfg, 1))

new[1].delta               # dimensionless overlap δ = Vᴴb̃ / B_T0
new[1].resonant_fraction_percent    # how much of this coil's own spectrum is resonant
new[1].spectrum            # b̃ itself, for the diagnostics below

Every normalization the quantity is quoted in travels on the CoilOverlap, so a caller never has to work out which one a bare number was in.

A revision can rename or split coils. Because the field is linear in the currents, a design current pattern is applied afterwards rather than by re-running the geometry, and a split is written as the combination it is:

sets = ForcingTerms.conductors(ForcingTerms.read_coil_dat("old_pair.dat"))   # one file, two coils
ovs  = EF.coil_overlaps(ctx, sets)                                          # each at 1 kA, say

EF.combine_overlaps(ovs, "old_pair_1" => 20.0, "old_pair_2" => -15.0)       # weights in kA

Weights multiply each spectrum as it was evaluated; they are not absolute currents, since the sets may have been swept at any current. Evaluating every conductor at 1 kA and weighting by the design current in kA is the idiom that makes them read as currents. An unmatched name raises rather than contributing zero, so a coil renamed between revisions cannot quietly drop out of a comparison.

For sensitivities to rigid motion as well, compute_coil_sensitivities takes the same context. It costs thirteen spectrum evaluations per coil set instead of one, so reach for coil_overlaps when only the overlaps are wanted.

Result names

Every result struct in ErrorFields (CoilSensitivities, SensitivityTable, CoilOverlap, MonteCarloResult, RiskResult, ToleranceScan, PhasingMap, EFCCoupling) names its fields by one grammar, and every field name is also its HDF5 dataset name under ErrorFields/:

<quantity>[_<instance>][_efc][_<statistic>][_<unit>]

quantity is the physical thing, spelled out: delta is the dimensionless dominant-mode overlap v·b̃ / B_T0, a complex phasor per coil set, and abs_delta its magnitude; threshold is the penetration threshold in the same units; locking_probability is the probability of locking (never abbreviated, and never called a risk in a name); resonant_fraction is 100·|v·b̃| / ‖b̃‖; field, spectrum, torque, current, shift_sensitivity and tilt_sensitivity are what they say.
instance says which geometry or value the quantity is evaluated for: _as_designed, the unperturbed design; _sampled, over the Monte Carlo misalignment samples; _worst_case, every term at its tolerance edge and in phase; _fit, a scaling-law value at the fitted exponents without the fit's scatter; _total, the coherent sum over coil sets (a per-coil vector carries no qualifier and documents its shape as [ncoil_set]). A result evaluated on whatever geometry the caller passed (CoilOverlap, EFCCoupling) carries no instance qualifier. The word nominal is not used in names.
_efc marks a value with error-field correction applied, as in the TOML keys efc_factor and efc_coils.
statistic: _mean, _pdf, _batches, _spread (range over the Monte Carlo batches), _samples, _bin_edges.
unit is present whenever the quantity has one: _percent, _t (tesla), _m, _mm, _deg, _per_m, _per_deg, _per_mm, _per_kat, _per_kat2. delta and probabilities given as fractions carry none; a probability stored as a percentage always carries _percent.

The order is fixed, so locking_probability_efc_batches_percent reads as "the locking probability, corrected, per batch, in percent", and abs_delta_total_as_designed as "the magnitude of the total overlap of the design as built".

API Reference

GeneralizedPerturbedEquilibrium.ErrorFields.MIN_NZETA_PER_PERIOD — Constant

Minimum toroidal boundary points per field period below which a post-hoc evaluation warns.

Helical windings alias on a coarse toroidal grid. A scan on multi-turn coils against a 128-point reference put 16 points per period about 1 % off, 20 about 0.4 %, 24 about 0.1 %, and 32 within one part in 1e5. A few tenths of a percent is the same size as the differences a coil-revision comparison is trying to resolve, so a grid inherited silently from a stored deck and coarser than ForcingTerms.NZETA_POINTS_PER_PERIOD is worth a warning; the default itself is converged.