API Reference¶
Top-Level Exports¶
import scpn_control
scpn_control.__version__ # "0.20.7"
scpn_control.FusionKernel # Grad-Shafranov equilibrium solver
scpn_control.RUST_BACKEND # True if Rust acceleration available
scpn_control.TokamakConfig # Preset tokamak geometries
scpn_control.StochasticPetriNet
scpn_control.FusionCompiler
scpn_control.CompiledNet
scpn_control.NeuroSymbolicController
scpn_control.kuramoto_sakaguchi_step
scpn_control.order_parameter
scpn_control.KnmSpec
scpn_control.build_knm_paper27
scpn_control.UPDESystem
scpn_control.LyapunovGuard
scpn_control.RealtimeMonitor
scpn_control.PhysicsDebugAssistant
SCPN — AER Control Observation¶
scpn_control.scpn.observation adapts asynchronous Address-Event
Representation spike streams into bounded controller features without changing
the existing ControlObservation consumers. The Python surface provides
SpikeEvent, SpikeBuffer, rate/temporal/ISI decoders, and
AERControlObservation.to_features(). The matching Rust implementation lives in
control_core::spike_buffer, with optional PyO3 bindings exposed as
scpn_control_rs.PySpikeBuffer and scpn_control_rs.aer_decode_*.
SpikeBuffer.admission_report(), monotonic_input, and
out_of_order_event_count provide bounded ingress evidence. Observations may
set require_monotonic=True to fail closed before decoding if upstream AER
timestamps violate the monotonic admission contract.
This is an ingress adapter and feature-decoding contract. It does not claim hardware AER signal integrity, target neuromorphic-device deployment, FPGA timing closure, or PCS admission without separate hardware evidence.
dataclass
¶
Single Address-Event Representation spike event.
Bounded FIFO ring buffer for AER spike events.
Overflow drops the oldest event and latches overflowed until clear
is called. This keeps memory bounded under bursty neuromorphic input.
dataclass
¶
AERControlObservation(timestamp_ns, spike_stream, decode_window_ns, decode_strategy='rate', n_features=64, feature_normalisation='unit', require_monotonic=False)
Decode AER events as per-feature event fractions over the active window.
Decode by first-spike timing; earlier spikes map closer to one.
Decode by mean inter-spike interval; faster repeated spikes map higher.
SCPN — Geometry-Neutral Replay Schema¶
scpn_control.scpn.geometry_neutral_replay publishes deterministic
geometry-neutral replay reports and schema admission helpers. The v1.1 schema
extends the v1 report with optional pulsed-shot context: UUID pulse IDs,
capacitor initial energy, trigger timestamp, recovered energy, sorted
shot-phase logs, FRC diagnostic scalars, and digest-bound AER admission
metadata. Existing v1 reports remain loadable under the v1.1 schema bundle.
Use the dedicated Geometry-Neutral Replay guide for field semantics, example payloads, and claim boundaries.
Generate a deterministic compact SCPN-control replay report.
Validate the compact replay report contract without external packages.
Load a bundled geometry-neutral replay JSON Schema document.
Return the bundled v1.1 replay schema after self-identification checks.
Validate a v1 replay report while the v1.1 schema bundle is installed.
build_aer_admission_metadata(*, admission_report, decode_strategy, decode_window_ns, n_features, feature_normalisation='unit', require_monotonic=False, feature_vector=None)
Build replay-safe AER admission metadata from a decoded observation.
Attach digest-bound AER admission metadata to a replay v1.1 report.
Persist a validated geometry-neutral replay report as stable JSON.
Load and validate a geometry-neutral replay report with duplicate-key checks.
Control — Pulsed Scenario Scheduler v2¶
scpn_control.control.pulsed_scenario_scheduler_v2 owns the reusable
pulsed-fusion lifecycle contract that MIF-CORE incubated as MIF-004. The
scheduler models the adjacent state ring:
idle -> ramp_up -> flat_top -> burn -> expansion -> dump -> recharge -> cool_down -> idle
The Python surface provides the control-plane API and audit-log contract. The
matching Rust kernel lives in control_control::pulsed_scenario for the
compiled hot-path lane. When the optional extension is built,
scpn_control_rs.PyPulsedScenarioScheduler exposes that Rust kernel directly
to Python without changing the pure-Python API. All surfaces use the same state
names, action names, guard thresholds, monotone timestamp checks, and
transition reasons.
scpn_control.control.pulsed_scenario_scheduler is retained as the
SCPN-MIF-CORE compatibility import and re-exports the v2 implementation. The
finite-state topology is also captured in lean/SCPNControl/PulsedFSM.lean,
including the liveness theorem
pulsed_fsm_eventually_returns_to_idle.
The scheduler is a generic pulsed-reactor controller primitive. It is not a facility-validated PCS implementation by itself. Hardware timing, actuator mapping, capacitor-bank plant dynamics, and measured-shot validation remain separate admission surfaces.
Compatibility surface for the CONTROL-owned pulsed scenario scheduler.
SCPN-MIF-CORE originally contracted the module path
scpn_control.control.pulsed_scenario_scheduler. The production
implementation lives in pulsed_scenario_scheduler_v2; this module keeps the
cross-repository import stable without duplicating scheduler logic.
dataclass
¶
Capacitor-bank telemetry consumed by lifecycle transition guards.
property
¶
Return bank voltage as a fraction of the declared maximum.
dataclass
¶
PulsedPlasmaTelemetry(coil_current_A, temperature_eV, phase_lock_error_rad, reference_error_m, fusion_power_W, radial_velocity_m_s)
Plasma telemetry consumed by lifecycle transition guards.
Bases: StrEnum
Command emitted for the active lifecycle state.
dataclass
¶
Command emitted by one pulsed-scenario scheduler step.
CONTROL-owned reusable eight-state scheduler for pulsed-fusion shots.
Perform a validated manual adjacent transition.
dataclass
¶
PulsedScenarioSpec(min_precharge_energy_J, ramp_current_A, phase_tolerance_rad, spatial_tolerance_m, burn_temperature_eV, min_fusion_power_W, expansion_velocity_m_s, dump_energy_floor_J, recharge_voltage_fraction, cooldown_temperature_eV, cooldown_current_A, min_burn_duration_s=0.0)
Guard thresholds for a reusable pulsed-fusion shot lifecycle.
Bases: StrEnum
Canonical pulsed-fusion lifecycle states.
Bases: StrEnum
Canonical pulsed-fusion lifecycle states.
Bases: StrEnum
Command emitted for the active lifecycle state.
dataclass
¶
PulsedScenarioSpec(min_precharge_energy_J, ramp_current_A, phase_tolerance_rad, spatial_tolerance_m, burn_temperature_eV, min_fusion_power_W, expansion_velocity_m_s, dump_energy_floor_J, recharge_voltage_fraction, cooldown_temperature_eV, cooldown_current_A, min_burn_duration_s=0.0)
Guard thresholds for a reusable pulsed-fusion shot lifecycle.
dataclass
¶
PulsedPlasmaTelemetry(coil_current_A, temperature_eV, phase_lock_error_rad, reference_error_m, fusion_power_W, radial_velocity_m_s)
Plasma telemetry consumed by lifecycle transition guards.
dataclass
¶
Capacitor-bank telemetry consumed by lifecycle transition guards.
property
¶
Return bank voltage as a fraction of the declared maximum.
dataclass
¶
Single lifecycle transition audit entry.
dataclass
¶
Command emitted by one pulsed-scenario scheduler step.
CONTROL-owned reusable eight-state scheduler for pulsed-fusion shots.
Perform a validated manual adjacent transition.
Control — Capacitor Bank State Model¶
scpn_control.control.capacitor_bank_state owns the CONTROL-side bounded
series-RLC capacitor-bank contract for pulsed-shot admission. It mirrors the
MIF-CORE MIF-005 capacitor-bank mathematics at the control boundary: damping
regime classification, closed-form free response, Crank-Nicolson stepping,
midpoint-sampled discharge waveforms, conservative feasibility guards, and
constant-power recharge projection.
The state equation is:
where C is capacitance in farads, L is inductance in henries, R is series
resistance in ohms, v_C is capacitor voltage in volts, i is series current
in amperes, and i_load is the prescribed external load current. The numerical
step uses a Crank-Nicolson update so natural-response validation can compare
against the analytical underdamped, critical, and overdamped solutions.
CapacitorBank.telemetry() adapts the state model to
PulsedScenarioScheduler v2 by emitting scheduler-compatible capacitor-bank
telemetry with absolute voltage magnitude, declared voltage limit, and stored
energy. The matching Rust kernel lives in
control_control::capacitor_bank. When the optional extension is built,
scpn_control_rs.PyCapacitorBankModel, scpn_control_rs.PyCapacitorBankSpec,
and scpn_control_rs.capacitor_bank_free_response() expose the compiled
surface directly to Python.
scpn_control.control.capacitor_bank is retained as the SCPN-MIF-CORE
compatibility import and re-exports the state-model implementation without
duplicating the RLC mathematics.
CapacitorBank.discharge() now reports an explicit total RLC energy ledger for
admission and replay: initial total stored energy, remaining total stored
energy, remaining capacitor electric energy, remaining inductor magnetic energy,
integrated ohmic loss, integrated prescribed-load extraction, absolute residual,
relative residual, and a boolean pass/fail flag. The residual contract is:
The ledger uses midpoint quantities from the same Crank-Nicolson step that
advances the state, so the residual is a numerical consistency check on the
CONTROL admission model rather than a facility hardware protection claim.
CapacitorBankState.energy_J remains scheduler-facing capacitor electric
energy only.
This model is a bounded control admission and scheduling primitive. It is not a validated facility capacitor-bank driver, insulation model, switch model, or hardware interlock implementation. Facility deployment still requires target hardware timing, protection relays, isolation evidence, and shot-matched validation artefacts.
Compatibility surface for the CONTROL-owned capacitor-bank model.
SCPN-MIF-CORE originally contracted the module path
scpn_control.control.capacitor_bank. The production implementation lives in
capacitor_bank_state; this module keeps the cross-repository import stable
without duplicating the RLC mathematics.
Mutable bounded series-RLC bank with analytical and numerical stepping.
Return scheduler-compatible bank telemetry using absolute voltage magnitude.
Reset bank state at t = 0 with bounded initial voltage.
Advance the bank one Crank-Nicolson step with optional load current.
Drive the bank with pulse using midpoint-sampled load current.
Run conservative pulse admissibility guards against current bank state.
Project constant-power recharge state after non-negative time t.
dataclass
¶
CapacitorBankSpec(capacitance_F, inductance_H, series_resistance_ohm, voltage_max_V, recharge_power_kW=0.0)
Physical parameters for a bounded series RLC capacitor bank.
property
¶
Return :math:\omega_0 = 1/\sqrt{LC}.
dataclass
¶
Instantaneous bank state in SI units.
dataclass
¶
EnergyReport(energy_delivered_J, energy_initial_J, energy_remaining_J, capacitor_energy_remaining_J, inductor_energy_remaining_J, resistive_loss_J, load_energy_J, energy_balance_residual_J, energy_balance_relative_error, energy_balance_passed, peak_voltage_V, peak_current_A, discharge_duration_s, rlc_regime)
Energy bookkeeping returned by a discharge simulation.
dataclass
¶
Prescribed load-current waveform for bounded discharge simulations.
Bases: StrEnum
Canonical damping regimes for the series RLC bank model.
Evaluate the closed-form homogeneous series-RLC response at time t.
Bases: StrEnum
Canonical damping regimes for the series RLC bank model.
dataclass
¶
CapacitorBankSpec(capacitance_F, inductance_H, series_resistance_ohm, voltage_max_V, recharge_power_kW=0.0)
Physical parameters for a bounded series RLC capacitor bank.
property
¶
Return :math:\omega_0 = 1/\sqrt{LC}.
dataclass
¶
Instantaneous bank state in SI units.
dataclass
¶
Prescribed load-current waveform for bounded discharge simulations.
dataclass
¶
EnergyReport(energy_delivered_J, energy_initial_J, energy_remaining_J, capacitor_energy_remaining_J, inductor_energy_remaining_J, resistive_loss_J, load_energy_J, energy_balance_residual_J, energy_balance_relative_error, energy_balance_passed, peak_voltage_V, peak_current_A, discharge_duration_s, rlc_regime)
Energy bookkeeping returned by a discharge simulation.
Evaluate the closed-form homogeneous series-RLC response at time t.
Mutable bounded series-RLC bank with analytical and numerical stepping.
Return scheduler-compatible bank telemetry using absolute voltage magnitude.
Reset bank state at t = 0 with bounded initial voltage.
Advance the bank one Crank-Nicolson step with optional load current.
Drive the bank with pulse using midpoint-sampled load current.
Run conservative pulse admissibility guards against current bank state.
Project constant-power recharge state after non-negative time t.
Control — Pulsed-Shot MPC Admission Adapter¶
scpn_control.control.fusion_sota_mpc.PulsedShotMPCAdapter wraps the
CONTROL-owned gradient MPC surface and admits its first action through the
pulsed-scenario scheduler and capacitor-bank feasibility guard. It is a
control-boundary adapter over existing MPC output, not a new equilibrium solver
or a duplicate solver lane.
The adapter applies three deterministic checks:
- Non-
burnscheduler states replace selected burn-action components with the configured safe action. burnstate evaluates aPulseSpecagainstCapacitorBank.feasibility()unless the caller disables that policy.- Every step records an explainable decision dictionary with scheduler state, capacitor feasibility text, constraint slack, MPC objective, whether a safe action was applied, and a digest-bound admission evidence payload.
The matching Rust kernel lives in control_control::mpc::MPController as
plan_pulsed(). When the optional PyO3 extension is rebuilt,
scpn_control_rs.PyMpcController.plan_pulsed() exposes the same admission
fields to Python, including evidence_schema_version, action_sha256,
safe_action_sha256, burn_action_mask_sha256, peak_current_A, and
admission_digest.
Use the dedicated Pulsed MPC Adapter guide for examples, runtime boundaries, and benchmark evidence commands.
dataclass
¶
PulsedShotMPCDecision(action, mpc_objective, constraint_slack, scheduler_state, bank_feasibility, reason, bank_feasible, safe_action_applied, burn_components_masked, peak_current_A, evidence_schema_version, action_sha256, safe_action_sha256, burn_action_mask_sha256, admission_digest)
Admitted pulsed-shot MPC action and its control-boundary rationale.
PulsedShotMPCAdapter(nmpc, scheduler, bank, *, burn_action_mask=None, safe_action=None, pulse_duration_s=0.001, pulse_waveform='half_sine', refuse_burn_when_uncharged=True)
Gate an existing MPC action through pulsed-shot lifecycle and bank guards.
Control — Multi-Shot Campaign Orchestrator¶
scpn_control.control.multi_shot_campaign runs repeated pulsed-shot admission
over the scheduler, capacitor-bank telemetry, and replay v1.1 metadata
contracts. It accepts explicit shot telemetry, records command and transition
logs, requires the canonical pulsed lifecycle by default, and emits a
digest-bound campaign report. When supplied, per-shot
pulsed_mpc_admission_digest values are validated as lowercase SHA-256
digests and preserved in the campaign report plus replay v1.1 extension fields.
The matching Rust kernel lives in control_control::multi_shot_campaign. When
the optional PyO3 extension is rebuilt,
scpn_control_rs.PyMultiShotCampaignOrchestrator.run_table() exposes the Rust
surface through table-shaped NumPy inputs, including optional per-shot
pulsed_mpc_admission_digests.
Use the dedicated Multi-Shot Campaign guide for examples, claim boundaries, and benchmark commands.
Release admission is handled by
validation.validate_multi_shot_campaign_evidence. The gate admits the
Python/PyO3 benchmark report and the Rust benchmark report only when both carry
the pulsed-MPC digest chain, benchmark context, stable SHA-256 payload seals,
and a complete Python, PyO3, and Rust surface set. The top-level
scpn-control validate command runs this gate by default before a release
evidence JSON can be admitted.
dataclass
¶
One timestamped plasma and bank telemetry sample for a pulsed shot.
dataclass
¶
CampaignShotPlan(shot_id, samples, initial_bank_voltage_V, initial_bank_current_A=0.0, pulsed_mpc_admission_digest=None)
Input contract for one shot in a multi-shot campaign.
dataclass
¶
dataclass
¶
MultiShotCampaignOrchestrator(campaign_id, scheduler_spec, bank_spec, *, require_complete_lifecycle=True)
Run deterministic multi-shot admission over scheduler and bank contracts.
Control — PREEMPT_RT Runtime Admission¶
scpn_control.core.runtime_admission emits a schema-versioned admission report
before native hardware-campaign execution. It binds the observed Linux kernel,
PREEMPT_RT evidence, current process affinity, requested execution cores,
scheduler policy, CPU governors, heartbeat configuration, and memory-lock limits
to the campaign summary. run-hardware-campaign --runtime-admission-policy
require fails closed before native execution if the host is not
production-admissible.
The optional PyO3 counterpart scpn_control_rs.runtime_admission_snapshot()
exposes native-side kernel/core parsing when the Rust extension is installed.
validation.validate_runtime_admission_evidence admits persisted
runtime-admission benchmark reports into the top-level release gate only when
the report carries command, CPU-affinity, host-load, isolation, latency,
SHA-256 payload, and production-claim-boundary evidence.
Use the dedicated PREEMPT_RT Runtime Admission guide for policy semantics and operator examples.
dataclass
¶
RuntimeAdmissionRequest(execution_backend='auto', pacing_mode='sleep', transport_backend='std', formal_mode='async_drop', tick_interval_s=0.001, core_snn=1, core_z3=2, core_net=3, core_hb=4, heartbeat_port=0, native_backend_available=False, require_preempt_rt=False, require_realtime_scheduler=False, require_performance_governor=False, require_heartbeat=False, min_memlock_bytes=DEFAULT_MIN_MEMLOCK_BYTES)
Operator-requested runtime contract for a native control campaign.
dataclass
¶
RuntimeAdmissionProbe(generated_ns, platform_system, kernel_release, is_linux, preempt_rt, realtime_sysfs, affinity, available_parallelism, scheduler_policy, scheduler_priority, governors, memlock_soft_bytes, memlock_hard_bytes, native_snapshot=None)
Observed host runtime properties used for admission.
Collect and evaluate current host runtime admission evidence.
Evaluate a runtime admission request against an observed probe.
Validation — Benchmark Regression Gates¶
validation/validate_benchmark_regression_gates.py admits persisted benchmark
evidence before the preflight gate accepts a regression baseline. The gate does
not run benchmarks and does not create timing evidence. It validates
validation/reports/benchmark_regression_gates.json against the referenced
latency reports, SHA-256 digests, metric paths, bounded thresholds, sample
counts, hardware-context metadata, and explicit non-HIL claim boundaries.
The gate is intended to catch stale, tampered, overclaimed, or regressed local benchmark evidence before release preflight continues. Hardware-in-the-loop, target-device, cloud-GPU, or plant real-time claims remain blocked until those specific benchmark artefacts are generated and admitted separately.
Core — Native C++ Solver Bridge¶
scpn_control.core.hpc_bridge exposes the optional native Grad-Shafranov
solver bridge. Native compilation is disabled unless
SCPN_ALLOW_NATIVE_BUILD=1 is set. When enabled, the bridge admits only the
package-local solver.cpp whose SHA-256 digest matches
solver_manifest.json, resolves the compiler to an absolute regular
executable, strips dynamic-loader injection variables from the build
environment, rejects symlinked solver inputs and output targets, compiles to a
temporary package-local file, and publishes the shared library atomically after
the compiler produced a regular file.
External runtime solver libraries remain blocked unless
SCPN_ALLOW_EXTERNAL_SOLVER_LIB=1 is set for a vetted absolute path. The
default path searches package-local solver locations only.
Physics Debug Assistance¶
scpn_control.physics_debug provides a local-first advisory assistant for
physics validation gaps. The default provider policy admits loopback endpoints
only; facility or external gateways must be explicitly allowlisted.
build_local_provider() supplies loopback profiles for common onsite gateway
protocols: chat-completions-compatible, Ollama-style chat, direct JSON, and
text-generation endpoints. Reports are schema-versioned advisory evidence with
secret redaction, falsifiable hypothesis checks, campaign risk controls, and
risk-bound prompt-injection neutralization for untrusted evidence text before
provider prompting. Prompt-guard findings are recorded in the tamper-evident
payload digest. build_guardrail_provider() adds an optional hallucination
guardrail gateway with a director-ai default profile and explicit alternate
profiles for lab-owned guardrail solutions. Guardrail block decisions fail
closed before report persistence; allow findings are bound into the report
digest together with the SHA-256 digest of the reviewed provider draft.
The guardrail request also binds the provider metadata, safety policy, and
guardrail policy digests so reviews cannot be replayed across another provider
or relaxed policy. High-severity guardrail findings must use block actions, and
risk controls must meet the configured guardrail policy before persistence.
They are not validated physics truth, controller-parameter promotion, or
facility safety approval.
run_provider_quorum() runs multiple providers in local-first order and admits
only hypotheses corroborated by the required provider count over the same gap
and evidence set.
PhysicsDebugSafetyPolicy binds the human-review requirement, caps advisory
confidence, and rejects provider text that attempts controller promotion,
actuation, review bypass, or approval claims.
dataclass
¶
ProviderPolicy(allowed_endpoint_prefixes=('http://127.0.0.1', 'http://localhost', 'https://127.0.0.1', 'https://localhost'), allow_remote_providers=False, max_prompt_chars=20000, max_response_chars=50000)
Security policy for model-provider access.
dataclass
¶
PhysicsDebugGuardrailPolicy(required=False, block_actions=('block',), allowed_decisions=('allow', 'block'), blocking_severities=('high', 'critical'), minimum_risk_controls=1)
Admission policy for optional physics-debug hallucination guardrails.
dataclass
¶
Evidence item supplied to the physics debugging assistant.
dataclass
¶
Physics validation gap to analyse.
dataclass
¶
PhysicsDebugSafetyPolicy(human_review_required=True, max_advisory_confidence=0.95, forbidden_action_phrases=DEFAULT_FORBIDDEN_ACTION_PHRASES)
Fail-closed admission policy for advisory physics-debug output.
dataclass
¶
HTTPChatProvider(provider_name, endpoint, model, protocol='chat-completions', headers=dict(), timeout_seconds=30.0, transport=None)
Allowlisted HTTP chat gateway for local or facility-approved providers.
dataclass
¶
PhysicsDebugGuardrailProvider(provider_name, endpoint, model, profile='director-ai', protocol='direct-json', headers=dict(), timeout_seconds=30.0, transport=None)
Allowlisted guardrail gateway for advisory physics-debug reports.
dataclass
¶
PhysicsDebugAssistant(policy=ProviderPolicy(), safety_policy=PhysicsDebugSafetyPolicy(), guardrail_policy=PhysicsDebugGuardrailPolicy())
Evidence-first assistant for physics gap analysis.
Ask an allowlisted provider for advisory hypotheses and campaigns.
build_local_provider(*, family, model, provider_name, host='127.0.0.1', port=None, path=None, transport=None, headers=None, timeout_seconds=30.0)
Build a loopback-only provider profile for onsite model gateways.
build_guardrail_provider(*, model, provider_name, profile='director-ai', host='127.0.0.1', port=None, path=None, transport=None, headers=None, timeout_seconds=30.0)
Build a loopback-only guardrail profile for physics debug reviews.
build_physics_debug_report(*, provider, evidence, gaps, provider_output, redactions=None, safety_policy=None, prompt_guard_findings=None, guardrail_policy=None, guardrail=None)
Build a schema-versioned advisory physics debugging report.
run_provider_quorum(*, evidence, gaps, providers, policy=None, safety_policy=None, guardrail_policy=None, guardrail_provider=None, min_providers=2, min_local_providers=1)
Run multiple providers and admit only corroborated advisory evidence.
Validate a tamper-evident advisory physics debugging report.
Validate a tamper-evident provider-quorum physics debug report.
Persist a validated advisory physics debugging report as JSON.
Persist a validated provider-quorum physics debug report as JSON.
Control — Federated Disruption Prediction¶
scpn_control.control.federated_disruption supports FedAvg and FedProx
training across named tokamak clients without centralising facility arrays.
create_facility_clients_from_arrays() is the production ingestion boundary
for per-facility X_train, y_train, X_test, and y_test arrays. It
enforces the shared 8-feature disruption contract and binary labels before a
client joins the federation.
DifferentialPrivacyConfig enables facility-update clipping, Gaussian noise,
and a serialisable privacy ledger. The shipped benchmark
validation/benchmark_federated_disruption.py publishes deterministic
synthetic multi-facility evidence in
validation/reports/federated_disruption_benchmark.json and
validation/reports/federated_disruption_benchmark.md. Those artefacts test
federation, heterogeneity, and differential privacy contracts; they do not
claim measured cross-facility validation without external shot databases.
dataclass
¶
Facility-level differential privacy for federated client updates.
dataclass
¶
PrivacyLedgerEntry(round_index, participating_clients, epsilon_spent, cumulative_epsilon, delta, max_update_norm, noise_multiplier, clipped_clients)
Per-round facility-level privacy accounting record.
dataclass
¶
FacilityBenchmarkSummary(aggregation, machines, n_rounds, mean_accuracy, mean_loss, per_machine_accuracy, privacy_epsilon, privacy_delta, evidence_kind)
Deterministic benchmark summary for a federated disruption campaign.
dataclass
¶
FederatedConfig(n_rounds=10, local_epochs=5, learning_rate=0.01, aggregation='fedavg', mu_proximal=0.01, min_clients=2, machines=(lambda: ['DIII-D', 'JET', 'KSTAR'])(), dp_config=None)
Configuration for federated disruption prediction training.
Local training client for a single tokamak.
SGD on local data, starting from global_weights.
When mu_proximal > 0, adds the FedProx penalty (mu/2)||w - w_global||^2 to the loss gradient.
Orchestrates federated training across machine clients.
FedAvg: weighted average of model weights by dataset size.
McMahan et al., "Communication-Efficient Learning of Deep Networks from Decentralized Data", AISTATS 2017.
FedProx aggregation with proximal regularisation.
The proximal term is applied during local training (not aggregation), so aggregation itself is weighted averaging — the difference from FedAvg is in the client-side gradient update. This method exists for API symmetry; the mu parameter documents the proximal weight used.
Full federated training loop.
Returns per-round metrics including per-client accuracy, loss, n_samples.
Return cumulative facility-level DP spend for the current server.
Create clients from per-facility arrays without centralising raw shots.
Each facility payload must contain X_train, y_train, X_test, and
y_test. The constructor enforces the shared 8-feature disruption
contract and binary label boundary before the data can enter a federation.
run_synthetic_multifacility_benchmark(*, machines=('DIII-D', 'JET', 'KSTAR', 'EAST'), n_rounds=4, local_epochs=3, aggregation='fedprox', dp_config=None, seed=20240531)
Run a deterministic synthetic multi-facility disruption benchmark.
This benchmark exercises the production federation, heterogeneity, and privacy-accounting contracts. It is not measured cross-facility validation.
Control — Quantum Disruption Bridge¶
scpn_control.control.quantum_disruption_bridge is a fail-closed facade for
optional quantum-enhanced disruption prediction. Quantum circuit and provider
ownership stays in scpn-quantum-control; SCPN-CONTROL owns the control
feature contract, lazy optional import boundary, bounded claim metadata, and
tamper-evident advisory reports. The bridge maps the CONTROL 8-feature
disruption vector to the ITER 11-feature contract only when missing ITER fields
are either supplied explicitly or declared as bounded centre defaults. Reports
are not facility validation, controller promotion, or publication-safe evidence
without external disruption databases and benchmark artefacts.
quantum_disruption_kernel_matrix() emits a bounded amplitude-encoding kernel
report with symmetry, diagonal, and [0, 1] admission checks. The callable
quantum owner path uses scpn_quantum_control.control.q_disruption_iter
lazily; when that optional dependency is unavailable the report fails closed
with status="quantum-unavailable" and no quantum score. Every bridge report
also records advisory admission evidence: CONTROL feature digests, ITER mapping
digests, default-use reasons, and the external evidence still required before
facility or publication claims are admissible. Bridge and kernel reports carry
schema-versioned advisory certificates that bind report kind, repository
ownership, claim boundary, downstream non-admission policy, and the content
digest before the outer payload digest is accepted. The facade also publishes a
machine-readable dependency contract that names the scpn-quantum-control
backend module, required classifier surface, Qiskit core dependencies, optional
provider dependency families, report schemas, feature contract, and
non-admission policy for future backend hardening. Generated bridge and kernel
reports embed that dependency contract and bind its digest into the advisory
certificate so archived reports cannot be replayed against a different quantum
backend contract. When the optional backend exposes its own
scpn_control_bridge_dependency_contract() callable, CONTROL compares it
against the expected contract, records backend-contract attestation evidence,
and fails closed before report creation if the backend advertises a conflicting
contract. Bridge reports also include schema-versioned advisory decision
evidence that records whether the score came from the quantum backend or the
classical fallback, applies deterministic low/elevated/high risk-band
thresholds, records backend-contract validation state, fixes the control action
to blocked, and binds the decision digest into the advisory certificate.
dataclass
¶
QuantumDisruptionBridgeConfig(allow_center_defaults=False, require_quantum_backend=False, seed=20240531, backend_profile='statevector', quantum_module=QUANTUM_MODULE, claim_status='bounded_model', source_mode='control-facade')
Configuration for the optional quantum disruption bridge.
dataclass
¶
QuantumFeatureMapping(raw_iter_features, normalized_iter_features, control_feature_names, iter_feature_names, defaults_used=tuple(), unmapped_control_features=('dBp_dt', 'n1_rms'), claim_status='bounded_model', publication_safe=False)
CONTROL-to-ITER feature mapping with provenance.
Map the CONTROL 8-feature disruption contract to the ITER 11-feature contract.
Return a bounded amplitude-encoding kernel report for CONTROL feature samples.
Return the CONTROL-to-QUANTUM disruption bridge dependency contract.
Run the optional quantum disruption bridge and return an advisory report.
Validate the CONTROL-to-QUANTUM disruption bridge dependency contract.
Validate a tamper-evident quantum disruption bridge report.
Validate a tamper-evident quantum disruption kernel report.
Core — Physics Solvers¶
FusionKernel¶
FusionKernel validates its JSON configuration before grid construction:
the root must be an object, duplicate JSON keys are rejected, dimensions and
grid resolution must be physical, physics.plasma_current_target must be
positive finite, and physics.vacuum_permeability must be positive finite
when supplied.
Non-linear Grad-Shafranov equilibrium solver.
Parameters¶
config_path : str | Path Path to a JSON configuration file describing the reactor geometry, coil set, physics parameters and solver settings.
Attributes¶
Psi : FloatArray Poloidal flux on the (Z, R) grid. J_phi : FloatArray Toroidal current density on the (Z, R) grid. B_R, B_Z : FloatArray Radial and vertical magnetic field components (set after solve).
Load reactor configuration from a JSON file.
Parameters¶
path : str | Path Filesystem path to the configuration JSON.
Build a free-boundary CoilSet from the active JSON configuration.
solve(*, boundary_variant=None, coils=None, preserve_initial_state=False, boundary_flux=None, max_outer_iter=20, tol=0.0001, optimize_shape=False, tikhonov_alpha=0.0001)
Solve using the requested boundary variant.
The default is taken from solver.boundary_variant in the active config.
Explicit entry point for the fixed-boundary variant.
Compute the vacuum poloidal flux from the external coil set.
Uses elliptic integrals (toroidal Green's function) for each coil.
Returns¶
FloatArray Vacuum flux Psi_vac on the (NZ, NR) grid.
staticmethod
¶
Evaluate a modified-tanh pedestal profile (vectorised).
Parameters¶
psi_norm : FloatArray
Normalised poloidal flux (0 at axis, 1 at separatrix).
params : dict[str, float]
Profile shape parameters with keys ped_top, ped_width,
ped_height, core_alpha.
Returns¶
FloatArray Profile value; zero outside the plasma region.
staticmethod
¶
Evaluate d(mtanh_profile)/dpsi_norm for H-mode Newton linearisation.
Compute the toroidal current density J_phi from the GS source.
Uses J_phi = R p'(psi) + FF'(psi) / (mu0 R) with either
L-mode (linear) or H-mode (mtanh) profiles, then renormalises to
match the target plasma current.
Parameters¶
Psi_axis : float Poloidal flux at the magnetic axis (O-point). Psi_boundary : float Poloidal flux at the separatrix (X-point or limiter).
Returns¶
FloatArray
Updated J_phi on the (NZ, NR) grid.
Run the full Picard-iteration equilibrium solver.
Iterates: topology analysis -> source update -> elliptic solve -> under-relaxation until the residual drops below the configured convergence threshold or the maximum iteration count is reached.
When solver.solver_method is "anderson", Anderson
acceleration is applied every few Picard steps to speed up
convergence.
Returns¶
dict[str, Any]
{"psi": FloatArray, "converged": bool, "iterations": int,
"residual": float, "residual_history": list[float],
"gs_residual": float, "gs_residual_best": float,
"gs_residual_history": list[float],
"wall_time_s": float, "solver_method": str}
Parameters¶
preserve_initial_state : bool, optional
Keep current interior self.Psi values and only enforce boundary
conditions before the iterative solve. Default is False.
boundary_flux : FloatArray | None, optional
Explicit boundary map to enforce during solve. Must have shape
(NZ, NR) when provided.
Raises¶
RuntimeError
If the solver produces NaN or Inf and solver.fail_on_diverge
is enabled in the active configuration.
optimize_coil_currents(coils, target_flux, tikhonov_alpha=0.0001, *, x_point_flux_target=None, divertor_flux_targets=None)
Find coil currents that best satisfy free-boundary target constraints.
Solves a bounded linear least-squares problem that can include:
- boundary-flux targets at
target_flux_points - X-point isoflux and null-field constraints
-
divertor strike-point isoflux constraints
min_I || A I - b ||^2 + alpha * ||I||^2 s.t. -I_max <= I <= I_max (per coil)
where A stacks the active constraint blocks.
Parameters¶
coils : CoilSet
Coil geometry and optional objective targets.
target_flux : FloatArray, shape (n_pts,)
Desired poloidal flux at target_flux_points. Can be empty when
only X-point and/or divertor constraints are active.
tikhonov_alpha : float
Regularisation strength to penalise large currents.
x_point_flux_target : float or None
Scalar target flux at coils.x_point_target.
divertor_flux_targets : FloatArray or None
Flux targets at coils.divertor_strike_points.
Returns¶
FloatArray, shape (n_coils,) Optimised coil currents [A].
solve_free_boundary(coils, max_outer_iter=20, tol=0.0001, optimize_shape=False, tikhonov_alpha=0.0001, objective_tolerances=None)
Experimental external-coil outer loop around the fixed-boundary GS solve.
Iterates between updating boundary flux from coils and solving the
internal GS equation. When optimize_shape=True and the coil set
has target_flux_points, an additional outer loop optimises the
coil currents to match the desired plasma boundary shape.
This helper should not be interpreted as a complete production-grade free-boundary solver. The current project standard for closing the free-boundary roadmap item is higher: shape control, X-point geometry, and divertor-configuration support must all be demonstrated.
Parameters¶
coils : CoilSet
External coil set.
max_outer_iter : int
Maximum outer-loop iterations for the experimental coil-coupled path.
tol : float
Convergence tolerance on max |delta psi|.
optimize_shape : bool
When True, run coil-current optimisation at each outer step.
tikhonov_alpha : float
Tikhonov regularisation for coil optimisation.
objective_tolerances : dict or None
Optional convergence gates for free-boundary target objectives.
Supported keys are shape_rms, shape_max_abs,
x_point_position, x_point_gradient, x_point_flux,
divertor_rms, and divertor_max_abs. When omitted, the
method falls back to free_boundary.objective_tolerances in the
config, if present.
Returns¶
dict
{"outer_iterations": int, "final_diff": float,
"coil_currents": NDArray}
phase_sync_step(theta, omega, *, dt=0.001, K=None, alpha=None, zeta=None, psi_driver=None, psi_mode=None, actuation_gain=None)
Reduced-order plasma sync kernel (phase reduction).
dθ_i/dt = ω_i + K·R·sin(ψ_r − θ_i − α) + ζ·sin(Ψ − θ_i)
Ψ is exogenous when psi_mode="external" (no dotΨ equation). This is the reviewer's ζ sin(Ψ−θ) injection for plasma sync stability.
phase_sync_step_lyapunov(theta, omega, *, n_steps=100, dt=0.001, K=None, zeta=None, psi_driver=None, psi_mode=None)
Multi-step phase sync with Lyapunov stability tracking.
Returns final state, R trajectory, V trajectory, and λ exponent. λ < 0 ⟹ stable convergence toward Ψ.
TokamakConfig¶
dataclass
¶
Tokamak equilibrium operating point.
Units: lengths [m], field [T], current [MA], density [1e19 m^-3], temperature [keV], power [MW].
classmethod
¶
Return the JSON Schema for serialized tokamak configurations.
TransportSolver¶
TransportSolver(config_path, *, nr=50, multi_ion=False, transport_model='gyro_bohm', external_gk_allow_gyrobohm_fallback=False, tglf_native_allow_gyrobohm_fallback=False, allow_constant_transport_fallback=False, allow_simplified_bootstrap_fallback=False, allow_legacy_approximations=False)
Bases: FusionKernel
1.5D Integrated Transport Code. Solves Heat and Particle diffusion equations on flux surfaces, coupled self-consistently with the 2D Grad-Shafranov equilibrium.
When multi_ion=True, the solver evolves separate D/T fuel densities,
He-ash transport with pumping (configurable tau_He), independent
electron temperature Te, coronal-equilibrium tungsten radiation
(Pütterich et al. 2010), and per-cell Bremsstrahlung.
property
¶
Relative energy conservation error from the last evolution step.
property
¶
Relative particle conservation error from the last evolution step.
Only meaningful in multi-ion mode. Returns 0.0 otherwise.
Configure Chang-Hinton neoclassical transport model.
When set, update_transport_model uses the Chang-Hinton formula instead of the constant chi_base = 0.5.
Backward-compatible Chang-Hinton profile helper.
Older parity tests call this no-arg method on a partially-initialized transport object. Keep the method as a thin adapter over the module function so those tests remain stable.
Models impurity influx from PWI erosion. Simple diffusion model: Source at edge, diffuses inward.
Calculate bootstrap current using full Sauter closure by default.
If neoclassical configuration is missing, this method fails closed unless
allow_simplified_bootstrap_fallback=True was explicitly enabled.
Gyro-Bohm + neoclassical transport model with EPED-like pedestal.
When neoclassical params are set, uses: - Chang-Hinton neoclassical chi as additive floor - Gyro-Bohm anomalous transport (calibrated c_gB) - EPED-like pedestal model for H-mode boundary condition
Fails closed when neoclassical parameters are not configured, unless
allow_constant_transport_fallback=True is explicitly enabled.
Advance Ti by one time step using Crank-Nicolson implicit diffusion.
The scheme is unconditionally stable, allowing dt up to ~1.0 s without NaN. The full equation solved is:
(T^{n+1} - T^n)/dt = 0.5*[L_h(T^{n+1}) + L_h(T^n)] + S - Sink
Parameters¶
dt : float
Time step [s].
P_aux : float
Auxiliary heating power [MW].
enforce_conservation : bool
When True, raise :class:PhysicsError if the per-step energy
conservation error exceeds 1%.
ped_te : PedestalProfile, optional
Pedestal profile model for electron temperature.
ped_ti : PedestalProfile, optional
Pedestal profile model for ion temperature.
Projects the 1D radial profiles back onto the 2D Grad-Shafranov grid, including neoclassical bootstrap current.
Compute the energy confinement time from stored energy.
τ_E = W_stored / P_loss, where W_stored = ∫ 3/2 n (Ti+Te) dV and the volume element is estimated from the 1D radial profiles using cylindrical approximation.
Parameters¶
P_loss_MW : float Total loss power [MW]. Must be > 0.
Returns¶
float Energy confinement time [s].
Run self-consistent GS <-> transport iteration.
This implements the standard integrated-modelling loop used by codes such as ASTRA and JINTRAC: evolve the 1D transport for n_inner steps, project profiles onto the 2D grid, re-solve the Grad-Shafranov equilibrium, and repeat until the poloidal-flux change drops below psi_tol.
Algorithm¶
- Run transport for n_inner steps (evolve Ti/Te/ne).
- Call :meth:
map_profiles_to_2dto updateJ_phion the 2D grid. - Re-solve the Grad-Shafranov equilibrium with the updated source.
- Check psi convergence:
||Psi_new - Psi_old|| / ||Psi_old|| < psi_tol. - Repeat until converged or n_outer iterations exhausted.
Parameters¶
P_aux : float Auxiliary heating power [MW]. n_inner : int Number of transport evolution steps per outer iteration. n_outer : int Maximum number of outer (GS re-solve) iterations. dt : float Transport time step [s]. psi_tol : float Relative psi convergence tolerance.
Returns¶
dict
{"T_avg": float, "T_core": float, "tau_e": float,
"n_outer_converged": int, "psi_residuals": list[float],
"Ti_profile": ndarray, "ne_profile": ndarray,
"converged": bool}
run_to_steady_state(P_aux, n_steps=500, dt=0.01, adaptive=False, tol=0.001, self_consistent=False, sc_n_inner=100, sc_n_outer=10, sc_psi_tol=0.001)
Run transport evolution until approximate steady state.
Parameters¶
P_aux : float
Auxiliary heating power [MW].
n_steps : int
Number of evolution steps.
dt : float
Time step [s] (initial value when adaptive=True).
adaptive : bool
Use Richardson-extrapolation adaptive time stepping.
tol : float
Error tolerance for adaptive stepping.
self_consistent : bool
When True, delegate to :meth:run_self_consistent which
iterates GS <-> transport to convergence. The remaining
sc_* parameters are forwarded.
sc_n_inner : int
Transport steps per outer GS iteration (self-consistent mode).
sc_n_outer : int
Maximum outer GS iterations (self-consistent mode).
sc_psi_tol : float
Relative psi convergence tolerance (self-consistent mode).
Returns¶
dict
{"T_avg": float, "T_core": float, "tau_e": float,
"n_steps": int, "Ti_profile": ndarray,
"ne_profile": ndarray}
When adaptive=True, also includes dt_final,
dt_history, error_history.
When self_consistent=True, returns the
:meth:run_self_consistent dict instead.
Scaling Laws¶
Evaluate the IPB98(y,2) confinement time scaling law.
Parameters¶
Ip : float Plasma current [MA]. BT : float Toroidal magnetic field [T]. ne19 : float Line-averaged electron density [10^19 m^-3]. Ploss : float Loss power [MW]. Must be > 0. R : float Major radius [m]. kappa : float Elongation (κ). epsilon : float Inverse aspect ratio a/R. M : float Effective ion mass [AMU]. Default 2.5 (D-T). coefficients : dict, optional Pre-loaded coefficient dict. If None, loaded from disk.
Returns¶
float Predicted thermal energy confinement time τ_E [s].
Raises¶
ValueError If any input is non-finite or non-positive.
Compute the H-factor (enhancement factor over scaling law).
Parameters¶
tau_actual : float Positive measured or simulated confinement time [s]. tau_predicted : float Positive IPB98(y,2) predicted confinement time [s].
Returns¶
float H98(y,2) = tau_actual / tau_predicted.
Raises¶
ValueError If either confinement time is non-finite or non-positive.
GEQDSK I/O¶
dataclass
¶
GEqdsk(description='', nw=0, nh=0, rdim=0.0, zdim=0.0, rcentr=0.0, rleft=0.0, zmid=0.0, rmaxis=0.0, zmaxis=0.0, simag=0.0, sibry=0.0, bcentr=0.0, current=0.0, fpol=(lambda: array([], dtype=float64))(), pres=(lambda: array([], dtype=float64))(), ffprime=(lambda: array([], dtype=float64))(), pprime=(lambda: array([], dtype=float64))(), qpsi=(lambda: array([], dtype=float64))(), psirz=(lambda: array([], dtype=float64))(), rbdry=(lambda: array([], dtype=float64))(), zbdry=(lambda: array([], dtype=float64))(), rlim=(lambda: array([], dtype=float64))(), zlim=(lambda: array([], dtype=float64))())
Container for all data in a G-EQDSK file.
Read a G-EQDSK file and return a :class:GEqdsk container.
Handles both whitespace-separated and fixed-width Fortran formats,
including the common case where values run together without spaces
(e.g. 2.385E+00-1.216E+01).
Parameters¶
path : str or Path Path to the G-EQDSK file.
Returns¶
GEqdsk Parsed equilibrium data.
Write a :class:GEqdsk to a G-EQDSK file.
Parameters¶
eq : GEqdsk Equilibrium data. path : str or Path Output file path.
Uncertainty Quantification¶
Monte Carlo uncertainty quantification for fusion performance.
Samples scaling-law coefficients from their Gaussian posteriors and propagates through the confinement and fusion power models.
Parameters¶
scenario : PlasmaScenario Plasma parameters (held fixed). n_samples : int Number of Monte Carlo samples (default 10,000). seed : int, optional Random seed for reproducibility.
Returns¶
UQResult — central estimates + error bars + percentiles.
quantify_full_chain(scenario, n_samples=5000, seed=None, chi_gB_sigma=0.3, pedestal_sigma=0.2, boundary_sigma=0.02)
Full-chain Monte Carlo uncertainty propagation: equilibrium -> transport -> fusion power -> gain.
Now uses correlated sampling for IPB98 coefficients.
dataclass
¶
UQClaimEvidence(schema_version, source, source_id, scenario_source, prior_source, propagation_chain, sensitivity_source, model_id, seed, n_samples, tau_E_s, P_fusion_MW, Q, tau_E_sigma, P_fusion_sigma, Q_sigma, tau_E_percentiles_ordered, P_fusion_percentiles_ordered, Q_percentiles_ordered, finite_outputs, fuel_ion_fraction, dP_dn_MW_per_1e19m3, dP_dT_MW_per_keV, tau_E_relative_error, P_fusion_relative_error, Q_relative_error, sigma_relative_error, relative_tolerance, sigma_relative_tolerance, calibrated_uq_claim_allowed, claim_status)
Serialisable provenance and reference-comparison evidence for UQ claims.
uq_claim_evidence(scenario, result, *, source, source_id, scenario_source, prior_source, propagation_chain, sensitivity_source, seed, model_id='bounded_full_chain_uq', reference_tau_E=None, reference_P_fusion=None, reference_Q=None, reference_tau_E_sigma=None, relative_tolerance=0.05, sigma_relative_tolerance=0.1)
Build fail-closed evidence for bounded or calibrated UQ claims.
Raise when UQ evidence is insufficient for a calibrated predictive claim.
Persist UQ claim evidence as deterministic JSON.
JAX-Accelerated Transport Primitives¶
Requires pip install "scpn-control[jax]". GPU execution automatic when jaxlib has CUDA/ROCm.
thomas_solve(a, b, c, d, *, use_jax=True, allow_numpy_fallback=False, allow_legacy_numpy_fallback=False)
Tridiagonal solve with automatic JAX/GPU dispatch.
Parameters¶
a : sub-diagonal, length n-1 b : main diagonal, length n c : super-diagonal, length n-1 d : right-hand side, length n use_jax : attempt JAX backend (falls back to NumPy if unavailable)
diffusion_rhs(T, chi, rho, drho, *, use_jax=True, allow_numpy_fallback=False, allow_legacy_numpy_fallback=False)
Cylindrical diffusion operator L_h(T) with JAX/GPU dispatch.
crank_nicolson_step(T, chi, source, rho, drho, dt, T_edge=0.1, *, use_jax=True, allow_numpy_fallback=False, allow_legacy_numpy_fallback=False)
Single Crank-Nicolson transport step with JAX/GPU dispatch.
Parameters¶
T : temperature profile, length n chi : diffusivity profile, length n source : net heating source, length n rho : radial grid, length n drho : grid spacing dt : timestep T_edge : edge boundary condition (Dirichlet), keV
Differentiable Transport Facade¶
Requires pip install "scpn-control[jax]" for gradient evaluation. The NumPy
path is deterministic for parity checks and non-JAX deployments, but
transport_loss_gradient() fails closed without JAX.
transport_parameter_gradients() extends the same traced Crank-Nicolson
contract to source schedules, returning JAX gradients for both transport
coefficients and additive heating, fuelling, or impurity-source inputs.
differentiable_transport_rollout() advances a bounded multi-step source
schedule with the same four-channel boundary contract.
transport_rollout_source_gradients() returns fail-closed JAX gradients for
that full source schedule so controller tuning can optimise time-distributed
heating, fuelling, and impurity-source inputs without finite differences. The
rollout gradient path keeps the loss inside the traced JAX graph and enables
JAX x64 before importing jax.numpy, so persisted dtype evidence is not
silently downgraded.
audit_transport_rollout_source_gradients() and
assert_transport_rollout_source_gradients_consistent() compare those rollout
gradients against sampled NumPy finite-difference perturbations.
audit_transport_parameter_gradients() and
assert_transport_parameter_gradients_consistent() compare those JAX gradients
against sampled independent finite-difference perturbations before
controller-tuning admission.
transport_coefficients_from_neural_closure() maps bounded neural transport
closure outputs into the four-channel coefficient order used by the facade:
electron heat, ion heat, electron particle diffusivity, and a declared impurity
diffusivity fraction.
gyrokinetic_transport_closure_profiles() wraps the reduced gyrokinetic
transport profile evaluator as bounded closure provenance, and
transport_coefficients_from_gyrokinetic_closure() maps that closure into the
same four-channel coefficient order without promoting the reduced GK model to
an externally validated transport claim.
transport_campaign_metadata() records backend, dtype, radial grid, timestep,
boundary conditions, closure provenance, gradient tolerance, and optional
equilibrium-grid shape for reproducible controller-tuning campaigns.
save_transport_campaign_metadata() and load_transport_campaign_metadata()
persist the same contract as schema-versioned JSON and fail closed on malformed
or physically inconsistent replay metadata.
assert_transport_campaign_metadata_replay() compares archived campaign
metadata with a candidate setup and raises on backend, grid, boundary, closure,
gradient-tolerance, or equilibrium-shape drift before controller tuning reruns.
transport_differentiability_evidence() and
assert_transport_differentiability_claim_admissible() add a tamper-evident
admission envelope over campaign metadata and gradient-audit results. That
envelope requires JAX backend evidence, passed sampled finite-difference
gradient audit, stable SHA-256 digests for the campaign and audit payloads, and
an optional link to the safety-critical controller proof artifact digest.
Admission revalidates finite non-negative audit losses and errors, tolerance
agreement with campaign metadata, unique in-domain sampled audit indices,
pass/fail consistency with maximum audit error, and ordered latency percentiles
before persisted controller-tuning evidence is accepted. Latency reports also
bind runtime provenance for later CPU/GPU comparison: Python version, platform,
machine class, JAX and jaxlib versions, default backend, visible JAX devices,
and x64 state.
transport_full_fidelity_readiness_evidence() and
assert_transport_full_fidelity_claim_ready() add the stricter promotion gate
for full-fidelity differentiable-transport claims. The gate binds campaign
metadata, one-step gradient-latency evidence, rollout source-gradient latency
evidence, audit digests, controller formal-proof digests, equilibrium-coupled
metadata, and an admitted external reference artifact. The repository benchmark
binds the canonical payload SHA-256 from validation/reports/scpn_z3_formal.json
when that bounded Z3 formal report is present and passing. Without all of those
inputs the claim remains bounded local differentiability evidence.
equilibrium_weighted_transport_rollout_tracking_loss() extends the optional
Grad-Shafranov flux-map weighting from one transport step to a full source
rollout. equilibrium_weighted_transport_rollout_source_gradient() returns
fail-closed JAX gradients with respect to both the source schedule and the
equilibrium flux map for controller-tuning studies.
differentiable_transport_step(profiles, chi, sources, rho, dt, edge_values, *, use_jax=True, allow_numpy_fallback=False, allow_legacy_numpy_fallback=False)
Advance four transport channels by one differentiable radial step.
Channel order is electron temperature, ion temperature, electron density, and impurity density. The radial coordinate is a strictly increasing, uniformly spaced normalised axis from core-side interior to edge. The core boundary uses the inherited zero-gradient condition, and each channel uses its supplied Dirichlet edge value.
transport_tracking_loss(profiles, chi, sources, target_profiles, rho, dt, edge_values, *, weights=None, use_jax=True, allow_numpy_fallback=False, allow_legacy_numpy_fallback=False)
Return weighted one-step transport tracking loss for controller tuning.
transport_loss_gradient(profiles, chi, sources, target_profiles, rho, dt, edge_values, *, weights=None)
Return the tracking loss and JAX gradient with respect to chi.
transport_parameter_gradients(profiles, chi, sources, target_profiles, rho, dt, edge_values, *, weights=None)
Return JAX gradients with respect to chi and source schedules.
This is the controller-tuning primitive for differentiable auxiliary
heating, fuelling, and impurity-source schedules. It keeps the same
four-channel Crank-Nicolson, source-term, core zero-gradient, and edge
Dirichlet contracts as :func:differentiable_transport_step; unlike
:func:transport_loss_gradient, it exposes gradients for both turbulent
transport coefficients and additive source terms.
dataclass
¶
JAX gradients of transport tracking loss for tunable transport inputs.
differentiable_transport_rollout(initial_profiles, chi, source_sequence, rho, dt, edge_values, *, use_jax=True, allow_numpy_fallback=False, allow_legacy_numpy_fallback=False)
Advance a time-series of four-channel transport source schedules.
The returned array has shape (n_steps, 4, n_rho). Transport
coefficients are held fixed over the rollout, while source_sequence
supplies differentiable additive heating, fuelling, and impurity-source
schedules at each step. This is a bounded controller-tuning primitive, not
an externally validated integrated-transport campaign.
transport_rollout_tracking_loss(initial_profiles, chi, source_sequence, target_history, rho, dt, edge_values, *, weights=None, use_jax=True, allow_numpy_fallback=False, allow_legacy_numpy_fallback=False)
Return weighted multi-step transport tracking loss.
transport_rollout_source_gradients(initial_profiles, chi, source_sequence, target_history, rho, dt, edge_values, *, weights=None)
Return JAX gradients for a multi-step transport source schedule.
dataclass
¶
JAX gradients of multi-step transport loss for source schedules.
audit_transport_rollout_source_gradients(initial_profiles, chi, source_sequence, target_history, rho, dt, edge_values, *, weights=None, epsilon=1e-05, tolerance=0.0005, sample_indices=None)
Compare JAX rollout source gradients with sampled finite differences.
assert_transport_rollout_source_gradients_consistent(initial_profiles, chi, source_sequence, target_history, rho, dt, edge_values, *, weights=None, epsilon=1e-05, tolerance=0.0005, sample_indices=None)
Return rollout source-gradient audit evidence or fail closed.
dataclass
¶
TransportRolloutGradientAudit(loss, epsilon, tolerance, checked_indices, source_max_abs_error, passed)
Finite-difference audit of multi-step source-rollout gradients.
audit_transport_parameter_gradients(profiles, chi, sources, target_profiles, rho, dt, edge_values, *, weights=None, epsilon=1e-05, tolerance=0.0005, sample_indices=None)
Audit JAX transport gradients against independent finite differences.
The audit is intentionally sampled rather than exhaustive so it can run in controller-tuning admission checks. It evaluates deterministic interior radial points for every channel by default and compares JAX gradients for both turbulent transport coefficients and additive source schedules against independently perturbed NumPy losses.
assert_transport_parameter_gradients_consistent(profiles, chi, sources, target_profiles, rho, dt, edge_values, *, weights=None, epsilon=1e-05, tolerance=0.0005, sample_indices=None)
Return the gradient audit or fail closed when consistency is violated.
dataclass
¶
TransportGradientAudit(loss, epsilon, tolerance, checked_indices, chi_max_abs_error, source_max_abs_error, passed)
Finite-difference audit of differentiable transport tuning gradients.
benchmark_transport_parameter_gradient_latency(profiles, chi, sources, target_profiles, rho, dt, edge_values, *, weights=None, epsilon=1e-05, tolerance=0.0005, sample_indices=None, warmup_runs=1, timed_runs=5)
Measure audited JAX gradient-admission latency for controller tuning.
The measured path is intentionally the fail-closed admission contract: JAX gradients for transport coefficients and source schedules plus sampled independent finite-difference audit. The report is local timing evidence, not a real-time control-loop guarantee.
dataclass
¶
TransportGradientLatencyReport(schema_version, backend, dtype, n_rho, channel_count, warmup_runs, timed_runs, p50_ms, p95_ms, max_ms, runtime_metadata, audit, claim_status)
Latency evidence for audited differentiable transport tuning gradients.
Persist differentiable transport gradient-latency evidence as JSON.
benchmark_transport_rollout_source_gradient_latency(initial_profiles, chi, source_sequence, target_history, rho, dt, edge_values, *, weights=None, epsilon=1e-05, tolerance=0.0005, sample_indices=None, warmup_runs=1, timed_runs=5)
Measure audited multi-step source-rollout gradient latency.
The measured path is the controller-admission contract for source schedules: JAX rollout gradients with a sampled independent NumPy finite-difference audit. The report is local timing evidence, not a real-time control-loop or externally validated transport claim.
dataclass
¶
TransportRolloutGradientLatencyReport(schema_version, backend, dtype, n_rho, n_steps, channel_count, warmup_runs, timed_runs, p50_ms, p95_ms, max_ms, runtime_metadata, audit, claim_status)
Latency evidence for audited multi-step source-rollout gradients.
Persist rollout source-gradient latency evidence as JSON.
transport_coefficients_from_neural_closure(closure, *, impurity_diffusivity_fraction=1.0, chi_floor=0.0)
Map a neural transport closure into four facade coefficient channels.
The facade channel order is electron temperature, ion temperature, electron density, and impurity density. Neural transport closures provide electron heat, ion heat, and particle diffusivity profiles; the electron-density channel uses the particle diffusivity, while the impurity-density channel uses a declared bounded fraction of the same particle diffusivity until species-resolved impurity transport is externally validated.
Wrap a reduced gyrokinetic profile closure for controller tuning.
model must provide the existing evaluate_profile(rho, profiles)
contract and return ion heat, electron heat, and particle diffusivity
profiles. This adapter validates only the CONTROL boundary and does not
promote the reduced GK model to an externally validated transport claim.
transport_coefficients_from_gyrokinetic_closure(closure, *, impurity_diffusivity_fraction=1.0, chi_floor=0.0)
Map a reduced gyrokinetic closure into four facade channels.
dataclass
¶
Reduced gyrokinetic closure profiles for differentiable transport.
dataclass
¶
TransportCampaignMetadata(backend, dtype, channel_order, n_rho, rho_min, rho_max, rho_spacing, dt, core_boundary, edge_boundary, edge_values, closure_source, closure_weights_checksum, gradient_tolerance, equilibrium_grid_shape)
Validated provenance for a differentiable transport tuning campaign.
transport_campaign_metadata(profiles, chi, sources, rho, dt, edge_values, *, backend, closure=None, gradient_tolerance=None, equilibrium_psi=None)
Return serialisable provenance for differentiable transport campaigns.
Persist transport campaign metadata as schema-versioned JSON.
Load and validate schema-versioned transport campaign metadata JSON.
assert_transport_campaign_metadata_replay(archived, profiles, chi, sources, rho, dt, edge_values, *, backend, closure=None, gradient_tolerance=None, equilibrium_psi=None)
Validate that a candidate transport setup matches archived metadata.
This guard is intended for replaying differentiable transport tuning campaigns. It fails closed on backend, dtype, grid, timestep, boundary, closure-provenance, gradient-tolerance, or equilibrium-shape drift before a controller rerun can silently compare against a different physics setup.
dataclass
¶
TransportDifferentiabilityEvidence(schema_version, backend, campaign_sha256, gradient_audit_sha256, gradient_tolerance, audit_kind, audit_passed, n_rho, channel_order, equilibrium_coupled, controller_formal_artifact_sha256, claim_status)
Tamper-evident admission evidence for differentiable transport gradients.
Build tamper-evident evidence for audited differentiable transport.
Fail closed unless differentiable transport evidence matches inputs.
dataclass
¶
TransportFullFidelityReadinessEvidence(schema_version, backend, campaign_sha256, gradient_latency_report_sha256, gradient_audit_sha256, rollout_latency_report_sha256, rollout_audit_sha256, external_reference_artifact_sha256, external_reference_admitted, controller_formal_artifact_sha256, n_rho, rollout_steps, channel_order, equilibrium_coupled, full_fidelity_claim_admissible, blocked_reasons, claim_status)
Fail-closed readiness evidence for full-fidelity transport claims.
transport_full_fidelity_readiness_evidence(metadata, gradient_report, *, rollout_report=None, external_reference_artifact_sha256=None, external_reference_admitted=False, controller_formal_artifact_sha256=None)
Build fail-closed readiness evidence for full-fidelity transport claims.
Local gradient audits and timing reports are necessary but not sufficient for a full-fidelity claim. This certificate requires a JAX campaign, equilibrium coupling, one-step and rollout audit reports, a controller proof digest, and an independently admitted external reference artefact.
assert_transport_full_fidelity_claim_ready(evidence, metadata, gradient_report, *, rollout_report=None)
Fail closed unless readiness evidence admits a full-fidelity claim.
Return positive mean-one radial weights from a Grad-Shafranov flux map.
equilibrium_weighted_transport_tracking_loss(profiles, chi, sources, target_profiles, rho, dt, edge_values, equilibrium_psi, *, weights=None, use_jax=True, allow_numpy_fallback=False, allow_legacy_numpy_fallback=False)
Return transport tracking loss weighted by GS-equilibrium flux geometry.
equilibrium_weighted_transport_loss_gradient(profiles, chi, sources, target_profiles, rho, dt, edge_values, equilibrium_psi, *, weights=None)
Return JAX gradients of equilibrium-weighted transport loss.
The returned gradients are with respect to the transport coefficients and the supplied equilibrium flux map. If the flux map was produced inside an outer JAX graph by the Grad-Shafranov solver, this loss is compatible with further chain-rule propagation through that equilibrium solve.
dataclass
¶
JAX gradient of equilibrium-weighted transport tracking loss.
equilibrium_weighted_transport_rollout_tracking_loss(initial_profiles, chi, source_sequence, target_history, rho, dt, edge_values, equilibrium_psi, *, weights=None, use_jax=True, allow_numpy_fallback=False, allow_legacy_numpy_fallback=False)
Return multi-step transport rollout loss weighted by GS flux geometry.
equilibrium_weighted_transport_rollout_source_gradient(initial_profiles, chi, source_sequence, target_history, rho, dt, edge_values, equilibrium_psi, *, weights=None)
Return JAX gradients of GS-weighted rollout loss.
The returned gradients are with respect to the full source schedule and the supplied equilibrium flux map. If the flux map was produced inside an outer JAX graph by the Grad-Shafranov solver, this loss is compatible with chain-rule propagation through that equilibrium solve.
dataclass
¶
EquilibriumWeightedTransportRolloutGradient(loss, source_gradient, equilibrium_gradient, radial_weights, final_profiles)
JAX gradient of equilibrium-weighted multi-step transport rollout loss.
Neural Equilibrium¶
NeuralEquilibriumAccelerator.pretrain_from_synthetic_equilibria() trains
JAX-compatible PCA plus MLP weights on bounded synthetic Solovev-like
equilibria for pretraining. The corresponding real EFIT fine-tuning entry point
fine_tune_from_efit_reconstructions() fails closed unless the persisted
P-EFIT or documented-public-reference artefact validator passes.
PCA + MLP surrogate for Grad-Shafranov equilibrium.
Can be trained from SPARC GEQDSK files (preferred) or from a FusionKernel config with coil perturbations.
Evaluate on test set. Returns dict with mse, max_error, gs_residual.
pretrain_from_synthetic_equilibria(n_samples=2048, *, seed=20240531, ridge_alpha=1e-06, save_path=None)
Pretrain the surrogate on bounded synthetic equilibria.
Uses PCA compression followed by a deterministic ridge-regression readout. The saved weights are compatible with the existing JAX inference path because the readout is represented as a one-layer MLP.
fine_tune_from_efit_reconstructions(geqdsk_paths, *, reference_artifact_root=None, require_reference_artifacts=True, n_perturbations=5, seed=42)
Fine-tune on real EFIT/P-EFIT artefacts only after evidence admission.
Train on real SPARC GEQDSK equilibria with perturbations.
For each GEQDSK file, generates n_perturbations by scaling p'/FF' profiles, yielding n_files × n_perturbations training pairs.
Input features (12-dim): [I_p, B_t, R_axis, Z_axis, pprime_scale, ffprime_scale, simag, sibry, kappa, delta_upper, delta_lower, q95]
Output: flattened ψ(R,Z) → PCA coefficients
Uses 70/15/15 train/val/test split with val-loss early stopping (patience=20) and combined MSE + lambda_gs * GS residual loss.
dataclass
¶
NeuralEquilibriumClaimEvidence(schema_version, model_id, source, source_id, weights_path, weights_sha256, reference_source, reference_dataset_id, reference_artifact_sha256, reference_equilibria_count, grid_shape, feature_names, n_components, explained_variance, synthetic_test_mse, synthetic_gs_residual, psi_rmse_Wb, pressure_rmse_Pa, q_profile_rmse, boundary_rmse_m, axis_position_error_m, psi_tolerance_Wb, pressure_tolerance_Pa, q_profile_tolerance, boundary_tolerance_m, axis_position_tolerance_m, facility_claim_allowed, claim_status)
Serialisable admission evidence for neural-equilibrium predictive claims.
Generate bounded Solovev-like equilibria for pretraining.
The generated targets are synthetic Grad-Shafranov-shaped flux maps for pretraining only. They are not matched EFIT or P-EFIT validation evidence.
neural_equilibrium_claim_evidence(pretraining, *, weights_path, source, source_id, model_id='neural_equilibrium_pca_mlp', reference_artifact_path=None)
Build fail-closed evidence for neural-equilibrium predictive claims.
Synthetic pretraining may support bounded pretraining and inference-plumbing statements only. Facility or predictive claims require a persisted P-EFIT/documented-reference artefact that validates against the same weight checksum and declares psi, pressure, q-profile, boundary, and magnetic-axis error tolerances.
Return evidence only when strict matched-reference admission passed.
Persist neural-equilibrium claim evidence as deterministic JSON.
pretrain_neural_equilibrium_synthetic(*, n_samples=2048, save_path=REPO_ROOT / 'weights' / 'neural_equilibrium_synthetic_pretrain.npz', grid_shape=(65, 65), n_components=20, seed=20240531)
Convenience entry point for bounded synthetic pretraining.
dataclass
¶
PretrainingResult(n_samples, n_components, explained_variance, train_mse, val_mse, test_mse, test_max_error, gs_residual, train_time_s, weights_path, evidence_kind, campaign)
Synthetic neural-equilibrium pretraining summary.
dataclass
¶
Metadata for bounded synthetic equilibrium pretraining data.
JAX-Accelerated Neural Equilibrium¶
Requires pip install "scpn-control[jax]". GPU and autodiff via jax.grad.
Predict psi(R,Z) from input features using JAX-accelerated pipeline.
Parameters¶
features : shape (n_features,) or (batch, n_features) mlp_weights : ((W0, W1, ...), (b0, b1, ...)) pca_params : (pca_mean, pca_components) norm_params : (input_mean, input_std) grid_shape : (nh, nw) output grid dimensions
Returns¶
psi : shape (nh, nw) or (batch, nh, nw)
Load neural equilibrium .npz weights into JAX-compatible structures.
Returns¶
(mlp_weights, pca_params, norm_params, grid_shape) mlp_weights: ((W0, W1, ...), (b0, b1, ...)) pca_params: (pca_mean, pca_components) norm_params: (input_mean, input_std) grid_shape: (nh, nw)
Neural Transport¶
cross_validate_neural_transport() benchmarks the active surrogate against the
analytic critical-gradient reference across fixed regime cases and a canonical
profile, so shipped weights can be checked against a deterministic local
baseline instead of only reporting synthetic training RMSE.
neural_transport_closure_profiles() packages profile transport coefficients
for controller and differentiable-transport coupling. It validates finite
strictly ordered profile inputs, fails closed when neural weights are required
but unavailable, and records whether coefficients came from loaded weights or
the analytic fallback.
NeuralTransportModel(weights_path=None, *, auto_discover=True, allow_weight_load_fallback=False, allow_legacy_weight_load_fallback=False)
Neural transport surrogate with analytic fallback.
On construction, attempts to load MLP weights from weights_path.
If loading fails (file missing, wrong format), the model
transparently falls back to :func:critical_gradient_model.
Parameters¶
weights_path : str or Path, optional
Path to a .npz file containing MLP weights. The file must
contain arrays: w1, b1, w2, b2, w3, b3, input_mean,
input_std, output_scale.
Examples¶
model = NeuralTransportModel() # fallback mode inp = TransportInputs(grad_ti=8.0) fluxes = model.predict(inp) fluxes.chi_i > 0 True model.is_neural False
Predict turbulent transport fluxes for given local parameters.
Uses the neural MLP if weights are loaded, otherwise falls back to the analytic critical-gradient model.
Parameters¶
inp : TransportInputs Local plasma parameters at a single radial point.
Returns¶
TransportFluxes Predicted heat and particle diffusivities.
Predict transport coefficients on the full radial profile.
Computes normalised gradients from the profile arrays via central finite differences, then evaluates the surrogate at each radial point. When the MLP is loaded the entire profile is evaluated in a single batched forward pass (no Python loop).
Parameters¶
rho : FloatArray
Normalised radius grid (0 to 1), shape (N,).
te, ti : FloatArray
Electron/ion temperature profiles [keV], shape (N,).
ne : FloatArray
Electron density profile [10^19 m^-3], shape (N,).
q_profile : FloatArray
Safety factor profile, shape (N,).
s_hat_profile : FloatArray
Magnetic shear profile, shape (N,).
r_major : float
Major radius [m] for gradient normalisation.
Returns¶
chi_e, chi_i, d_e : FloatArray
Transport coefficient profiles, each shape (N,).
dataclass
¶
Profile transport closure with provenance for controller coupling.
Parameters¶
chi_e, chi_i, d_e : FloatArray
Electron heat, ion heat, and particle diffusivity profiles.
channel_weights : FloatArray
Per-radius normalised contribution weights for chi_e, chi_i,
and d_e with shape (3, N).
source : str
"neural" when loaded weights are active, otherwise
"analytic_fallback".
weights_checksum : str or None
Loaded neural-weight checksum, or None for analytic fallback.
dataclass
¶
NeuralTransportClaimEvidence(schema_version, model_id, source, source_id, surrogate_mode, weights_path, weights_sha256, internal_weights_checksum, feature_schema, reference_source, reference_dataset_id, reference_artifact_sha256, reference_sample_count, local_case_count, local_max_abs_error, local_channel_agreement, local_per_channel_relative_rmse, local_profile_per_channel_relative_rmse, chi_i_rmse_m2_s, chi_e_rmse_m2_s, d_e_rmse_m2_s, chi_i_relative_mae, unstable_branch_accuracy, chi_i_rmse_tolerance_m2_s, chi_e_rmse_tolerance_m2_s, d_e_rmse_tolerance_m2_s, chi_i_relative_mae_tolerance, unstable_branch_accuracy_min, quantitative_claim_allowed, claim_status)
Serialisable admission evidence for neural-transport surrogate claims.
neural_transport_closure_profiles(rho, te, ti, ne, q_profile, s_hat_profile, *, model=None, r_major=6.2, require_neural=False, allow_fallback=False, allow_legacy_fallback=False)
Return bounded neural-transport closure profiles with provenance.
The function packages :meth:NeuralTransportModel.predict_profile for
controller and differentiable-transport coupling. It validates the radial
grid and plasma profiles before evaluation, fails closed when neural weights
are required but unavailable, and records whether the result came from the
neural surrogate or the analytic critical-gradient fallback.
Compare the current transport surrogate to the analytic reference.
Returns pointwise and profile-level error metrics against the
:func:critical_gradient_model benchmark. This is a reference check, not a
claim of agreement with full QLKNN-10D or first-principles turbulence
solvers.
neural_transport_claim_evidence(validation_result, *, source, source_id, model_id='neural_transport_qlknn_facade', weights_path=None, reference_artifact_path=None)
Build fail-closed evidence for neural-transport quantitative claims.
Return evidence only when strict matched-reference admission passed.
Persist neural-transport claim evidence as deterministic JSON.
MHD Stability¶
Run all available MHD stability criteria and return a summary.
Mercier, ballooning, and Kruskal-Shafranov are always evaluated.
Troyon requires beta_t, Ip_MA, a, and B0.
NTM requires j_bs, j_total, and a.
Parameters¶
qp : QProfile — pre-computed safety-factor profile beta_t : float, optional — total toroidal beta (dimensionless) Ip_MA : float, optional — plasma current [MA] a : float, optional — minor radius [m] B0 : float, optional — toroidal field on axis [T] j_bs : array, optional — bootstrap current density [A/m^2] j_total : array, optional — total current density [A/m^2]
Returns¶
StabilitySummary
IMAS Adapter¶
dataclass
¶
Minimal subset of IMAS equilibrium IDS for scpn-control interop.
Fields match equilibrium.time_slice[0].profiles_2d[0].
Extract EquilibriumIDS from a solved FusionKernel instance.
HPC Bridge¶
HPCBridge loads compiled Grad-Shafranov solver libraries only from absolute
dynamic-library paths. Package-local solver libraries are trusted by default.
External paths provided through SCPN_SOLVER_LIB require the additional
operator gate SCPN_ALLOW_EXTERNAL_SOLVER_LIB=1; without that gate the bridge
fails closed before calling the dynamic loader.
Interface between Python and the compiled C++ Grad-Shafranov solver.
Loads the shared library (libscpn_solver.so / scpn_solver.dll)
at construction time. If the library is not found the bridge
gracefully degrades — :meth:is_available returns False and the
caller falls back to Python.
Parameters¶
lib_path : str, optional
Explicit path to the shared library. When None (default) the
bridge searches trusted package-local locations only. SCPN_SOLVER_LIB
must be an absolute dynamic-library path and is accepted only when it
resolves under the package-local solver directories, unless the operator
also sets SCPN_ALLOW_EXTERNAL_SOLVER_LIB=1 for a vetted external
library.
Create the C++ solver instance for the given grid dimensions.
Set a fixed Dirichlet boundary value for psi edges, if supported.
Run the C++ solver for iterations sweeps.
Returns None if the library is not loaded (caller should fall back to a Python solver).
Run the C++ solver and write results into psi_out in-place.
Run solver until convergence, if native API is available.
Returns (psi, iterations_used, final_delta). If the library is
unavailable or uninitialized, returns None.
Run convergence API and write results into psi_out in-place.
Gyrokinetic Transport (v0.16.0)¶
Drop-in replacement for Gyro-Bohm transport scaling.
Public evaluation points use rho in [0, 1] and positive finite local
geometry, temperature, density, safety-factor, and charge inputs.
GK Solver Interface (v0.17.0)¶
Bases: ABC
Abstract base for gyrokinetic solvers.
Concrete subclasses must implement prepare_input, run,
and is_available.
dataclass
¶
GKLocalParams(R_L_Ti, R_L_Te, R_L_ne, q, s_hat, alpha_MHD=0.0, Te_Ti=1.0, Z_eff=1.5, nu_star=0.1, beta_e=0.01, epsilon=0.1, kappa=1.0, delta=0.0, rho=0.5, R0=6.2, a=2.0, B0=5.3, n_e=10.0, T_e_keV=8.0, T_i_keV=8.0)
Local plasma parameters at a single flux surface.
The first 11 fields match the GyrokineticsParams convention already
used by gyrokinetic_transport.py. The remaining fields add Miller
geometry and dimensional quantities needed by external GK codes.
dataclass
¶
GKOutput(chi_i, chi_e, D_e, D_i=0.0, gamma=(lambda: empty(0))(), omega_r=(lambda: empty(0))(), k_y=(lambda: empty(0))(), dominant_mode='stable', converged=True, electromagnetic=False)
Gyrokinetic solver output for a single flux surface.
Fluxes are in physical units [m^2/s]. Growth rate arrays are normalised to c_s / a.
Native Linear GK Solver (v0.17.0)¶
solve_linear_gk(species_list=None, geom=None, vgrid=None, R0=2.78, a=1.0, B0=2.0, q=1.4, s_hat=0.78, Z_eff=1.0, nu_star=0.01, n_ky_ion=16, n_ky_etg=0, n_theta=64, n_period=2, electromagnetic=False, beta_e=0.0, alpha_MHD=0.0)
Full k_y scan of the linear GK eigenvalue solver.
Parameters¶
species_list : list of GKSpecies Plasma species. Default: [deuterium, adiabatic electron]. geom : MillerGeometry Flux-tube geometry. Default: circular with given (R0, a, q, s_hat). vgrid : VelocityGrid Velocity-space grid. Default: (16, 24). n_ky_ion, n_ky_etg : int Number of k_y points on ion and electron scales. electromagnetic : bool Enable finite-beta EM corrections (A_∥, delta-B_∥). Default False. beta_e : float Electron beta = 2 mu_0 n_e T_e / B_0^2. alpha_MHD : float Normalised pressure gradient alpha = -q^2 R dp/dr / (B^2/2mu_0).
quasilinear_fluxes_from_spectrum(result, ion, R0=2.78, a=1.0, B0=2.0, saturation='mixing_length', electron=None)
Convert linear GK spectrum to physical transport fluxes.
Parameters¶
result : LinearGKResult Output from solve_linear_gk(). ion : GKSpecies Main ion species (for gyro-Bohm normalisation). R0, a, B0 : float Geometry and field for dimensional conversion. saturation : str Saturation model: "mixing_length" (default). electron : GKSpecies | None Electron species for Te/Ti ratio. Defaults to ion temperature (Te/Ti=1).
GK Hybrid Validation (v0.17.0)¶
OODDetector(*, mahalanobis_threshold=4.0, soft_sigma_threshold=2.0, ensemble_disagreement_threshold=0.3, training_mean=None, training_cov_inv=None)
Three-method OOD detector for surrogate transport models.
Parameters¶
mahalanobis_threshold : float Mahalanobis distance above which an input is flagged OOD. Default 4.0 corresponds to ~99.99% chi-squared quantile for 10D. soft_sigma_threshold : float Number of training-set standard deviations for soft range check. ensemble_disagreement_threshold : float Relative std across ensemble predictions above which to flag OOD.
Spot-check scheduler for hybrid surrogate+GK transport.
Decide whether to run GK spot-checks this step.
Parameters¶
rho : array Radial grid. chi_i : array Current ion diffusivity profile (from surrogate). ood_results : list of OODResult or None Per-surface OOD detector results (adaptive mode).
Returns¶
SpotCheckRequest or None None if no validation needed this step.
Apply corrections to surrogate transport based on GK spot-checks.
Ballooning Solver (v0.16.0)¶
Second-order ODE for ideal MHD ballooning stability in the s-alpha model.
Performs ballooning stability analysis given a QProfile.
Extracts (s, alpha) at each radial point and returns per-radius stability margin. Positive margin means stable. margin = alpha_crit(s) - alpha_actual
Binary search for alpha_crit at fixed shear s. Returns the critical alpha (first stability boundary).
Current Diffusion (v0.16.0)¶
Implicit Crank-Nicolson solver for poloidal flux diffusion.
Governing equation (1D in ρ, cylindrical limit): ∂ψ/∂t = (η_∥/μ₀) · [∂²ψ/∂ρ² + (1/ρ)·∂ψ/∂ρ] + R₀·η_∥·j_source Reference: Jardin (2010), Ch. 8, Eq. 8.8.
η_∥ from neoclassical_resistivity() — Sauter 1999/2002 + Hirshman 1981. Crank-Nicolson discretisation — θ=½ implicit-explicit split, second-order accurate in time and space.
Advance ψ by one timestep dt [s].
Source term: j_source = j_bs + j_cd + j_ext Boundary conditions: axis (ρ=0) — regularity (L'Hôpital); edge (ρ=1) — Dirichlet ψ=const.
Current Drive (v0.16.0)¶
Electron Cyclotron Current Drive.
Default η_cd: Fisch & Boozer 1980, PRL 45, 720 (theory); range 0.02–0.05 A/W depending on launch angle and T_e. For angle-resolved efficiency, use eccd_efficiency().
Neutral Beam Injection.
Current-drive efficiency formula: Ehst & Karney 1991, Nucl. Fusion 31, 1933. Slowing-down time: Stix 1972, Plasma Physics 14, 367, Eq. 6.
Beam heating power per unit rho [W] using a grid-normalised finite-width deposition kernel.
Beam-driven current density [A/m^2].
Fast-ion density from steady-state balance: n_fast = P_heat · τ_s / E_beam τ_s from Stix 1972, Eq. 6 (nbi_slowing_down_time). j_cd = e · n_fast · v_∥ / Z_beam Ehst & Karney 1991, Nucl. Fusion 31, 1933.
Note: Ti_keV enters via Z_eff-dependent collisionality; here held constant.
Superposition of multiple current-drive sources.
I_cd = ∫ j_cd · 2π r dr [A]
r = ρ · a, dr = dρ · a → dA = 2π ρ a² dρ
NTM Dynamics (v0.16.0)¶
ECCD-based NTM controller: triggers at onset, deactivates at target.
Sawtooth Model (v0.16.0)¶
Tracks and triggers sawtooth crashes.
Apply Kadomtsev full-reconnection crash.
Returns (T_new, n_new, q_new, rho_1, rho_mix).
Helical flux proxy: dψ/dρ = ρ(1/q − 1), integrated from axis. rho_mix is where ψ returns to zero outside the q=1 surface. Inside rho_mix, T, n are volume-averaged and q is reset to 1.01.
Kadomtsev 1975, Sov. J. Plasma Phys. 1, 389. Porcelli et al. 1996, PPCF 38, 2163, §2.
SOL Model (v0.16.0)¶
Two-point model for the Scrape-Off Layer.
Parallel heat conduction: q_∥ = κ_0 T_u^{7/2} / (7 L_∥). Stangeby 2000, Eq. 5.69.
Solve two-point model for upstream and target conditions.
Upstream temperature from conduction integral (Stangeby 2000, Eq. 5.69): T_u = (3.5 · L_∥ · q_∥ / κ_0)^{2/7}
Target heat flux: q_target = P_SOL / (4π R λ_q f_exp). Stangeby 2000, Ch. 5.
Integrated Scenario (v0.16.0)¶
Time-dependent tokamak shot simulator.
Operator splitting sequence per step() (Strang 1968): a) half-step transport diffusion b) full-step sources c) full-step MHD events d) half-step transport diffusion
State variables carried on 1D rho grid: Te, Ti, ne, psi (→ q, j).
Advance one time step using Strang operator splitting.
Strang 1968, SIAM J. Numer. Anal. 5, 506: a) half-step transport b) full-step sources c) full-step MHD events d) half-step transport
audit_scenario_coupling(states, config, *, scenario_name='integrated_scenario', current_relative_tolerance=1e-09, energy_relative_step_limit=5.0)
Audit a scenario replay for finite, bounded, traceable module coupling.
The audit is deliberately a bounded controller-facing contract. It proves replay consistency, finite state exchange, declared timestep semantics, and conservation-bound diagnostics for a completed trajectory. It is not a measured-discharge or external integrated-modelling validation claim.
Persist a scenario-coupling audit report with deterministic JSON layout.
SCPN — Petri Net Compiler¶
StochasticPetriNet¶
Stochastic Petri Net with sparse matrix representation.
Usage::
net = StochasticPetriNet()
net.add_place("P_red", initial_tokens=1.0)
net.add_place("P_green", initial_tokens=0.0)
net.add_transition("T_r2g", threshold=0.5)
net.add_arc("P_red", "T_r2g", weight=1.0)
net.add_arc("T_r2g", "P_green", weight=1.0)
net.compile()
Add a place (state variable) with token density in [0, 1].
Add a transition (logic gate) with a firing threshold.
Add a directed arc between a Place and a Transition (either direction).
Valid arcs
Place -> Transition (input arc, stored in W_in) Transition -> Place (output arc, stored in W_out)
Parameters¶
inhibitor : bool, default False
If True, arc is encoded as a negative weight inhibitor input arc.
Only valid for Place -> Transition edges.
Raises ValueError for same-kind connections or unknown nodes.
Build sparse W_in and W_out matrices from the arc list.
Parameters¶
validate_topology : bool, default False
Compute and store topology diagnostics in
:attr:last_validation_report.
strict_validation : bool, default False
If True, raise ValueError when topology diagnostics contain
dead nodes or unseeded place cycles.
allow_inhibitor : bool, default False
Allow negative Place -> Transition weights for inhibitor arcs.
If False, compile fails when inhibitor arcs are present.
Return topology diagnostics without mutating compiled matrices.
Diagnostics:
- dead_places: places with zero total degree.
- dead_transitions: transitions with zero total degree.
- unseeded_place_cycles: place-level SCC cycles with no initial
token mass.
- input_weight_overflow_transitions: transitions whose positive
input arc-weight sum exceeds 1.0.
Return (n_places,) float64 vector of initial token densities.
Formal Verification¶
Bounded formal verifier for compiled StochasticPetriNet objects.
The current production backend is exact explicit-state reachability over
rational markings. backend='z3' is accepted only when the optional solver
package is installed; the same finite transition relation is then eligible
for SMT-backed proof extensions without changing public contracts.
Enumerate all markings reachable within max_depth firings.
Prove that all reachable markings satisfy per-place bounds.
Prove algebraic place invariants and bounded reachable consistency.
Prove bounded transition liveness by checking every transition fires on some path.
Verify bounded temporal safety/liveness specifications.
Verify bounded CTL formulas over the compiled Petri-net relation.
Verify bounded LTL formulas over all finite firing paths.
Run reachability, marking-bound, liveness, and temporal proof obligations.
dataclass
¶
Algebraic P-invariant over named places.
The verifier proves that sum(weights[p] * marking[p]) is unchanged by
every transition and then checks the bounded reachable graph against the
same invariant value. This is a structural proof obligation, not a sampled
simulation assertion.
dataclass
¶
SafetyCertificatePolicy(name, min_depth=0, require_artifact_binding=False, require_ctl=False, require_ltl=False, required_checked_specs=())
Admission policy for generated formal safety certificates.
classmethod
¶
Create a fail-closed policy requiring artifact binding plus CTL and LTL evidence.
dataclass
¶
SafetyCertificateBundlePolicy(name, min_certificates=1, required_policy_name=None, require_same_artifact=False, require_same_backend=False)
Admission policy for a bundle of independent safety certificates.
dataclass
¶
Bounded CTL formula over a compiled SCPN transition system.
Supported operators are the certification-relevant bounded subset that can
produce exact counterexamples over the finite firing-depth transition
relation: AG safety, EF reachability, and AG_EF recurrence.
classmethod
¶
Create AG marking-bound safety: all paths always satisfy bounds.
classmethod
¶
Create EF transition reachability: some path eventually fires.
classmethod
¶
Create AG mutual-exclusion safety for two marked places.
classmethod
¶
Create AG EF recurrence: every reachable state can recover marking.
dataclass
¶
Bounded LTL formula over all finite firing paths of a compiled SCPN.
dataclass
¶
Temporal safety spec: all reachable markings keep named places in bounds.
dataclass
¶
Bounded recurrence spec: every reachable path can reach a marked place.
dataclass
¶
Bounded liveness spec: a transition fires on at least one reachable path.
dataclass
¶
Bounded response spec: a transition firing is followed by a marked place.
dataclass
¶
Temporal safety spec: two places are never simultaneously marked above threshold.
build_safety_certificate_payload(report, *, ctl_report=None, ltl_report=None, artifact_sha256=None, issuer='scpn-control', policy=None)
Build a schema-versioned tamper-evident formal safety certificate.
Build a schema-versioned tamper-evident safety certificate bundle.
Build a schema-versioned reference to a persisted certificate bundle.
generate_safety_certificate(net, *, max_depth, marking_bounds, json_path, markdown_path, temporal_specs=None, ctl_specs=None, ltl_specs=None, artifact_path=None, artifact_sha256=None, issuer='scpn-control', backend='explicit-state', policy=None)
Run proof obligations and persist a bound formal safety certificate.
The workflow resolves one verifier backend, runs base Petri-net safety, liveness, optional CTL, and optional LTL obligations through that same transition relation, binds the certificate to optional controller artifact bytes, and writes JSON/Markdown evidence in one call.
Validate a schema-versioned formal safety certificate payload.
Validate a schema-versioned formal safety certificate bundle.
Validate a persisted certificate-bundle artifact reference.
Load and validate the bundle referenced by a certificate-bundle artifact.
write_safety_certificate(report, *, json_path, markdown_path, ctl_report=None, ltl_report=None, artifact_sha256=None, issuer='scpn-control', policy=None)
Persist a formal safety certificate as JSON and Markdown artifacts.
Persist a formal safety certificate bundle as JSON and Markdown artifacts.
Bounded SMT model checker for compiled StochasticPetriNet graphs.
Prove marking bounds by asking Z3 for a reachable counterexample.
Verify supported bounded temporal specifications with Z3.
Verify bounded CTL formulas with the Z3 transition relation.
Verify bounded LTL formulas with the Z3 transition relation.
verify_trigger_response_latency(*, trigger_transition, response_transition, max_latency_ns, tick_period_ns=10.0, max_depth, weight_bounds=None)
Verify a bounded trigger->response latency contract in transition steps.
The bound is interpreted as:
at most ceil(max_latency_ns / tick_period_ns) non-idle transition
steps from a trigger fire to the first subsequent response fire. Idle
stalls are intentionally discharged by the separate no-stall SymbiYosys
contract rather than conflated with trigger-response latency.
build_symbiyosys_contract(*, max_depth, weight_bounds=None, trigger_transition=None, response_transition=None, max_latency_ns=50.0, tick_period_ns=10.0, no_stall_window_ns=50.0)
Build a bounded symbolic contract for SymbiYosys/SMTBMC review.
The emitted contract uses assume bounds for hot-swappable transition
weights and assert obligations for bounded trigger->response latency and
bounded non-stall progression.
verify_z3_formal_contracts(net, *, max_depth, marking_bounds, temporal_specs=None, weight_bounds=None)
Run bounded Z3 safety and temporal proof obligations.
Return the schema-versioned JSON payload for a bounded Z3 report.
Return a schema-versioned blocked report for unavailable SMT evidence.
Validate a schema-versioned Z3 formal evidence report.
Load and validate a duplicate-key-safe Z3 formal evidence report.
Persist a bounded SMT report as JSON and Markdown evidence artifacts.
FusionCompiler¶
FusionCompiler(bitstream_length=1024, seed=42, *, lif_tau_mem=1000000.0, lif_noise_std=0.0, lif_dt=1.0, lif_resistance=1.0, lif_refractory_period=0)
Compiles a StochasticPetriNet into a CompiledNet.
Parameters¶
bitstream_length : Number of bits per stochastic stream (default 1024). seed : Base RNG seed for reproducibility.
classmethod
¶
with_reactor_lif_defaults(bitstream_length=1024, seed=42, *, lif_dt=1.0, lif_resistance=1.0, lif_refractory_period=0)
Create a compiler with reactor-like LIF defaults.
Keeps the legacy constructor defaults untouched for compatibility while exposing a realistic preset for new control experiments.
staticmethod
¶
traceable_runtime_kwargs(*, runtime_backend='auto', allow_runtime_backend_fallback=False, allow_legacy_runtime_backend_fallback=False)
Recommended controller kwargs for traceable runtime loops.
compile(net, firing_mode='binary', firing_margin=0.05, *, allow_inhibitor=False, validate_topology=False, strict_topology=False)
Compile the Petri Net into sc_neurocore artifacts.
Parameters¶
net : compiled StochasticPetriNet.
firing_mode : "binary" (default) or "fractional".
firing_margin : margin for fractional firing (ignored in binary mode).
allow_inhibitor : enable inhibitor arc compilation.
validate_topology : run topology diagnostics during compile.
strict_topology : raise if topology diagnostics detect issues.
Steps
- Extract dense W_in (nT x nP) and W_out (nP x nT).
- Create one LIF neuron per transition (pure threshold comparator).
- Pre-encode weight matrices as packed uint64 bitstreams.
- Return
CompiledNetwith all artifacts.
CompiledNet¶
dataclass
¶
CompiledNet(n_places, n_transitions, place_names, transition_names, W_in, W_out, W_in_packed=None, W_out_packed=None, neurons=list(), bitstream_length=1024, thresholds=(lambda: array([], dtype=float64))(), transition_delay_ticks=(lambda: array([], dtype=int64))(), initial_marking=(lambda: array([], dtype=float64))(), seed=42, firing_mode='binary', firing_margin=0.05, lif_tau_mem=1000000.0, lif_noise_std=0.0, lif_dt=1.0, lif_resistance=1.0, lif_refractory_period=0)
Compiled Petri Net ready for sc_neurocore execution.
Holds both the dense float matrices (for validation / fallback) and pre-packed uint64 weight bitstreams (for the stochastic path).
Run LIF threshold detection on all transitions.
Binary mode: f_t = 1 if current >= threshold else 0
Fractional mode: f_t = clip((current - threshold) / margin, 0, 1)
Parameters¶
currents : (n_transitions,) float64 — weighted-sum activations.
Returns¶
fired : (n_transitions,) float64 vector. Binary mode → values in {0.0, 1.0}. Fractional mode → values in [0.0, 1.0].
Build an Artifact from compiled state + user-provided config.
Parameters¶
name : artifact name.
dt_control_s : control tick period (s).
readout_config : dict with actions, gains, abs_max,
slew_per_s lists. Required for a complete artifact.
injection_config : list of place-injection dicts.
NeuroSymbolicController¶
NeuroSymbolicController rejects nonzero sc_bitflip_rate unless
allow_fault_injection=True is supplied explicitly and the process environment
sets SCPN_ALLOW_CONTROLLER_FAULT_INJECTION=1. Bit-flip mutation is a
double-gated fault-injection test mode, not a production control default.
Controller JSONL logging requires an explicit log_root whenever log_path is
provided. Relative and absolute log paths must resolve under that root and use a
.jsonl suffix before any file is opened. Log appends use a constrained append
helper that rejects symlink targets where the platform exposes no-follow open
semantics.
NeuroSymbolicController(artifact, seed_base, targets, scales, sc_n_passes=8, sc_bitflip_rate=0.0, allow_fault_injection=False, sc_binary_margin=None, sc_antithetic=True, enable_oracle_diagnostics=True, feature_axes=None, runtime_profile='adaptive', runtime_backend='auto', allow_runtime_backend_fallback=False, allow_legacy_runtime_backend_fallback=False, rust_backend_min_problem_size=1, sc_antithetic_chunk_size=2048)
Reference controller with oracle float and stochastic paths.
Parameters¶
artifact : loaded .scpnctl.json artifact.
seed_base : 64-bit base seed for deterministic stochastic execution.
targets : control setpoint targets.
scales : normalisation scales.
allow_fault_injection : explicit opt-in for stochastic bit-flip fault
injection. Nonzero sc_bitflip_rate is rejected unless this is true.
Execute one control tick.
Steps
extract_features(obs)→ 4 unipolar features_inject_places(features)_oracle_step()— float path (optional)_sc_step(k)— deterministic stochastic path_decode_actions()— gain × differencing, slew + abs clamp- Optional JSONL logging
Execute one control tick from a fixed-order observation vector.
The vector order is self._axis_obs_keys. This avoids per-step key
lookups/dict allocation in tight control loops.
Contracts¶
Bases: TypedDict
Plant observation at a single control tick.
Bases: TypedDict
Actuator command output for a single control tick.
Keys are dynamic and depend on the Petri Net readout configuration.
dataclass
¶
Setpoint targets for the control loop.
Map observation → unipolar [0, 1] feature sources.
By default returns keys x_R_pos, x_R_neg, x_Z_pos, x_Z_neg.
Custom feature mappings can be supplied via feature_axes for arbitrary
observation dictionaries.
Decode marking → actuator commands with gain, slew-rate, and abs clamp.
Parameters¶
marking : current marking vector (len >= max place index). actions_spec : per-action pos/neg place definitions. gains : per-action gain multiplier. abs_max : per-action absolute saturation. slew_per_s : per-action max change rate (units/s). dt : control tick period (s). prev : previous action outputs (same length as actions_spec).
Returns¶
dict mapping action name → clamped value. Also mutates prev in-place.
Artifacts¶
Safety-critical controller admission must call load_artifact(...,
require_formal_verification=True) or validate_safety_critical_artifact().
That gate rejects missing, blocked, failed, malformed, or unbounded proof
evidence and accepts only hash-addressed bounded proof manifests tied to the
compiled controller artifact. The proof manifest must include the canonical
artifact payload SHA-256, report SHA-256, bounded proof depth, checked
specification names, backend/solver metadata, and a safe relative report URI.
When callers provide formal_report_root, the loader resolves the report URI
under that root and verifies the report bytes against the declared SHA-256.
Z3 and Lean report files reached through an artifact manifest are loaded through
duplicate-key-safe and schema-strict public report loaders.
Z3 reports are additionally schema-versioned as
scpn-control.z3-formal-report.v1, carry a canonical payload SHA-256 over the
proof payload, reject unknown top-level and proof-section fields, schema-check
serialized counterexample records, enforce solver-status/holds/counterexample
consistency, reject counterexamples on unknown solver sections, and must
match the manifest status, solver, proof depth, and
checked specification list before a safety-critical artifact is admitted.
Each Z3 proof section must also carry unique non-empty checked_specs; duplicate
section obligations are rejected before top-level report/spec matching.
Blocked Z3 reports are not proof evidence: they must use the unavailable solver
label, zero proof depth, and only the z3_solver_available checked
specification. Pass/fail Z3 reports must identify z3-solver and must not
use the unavailable-solver label; a missing z3-solver dependency is
admitted only as blocked SMT evidence.
Lean 4 reports are admitted only through the bounded lean4 manifest path:
the manifest must bind a solver string that identifies Lean and includes the
declared Lean version, Lake file SHA-256, proof-source SHA-256, theorem names,
theorem modules, proved contracts, linked production module paths, safety-case
identifiers, checked specification list, report SHA-256, and compiled artifact
SHA-256. The Lean report schema is
scpn-control.lean4-formal-report.v1; when a report root is supplied, every
manifest field above must match the report before admission. The current
required Lean proof-contract surface covers PID actuator-saturation preservation
and SNN/neuro-symbolic marking-bound preservation, and reports cannot list
unsupported proved contracts to imply wider formal coverage. This is an evidence
admission contract; it does not claim certification unless the referenced
machine-checked proof artefacts are present and verified.
The admitted PID/SNN surface is also exact-linked: theorem modules, theorem
names, linked production module paths, and safety-case identifiers must remain
inside the expected PID/SNN namespaces and module paths instead of padding a
valid report with unrelated evidence.
dataclass
¶
Full SCPN controller artifact (.scpnctl.json).
dataclass
¶
FormalVerificationEvidence(required, status, backend, solver, max_depth, checked_specs, artifact_sha256, report_sha256, claim_boundary, report_uri=None, generated_utc=None, counterexample_path=None, counterexample_property=None, lean_version=None, lakefile_sha256=None, proof_source_sha256=None, theorem_names=None, theorem_modules=None, proved_contracts=None, module_paths=None, safety_case_ids=None, proof_assumptions=None, assumption_sha256=None)
Hashable bounded-proof admission evidence for controller artifacts.
Compute a canonical SHA-256 digest over artifact payload sections.
Serialize an Artifact to indented JSON.
Parse a .scpnctl.json file into an Artifact dataclass.
Fail closed unless a controller artifact carries passing bounded-proof evidence.
dataclass
¶
LeanFormalVerificationReport(status, solver, lean_version, checked_specs, artifact_sha256, proof_source_sha256, lakefile_sha256, theorem_names, theorem_modules, proved_contracts, module_paths, safety_case_ids, claim_boundary, proof_assumptions)
Machine-checkable Lean 4 report metadata for bounded controller claims.
Build a canonical Lean formal-verification report payload.
Validate a Lean formal-verification report payload.
Write a canonical Lean formal-verification report and return the payload.
Load and validate a Lean formal-verification report payload.
load_lean_formal_report() and the validation executable
validation/validate_scpn_lean_formal.py validate Lean report JSON with
duplicate-key rejection and can additionally call
load_artifact(..., require_formal_verification=True) against a supplied
artifact and report root. This is the release-gate path for admitting Lean proof
evidence without running long proof jobs inside the Python test suite.
The Lean report and artifact manifest validators also enforce namespace
coverage for required controller contracts: PID actuator-saturation evidence
must include a ScpnControl.PID theorem module and theorem name, while
SNN/neuro-symbolic marking-bound evidence must include a ScpnControl.SNN
theorem module and theorem name.
Lean evidence must also declare explicit bounded proof_assumptions and a
canonical assumption_sha256; report and artifact admission reject unbounded
assumptions, certification overclaims, malformed assumption digests, or
manifest/report assumption mismatches.
Admission also rejects non-Lean solver declarations, Lean solver strings that
do not include the reported lean_version, and any proved_contracts outside
the currently admitted PID/SNN proof surface.
Reports also fail closed when theorem_names, theorem_modules,
module_paths, or safety_case_ids include unrelated entries outside the
admitted PID/SNN proof boundary.
The same checks run on the .scpnctl artifact manifest itself, even without a
report root, so safety-critical artifact loading cannot admit a stale or padded
Lean manifest before report-byte comparison is available.
Both Lean report payloads and artifact formal_verification manifests are
closed schemas: unknown proof fields are rejected rather than ignored.
External Lean reports must also carry the canonical payload_sha256 self-digest;
reports that omit it are rejected before safety-critical artifact admission.
Phase — Paper 27 Dynamics¶
Kuramoto-Sakaguchi Step¶
kuramoto_sakaguchi_step(theta, omega, *, dt, K, alpha=0.0, zeta=0.0, psi_driver=None, psi_mode='external', wrap=True)
Single Euler step of mean-field Kuramoto-Sakaguchi + global driver.
dθ_i/dt = ω_i + K·R·sin(ψ_r − θ_i − α) + ζ·sin(Ψ − θ_i)
Uses Rust backend when available (rayon-parallelised, sub-ms for N>1000).
Kuramoto order parameter R·exp(i·ψ_r) =
Returns (R, ψ_r).
Lyapunov candidate V(t) = (1/N) Σ (1 − cos(θ_i − Ψ)).
V=0 at perfect sync (all θ_i = Ψ), V=2 at maximal desync. Range: [0, 2]. Mirror of control-math/kuramoto.rs::lyapunov_v.
λ = (1/T) · ln(V_final / V_initial). λ < 0 ⟹ stable.
dataclass
¶
Resolve the global field phase Ψ.
mode="external" : Ψ supplied by caller (intention/carrier, no dotΨ).
mode="mean_field" : Ψ = arg(
dataclass
¶
KuramotoRuntimeEvidence(schema_version, generated_utc, target_id, oscillator_count, deployment_target_oscillators, dt_s, K, alpha, zeta, psi_mode, wrap, input_sha256, python_reference_sha256, rust_available, rust_parity_checked, rust_theta_max_abs_error_rad, rust_order_parameter_abs_error, parity_tolerance, parity_passed, timestep_refinement_checked, timestep_refinement_error_rad, timestep_refinement_tolerance, timestep_refinement_passed, deployment_claim_allowed, claim_status, payload_sha256)
Tamper-evident phase-runtime parity and timestep-refinement evidence.
kuramoto_runtime_evidence(theta, omega, *, dt, K, alpha=0.0, zeta=0.0, psi_driver=None, psi_mode='external', wrap=True, target_id='local-python-runtime', deployment_target_oscillators=None, parity_tolerance=1e-10, timestep_refinement_tolerance=0.005, generated_utc=None, deployment_claim_allowed=False)
Build bounded phase-runtime parity and timestep-refinement evidence.
Fail closed unless phase-runtime evidence supports deployment claims.
Persist Kuramoto runtime evidence as sorted JSON.
Load Kuramoto runtime evidence with duplicate-key and digest admission.
Knm Coupling Matrix¶
dataclass
¶
Paper 27 coupling specification.
K : (L, L) coupling matrix. K[n, m] = source n -> target m. alpha : (L, L) Sakaguchi phase-lag (optional). zeta : (L,) per-layer global-driver gain ζ_m (optional).
Build the canonical Paper 27 Knm with exponential distance decay.
K[i,j] = K_base · exp(−K_alpha · |i − j|), diag(K) kept for intra-layer sync (unlike the inter-oscillator Knm which zeros diag).
K_base=0.45 and K_alpha=0.3 from Paper 27 §3.2, Eq. 12. Calibration anchors from Paper 27, Table 2. Cross-hierarchy boosts from Paper 27 §4.3.
UPDE Multi-Layer Solver¶
dataclass
¶
Multi-layer UPDE driven by a KnmSpec.
step(theta_layers, omega_layers, *, psi_driver=None, actuation_gain=1.0, pac_gamma=0.0, K_override=None)
Advance all L layers by one Euler step.
Parameters¶
theta_layers : sequence of 1D arrays Phase vectors per layer. omega_layers : sequence of 1D arrays Natural frequencies per layer. psi_driver : float or None External global field phase Ψ (required if psi_mode="external"). actuation_gain : float Multiplicative gain on all coupling terms. pac_gamma : float PAC-like gating: boost inter-layer coupling by (1 + pac_gamma·(1 − R_source)). K_override : array or None Per-tick replacement for spec.K (adaptive coupling).
run(n_steps, theta_layers, omega_layers, *, psi_driver=None, actuation_gain=1.0, pac_gamma=0.0, K_override=None)
Run n_steps and return trajectory of per-layer R and global R.
run_lyapunov(n_steps, theta_layers, omega_layers, *, psi_driver=None, actuation_gain=1.0, pac_gamma=0.0, K_override=None)
Run n_steps with Lyapunov tracking.
Returns R histories, V histories, and per-layer + global λ. λ < 0 ⟹ stable convergence toward Ψ.
Lyapunov Guard¶
Sliding-window Lyapunov stability monitor.
Parameters¶
window : int Number of V(t) samples in the sliding window. dt : float Timestep between samples (for λ computation). lambda_threshold : float λ above this value counts as a violation. Default 0 (any growth). max_violations : int Consecutive violations before guard refuses.
Realtime Monitor¶
dataclass
¶
RealtimeMonitor(upde, guard, theta_layers, omega_layers, psi_driver=0.0, pac_gamma=0.0, adaptive_engine=None)
Tick-by-tick UPDE monitor with LyapunovGuard.
classmethod
¶
from_paper27(L=16, N_per=50, *, dt=0.001, zeta_uniform=0.5, psi_driver=0.0, pac_gamma=0.0, guard_window=50, guard_max_violations=3, seed=42)
Build from Paper 27 defaults.
classmethod
¶
from_plasma(L=8, N_per=50, *, mode='baseline', dt=0.001, zeta_uniform=0.0, psi_driver=0.0, pac_gamma=0.0, guard_window=50, guard_max_violations=3, adaptive_engine=None, seed=42)
Build from plasma-native Knm defaults with optional adaptive engine.
Advance one UPDE step and return dashboard snapshot.
Export recorded trajectory to HDF5.
Datasets: R_global, R_layer, V_global, V_layer, lambda_exp, guard_approved, latency_us, Psi_global. Attributes: L, N_per, psi_driver, pac_gamma, n_ticks.
Requires h5py (pip install h5py).
dataclass
¶
TrajectoryRecorder(R_global=list(), R_layer=list(), V_global=list(), V_layer=list(), lambda_exp=list(), guard_approved=list(), latency_us=list(), Psi_global=list())
Accumulates tick snapshots for batch export.
Adaptive Knm Engine¶
Diagnostic-driven online adaptation of the Knm coupling matrix.
Holds a baseline K from a KnmSpec and produces K_adapted each tick.
dataclass
¶
AdaptiveKnmConfig(beta_scale=0.3, beta_max_boost=0.5, risk_pairs=((2, 5), (3, 5), (2, 4)), risk_gain=0.4, coherence_Kp=0.15, coherence_Ki=0.02, coherence_R_target=0.6, coherence_max_boost=0.3, max_delta_per_tick=0.02, revert_on_guard_refusal=True)
Tuning knobs for each adaptation channel.
dataclass
¶
DiagnosticSnapshot(R_layer, V_layer, lambda_exp, beta_n, q95, disruption_risk, mirnov_rms, guard_approved)
Per-tick plasma diagnostic bundle.
Plasma Knm¶
build_knm_plasma(mode='baseline', L=8, K_base=0.3, zeta_uniform=0.0, custom_overrides=None, layer_names=None)
Build a plasma-native Knm coupling matrix.
Parameters¶
mode : str Instability scenario bias. One of: baseline, elm, ntm, sawtooth, hybrid. L : int Number of plasma layers (default 8). K_base : float Base coupling strength for exponential decay backbone. zeta_uniform : float Global driver gain ζ applied uniformly to all layers. custom_overrides : dict, optional Explicit (i, j) → value overrides applied last. Automatically symmetrised. layer_names : sequence of str, optional Layer labels. Defaults to PLASMA_LAYER_NAMES[:L] or PLASMA_LAYER_NAMES_16[:L].
Returns¶
KnmSpec Ready for UPDESystem consumption.
WebSocket Stream¶
PhaseStreamServer binds to loopback by default and requires authenticated
clients by default. Operators must supply an API key or explicitly disable
client authentication for local development. Non-loopback binds require an API
key, command frames are capped by max_payload_bytes, accepted commands are
rate-limited with token buckets per connection and per network peer, and
production remote exposure should enable TLS with require_tls=True.
Authentication, origin, payload, rate-limit, capacity, and command-authority
rejections emit structured security audit log events without logging tokens.
Query-string token authentication and plaintext
non-loopback binds are disabled by default and require explicit operator
opt-ins for constrained development or isolated lab environments. Browser
clients that send an Origin header are rejected unless the origin is
allowlisted, and deployments may restrict command authority with
allowed_actions.
websocket_runtime_evidence() emits a tamper-evident deployment artifact that
binds WebSocket configuration, authenticated sessions, accepted commands,
successful broadcasts, audit counters, TLS enforcement, payload caps, and
backpressure state without storing API-key material. Qualified facility
admission requires client authentication, configured TLS enforcement, observed
commands, observed broadcasts, no query-token authentication, no insecure remote
binding, and zero backpressure disconnects.
dataclass
¶
PhaseStreamServer(monitor, tick_interval_s=0.001, api_key=None, command_rate_limit=20, command_rate_window_s=1.0, max_payload_bytes=65536, client_send_timeout_s=0.25, max_clients=128, max_client_write_buffer_bytes=262144, require_client_auth=True, allow_query_token_auth=False, require_tls=False, allow_insecure_remote=False, allowed_origins=(), allowed_actions=_DEFAULT_ALLOWED_ACTIONS)
dataclass
¶
WebSocketRuntimeEvidence(schema_version, generated_utc, deployment_id, bind_host, uses_tls, require_client_auth, api_key_configured, api_key_min_bytes, require_tls, allow_query_token_auth, allow_insecure_remote, command_rate_limit, command_rate_window_s, max_payload_bytes, client_send_timeout_s, max_clients, max_client_write_buffer_bytes, allowed_actions, allowed_origins, config_sha256, auth_successes, command_frames, broadcast_frames, peak_connected_clients, authentication_failures, origin_rejections, payload_rejections, rate_limit_rejections, capacity_rejections, action_rejections, backpressure_disconnects, facility_claim_allowed, claim_status, payload_sha256)
Tamper-evident WebSocket runtime security and backpressure evidence.
websocket_runtime_evidence(server, *, deployment_id, bind_host='127.0.0.1', uses_tls=False, generated_utc=None, counters=None, facility_claim_allowed=False)
Build tamper-evident runtime evidence for a phase WebSocket deployment.
Fail closed unless WebSocket evidence supports a facility runtime claim.
Persist WebSocket runtime evidence as sorted JSON without secret material.
Load WebSocket runtime evidence with duplicate-key and digest admission.
Control — Controllers¶
H-infinity (Riccati DARE)¶
HInfinityController(A, B1, B2, C1, C2, gamma=None, D12=None, D21=None, enforce_robust_feasibility=False)
Riccati-based H-infinity controller for tokamak vertical stability.
Synthesis: Doyle et al. 1989, IEEE TAC 34, 831 (state-space ARE formulation). Riccati sign convention: Zhou, Doyle & Glover 1996, Eq. 14.18. γ-iteration: Green & Limebeer 1995, Ch. 13 (binary search on ||T_{zw}||_∞). Weighting selection: Ariola & Pironti 2008, Ch. 5 (tokamak vertical loop).
Parameters¶
A : array_like, shape (n, n) Plant state matrix. B1 : array_like, shape (n, p) Disturbance input matrix. B2 : array_like, shape (n, m) Control input matrix. C1 : array_like, shape (q, n) Performance output matrix. C2 : array_like, shape (l, n) Measurement output matrix. gamma : float, optional H-infinity attenuation level. If None, found by bisection. D12 : array_like, optional Feedthrough from control to performance. Default: identity-like. D21 : array_like, optional Feedthrough from disturbance to measurement. Default: identity-like. enforce_robust_feasibility : bool, optional If True, raise ValueError unless rho(XY) < gamma^2 after synthesis.
property
¶
True iff all eigenvalues of A + B2 F lie in the open left-half plane.
property
¶
Eigenvalue-based stability margin in dB.
Not the classical Bode gain margin. Measures the ratio of closed-loop to open-loop dominant eigenvalue real parts. For frequency-domain gain margin, use a Bode analysis on the loop transfer function.
Frobenius norms of the two H-infinity ARE residuals.
Residual form: X A + A^T X + Q - X B R^{-1} B^T X (Zhou, Doyle & Glover 1996, Eq. 14.18). Small residuals (< 1e-6 relative) confirm solver accuracy.
Return gamma^2 - rho(XY); positive values satisfy the strict condition.
Strict condition: rho(XY) < γ^2 — Doyle et al. 1989, §IV, Theorem 1.
H-infinity controller for tokamak vertical stability.
Plant: second-order vertical instability model from Ariola & Pironti 2008, Ch. 5, linearised about an unstable equilibrium:
A = [[0, 1 ],
[γ_v², -damping ]]
where γ_v is the vertical growth rate. The double integrator structure captures the leading-order Shafranov instability at low plasma current.
Parameters¶
gamma_growth : float Vertical instability growth rate [1/s]. ITER-like: ~100/s. SPARC: ~1000/s (Ariola & Pironti 2008, Ch. 5). damping : float Passive damping coefficient [1/s]. Resistive-wall contribution; 10.0 is representative for ITER-like wall proximity (Ariola & Pironti 2008). enforce_robust_feasibility : bool, optional If True, require rho(XY) < gamma^2 and raise on infeasible synthesis.
Returns¶
HInfinityController Riccati-synthesized robust controller.
Model Predictive Control¶
Linearized surrogate model around current operating point.
ModelPredictiveController(surrogate, target_state, *, prediction_horizon=PREDICTION_HORIZON, learning_rate=0.5, iterations=20, action_limit=2.0, action_regularization=0.1)
Gradient-based MPC planner over surrogate dynamics.
Optimal Control¶
OptimalController(config_file, *, kernel_factory=FusionKernel, verbose=True, correction_limit=5.0, coil_current_limits=(-40.0, 40.0), current_target_limits=(5.0, 16.0))
Digital Twin¶
run_digital_twin() now supports persistent sensor calibration bias and drift
in addition to dropout and white-noise corruption, and it can now stress the
command path with deterministic actuator bias, drift, first-order lag, and
rate limiting. The returned summary exposes both commanded and applied actions
plus actuator-lag telemetry so replay tests can see what the plant actually
received. Density and effective-charge knobs are explicit model-update
parameters.
digital_twin_online_update adds fail-closed TRANSP/TSC simulator artifact
metadata validation and deterministic Bayesian optimisation over bounded twin
parameters. The shipped benchmark is synthetic online-update evidence only;
external simulator replay claims require validated artifact metadata and the
strict digital-twin reference gate.
digital_twin_update_evidence() and
assert_digital_twin_update_claim_admissible() bind a bounded Bayesian update
to TRANSP and TSC simulator metadata digests, observation and prior digests,
result digest, baseline-improvement evidence, and an optional
safety-critical controller proof artifact digest. Admission also revalidates
source binding, finite non-negative loss history, minimum-loss consistency,
best-parameter bounds, strict integer campaign settings, and simulator unit
coverage for every observation target.
run_digital_twin(time_steps=TIME_STEPS, seed=42, save_plot=True, output_path='Tokamak_Digital_Twin.png', verbose=True, gyro_surrogate=None, chaos_monkey=False, sensor_dropout_prob=0.0, sensor_noise_std=0.0, sensor_bias_std=0.0, sensor_drift_std=0.0, sensor_thermal_lag_tau=0.05, dt=0.001, actuator_bias=0.0, actuator_drift_std=0.0, actuator_tau_steps=0.0, actuator_rate_limit=2.0, n_e=1e+20, Z_eff=1.5, rng=None)
Run deterministic digital-twin control simulation.
Returns a summary dict so callers can use the simulation without relying on console text or plot artifacts.
Validate TRANSP/TSC-style simulator artefact metadata before coupling.
Run deterministic Bayesian optimisation for online twin calibration.
dataclass
¶
DigitalTwinUpdateEvidence(schema_version, simulator_codes, external_artifacts_sha256, observation_sha256, priors_sha256, result_sha256, best_loss, baseline_loss, improved_over_baseline, evaluated_points, controller_formal_artifact_sha256, claim_status)
Tamper-evident admission evidence for bounded online twin updates.
digital_twin_update_evidence(observation, priors, result, external_artifacts, *, controller_formal_artifact_sha256=None)
Build tamper-evident evidence for bounded online digital-twin updates.
assert_digital_twin_update_claim_admissible(evidence, observation, priors, result, external_artifacts, *, require_simulators=SUPPORTED_EXTERNAL_SIMULATORS)
Fail closed unless online twin-update evidence matches replay inputs.
Deterministic local benchmark for the online-update contract.
Controller Safety-Case Evidence¶
control.safety_case defines the bounded safety-case workflow contract that
links a passing safety-critical controller proof manifest, audited JAX
differentiable-transport evidence, and TRANSP/TSC-backed bounded digital-twin
online-update evidence. The bundle is tamper-evident and fails closed unless
all evidence items bind to the same canonical controller artifact SHA-256. This
is a repository safety-package admission boundary, not an external
certification claim. save_controller_safety_case_evidence() persists the
bundle with a manifest integrity digest, and
load_controller_safety_case_evidence() rejects schema drift, malformed fields,
or edited evidence payloads before replay admission.
evaluate_controller_safety_case_readiness() separates linked bounded evidence
from promotion readiness: external physics validation, target-hardware timing
evidence, qualified HIL replay evidence, qualified CODAC/EPICS runtime
evidence, qualified WebSocket runtime evidence, qualified HDL export evidence,
and independent safety-review digests are all required before
assert_controller_safety_case_readiness_admissible() accepts the package.
ReadinessArtifactEvidence and
evaluate_controller_safety_case_readiness_from_artifacts() provide the normal
promotion path: each required readiness input must be a typed artifact with a
known kind, SHA-256 digest, safe relative artifact URI, producer, and generation
timestamp before it can satisfy the promotion gate. The evaluator also requires
an explicit artifact_root: each URI must resolve below that root and match the
declared bytes. target_hardware_timing artifacts must additionally pass the
schema-versioned E2E latency evidence validator with qualified target hardware
and the configured p95 limit. hil_replay_evidence artifacts must pass the
schema-versioned HIL replay admission loader with qualified target hardware, so
local workstation replay cannot satisfy deployment promotion readiness.
codac_runtime_evidence artifacts must pass the schema-versioned CODAC runtime
admission loader with qualified facility-claim status, deadline-clean cycle
evidence, exercised interlock blocking, zero backpressure events, and
hash-bound EPICS/OPC-UA exports.
websocket_runtime_evidence artifacts must pass the schema-versioned WebSocket
runtime admission loader with authenticated command evidence, TLS enforcement,
token-bucket and payload-cap configuration, successful broadcast counters, and
zero backpressure disconnects.
hdl_export_evidence artifacts must pass the schema-versioned FPGA export
admission loader with controller-artifact binding, generated project file
digests, synthesis-report digest binding, and non-negative timing slack.
save_controller_safety_case_readiness() and
load_controller_safety_case_readiness() persist that readiness decision with
the same schema-versioned integrity-digest semantics as the safety-case bundle.
dataclass
¶
ControllerSafetyCaseEvidence(schema_version, controller_artifact_sha256, formal_report_sha256, formal_backend, formal_max_depth, transport_evidence_sha256, digital_twin_evidence_sha256, claim_status)
Tamper-evident bounded controller safety-case evidence bundle.
dataclass
¶
SafetyCaseReadinessEvidence(schema_version, safety_case_sha256, status, external_physics_validation_sha256, target_hardware_timing_sha256, hil_replay_evidence_sha256, hdl_export_evidence_sha256, codac_runtime_evidence_sha256, websocket_runtime_evidence_sha256, independent_safety_review_sha256, blocking_reasons, claim_status)
Promotion-readiness gate for a bounded controller safety-case bundle.
dataclass
¶
Typed artifact evidence required for safety-case promotion readiness.
Build a bounded safety-case evidence bundle for one controller artifact.
assert_controller_safety_case_admissible(evidence, controller_artifact, transport_evidence, digital_twin_evidence)
Fail closed unless the safety-case evidence matches all supplied inputs.
Persist controller safety-case evidence with an integrity digest.
Load controller safety-case evidence and verify manifest integrity.
evaluate_controller_safety_case_readiness(safety_case, *, external_physics_validation_sha256=None, target_hardware_timing_sha256=None, hil_replay_evidence_sha256=None, hdl_export_evidence_sha256=None, codac_runtime_evidence_sha256=None, websocket_runtime_evidence_sha256=None, independent_safety_review_sha256=None)
Evaluate whether a bounded safety-case bundle is promotion-ready.
The linked internal evidence chain is necessary but not sufficient for promotion readiness. This gate requires external physics validation, target-hardware timing evidence, HIL replay evidence, HDL export evidence, CODAC/EPICS runtime evidence, WebSocket runtime evidence, and an independent safety review digest.
evaluate_controller_safety_case_readiness_from_artifacts(safety_case, artifacts, *, artifact_root, max_target_hardware_e2e_p95_us=1000.0)
Evaluate promotion readiness from typed external evidence artifacts.
Fail closed unless readiness evidence matches the safety case and is complete.
Persist controller safety-case readiness with an integrity digest.
Load controller safety-case readiness and verify manifest integrity.
Flight Simulator¶
IsoFluxController(config_file, kernel_factory=FusionKernel, verbose=True, actuator_tau_s=0.06, heating_actuator_tau_s=None, actuator_current_delta_limit=1000000000.0, heating_beta_max=5.0, control_dt_s=0.05)
Simulates the Plasma Control System (PCS). Uses PID loops to adjust Coil Currents to maintain plasma shape.
run_flight_sim(config_file=None, shot_duration=SHOT_DURATION, seed=42, save_plot=True, output_path='Tokamak_Flight_Report.png', verbose=True, actuator_tau_s=0.06, heating_actuator_tau_s=None, actuator_current_delta_limit=1000000000.0, heating_beta_max=5.0, control_dt_s=0.05, kernel_factory=FusionKernel)
Run deterministic tokamak flight-sim control loop and return summary.
Free-Boundary Tracking¶
Experimental closed-loop free-boundary tracking that keeps the full
FusionKernel in the loop and re-identifies the local coil-response map from
repeated solves. Safe-current fallback targets can be supplied through the
free_boundary_tracking.fallback_currents config block when supervisor
rejection should ramp the coils toward a predefined safe state. Persistent
objective residuals can also be accumulated with the config-driven
free_boundary_tracking.observer_gain and observer_max_abs settings. When
free-boundary objective tolerances are configured, the controller also uses
them directly in its correction and accept/reject logic so tighter X-point or
divertor targets take precedence over looser shape goals, and it refuses trial
steps that would push an already-satisfied objective back outside tolerance. If
the identified coil-response map loses authority entirely, the controller marks
that degraded state explicitly and drops into the safe-state recovery path
instead of silently accepting a zero-action step. Residuals already inside the
configured tolerances are also treated as deadband, so the controller stops
chattering the coils once the protected objectives are met. Coil allocation is
also headroom-aware, so the regularized solve prefers actuators that still have
current authority instead of leaning equally on a nearly saturated coil.
Deterministic objective-space sensor bias and per-step drift can be injected
through free_boundary_tracking.measurement_bias and
measurement_drift_per_step, and known calibration corrections can be applied
with measurement_correction_bias and measurement_correction_drift_per_step.
The run summary reports both measured and hidden true objective errors so
calibration faults cannot masquerade as control success in acceptance tests.
from scpn_control.control.free_boundary_tracking import run_free_boundary_tracking
summary = run_free_boundary_tracking(
"iter_config.json",
shot_steps=5,
gain=0.8,
verbose=False,
coil_slew_limits=2.5e5,
supervisor_limits={"x_point_position": 0.15, "max_abs_actuator_lag": 1.0e5},
hold_steps_after_reject=2,
)
print(summary["shape_rms"], summary["objective_converged"], summary["supervisor_intervention_count"])
FreeBoundaryTrackingController(config_file, *, kernel_factory=FusionKernel, verbose=True, identification_perturbation=0.25, correction_limit=2.0, response_regularization=0.001, response_refresh_steps=1, solve_max_outer_iter=10, solve_tol=0.0001, objective_tolerances=None, control_dt_s=None, coil_actuator_tau_s=None, coil_slew_limits=None, supervisor_limits=None, hold_steps_after_reject=None, state_estimator=None)
Direct free-boundary controller using local coil-response identification.
This path keeps the full Grad-Shafranov kernel in the loop instead of replacing it with a reduced-order plant.
Examples¶
from scpn_control.control.free_boundary_tracking import run_free_boundary_tracking summary = run_free_boundary_tracking( ... "iter_config.json", ... shot_steps=3, ... gain=0.8, ... verbose=False, ... ) summary["boundary_variant"] 'free_boundary'
run_free_boundary_tracking(config_file=None, *, shot_steps=10, gain=1.0, verbose=True, kernel_factory=FusionKernel, objective_tolerances=None, control_dt_s=None, coil_actuator_tau_s=None, coil_slew_limits=None, supervisor_limits=None, hold_steps_after_reject=None, disturbance_callback=None, stop_on_convergence=False)
Run deterministic free-boundary tracking over the configured objectives.
Examples¶
summary = run_free_boundary_tracking("iter_config.json", shot_steps=2, gain=0.7, verbose=False) bool(summary["objective_convergence_active"]) True
Disruption Predictor¶
predict_disruption_risk_safe() still returns a bounded scalar risk, but its
metadata now includes deterministic sigma-point uncertainty summaries
(risk_p05, risk_p50, risk_p95, risk_std, risk_interval) for both
fallback and checkpoint inference paths.
Lightweight deterministic disruption risk estimator returning a value in [0, 1].
Supplements the Transformer pathway by explicitly consuming toroidal asymmetry observables from 3D diagnostics (n=1,2,3 mode amplitudes).
Feature weights tuned on synthetic DIII-D/JET validation shots; see validation/reports/disruption_replay_pipeline_benchmark.md. Logit bias: sigmoid(−4.0) ≈ 0.018, giving low base risk on zero features.
predict_disruption_risk_safe(signal, toroidal_observables=None, *, model_path=None, seq_len=DEFAULT_SEQ_LEN, train_if_missing=False, mc_samples=10, allow_legacy_fallback=False, allow_inference_fallback=False, allow_legacy_inference_fallback=False)
Predict disruption risk with MC dropout uncertainty if model is available.
Returns¶
risk, metadata
risk is the MC mean risk in [0, 1].
metadata includes risk_std (epistemic uncertainty).
Disruption Contracts¶
Lightweight deterministic disruption risk estimator returning a value in [0, 1].
Supplements the Transformer pathway by explicitly consuming toroidal asymmetry observables from 3D diagnostics (n=1,2,3 mode amplitudes).
Feature weights tuned on synthetic DIII-D/JET validation shots; see validation/reports/disruption_replay_pipeline_benchmark.md. Logit bias: sigmoid(−4.0) ≈ 0.018, giving low base risk on zero features.
SPI Mitigation¶
Reduced SPI mitigation model for thermal/current quench campaigns.
SPI concept: Commaux et al. 2010, Nucl. Fusion 50, 112001.
staticmethod
¶
Estimate Z_eff for a noble-gas SPI cocktail.
Radiative efficiency weighting reflects that higher-Z species dominate coronal radiation. Lehnen et al. 2015, J. Nucl. Mater. 463, 39.
staticmethod
¶
Current quench time scale.
Scaling: τ_CQ ∝ T_e^0.25 / Z_eff, consistent with resistive decay. τ_TQ ~ 1–3 ms, τ_CQ ~ 20–150 ms for ITER. Hender et al. 2007, Nucl. Fusion 47, S128, Sec. 3.
run_spi_mitigation(*, plasma_energy_mj=300.0, plasma_current_ma=15.0, neon_quantity_mol=0.1, argon_quantity_mol=0.0, xenon_quantity_mol=0.0, duration_s=0.05, dt_s=1e-05, save_plot=True, output_path='SPI_Mitigation_Result.png', verbose=True)
Run SPI mitigation simulation and return deterministic summary metrics.
Fusion Control Room¶
run_control_room(sim_duration=SIM_DURATION, *, seed=42, save_animation=True, save_report=True, output_gif='SCPN_Fusion_Control_Room.gif', output_report='SCPN_Fusion_Status_Report.png', verbose=True, kernel_factory=None, config_file=None, allow_analytic_fallback=False, allow_legacy_analytic_fallback=False, allow_runtime_solve_fallback=False, allow_legacy_runtime_solve_fallback=False, allow_kernel_psi_fallback=False, allow_legacy_kernel_psi_fallback=False, allow_coil_update_fallback=False, allow_legacy_coil_update_fallback=False, allow_output_render_fallback=False, allow_legacy_output_render_fallback=False)
Run the control-room loop and return deterministic summary metrics.
Reduced Grad-Shafranov geometry model with optional kernel-Psi ingestion.
Gymnasium Environment¶
TokamakEnv(dt=0.001, max_steps=500, T_target=20.0, noise_std=0.01, seed=42, n_e_20=1.0, V_plasma=830.0)
Minimal Gymnasium-compatible tokamak control environment.
Implements a bounded reduced-order plasma response model for energy and current evolution based on 0D lumped-parameter equations. Ref: Wesson, J. (2011). Tokamaks. 4th Edition, Chapter 1.
Follows the gymnasium.Env interface (reset/step/render) without
requiring gymnasium as a hard dependency. If gymnasium is installed,
this class can be registered via gymnasium.register().
Parameters¶
dt : float Timestep per step call [s]. max_steps : int Episode length. T_target : float Target axis temperature [keV]. noise_std : float Observation noise standard deviation. n_e_20 : float Line-averaged electron density [10^20 m^-3]. Default 1.0. V_plasma : float Plasma volume [m^3]. Default 830.0 (ITER).
Analytic Solver¶
AnalyticEquilibriumSolver(config_path, *, kernel_factory=FusionKernel, verbose=True)
Analytic vertical-field target and least-norm coil-current solve.
Bio-Holonomic Controller¶
Biological Feedback Controller mapping L4/L5 states to clinical actions.
Rather than mitigating plasma disruptions, this controller mitigates biological decoherence by tracking autonomic tone and triggering resonant acoustic interventions.
Digital Twin Ingest¶
In-memory realtime ingest + SNN planning hook.
Director Interface¶
DirectorInterface(config_path, *, allow_fallback=False, allow_legacy_fallback=False, director=None, controller_factory=NeuroCyberneticController, entropy_threshold=0.3, history_window=10)
Interfaces the 'Director' (Layer 16: Coherence Oversight) with the Fusion Reactor.
Role: The Director does NOT control the coils (Layer 2 does that). The Director controls the Controller. It sets the strategy and monitors for "Backfire".
Mechanism: 1. Sample System State (Physics + Neural Activity). 2. Format as a "Prompt" for the Director. 3. Director calculates Entropy/Risk. 4. If Safe: Director updates Target Parameters. 5. If Unsafe: Director triggers corrective action.
Translate physical telemetry into a semantic prompt for the Director.
Fueling Mode Controller¶
Hybrid Petri-to-SNN + PI fueling controller.
Halo RE Physics¶
HaloCurrentModel(plasma_current_ma=15.0, minor_radius_m=2.0, major_radius_m=6.2, wall_resistivity_ohm_m=7e-07, wall_thickness_m=0.06, tpf=2.0, contact_fraction=0.3)
Fitzpatrick-style L/R circuit halo current model.
Parameters¶
plasma_current_ma : float Pre-disruption plasma current (MA). minor_radius_m : float Plasma minor radius (m). major_radius_m : float Plasma major radius (m). wall_resistivity_ohm_m : float Wall resistivity (Ohm·m). Default: stainless steel ~7e-7. wall_thickness_m : float Wall thickness (m). tpf : float Toroidal peaking factor (1.0 = uniform, up to ~2.5 in severe VDEs). contact_fraction : float Fraction of plasma cross-section in wall contact (0–1).
dataclass
¶
DisruptionMitigationClaimEvidence(schema_version, model_id, source, source_id, ensemble_runs, ensemble_seed, prevention_rate, mean_halo_peak_ma, p95_halo_peak_ma, mean_re_peak_ma, p95_re_peak_ma, mean_tpf_product, passes_iter_limits, reference_source, reference_dataset_id, reference_artifact_sha256, reference_case_count, risk_after_abs_error, detection_lead_time_abs_error_ms, halo_current_relative_error, runaway_beam_relative_error, tbr_abs_error, risk_after_abs_tolerance, detection_lead_time_abs_tolerance_ms, halo_current_relative_tolerance, runaway_beam_relative_tolerance, tbr_abs_tolerance, mitigation_claim_allowed, claim_status)
Serialisable admission evidence for halo/runaway mitigation claims.
disruption_mitigation_claim_evidence(report, *, source, source_id, ensemble_seed, model_id='halo_runaway_disruption_mitigation', reference_artifact_path=None)
Build fail-closed evidence for halo/runaway mitigation claims.
Return evidence only when strict matched-reference admission passed.
Persist disruption mitigation claim evidence as deterministic JSON.
HIL Test Harness¶
Real-time 1 kHz control loop with measured timing.
Executes a user-supplied control callback at a target rate (default 1 kHz)
and measures actual loop timing using time.perf_counter_ns.
Parameters¶
target_rate_hz : float Target control loop rate (Hz). Default: 1000 (1 kHz). sensor : SensorInterface or None Sensor/actuator interface. Created with defaults if None.
Execute the control loop and measure timing.
Parameters¶
iterations : int Number of control loop iterations. plant_fn : callable or None Plant model: state_new = plant_fn(state, command). Default: simple integrator (state += command * dt). initial_state : float Initial plant state. setpoint : float Control target.
dataclass
¶
HILBenchmarkResult(control_metrics, sensor_latency_us, controller_latency_us, actuator_latency_us, total_loop_latency_us, passes_sub_ms, passes_1khz, fpga_register_map)
Benchmark result for HIL sub-ms timing validation.
hil_replay_evidence(metrics, *, controller_id, target_hardware_id='local-unqualified-host', target_hardware_class='local-userspace-replay', rt_kernel='generic-userspace', clock_source='time.perf_counter_ns', interlock_events=0, backpressure_events=0, deployment_claim_allowed=False, generated_at=None)
Build tamper-evident HIL replay evidence with claim admission metadata.
Local workstation timing is useful regression evidence, but it is not
production deployment evidence. Setting deployment_claim_allowed=True
therefore fail-closes unless qualified target hardware metadata, zero
safety events, zero overruns, and target-rate latency bounds are present.
Validate a HIL replay evidence payload and return an admitted copy.
Persist an admitted HIL replay evidence payload as canonical JSON.
Load and admit a HIL replay evidence payload from JSON.
JAX Traceable Runtime¶
Requires pip install "scpn-control[jax]".
dataclass
¶
Configuration for reduced traceable first-order actuator dynamics.
LIF+NEF SNN Controller¶
SNN controller for tokamak position control.
Pure-NumPy LIF + NEF implementation. Per channel: error ensemble decodes gain·x into control ensemble, which decodes identity to output.
Parameters¶
config : NengoSNNConfig or None Controller configuration. Uses defaults if None.
Neuro-Cybernetic Controller¶
TORAX Hybrid Loop¶
Run deterministic NSTX-U-like realtime hybrid control campaign.
Advanced SOC Learning¶
run_advanced_learning_sim(size=L, time_steps=TIME_STEPS, seed=42, *, epsilon=EPSILON, noise_probability=0.01, shear_step=0.05, shear_bounds=(0.0, 1.0), save_plot=True, output_path='Advanced_SOC_Learning.png', verbose=True)
Run deterministic SOC+Q-learning control simulation and return summary metrics.
NMPC Controller (v0.16.0)¶
NonlinearMPC validates the NMPC quadratic program contract before
optimization: Q, R, and optional terminal P must be finite symmetric
positive-definite matrices with tokamak state/input dimensions; state, input,
and slew bounds must be finite and ordered; and plant-model evaluations must
return finite state vectors. Invalid math contracts fail closed instead of
propagating undefined SQP or PGD iterates.
The public compute_cost() evaluator includes the finite-horizon terminal
penalty, using configured P when supplied and the controller's conservative
fallback terminal weight otherwise.
Production plant models may provide an analytic linearization_model(x, u)
contract returning finite (6, 6) state and (6, 3) input Jacobians. The
controller validates those matrices before use and records
last_linearization_source == "analytic". If no analytic provider is supplied,
the controller can use linearization_backend="jax" for JAX-traceable plant
models and records last_linearization_source == "jax"; otherwise it falls
back to bounded central finite differences and records
last_linearization_source == "finite_difference". Quadratic weights use a
strict symmetry gate before positive-definite projection so near-zero
off-diagonal asymmetry cannot pass as a valid cost matrix.
DARE-derived terminal matrices are accepted only when finite, symmetric, and
positive definite; invalid solver output falls back to the conservative terminal
weight.
Explicit terminal state sets are configured with paired terminal_x_min and
terminal_x_max vectors. These bounds must lie inside the configured physics
state envelope and currently require qp_backend="scipy", qp_backend="osqp",
qp_backend="casadi", or qp_backend="acados" so the coupled terminal-state
inequality is enforced inside the constrained QP solve rather than checked
after the fact. casadi is a repository-local optional dependency path.
The acados backend is a full optional OCP interface: deployments may inject a
pre-built acados OCP/solver factory, or provide symbolic_dynamics_model(ca, x,
u) so the controller builds a discrete augmented-state acados model. The
augmented state stores the previous actuator vector, making |Δu| <= du_max
a native acados path constraint instead of a post-solve clamp. The default
builder configures SQP, partial-condensing HPIPM, exact Hessian mode, linear
least-squares stage/terminal costs, state/input bounds, terminal state bounds,
warm starts, fail-closed solver-status handling, and a runtime plant-consistency
gate. The returned acados state trajectory must start from the commanded state,
remain inside configured state bounds, satisfy any terminal admissible set, and
match plant_model transitions within acados_dynamics_residual_tol before the
first actuator command is admitted.
The previous input supplied to step() must already satisfy actuator bounds so
the slew-rate projection cannot propagate an unsafe actuator state.
The accepted horizon=1 case is handled as a valid one-step receding-horizon
controller and warm-starts from the bounded previous input.
Each QP solve records last_qp_iterations and last_qp_converged, making
projection-tolerance convergence distinguishable from iteration-budget
exhaustion.
The projected-gradient QP iteration budget is configured by qp_max_iter
instead of being an unobservable hard-coded loop bound.
Linearization perturbations are clipped to the configured state/input domain:
interior points use central differences, while boundary points use one-sided
finite differences.
tune_transport_coefficients_for_tracking() connects NMPC controller tuning to
the differentiable transport facade. It updates four-channel transport
coefficients from the JAX gradient of the transport tracking loss, applies
non-negative coefficient bounds and fractional update caps, and fails closed
when JAX gradients are unavailable. By default, coefficient tuning also runs the
differentiable-transport finite-difference gradient audit before admission and
stores the audit result beside the validated transport campaign metadata for
backend, dtype, radial grid, boundary conditions, closure provenance, and
gradient tolerance.
tune_neural_transport_closure_for_tracking() initialises the same tuning path
from a bounded neural transport closure, preserving the differentiable facade's
four-channel coefficient contract, the explicit JAX-gradient requirement, and
the default gradient-audit admission gate.
tune_transport_sources_for_tracking() applies the audited JAX gradient path to
additive heating, fuelling, and impurity-source schedules. Source lower and
upper bounds are explicit because replay studies may include physically valid
sink terms, and every accepted update carries campaign metadata plus the
gradient-audit result.
tune_transport_source_rollout_for_tracking() extends that admission boundary
from a single transport step to a complete (n_steps, 4, n_rho) source
schedule. It uses JAX for the multi-step rollout gradient, requires a sampled
NumPy finite-difference audit by default, clips per-entry source updates when
configured, and records bounded campaign metadata before the schedule can enter
NMPC tuning.
NonlinearMPC(plant_model, config, linearization_model=None, symbolic_dynamics_model=None, acados_ocp_factory=None, acados_solver_factory=None)
SQP-based NMPC with validated plant linearization contracts.
Each SQP outer iteration linearizes f around the nominal trajectory with an optional analytic Jacobian provider. When no provider is configured, the controller falls back to bounded finite differences. The condensed QP is solved by either SciPy SLSQP or curvature-scaled projected gradient.
Release external solver resources without suppressing control-loop faults.
Evaluate the NMPC cost J over a trajectory.
J = Σ_{k=0}^{N-1} ‖x_k − x_ref‖²_Q + ‖u_k‖²_R Rawlings, Mayne & Diehl 2017, Ch. 1, Eq. (1.2).
Compute optimal first control action via SQP.
Warm-started from the previous solution shifted by one step.
dataclass
¶
Result of a gradient-based transport-coefficient tuning step.
dataclass
¶
TransportSourceScheduleTuningResult(loss, gradient, updated_sources, step_norm, metadata, gradient_audit)
Result of a gradient-based transport source-schedule tuning step.
dataclass
¶
TransportSourceRolloutGradientAudit(loss, epsilon, tolerance, checked_indices, source_max_abs_error, passed)
Finite-difference audit for multi-step source-schedule gradients.
dataclass
¶
TransportSourceRolloutTuningResult(loss, gradient, updated_sources, final_profiles, step_norm, metadata, gradient_audit)
Result of a gradient-based multi-step transport source rollout update.
tune_transport_coefficients_for_tracking(profiles, chi, sources, target_profiles, rho, dt, edge_values, *, weights=None, learning_rate, chi_min=0.0, max_fractional_update=0.1, gradient_tolerance=None, require_gradient_audit=True, gradient_audit_epsilon=1e-05, gradient_audit_tolerance=0.0005, gradient_audit_sample_indices=None, equilibrium_psi=None, _closure_for_metadata=None)
Tune transport coefficients for NMPC tracking through JAX autodiff.
The gradient is taken with respect to the four-channel transport
coefficient profile used by
scpn_control.core.differentiable_transport.transport_loss_gradient.
This function intentionally has no finite-difference fallback: coefficient
tuning is exposed to NMPC only when the differentiable JAX path is present.
tune_transport_sources_for_tracking(profiles, chi, sources, target_profiles, rho, dt, edge_values, *, weights=None, learning_rate, source_min=None, source_max=None, max_absolute_update=None, gradient_tolerance=None, require_gradient_audit=True, gradient_audit_epsilon=1e-05, gradient_audit_tolerance=0.0005, gradient_audit_sample_indices=None, equilibrium_psi=None, _closure_for_metadata=None)
Tune additive transport source schedules through JAX autodiff.
This is the NMPC admission path for heating, fuelling, and impurity-source schedules. It uses the same differentiable transport loss as coefficient tuning, but applies the update to the source array and keeps source bounds explicit because sinks or negative feedback terms may be physically valid in reduced replay studies.
tune_transport_source_rollout_for_tracking(initial_profiles, chi, source_sequence, target_history, rho, dt, edge_values, *, weights=None, learning_rate, source_min=None, source_max=None, max_absolute_update=None, gradient_tolerance=None, require_gradient_audit=True, gradient_audit_epsilon=1e-05, gradient_audit_tolerance=0.0005, gradient_audit_sample_indices=None, equilibrium_psi=None, _closure_for_metadata=None)
Tune a full NMPC transport source rollout through JAX autodiff.
The update acts on the complete (n_steps, 4, n_rho) heating, fuelling,
and impurity-source schedule. A sampled NumPy finite-difference audit is
required by default so the controller does not admit unaudited JAX
gradients for multi-step source optimisation.
tune_neural_transport_closure_for_tracking(profiles, closure, sources, target_profiles, rho, dt, edge_values, *, weights=None, learning_rate, impurity_diffusivity_fraction=1.0, chi_min=0.0, max_fractional_update=0.1, gradient_tolerance=None, require_gradient_audit=True, gradient_audit_epsilon=1e-05, gradient_audit_tolerance=0.0005, gradient_audit_sample_indices=None, equilibrium_psi=None)
Tune NMPC transport coefficients initialised from a neural closure.
Mu-Synthesis (v0.16.0)¶
Structured robust controller with bounded static mu-analysis evidence.
Theory: Doyle 1982 (μ definition); Balas et al. 1993 (DK algorithm); Skogestad & Postlethwaite 2005, §8.5 (convergence and practical use).
Physical uncertainty model for tokamak control follows Ariola & Pironti 2008, Ch. 7: - plasma_position real_scalar ±2 cm - plasma_current real_scalar ±3 % - plasma_shape full ±5 %
D-scaling upper bound on μ(M).
Computes min_D σ̄(D M D^{-1}) where D is block-diagonal with positive real scalars matching delta_structure. This is always ≥ μ(M) and equals μ(M) for complex full blocks (Doyle 1982, IEE Proc. D 129, 242).
A finite-difference descent on log(D) is used to fit the static D-scaling.
Parameters¶
M : np.ndarray, shape (n, n) Closed-loop transfer matrix evaluated at a single frequency. delta_structure : list of (size, block_type) Block sizes and types from StructuredUncertainty.build_Delta_structure().
Returns¶
float Upper bound μ̄ ≥ μ(M).
dataclass
¶
MuSynthesisClaimEvidence(schema_version, source, source_id, model_id, state_dimension, control_dimension, output_dimension, uncertainty_block_count, uncertainty_total_size, max_uncertainty_bound, block_structure, mu_peak_upper_bound, robustness_margin, controller_gain_frobenius_norm, d_scalings, closed_loop_spectral_abscissa, static_dc_analysis_only, reference_source, reference_dataset_id, reference_artifact_sha256, reference_case_count, mu_upper_bound_relative_error, robustness_margin_abs_error, controller_gain_relative_error, d_scaling_relative_error, closed_loop_spectral_abscissa_abs_error, mu_upper_bound_relative_tolerance, robustness_margin_abs_tolerance, controller_gain_relative_tolerance, d_scaling_relative_tolerance, closed_loop_spectral_abscissa_abs_tolerance, validated_claim_allowed, claim_status, payload_sha256='')
Serialisable evidence for bounded or externally validated μ-analysis claims.
mu_synthesis_claim_evidence(controller, *, source, source_id, model_id='bounded_static_mu_synthesis', reference_artifact=None, mu_upper_bound_relative_tolerance=0.05, robustness_margin_abs_tolerance=0.05, controller_gain_relative_tolerance=0.1, d_scaling_relative_tolerance=0.1, closed_loop_spectral_abscissa_abs_tolerance=0.05)
Build fail-closed evidence for bounded static μ-analysis claims.
Raise when μ-analysis evidence is insufficient for validated claims.
Persist μ-analysis claim evidence as deterministic JSON.
Load μ-analysis claim evidence with duplicate-key and digest admission.
Real-Time EFIT (v0.16.0)¶
dataclass
¶
EFITLiteClaimEvidence(schema_version, source, source_id, diagnostic_source, model_id, grid_shape, n_flux_loops, n_b_probes, rogowski_radius_m, chi_squared, n_iterations, wall_time_ms, ip_reconstructed_A, q95, beta_pol, li, psi_relative_error, ip_relative_error, q95_abs_error, beta_pol_abs_error, li_abs_error, psi_relative_tolerance, ip_relative_tolerance, q95_abs_tolerance, beta_pol_abs_tolerance, li_abs_tolerance, facility_claim_allowed, claim_status)
Serialisable admission evidence for EFIT-lite reconstruction claims.
efit_lite_claim_evidence(result, diagnostics, *, source, source_id, diagnostic_source, model_id='bounded_efit_lite', reference_psi=None, reference_shape=None, psi_relative_tolerance=0.05, ip_relative_tolerance=0.02, q95_abs_tolerance=0.1, beta_pol_abs_tolerance=0.1, li_abs_tolerance=0.1)
Build fail-closed evidence for EFIT-lite claim admission.
Return evidence or fail closed before an EFIT-lite facility claim.
Persist EFIT-lite claim evidence as deterministic JSON.
Gain-Scheduled Controller (v0.16.0)¶
Multi-regime PID controller with bumpless gain interpolation.
Scheduling approach: Rugh & Shamma 2000, Automatica 36, 1401, §3. Bumpless transfer via linear interpolation over _TAU_SWITCH seconds avoids step transients at regime boundaries (Walker et al. 2006, §3.2).
LPV interpretation: gains are piecewise-affine in the scheduling vector (I_p, β_N, l_i) per Packard 1994, Systems & Control Letters 22, 79.
Compute PID output with bumpless gain interpolation.
On regime switch: α = (t - t_switch) / τ_switch ∈ [0,1]. Gains interpolated linearly: K(α) = (1-α) K_old + α K_new. Walker et al. 2006, §3.2, Eq. (4).
Shape Controller (v0.16.0)¶
Real-time shape controller using Tikhonov-regularized pseudoinverse.
Gain matrix
K = (JᵀWJ + λI)⁻¹ JᵀW
Ferron et al. 1998, Nucl. Fusion 38, 1055: ISOFLUX control law maps ψ differences at boundary points to coil current corrections.
Safe RL Controller (v0.16.0)¶
PPO augmented with Lagrangian multipliers for constraint satisfaction.
Objective (Tessler et al. 2018, Algorithm 1): r_aug = r - Σ_i λ_i c_i
Multiplier update (dual gradient ascent): λ_i ← max(0, λ_i + lr (C_i - d_i))
Safety guarantee: at convergence J_{C_i}(π*) ≤ d_i for all i (Achiam et al. 2017, Theorem 1, assuming feasibility).
Degrave et al. 2022, Nature 602, 414: penalty-based safety shaping used in DeepMind tokamak controller for disruption avoidance.
Dual gradient ascent step.
λ_i ← max(0, λ_i + lr (C_i - d_i)) Achiam et al. 2017, §5, Lagrangian update.
Stub training loop demonstrating λ update mechanics.
A production implementation replaces the random-action rollout with a PPO policy gradient update on the augmented reward.
Sliding-Mode Vertical (v0.16.0)¶
Vertical-position (VS) controller for an elongated tokamak.
Vertical instability growth rate for elongated plasmas
γ ≈ (κ - 1) / τ_wall
where κ is elongation, τ_wall is the resistive wall time. Lazarus et al. 1990, Nucl. Fusion 30, 111, §2.
Real-time VS implementation at DIII-D: Humphreys et al. 2009, Nucl. Fusion 49, 115003.
Restoring force on the plasma column (rigid-body model): F = -n μ₀ I_p² / (4π R₀) · Z ≡ -K_vs · Z where n is the field-index (n < 0 for unstable equilibria), μ₀ = 4π×10⁻⁷ H/m (SI), I_p in A, R₀ in m. (Wesson 2004, "Tokamaks", Oxford, §3.7)
property
¶
Vertical restoring force coefficient [N/m].
K_vs = -n μ₀ I_p² / (4π R₀) Wesson 2004, §3.7, Eq. (3.7.4). Negative n_index (n < 0) makes K_vs > 0 (unstable).
Compute coil correction voltage.
Plant: m_eff Z̈ = -K_vs Z + F_coil (linearised, Humphreys 2009). e = Z_meas - Z_ref; SMC drives e → 0.
Scenario Scheduler (v0.16.0)¶
Offline trajectory design.
Fault-Tolerant Control (v0.16.0)¶
Control-allocation reconfiguration after actuator or sensor faults.
Weighted pseudo-inverse gain: u = (J^T W J + λI)^{-1} J^T W v. Reference: Bodson 2002, J. Guidance 25, 307, Eq. 12.
Faulted coil columns are zeroed in J before the gain is recomputed, matching the reconfiguration procedure in Ambrosino et al. 2008, Fusion Eng. Des. 83, 1485 for the ITER VS system.
True if remaining actuators span the minimum required target space.
MIN_REQUIRED_RANK = 2 covers control of Ip and vertical position, the two safety-critical outputs identified by Ambrosino et al. 2008.
RZIp Model (v0.16.0)¶
dataclass
¶
RZIPCalibrationEvidence(schema_version, source, source_id, model_id, vertical_inertia_kg, wall_time_constant_s, growth_rate_s_inv, growth_time_ms, reference_growth_rate_s_inv, growth_rate_relative_error, growth_rate_relative_tolerance, facility_claim_allowed, claim_status, evidence_payload_sha256)
Serialisable calibration/admission evidence for a bounded RZIP plant.
rzip_calibration_evidence(rzip, *, source, source_id, wall_time_constant_s, model_id='bounded_rzip', reference_growth_rate_s_inv=None, growth_rate_relative_tolerance=0.2)
Build fail-closed calibration evidence for RZIP claim admission.
local_regression_reference evidence is useful for deterministic
regression reports but never permits facility claims. Facility claims require
documented public, external-code, or measured-discharge reference sources
and must satisfy the declared growth-rate tolerance.
Return evidence or fail closed before a RZIP facility-control claim.
Persist RZIP calibration evidence as deterministic JSON.
RWM Feedback (v0.16.0)¶
Sensor-coil PD feedback for RWM stabilization.
The closed-loop effective growth rate combines wall physics (including rotation) with proportional feedback suppression:
γ_eff = γ_total − G_p M_coil γ_total / (1 + γ_total τ_ctrl)
Garofalo et al. 2002, Phys. Plasmas 9, 1997, Sec. III.
Compute feedback coil currents from sensor measurements.
I = G_p B_r + G_d dB_r/dt
Closed-loop growth rate.
γ_eff = γ_total − G_p M_coil γ_total / (1 + γ_total τ_ctrl)
Garofalo et al. 2002, Phys. Plasmas 9, 1997, Eq. (5). Uses γ_total from RWMPhysics (includes rotation if Ω_φ ≠ 0).
dataclass
¶
RWMClaimEvidence(schema_version, source, source_id, model_id, beta_n, beta_n_nowall, beta_n_wall, tau_wall_s, tau_eff_s, omega_phi_rad_s, wall_radius_m, plasma_radius_m, n_sensors, n_coils, proportional_gain, derivative_gain, controller_latency_s, coil_coupling, open_loop_growth_rate_s_inv, closed_loop_growth_rate_s_inv, required_proportional_gain, reference_closed_loop_growth_rate_s_inv, closed_loop_growth_rate_abs_error, closed_loop_growth_rate_abs_tolerance, facility_claim_allowed, claim_status)
Serialisable admission evidence for RWM feedback claims.
rwm_claim_evidence(rwm, controller, *, source, source_id, model_id='bounded_rwm_feedback', reference_closed_loop_growth_rate_s_inv=None, closed_loop_growth_rate_abs_tolerance=5.0)
Build fail-closed evidence for RWM feedback claim admission.
Local regression evidence is admissible only for bounded repository claims. Facility claims require documented public, external MHD, or measured-shot references plus a finite closed-loop growth-rate comparison within the declared tolerance.
Return evidence or fail closed before a RWM facility-control claim.
Persist RWM claim evidence as deterministic JSON.
Complete Module Index¶
This index keeps the published API reference aligned with every tracked Python module under src/scpn_control/. Domain pages above highlight primary entry points; this section exposes the remaining module surfaces through mkdocstrings so public signatures and docstrings stay visible as the codebase grows.
Top-Level CLI¶
Cli¶
CLI for scpn-control.
Usage::
scpn-control demo --scenario combined --steps 1000
scpn-control benchmark --n-bench 5000
scpn-control validate
scpn-control validate-manifest real_manifest.json --json-out
scpn-control validate-data-manifests --json-out
scpn-control validate-release-evidence artifacts/release_evidence_report.json --json-out
scpn-control validate-physics-traceability --json-out
scpn-control validate-gk-geometry-reference --json-out
scpn-control validate-gk-species-reference --json-out
scpn-control validate-jax-gk-parity --require-parity-artifacts --json-out
scpn-control validate-gk-ood-calibration --require-campaign-artifacts --json-out
scpn-control validate-gk-interface-artifacts --require-interface-artifacts --json-out
scpn-control validate-blob-transport-reference --require-reference-artifacts --json-out
scpn-control validate-elm-reference --require-reference-artifacts --json-out
scpn-control validate-eped-reference --require-reference-artifacts --json-out
scpn-control validate-marfe-reference --require-reference-artifacts --json-out
scpn-control validate-ntm-reference --require-reference-artifacts --json-out
scpn-control validate-neural-equilibrium-reference --require-reference-artifacts --json-out
scpn-control validate-neural-transport-reference --require-reference-artifacts --json-out
scpn-control validate-neural-turbulence-reference --require-reference-artifacts --json-out
scpn-control validate-orbit-reference --require-reference-artifacts --json-out
scpn-control validate-uncertainty-reference --require-reference-artifacts --json-out
scpn-control validate-vmec-reference --require-reference-artifacts --json-out
scpn-control validate-rzip-reference --require-reference-artifacts --json-out
scpn-control validate-density-reference --require-reference-artifacts --json-out
scpn-control validate-burn-reference --require-reference-artifacts --json-out
scpn-control validate-volt-second-reference --require-reference-artifacts --json-out
scpn-control validate-current-drive-reference --require-reference-artifacts --json-out
scpn-control validate-mu-synthesis-reference --require-reference-artifacts --json-out
scpn-control validate-free-boundary-reference --require-reference-artifacts --json-out
scpn-control validate-disruption-reference --require-reference-artifacts --json-out
scpn-control validate-digital-twin-reference --require-reference-artifacts --json-out
scpn-control validate-soc-reference --require-reference-artifacts --json-out
scpn-control acquire-mdsplus-shot --spec-json validation/reference_data/diiid/acquisition_specs/shot_163303_mdsplus.json
scpn-control live --port 8765 --zeta 0.5
scpn-control hil-test --shots-dir validation/reference_data/diiid/disruption_shots
scpn-control run-hardware-campaign --neurons 64 --core-snn 1 --core-z3 2 --core-net 3 --core-hb 4
Timing benchmark: PID vs SNN control step latency.
A warm-up phase is executed but excluded from timing to isolate steady-state micro-benchmark behavior from cold-start effects (JIT, allocator, cache, init).
run_hardware_campaign(neurons, seed, state_init_r, state_init_z, plant_gain, target_r, target_z, beta_n_limit, itpa_cgb, itpa_path, kuramoto_weights, transport_endpoint, transport_port, transport_ttl, transport_max_queue, transport_backend, execution_backend, heartbeat_port, heartbeat_timeout_ms, core_snn, core_z3, core_net, core_hb, steps, runtime_s, tick_interval_s, pacing_mode, runtime_admission_policy, max_publish_failures, formal_mode, formal_stride, formal_depth, formal_max_marking, formal_channel_capacity, require_native, json_out)
Launch the Rust-native execution plane from Python control layer.
validate(json_out, data_manifest_root, no_data_manifests, no_verify_artifacts, jax_gk_parity_root, no_jax_gk_parity, physics_traceability_registry, no_physics_traceability, multi_shot_campaign_python_report, multi_shot_campaign_rust_report, multi_shot_min_digest_count, no_multi_shot_campaign_evidence, runtime_admission_report, no_runtime_admission_evidence, native_formal_certificate_report, native_formal_max_aot_p99_cycle_us, no_native_formal_certificate)
Run import hygiene, provenance, parity, traceability, and formal evidence validation.
Validate top-level release evidence report admission.
Validate real-shot or synthetic-shot data manifest provenance.
validate_data_manifests_command(root, output_json, no_verify_artifacts, require_real_acquisition, json_out)
Validate repository data manifests, local artefacts, and acquisition specs.
Validate physics fidelity traceability and bounded-claim contracts.
Validate real external-code evidence for linear GK agreement.
Validate Miller geometry output against immutable reference cases.
Validate GK species normalisation and collision reference cases.
validate_jax_gk_parity_command(artifact_root, require_parity_artifacts, require_cases, require_backends, output_json, json_out)
Validate persisted JAX/native GK parity artifacts.
validate_gk_ood_calibration_command(artifact_root, require_campaign_artifacts, output_json, json_out)
Validate persisted GK OOD calibration campaign artifacts.
validate_gk_interface_artifacts_command(artifact_root, require_interface_artifacts, output_json, json_out)
Validate persisted external GK interface parser artifacts.
validate_blob_transport_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted blob-transport reference artifacts.
Validate persisted ELM reference artifacts.
Validate persisted EPED reference artifacts.
Validate persisted MARFE reference artifacts.
Validate persisted NTM reference artifacts.
validate_neural_equilibrium_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted neural-equilibrium reference artifacts.
validate_neural_transport_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted neural-transport reference artifacts.
validate_neural_turbulence_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted neural-turbulence reference artifacts.
Validate persisted orbit-following reference artifacts.
validate_uncertainty_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted uncertainty-quantification reference artifacts.
Validate persisted VMEC-lite reference artifacts.
Validate persisted RZIP vertical-stability reference artifacts.
validate_current_drive_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted current-drive reference artifacts.
validate_mu_synthesis_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted mu-synthesis reference artifacts.
validate_volt_second_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted volt-second reference artifacts.
Validate persisted DT burn-control reference artifacts.
validate_density_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted density-control reference artifacts.
validate_free_boundary_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted free-boundary tracking reference artifacts.
validate_disruption_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted disruption-mitigation reference artifacts.
validate_digital_twin_reference_command(artifact_root, require_reference_artifacts, output_json, json_out)
Validate persisted tokamak digital-twin reference artifacts.
Validate persisted SOC turbulence-learning reference artifacts.
acquire_mdsplus_shot_command(spec_json, tree, shot, signals_json, output_npz, manifest_json, source_uri, access_policy, licence, retrieved_at, json_out)
Acquire selected MDSplus shot signals and write a validated manifest.
Run full RMSE validation dashboard against reference data.
live(port, host, layers, n_per, zeta, psi, tick_interval, api_key, command_rate_limit, command_rate_window_s, max_payload_bytes, max_client_write_buffer_bytes, allow_unauthenticated_clients, allow_query_token_auth, require_tls, allow_insecure_remote, allowed_origin, allowed_action, tls_cert, tls_key)
Start real-time WebSocket phase sync server.
Streams Kuramoto R/V/lambda tick snapshots over ws://
Hardware-in-the-loop test campaign against reference shots.
Control Modules¶
Burn Controller¶
D-T burn-control, alpha-heating, Lawson, and auxiliary-heating utilities.
dataclass
¶
BurnControlClaimEvidence(schema_version, source, source_id, model_id, R0_m, a_m, kappa, profile_points, ne_volume_average_1e20_m3, Te_volume_average_keV, Ti_volume_average_keV, P_alpha_MW, P_aux_MW, Q, lawson_triple_product, lawson_margin, burn_fraction, reactivity_exponent, thermally_stable, controller_Q_target, controller_T_target_keV, controller_P_aux_max_MW, controller_command_MW, reference_source, reference_dataset_id, reference_artifact_sha256, reference_case_count, P_alpha_relative_error, Q_abs_error, lawson_margin_abs_error, burn_fraction_relative_error, reactivity_exponent_abs_error, P_alpha_relative_tolerance, Q_abs_tolerance, lawson_margin_abs_tolerance, burn_fraction_relative_tolerance, reactivity_exponent_abs_tolerance, reactor_claim_allowed, claim_status)
Serialisable evidence for bounded or validated burn-control claims.
Alpha heating power from D-T fusion.
P_α = n_D n_T <σv> E_α × V. ITER Physics Basis 1999, Nucl. Fusion 39, 2137. Reactivity <σv>_DT: Bosch & Hale 1992, Nucl. Fusion 32, 611.
Alpha power density [MW/m^3] for 50:50 D-T mixture.
p_α = n_D n_T <σv>(T_i) E_α. ITER Physics Basis 1999, Nucl. Fusion 39, 2137, Eq. (2.2.1).
P_alpha [MW] integrated over plasma volume.
dV = 4π² R₀ a² κ ρ dρ (torus shell element in normalised radius).
Fusion energy gain Q = P_fus / P_aux = 5 P_α / P_aux.
P_fus = 5 P_α because α carries 1/5 of total DT energy (E_α/E_fus = 3.52/17.6). Delayed feedback for alpha-heating runaway: Mitarai & Muraoka 1999, Nucl. Fusion 39, 725.
Thermal burn stability based on reactivity exponent.
Delayed alpha-heating feedback control: Mitarai & Muraoka 1999, Nucl. Fusion 39, 725.
d(ln <σv>) / d(ln T), evaluated via finite difference.
Stability requires this exponent < 2. Mitarai & Muraoka 1999, Nucl. Fusion 39, 725, Eq. (5).
True if d(ln <σv>)/d(ln T) < 2 (Mitarai & Muraoka 1999).
T where d(ln <σv>)/d(ln T) = 2 (bisection, 5–30 keV).
PI burn controller targeting Q and T_i.
Delayed feedback prevents alpha-heating runaway. Mitarai & Muraoka 1999, Nucl. Fusion 39, 725. P_aux_max = 73 MW matches ITER heating system capacity. ITER Physics Basis 1999, Nucl. Fusion 39, 2137, Table I.
Scan T to find P_α(T) + P_aux = P_loss(T) intersections.
P_loss = 3 n T V / τ_E (energy balance, 0-D model). Lawson criterion: n τ_E T > 3×10^21 m^-3 s keV for ignition. Lawson 1957, Proc. Phys. Soc. B 70, 6.
burn_control_claim_evidence(alpha_heating, controller, *, rho, ne_20, Te_keV, Ti_keV, tau_E_s, P_aux_MW, source, source_id, model_id='bounded_dt_burn_control', reference_artifact=None, P_alpha_relative_tolerance=0.05, Q_abs_tolerance=0.5, lawson_margin_abs_tolerance=0.2, burn_fraction_relative_tolerance=0.1, reactivity_exponent_abs_tolerance=0.25)
Build fail-closed DT burn-control evidence from explicit profiles.
Raise when burn-control evidence is insufficient for reactor-control claims.
Persist burn-control claim evidence as deterministic JSON.
Return n τ_E T [m^-3 s keV].
Ignition requires n τ_E T > 3×10^21 m^-3 s keV. Lawson 1957, Proc. Phys. Soc. B 70, 6.
Approximate DT burn fraction.
f_b ≈ a² n_DT <σv> / (4 v_th). Wesson 2011, "Tokamaks" 4th ed., Eq. 1.7.3.
Codac Interface¶
ITER CODAC/EPICS integration prototype.
Generates EPICS .db and OPC-UA XML nodeset files for binding a NeuroSymbolicController to the ITER Plant Instrumentation & Control (I&C) infrastructure. Does NOT require pyepics or any EPICS runtime.
References
ITER CODAC Handbook v7.0, §3.2 (PV naming conventions) IEC 62541 (OPC Unified Architecture)
dataclass
¶
CODACConfig(pv_prefix='ITER-SCPN', cycle_hz=1000.0, timeout_ms=1.5, heartbeat_pv='ITER-SCPN:HEARTBEAT', interlock_pvs=('ITER-CIS:INTERLOCK:VDE', 'ITER-CIS:INTERLOCK:HALO', 'ITER-CIS:INTERLOCK:DISRUPTION'))
CODAC plant-system configuration.
dataclass
¶
Single EPICS process variable definition.
dataclass
¶
CODACRuntimeEvidence(schema_version, generated_utc, controller_id, plant_system, pv_prefix, cycle_hz, cycle_budget_us, timeout_budget_us, deadline_us, observed_cycle_p50_us, observed_cycle_p95_us, observed_cycle_p99_us, observed_cycle_max_us, input_channel_count, output_channel_count, interlock_pv_count, interlock_checks, interlock_blocks, backpressure_events, epics_db_sha256, opcua_nodeset_sha256, facility_claim_allowed, claim_status, payload_sha256)
Tamper-evident CODAC/EPICS runtime admission evidence.
Binds a NeuroSymbolicController to ITER CODAC I&C channels.
Convert raw EPICS PV dict to ControlObservation for the controller.
Return True if any hard limit is violated (actuation must be blocked).
Write EPICS .db file with record() entries for all channels.
Return OPC-UA XML nodeset text for ITER SDN integration.
Write OPC-UA XML nodeset for ITER SDN integration.
codac_runtime_evidence(interface, *, controller_id, observed_cycle_us, interlock_checks, interlock_blocks, backpressure_events=0, plant_system='ITER-CODAC-EPICS', generated_utc=None, facility_claim_allowed=False)
Build tamper-evident runtime evidence for a CODAC/EPICS boundary.
Fail closed unless CODAC runtime evidence can support a facility claim.
Persist CODAC runtime evidence as canonical sorted JSON.
Load CODAC runtime evidence with duplicate-key and digest admission.
dataclass
¶
CODACRuntimeEvidence(schema_version, generated_utc, controller_id, plant_system, pv_prefix, cycle_hz, cycle_budget_us, timeout_budget_us, deadline_us, observed_cycle_p50_us, observed_cycle_p95_us, observed_cycle_p99_us, observed_cycle_max_us, input_channel_count, output_channel_count, interlock_pv_count, interlock_checks, interlock_blocks, backpressure_events, epics_db_sha256, opcua_nodeset_sha256, facility_claim_allowed, claim_status, payload_sha256)
Tamper-evident CODAC/EPICS runtime admission evidence.
codac_runtime_evidence(interface, *, controller_id, observed_cycle_us, interlock_checks, interlock_blocks, backpressure_events=0, plant_system='ITER-CODAC-EPICS', generated_utc=None, facility_claim_allowed=False)
Build tamper-evident runtime evidence for a CODAC/EPICS boundary.
Fail closed unless CODAC runtime evidence can support a facility claim.
Persist CODAC runtime evidence as canonical sorted JSON.
Load CODAC runtime evidence with duplicate-key and digest admission.
Controller Tuning¶
Density Controller¶
Density-control plant, estimator, fuelling actuators, and optimisation utilities.
Particles/s. Edge-localised source from gas injection.
Pacher et al. 2007, Nucl. Fusion 47, 469: ITER gas-injection modelling places the effective source within ~3% of the minor radius from the wall.
pellet_source(speed_ms, radius_mm, launch_angle_deg=0.0, *, ne_profile=None, Te_eV_profile=None, B0_T=5.3, injection_side='HFS')
NGS trajectory deposition profile from a single pellet.
The deposition path is delegated to :class:PelletTrajectory, which
integrates Parks-Turnbull neutral-gas-shielding ablation along the
radial pellet trajectory and applies the repository drift correction.
Core-peaked particle source from neutral beam injection.
Recycling_coeff = 0.97 is the standard ITER assumption for a metal wall. ITER Physics Basis 1999, Nucl. Fusion 39, 2175, §4.2.
dataclass
¶
DensityControlClaimEvidence(schema_version, source, source_id, geometry_source, transport_source, actuator_source, diagnostic_source, model_id, n_rho, R0_m, a_m, dt_requested_s, dt_cfl_s, cfl_limited, greenwald_fraction, below_iter_greenwald_margin, gas_puff_rate, pellet_frequency, pellet_speed, cryopump_speed, total_source_particles_per_s, particle_inventory_delta, greenwald_fraction_abs_error, inventory_delta_relative_error, greenwald_fraction_abs_tolerance, inventory_delta_relative_tolerance, facility_density_claim_allowed, claim_status)
Serialisable provenance and reference-comparison evidence for density-control claims.
PI density controller with Greenwald limit enforcement.
Greenwald limit: n_GW = I_p / (π a²) [10^20 m^-3] Greenwald 2002, PPCF 44, R27, Eq. 1.
ITER operational margin: n/n_GW < 0.85. ITER Physics Basis 1999, Nucl. Fusion 39, 2175, §2.3.
staticmethod
¶
n_GW = I_p / (π a²) [10^20 m^-3], converted to [m^-3].
Greenwald 2002, PPCF 44, R27, Eq. 1.
Volume-averaged n / n_GW.
Greenwald 2002, PPCF 44, R27, Eq. 1. ITER safe operating limit: fraction < 0.85. ITER Physics Basis 1999, Nucl. Fusion 39, 2175, §2.3.
True if volume-averaged density is within the ITER safety margin.
ITER Physics Basis 1999, Nucl. Fusion 39, 2175, §2.3: n/n_GW < 0.85.
Abel-transform projection matrix for interferometry chords.
density_control_claim_evidence(model, controller, *, source, source_id, geometry_source, transport_source, actuator_source, diagnostic_source, ne_before, ne_after, sources, command, dt_requested_s, model_id='bounded_density_control', reference_greenwald_fraction=None, reference_inventory_delta=None, greenwald_fraction_abs_tolerance=0.02, inventory_delta_relative_tolerance=0.05)
Build fail-closed density-control evidence for bounded or facility claims.
Raise when density-control evidence is insufficient for a facility claim.
Persist density-control claim evidence as deterministic JSON.
dataclass
¶
DensityControlClaimEvidence(schema_version, source, source_id, geometry_source, transport_source, actuator_source, diagnostic_source, model_id, n_rho, R0_m, a_m, dt_requested_s, dt_cfl_s, cfl_limited, greenwald_fraction, below_iter_greenwald_margin, gas_puff_rate, pellet_frequency, pellet_speed, cryopump_speed, total_source_particles_per_s, particle_inventory_delta, greenwald_fraction_abs_error, inventory_delta_relative_error, greenwald_fraction_abs_tolerance, inventory_delta_relative_tolerance, facility_density_claim_allowed, claim_status)
Serialisable provenance and reference-comparison evidence for density-control claims.
density_control_claim_evidence(model, controller, *, source, source_id, geometry_source, transport_source, actuator_source, diagnostic_source, ne_before, ne_after, sources, command, dt_requested_s, model_id='bounded_density_control', reference_greenwald_fraction=None, reference_inventory_delta=None, greenwald_fraction_abs_tolerance=0.02, inventory_delta_relative_tolerance=0.05)
Build fail-closed density-control evidence for bounded or facility claims.
Raise when density-control evidence is insufficient for a facility claim.
Persist density-control claim evidence as deterministic JSON.
Detachment Controller¶
Divertor detachment state estimation and impurity-seeding control utilities.
Peak radiation temperature [eV] for each impurity species.
N₂: ~10 eV, Ne: ~30 eV, Ar: ~100 eV. Kallenbach et al. 2015, Nucl. Fusion 55, 053026, §3.
Radiation front location: 0 = target, 1 = X-point.
Higher seeding rate moves the front toward the X-point. Higher P_SOL pushes it back toward the target. Lipschultz et al. 1999, PPCF 41, A585: empirical front-position scaling.
DOD = Γ_t,attached / Γ_t,actual via T_target rollover.
DOD = 1 when attached (T_t > 5 eV). Below the onset (Stangeby 2000, Ch. 16) DOD rises as T_t falls, reflecting reduced ion flux from volumetric recombination.
PI controller driving divertor detachment via impurity seeding.
Target: T_div ≈ 3 eV (below the 5 eV onset, above the MARFE threshold). Stangeby 2000, Ch. 16: T_div < 5 eV criterion for detachment onset. Lipschultz et al. 1999, PPCF 41, A585: thermal bifurcation stability.
Nitrogen seeding rates up to 10^22 molecules/s for ITER. Kallenbach et al. 2015, Nucl. Fusion 55, 053026.
Steady-state scan across seeding rates to locate the thermal bifurcation.
Thermal bifurcation (S-curve) in T_target vs seeding_rate described by: Lipschultz et al. 1999, PPCF 41, A585 (Alcator C-Mod data + model). Two-point model (Stangeby 2000, Eq. 5.69) provides the upstream–target link.
Seeding rate where ion flux Γ ~ n_t √T_t peaks (flux rollover).
Rollover marks the detachment onset on the bifurcation S-curve. Stangeby 2000, Ch. 16; Lipschultz et al. 1999, PPCF 41, A585.
Parallel heat flux from the Spitzer conduction model.
q_∥ = κ₀ T_u^(7/2) / (7 L_∥) [W m^-2]
Stangeby 2000, "The Plasma Boundary of Magnetic Fusion Devices", Eq. 5.69. κ₀ = 2390 W m^-1 eV^(-7/2) (electron Spitzer conductivity, Stangeby Eq. 5.67).
Federated Disruption¶
Federated learning framework for cross-machine disruption prediction.
Trains a shared MLP disruption classifier across heterogeneous tokamak datasets (DIII-D, JET, KSTAR, EAST, SPARC) without centralising raw data. Supports FedAvg (McMahan et al., AISTATS 2017) and FedProx (Li et al., MLSys 2020) aggregation strategies.
Disruption features: Ip, beta_N, q95, n/n_GW, li, dBp/dt, locked_mode_amplitude, n1_rms — 8-dimensional input space whose distributions differ across machines (JET: higher Ip; DIII-D: more shaping; KSTAR: longer pulses).
dataclass
¶
Facility-level differential privacy for federated client updates.
dataclass
¶
PrivacyLedgerEntry(round_index, participating_clients, epsilon_spent, cumulative_epsilon, delta, max_update_norm, noise_multiplier, clipped_clients)
Per-round facility-level privacy accounting record.
dataclass
¶
FacilityBenchmarkSummary(aggregation, machines, n_rounds, mean_accuracy, mean_loss, per_machine_accuracy, privacy_epsilon, privacy_delta, evidence_kind)
Deterministic benchmark summary for a federated disruption campaign.
dataclass
¶
FederatedConfig(n_rounds=10, local_epochs=5, learning_rate=0.01, aggregation='fedavg', mu_proximal=0.01, min_clients=2, machines=(lambda: ['DIII-D', 'JET', 'KSTAR'])(), dp_config=None)
Configuration for federated disruption prediction training.
Local training client for a single tokamak.
SGD on local data, starting from global_weights.
When mu_proximal > 0, adds the FedProx penalty (mu/2)||w - w_global||^2 to the loss gradient.
Orchestrates federated training across machine clients.
FedAvg: weighted average of model weights by dataset size.
McMahan et al., "Communication-Efficient Learning of Deep Networks from Decentralized Data", AISTATS 2017.
FedProx aggregation with proximal regularisation.
The proximal term is applied during local training (not aggregation), so aggregation itself is weighted averaging — the difference from FedAvg is in the client-side gradient update. This method exists for API symmetry; the mu parameter documents the proximal weight used.
Full federated training loop.
Returns per-round metrics including per-client accuracy, loss, n_samples.
Return cumulative facility-level DP spend for the current server.
Return the conservative Gaussian-mechanism epsilon for one round.
Linearly compose the conservative per-round Gaussian epsilon budget.
Clip per-parameter gradient norms and add Gaussian noise (DP-SGD).
Reference: Abadi et al., "Deep Learning with Differential Privacy", CCS 2016.
Create MachineClient instances with synthetic disruption data.
Parameters¶
machine_configs : list of dicts Each dict must have "machine" (str) and optionally "n_train" (int, default 200), "n_test" (int, default 50), "disruption_fraction" (float, default 0.4), "learning_rate" (float, default 0.01). seed : int Base RNG seed; each machine gets seed + index.
Create clients from per-facility arrays without centralising raw shots.
Each facility payload must contain X_train, y_train, X_test, and
y_test. The constructor enforces the shared 8-feature disruption
contract and binary label boundary before the data can enter a federation.
run_synthetic_multifacility_benchmark(*, machines=('DIII-D', 'JET', 'KSTAR', 'EAST'), n_rounds=4, local_epochs=3, aggregation='fedprox', dp_config=None, seed=20240531)
Run a deterministic synthetic multi-facility disruption benchmark.
This benchmark exercises the production federation, heterogeneity, and privacy-accounting contracts. It is not measured cross-facility validation.
State Estimator¶
Extended Kalman Filter (EKF) for plasma state estimation.
EKF predict-update cycle
x_{k+1|k} = f(x_k, u_k) P_{k+1|k} = F P_k Fᵀ + Q y_k = z_k − H x_{k|k−1} S_k = H P_{k|k−1} Hᵀ + R K_k = P_{k|k−1} Hᵀ S_k⁻¹ x_{k|k} = x_{k|k−1} + K_k y_k P_{k|k} = (I − K_k H) P_{k|k−1}
Simon 2006, "Optimal State Estimation", Ch. 13.
Magnetic diagnostic-based tokamak state estimation: Lister et al. 1997, Nucl. Fusion 37, 1633 — real-time equilibrium reconstruction from magnetic measurements on TCV.
Observer-based current profile estimation: Moreau et al. 2008, Nucl. Fusion 48, 106001 — model-based estimation of the safety-factor profile q(ρ) from partial magnetic and kinetic data.
EKF for 6D plasma state estimation.
State vector: [R, Z, vR, vZ, Ip, Te_core] Measurement vector: [R, Z, Ip, Te_core]
Process model: constant-velocity for (R, Z); Ip and Te_core modelled as random walks (zero-input prediction). The linearization is exact because the process model is affine, so this degenerates to a standard Kalman filter for the current implementation.
Lister et al. 1997, Nucl. Fusion 37, 1633: magnetic flux and field measurements provide [R, Z, Ip] observability in real time. Moreau et al. 2008, Nucl. Fusion 48, 106001: Te_core estimated via ECE + model observer.
Parameters¶
x0 : np.ndarray Initial state estimate (6D). P0 : np.ndarray Initial error covariance (6×6). Q : np.ndarray Process noise covariance (6×6). R_cov : np.ndarray Measurement noise covariance (4×4).
Advance state estimate and covariance.
x_{k+1|k} = F x_k (constant-velocity kinematics + random walk) P_{k+1|k} = F P_k Fᵀ + Q dt
Simon 2006, Ch. 13, Eqs. (13.7)–(13.8).
Parameters¶
dt : float Time step [s]. u : np.ndarray, optional Control input (unused in constant-velocity model).
Returns¶
np.ndarray — Predicted state.
Correct state estimate with a new measurement.
y = z − H x_{k|k−1} S = H P_{k|k−1} Hᵀ + R K = P_{k|k−1} Hᵀ S⁻¹ x = x_{k|k−1} + K y P = (I − K H) P_{k|k−1}
Simon 2006, Ch. 13, Eqs. (13.9)–(13.13).
Parameters¶
z : np.ndarray Measured values [R, Z, Ip, Te_core] (4D).
Returns¶
np.ndarray — Updated state.
Volt Second Manager¶
Volt-second flux budgeting, consumption monitoring, and scenario feasibility utilities.
dataclass
¶
VoltSecondClaimEvidence(schema_version, source, source_id, model_id, Phi_CS_Vs, L_plasma_H, R_plasma_Ohm, Ip_MA, I_bs_MA, ramp_duration_s, flat_duration_s, ramp_down_duration_s, ramp_flux_Vs, flat_top_flux_Vs, ramp_down_flux_Vs, total_flux_Vs, margin_Vs, within_budget, ejima_startup_flux_Vs, max_flattop_duration_s, bootstrap_source, reference_source, reference_dataset_id, reference_artifact_sha256, reference_case_count, total_flux_relative_error, flat_top_duration_relative_error, ejima_flux_relative_error, bootstrap_current_abs_error_MA, margin_abs_error_Vs, total_flux_relative_tolerance, flat_top_duration_relative_tolerance, ejima_flux_relative_tolerance, bootstrap_current_abs_tolerance_MA, margin_abs_tolerance_Vs, facility_claim_allowed, claim_status)
Serialisable evidence for bounded or facility volt-second claims.
Volt-second budget tracker.
Implements the balance
∫ V_loop dt = L_p dI_p + R_p I_p dt
Reference: Wesson 2011, Tokamaks 4th ed., Eq. 3.7.4.
L_p · I_p — inductive volt-second consumption.
Wesson 2011, Tokamaks 4th ed., Eq. 3.7.4 (first term).
∫ R_p I_p dt — resistive volt-second consumption during ramp.
Wesson 2011, Tokamaks 4th ed., Eq. 3.7.4 (second term).
Startup flux via Ejima coefficient: ΔΨ = C_Ejima · μ₀ · R₀ · I_p.
Ejima et al. 1982, Nucl. Fusion 22, 1313, Eq. 2. C_EJIMA = 0.4 is the ITER design value.
τ_flat = (Ψ_avail − Ψ_startup) / (R_p I_p).
ITER Physics Basis 1999, Nucl. Fusion 39, 2137, §3.
staticmethod
¶
Simplified bootstrap current proxy: I_bs ~ ε^{1/2} · ∫ dp/dr dr.
Full neoclassical expression in Wesson 2011, Ch. 4.9. ε = a/R₀ is the inverse aspect ratio.
Integrate V_loop dt to track consumed volt-seconds.
Wesson 2011, Tokamaks 4th ed., Eq. 3.7.4 — V_loop drives both inductive and resistive flux consumption.
Decompose total flux consumption into ramp / flat-top / ramp-down.
Ramp: L_p I_p (inductive) + R_p · 0.5 I_p · t_ramp (resistive at mean current). Flat-top: R_p · (I_p − I_bs) · t_flat. Ramp-down: resistive loss minus partial inductive recovery. Reference: ITER Physics Basis 1999, Nucl. Fusion 39, 2137, §3.
volt_second_claim_evidence(budget, report, *, Ip_MA, I_bs_MA, ramp_duration_s, flat_duration_s, ramp_down_duration_s, R0_m, ramp_flux_for_flattop_Vs, source, source_id, bootstrap_source='repository bootstrap-current proxy', model_id='bounded_volt_second_manager', reference_artifact=None, total_flux_relative_tolerance=0.03, flat_top_duration_relative_tolerance=0.05, ejima_flux_relative_tolerance=0.05, bootstrap_current_abs_tolerance_MA=0.25, margin_abs_tolerance_Vs=0.5)
Build fail-closed evidence for scenario volt-second claims.
Raise when volt-second evidence is insufficient for facility scenario claims.
Persist volt-second claim evidence as deterministic JSON.
Core Support and Physics Modules¶
Rust Compatibility¶
Backward compatibility layer: imports from Rust (scpn_control_rs) if available, falls back to pure-Python implementations.
Usage
from scpn_control.core._rust_compat import FusionKernel, RUST_BACKEND
Drop-in wrapper around Rust PyFusionKernel that mirrors the Python FusionKernel attribute interface (.Psi, .R, .Z, .RR, .ZZ, .cfg, etc.).
Delegates equilibrium solve to Rust for ~20x speedup while keeping all attribute accesses compatible with downstream code.
Interpolate psi at a single (R, Z) point via Rust bilinear interpolation.
Compute (B_R, B_Z) from current psi, delegating to Rust when available.
Locate the null point (B=0) using local minimization. Matches Python FusionKernel.find_x_point() interface.
Compute vacuum field via Python FusionKernel (not yet in PyO3).
D-T fusion thermodynamics from current equilibrium (Rust backend).
Python wrapper for Rust SpikingControllerPool (LIF neuron population).
Falls back with ImportError if the Rust extension is not compiled.
Parameters¶
n_neurons : int Number of LIF neurons per sub-population (positive/negative). gain : float Output scaling factor. window_size : int Sliding window length for rate-code averaging.
Python wrapper for Rust NeuroCyberneticController (dual R+Z SNN pools).
Falls back with ImportError if the Rust extension is not compiled.
Parameters¶
target_r : float Target major-radius position [m]. target_z : float Target vertical position [m].
Rust PID controller (kp, ki, kd gains with finite-input validation).
Parameters¶
kp, ki, kd : float Proportional / integral / derivative gains.
Rust iso-flux controller (decoupled R + Z PID).
Parameters¶
target_r, target_z : float Target R, Z position [m].
Rust H-infinity observer-based controller for vertical stability.
Uses LQR-approximated gains for the 2-state VDE plant (full DARE solver pending ndarray-linalg).
Parameters¶
gamma_growth : float Unstable growth rate [1/s]. damping : float Passive damping coefficient. gamma : float H-infinity performance level. u_max : float Actuator saturation limit [A]. dt : float Nominal timestep [s].
SPI disruption-mitigation simulator. Uses Rust backend when available, falls back to pure-Python implementation matching Rust constants.
Parameters¶
w_th_mj : float Initial stored thermal energy [MJ]. ip_ma : float Initial plasma current [MA]. te_kev : float Initial electron temperature [keV].
Python tearing-mode simulation fallback.
Bosch-Hale D-T reaction rate [m³/s] at temperature t_kev [keV].
rust_multigrid_vcycle(source, psi_bc, r_min, r_max, z_min, z_max, nr, nz, tol=1e-06, max_cycles=500)
Multigrid V-cycle GS solver. Uses Rust backend when available, falls back to FusionKernel's Python multigrid.
Returns¶
tuple of (psi, residual, n_cycles, converged)
SVD-based coil current correction. Uses Rust backend when available, falls back to NumPy SVD pseudoinverse.
Parameters¶
response_matrix : ndarray, shape (m, n) Plant response Jacobian (typically 2 x n_coils). error : ndarray, shape (m,) Position error vector [R_err, Z_err]. gain : float Correction gain factor.
Returns¶
ndarray, shape (n,) Coil current deltas.
Validators¶
Shared input validators for public API boundaries.
require_bounded_float(name, value, *, low=-inf, high=inf, low_exclusive=False, high_exclusive=False)
Validate scalar is finite and within [low, high] (or exclusive bounds).
Subsumes require_positive_float, require_non_negative_float, require_fraction for arbitrary bound combinations.
Validate array is finite with optional shape/ndim constraint.
Alfven Eigenmodes¶
Alfven-eigenmode gap, drive, damping, and stability-screening utilities.
Shear Alfvén continuum and TAE gap structure.
Cheng & Chance 1986, Phys. Fluids 29, 3695 — gap condition q = (2m+1)/(2n).
Toroidal Alfvén eigenmode frequency and damping.
Cheng & Chance 1986, Phys. Fluids 29, 3695.
Electron Landau damping rate for a TAE.
γ_e = −(π/4)·(ω / |k_∥ v_the|) · ω · exp(−(ω / k_∥ v_the)²)
where v_the = sqrt(2 T_e / m_e), T_e in SI k_∥ = (n − m/q) / (q R₀) (parallel wave vector)
Returns |γ_e| (positive damping magnitude).
Rosenbluth & Rutherford 1975, Phys. Rev. Lett. 34, 1428, Eq. 9.
Fast-particle (alpha/beam-ion) drive for Alfvén eigenmodes.
Fu & Van Dam 1989, Phys. Fluids B 1, 1949. Heidbrink 2008, Phys. Plasmas 15, 055501 — review.
Analytic fast-ion resonance function.
Primary resonance (v_f ≃ v_A): F_main = (v_f/v_A)³ · exp(−(v_f/v_A − 1)² / (2·0.2²))
Sideband resonance (v_f ≃ v_A/3): F_sb = 0.15·(v_f/v_A)² · exp(−(v_f/v_A − 1/3)² / (2·0.1²))
Fu & Van Dam 1989, Phys. Fluids B 1, 1949, Eq. 28. σ_main = 0.2, σ_side = 0.1, A_side = 0.15 — ibid., Table 1.
γ_fast / ω ≃ β_f · q² · F(v_f/v_A) Fu & Van Dam 1989, Phys. Fluids B 1, 1949, Eq. 15. Heidbrink 2008, Phys. Plasmas 15, 055501, Eq. 10.
Combined TAE stability: electron Landau damping vs. fast-particle drive.
Rosenbluth & Rutherford 1975, Phys. Rev. Lett. 34, 1428. Fu & Van Dam 1989, Phys. Fluids B 1, 1949. Heidbrink 2008, Phys. Plasmas 15, 055501.
Reversed-Shear Alfvén Eigenmode (RSAE) frequency.
ω²_RSAE = ω²_BAE + (v_A / (2 q_min R₀))² · (n − m/q_min)²
where the beta-induced acoustic gap (BAE) frequency is ω_BAE = v_thi · sqrt(7/4 + T_e/T_i) / R₀, v_thi = sqrt(T_i / m_i).
Sharapov et al. 2001, Phys. Lett. A 289, 127, Eq. 5. Geodesic acoustic coupling: Zonca & Chen 1996, Plasma Phys. 38, 2011.
Blob Transport¶
SOL filament / blob transport.
References¶
D'Ippolito, D. A., Myra, J. R. & Zweben, S. J., Phys. Plasmas 18 (2011) 060501. Comprehensive review of convective SOL transport by blobs/filaments. Krasheninnikov, S. I., Phys. Lett. A 283 (2001) 368. Blob velocity scaling: v_b = 2T_e/(e B R δ_b) × (c_s/Ω_i) [Eq. 5]. Myra, J. R. et al., Phys. Plasmas 13 (2006) 112502. Blob size scaling: δ_b ≈ (ρ_s L_∥)^(2/5) R^(1/5) — critical size at sheath-inertial transition. Garcia, O. E. et al., Phys. Plasmas 19 (2012) 092306. PDF of SOL density fluctuations follows a gamma distribution; log-normal approximation valid for large intermittency parameter.
Blob velocity and size in sheath-connected and inertial regimes.
D'Ippolito et al. 2011, Phys. Plasmas 18, 060501 — two-regime picture: • sheath regime (small blobs): v_b ∝ δ_b^{-1/2} • inertial regime (large blobs): v_b ∝ δ_b^{+1/2}
Blob critical radius δ_b* at sheath–inertial transition [m].
Myra et al. 2006, Phys. Plasmas 13, 112502, Eq. 12: δ_b* ≈ 2 ρ_s (L_∥ / (R ρ_s))^(1/5)
Radial blob velocity in the sheath-connected regime [m/s].
Krasheninnikov 2001, Phys. Lett. A 283, 368, Eq. 5: v_b = 2 T_e / (e B R δ_b) × (c_s / Ω_i) ≡ 2 c_s ρ_s / (R √(δ_b / R)) [in the limit δ_b ≪ R]
Radial blob velocity in the inertia-limited regime [m/s].
D'Ippolito et al. 2011, Phys. Plasmas 18, 060501, Eq. 4: v_b ∝ (δ_b / R)^{1/2} c_s
Select sheath or inertial velocity based on δ_b vs δ_b*.
D'Ippolito et al. 2011, Phys. Plasmas 18, 060501, Fig. 2.
Statistical ensemble of blobs with gamma-distributed amplitudes.
Garcia et al. 2012, Phys. Plasmas 19, 092306: SOL density PDF is a gamma distribution; log-normal is the large-intermittency limit. Exponential waiting times produce a Poisson blob arrival process.
staticmethod
¶
SOL density profile with blob-enhanced transport.
Without blobs: n = n₀ exp(−r / λ_n). Blob transport broadens the profile via an effective λ. D'Ippolito et al. 2011, Phys. Plasmas 18, 060501, Sec. IV.
Threshold-crossing blob detector with conditional averaging.
Garcia et al. 2012, Phys. Plasmas 19, 092306 — standard method for identifying intermittent transport events in probe signals.
Checkpoint¶
Utilities for saving and restoring simulation and campaign states.
Enables resuming long-running stress campaigns and persisting solver hot-starts across sessions.
Disruption Sequence¶
Thermal-quench, current-quench, runaway-electron, and halo-current disruption sequence model.
Elm Model¶
ELM crash, recovery, suppression, and pedestal-coupling model utilities.
Stability against peeling-ballooning modes.
ELM onset: α > α_crit (ballooning) AND j_edge > j_peel (peeling). Snyder et al. 2002, Phys. Plasmas 9, 2037, Eq. 8.
Critical edge current density j_peel.
Geometry- and mode-aware peeling threshold proxy:
j_peel,crit ∝ (R0/a) * F(κ, δ) / (q95 * sqrt(n_mode))
where F(κ, δ) captures edge-shaping stabilisation and n_mode
captures reduced stability margin for higher-n peeling harmonics.
This keeps the Snyder-type coupled PB scaling structure while
avoiding a single-parameter 1/q95 closure.
Critical normalised pressure gradient α_crit.
Shaping factor (1 + κ²(1 + 2δ²)) from Sauter et al. 1999, Phys. Plasmas 6, 2834.
Combined peeling-ballooning criterion (elliptical boundary). Snyder et al. 2002, Phys. Plasmas 9, 2037, Eq. 8.
Distance to boundary (positive = stable).
Applies a Type I ELM crash to pedestal profiles.
ΔW_ELM / W_ped ≈ 0.04–0.15 (Type I). Loarte et al. 2003, PPCF 45, 1549, Fig. 12.
ΔW_ELM = f × W_ped; T and n drop by √(1 − f) (W ∝ n T). Loarte et al. 2003, PPCF 45, 1549, Fig. 12.
Bases: ELMCrashModel
Type III ELM: smaller energy loss, higher frequency, triggered at lower pedestal pressure. Zohm 1996, PPCF 38, 105.
Resonant Magnetic Perturbation effects on ELMs.
σ_Chir = Σ w_mn / dr_mn; w_mn = 4√(R0 q' δB_r / (n B0)).
f_ELM ∝ (P_loss − P_LH) / ΔW_ELM ≈ P_SOL / (f × W_ped). Leonard et al. 1999, J. Nucl. Mater. 266-269, 109.
Eped Pedestal¶
EPED-style pedestal prediction and validation-point utilities.
Bases: BaseModel
Pydantic v2 schema for EPED pedestal operating-point configuration.
dataclass
¶
EPEDConfig(R0, a, B0, kappa, delta, Ip_MA, ne_ped_19, B_pol_ped, C_KBM=C_KBM_DEFAULT, n_mode_min=5, n_mode_max=30, nu_star_e=0.0)
classmethod
¶
Return the JSON Schema for serialized EPED pedestal configurations.
Generates mtanh profiles from EPED predictions.
Produce (Te, ne) via mtanh.
rho_sym = 1 − Δ_ped/2 (centre of gradient region)
Wrapper for integrated_transport_solver compatibility.
Synthetic test data for EPED model validation. Values are approximate and do not represent real experimental measurements.
Self-consistent (p_ped, Δ_ped) from PB + KBM constraints. [S09]
Iteration: guess Δ → α_crit(Δ) → p_ped → β_p,ped → Δ_KBM. Collisionality correction applied after convergence. [S11]
GK CGYRO¶
CGYRO external solver interface.
Reference: Candy & Waltz, J. Comp. Phys. 186 (2003) 545.
GK GENE¶
GENE (Gyrokinetic Electromagnetic Numerical Experiment) external solver.
Generates parameters_gene namelist, executes via subprocess, parses
nrg_xxxx eigenvalue files.
Reference: Jenko et al., Phys. Plasmas 7 (2000) 1904.
GK Geometry¶
Miller parameterisation of local magnetic equilibrium geometry for flux-tube gyrokinetic calculations.
Computes metric coefficients, field-line curvature, and the Jacobian on a ballooning-angle grid from (R0, a, kappa, delta, q, s_hat, ...).
Reference: Miller et al., Phys. Plasmas 5 (1998) 973.
dataclass
¶
MillerGeometry(theta, R, Z, B_mag, jacobian, g_rr, g_rt, g_tt, metric_determinant, kappa_n, kappa_g, b_dot_grad_theta)
Flux-tube geometry on a ballooning-angle grid.
All arrays have shape (n_theta,).
miller_geometry(R0, a, rho, kappa=1.0, delta=0.0, s_kappa=0.0, s_delta=0.0, q=1.4, s_hat=0.78, alpha_MHD=0.0, dR_dr=0.0, B0=5.3, n_theta=64, n_period=2)
Compute Miller geometry on a ballooning-angle grid.
Parameters¶
R0 : float Major radius [m]. a : float Minor radius [m]. rho : float Normalised flux coordinate (r/a). kappa, delta : float Elongation and triangularity. s_kappa, s_delta : float Shear of elongation/triangularity: (r/kappa)(dkappa/dr), etc. q, s_hat : float Safety factor and magnetic shear s = (r/q)(dq/dr). alpha_MHD : float Shafranov shift parameter alpha = -q^2 R0 (dp/dr) / (B0^2 / 2mu0). dR_dr : float Shafranov shift gradient dR_axis/dr (typically negative). B0 : float Toroidal field at R0 [T]. n_theta : int Grid points per 2*pi period. n_period : int Number of poloidal periods (ballooning copies).
Circular cross-section limit (kappa=1, delta=0).
Useful for verification against analytic results and the Cyclone Base Case (Dimits et al. 2000).
GK GS2¶
GS2 external solver interface.
Reference: Kotschenreuther et al., Comp. Phys. Comm. 88 (1995) 128.
GK Nonlinear¶
Nonlinear δf gyrokinetic solver in flux-tube geometry.
Solves the gyrokinetic Vlasov equation for the perturbed distribution function δf(k_x, k_y, θ, v_∥, μ) with E×B nonlinearity computed via dealiased 2D FFT (Orszag 1971 2/3 rule).
Physics
- Quasineutrality field solve (adiabatic electrons)
- E×B advection: dealiased Arakawa bracket via FFT
- Parallel streaming: 4th-order compact finite differences
- Magnetic curvature/grad-B drift
- Simplified Sugama collision operator
- 4th-order hyperdiffusion for numerical stability
- RK4 time stepping with CFL-adaptive dt
References
- Dimits et al., Phys. Plasmas 7 (2000) 969 — CBC benchmark
- Orszag, J. Atmos. Sci. 28 (1971) 1074 — dealiasing
- Rosenbluth & Hinton, Phys. Rev. Lett. 80 (1998) 724 — zonal flows
- Sugama & Watanabe, Phys. Plasmas 13 (2006) 012501 — collisions
dataclass
¶
NonlinearGKConfig(n_kx=16, n_ky=16, n_theta=64, n_vpar=16, n_mu=8, n_species=2, dt=0.05, n_steps=5000, save_interval=100, Lx=80.0, Ly=62.83, vpar_max=3.0, mu_max=9.0, dealiasing='2/3', hyper_order=4, hyper_coeff=0.1, cfl_factor=0.5, cfl_adapt=True, collisions=True, nu_collision=0.01, collision_model='krook', nonlinear=True, kinetic_electrons=False, mass_ratio_me_mi=1.0 / 400.0, implicit_electrons=False, electromagnetic=False, beta_e=0.01, d_e_sq=0.0, R0=2.78, a=1.0, B0=2.0, q=1.4, s_hat=0.78, R_L_Ti=6.9, R_L_Te=6.9, R_L_ne=2.2)
Grid and physics parameters for the nonlinear solver.
dataclass
¶
Full 5D+1 state of the nonlinear solver.
dataclass
¶
NonlinearGKResult(chi_i, chi_e, chi_i_gB=0.0, Q_i_t=(lambda: empty(0))(), Q_e_t=(lambda: empty(0))(), phi_rms_t=(lambda: empty(0))(), zonal_rms_t=(lambda: empty(0))(), time=(lambda: empty(0))(), converged=False, final_state=None)
Time-averaged transport and diagnostics.
Nonlinear δf gyrokinetic solver.
Solve quasineutrality for φ(k_x, k_y, θ).
Adiabatic: [(1-Γ₀_i) + 1] φ = ∫ J₀_i h_i dv Kinetic e: (1-Γ₀_i + 1-Γ₀_e) φ = ∫ J₀_i h_i dv - ∫ J₀_e h_e dv
Ampere's law for A_∥: k_perp² A_∥ = β_e Σ_s q_s ∫ v_∥ J₀ h_s dv.
Returns A_par(n_kx, n_ky, n_theta). Zero when electromagnetic=False.
Poisson bracket {φ, f} = ∂φ/∂x ∂f/∂y - ∂φ/∂y ∂f/∂x.
Computed via 2D FFT with 2/3 dealiasing (Orszag 1971). phi: (nkx, nky, nθ) f_s: (nkx, nky, nθ, nvpar, nμ) Returns: (nkx, nky, nθ, nvpar, nμ)
4th-order FD for v_∥ b·∇θ ∂f/∂θ with ballooning connection BC.
At θ boundaries, kx shifts by ±s_hat×ky per poloidal turn. f_s: (nkx, nky, nθ, nvpar, nμ)
Curvature and grad-B drift contribution.
ω_D = k_y × 2(E/T) × [κ_n(1-λB/B₀) + κ_g√(1-λB/B₀)]
Background gradient drive: -ik_y ω_* × (φ - v_∥ A_∥) × F_M.
EM contribution: the effective potential is φ - v_∥ A_∥/c, which adds the magnetic flutter drive for KBM/MTM at finite β.
Ion and electron heat flux in gyro-Bohm units.
Q_i = Re[Σ_{ky>0} ik_y conj(φ) p_i] — radial E×B flux of pressure. The ik_y factor gives the radial component of v_E = ∇φ×b̂/B.
Single-mode initial condition for linear growth rate recovery.
GK Online Learner¶
OnlineLearner admits finite nonnegative transport targets only when the
caller-supplied OOD score is inside the configured threshold. Retraining uses a
validation holdout, rolls back on non-improvement, and can persist an auditable
JSON report containing every accepted or rejected update decision.
Online learning layer for surrogate transport model improvement.
Accumulates GK spot-check results as training data and periodically fine-tunes the QLKNN surrogate. Includes validation holdout and automatic rollback if performance degrades.
dataclass
¶
Single (input, target) pair from a GK spot-check.
dataclass
¶
LearnerConfig(buffer_size=100, validation_fraction=0.2, n_epochs=10, learning_rate=0.0001, max_generations=50, max_ood_score=3.0, min_validation_improvement=0.0)
Online learner parameters.
dataclass
¶
RetrainDecision(generation, accepted, val_loss, previous_best_val_loss, train_samples, validation_samples, buffer_size, reason, max_ood_score, learning_rate, n_epochs)
Auditable retraining admission decision.
GK QuaLiKiz¶
QuaLiKiz quasilinear gyrokinetic solver interface.
Unlike TGLF/GENE/GS2/CGYRO, QuaLiKiz has a Python API (qlknn-hyper). Falls back to subprocess if the Python module is unavailable.
Reference: Bourdelle et al., Phys. Plasmas 14 (2007) 112501.
GK Species¶
Species definitions and collision operator for the linear GK solver.
Velocity-space discretisation uses (energy, lambda) coordinates with Gauss-Legendre quadrature. The collision operator implements a bounded Sugama-style test-particle operator with separate pitch-angle deflection and thermal energy-relaxation coefficients.
Limitations
- Test-particle collision coefficients only (no field-particle response)
- No inter-species collisions (no electron-ion drag)
- Not the full linearised Fokker-Planck operator This bounded operator is adequate for local ITG/TEM regression ordering but not for quantitative collisional damping rates.
References
- Sugama & Watanabe, Phys. Plasmas 13 (2006) 012501
- Abel et al., Rep. Prog. Phys. 76 (2013) 116201, §4.4
dataclass
¶
Single plasma species for the GK solver.
Parameters¶
mass_amu : float Mass in atomic mass units (2.0 for deuterium). charge_e : float Charge in units of e (+1 for ions, -1 for electrons). temperature_keV : float Temperature [keV]. density_19 : float Density [10^19 m^-3]. R_L_T : float Normalised temperature gradient R/L_T. R_L_n : float Normalised density gradient R/L_n. is_adiabatic : bool If True, use adiabatic response instead of kinetic.
dataclass
¶
Energy-lambda velocity-space grid with Gauss-Legendre weights.
Energy E normalised to T_s. Lambda = mu B0 / E ∈ [0, 1].
dataclass
¶
Normalised species diamagnetic frequencies for GK drive terms.
Frequencies are normalised to v_th / R. Positive ion temperature or
density gradients therefore propagate in the ion diamagnetic direction,
while electron frequencies carry the opposite sign through charge_e.
Return density, temperature, and pressure diamagnetic frequencies.
The returned values are dimensionless frequencies normalised to v_th/R:
omega_*n = -sign(q_s) k_y rho_s R/L_n
omega_*T = -sign(q_s) k_y rho_s R/L_T
omega_*p = omega_*n + omega_*T
This is a local-drive bookkeeping utility for bounded GK regression and controller-adapter contracts; it does not introduce an external-code validation claim.
J_0(x) via scipy or polynomial approximation.
For the GK solver, x = k_perp * rho_s * sqrt(2 * lambda * E). At small x (long wavelength limit), J_0 → 1.
Compute deflection and energy-diffusion collision frequencies.
Returns (nu_D, nu_E) normalised to v_th / R.
Bounded Sugama-style test-particle coefficients
nu_D is the pitch-angle deflection rate. nu_E is a thermal energy-relaxation rate with the expected electron-ion mass-ratio suppression for heavy species and field temperature dependence. Field-particle and inter-species momentum conservation terms are intentionally excluded by the traceability contract.
Pitch-angle scattering operator matrix L_pitch on lambda grid.
L_pitch f = (1/2) d/d(xi) (1-xi^2) df/d(xi) where xi = v_∥ / v = sqrt(1 - lambda * B/B0).
Returns shape (n_lambda, n_lambda).
GK TGLF¶
TGLF (Trapped Gyro-Landau Fluid) external solver interface.
Generates TGLF input namelists, executes the tglf binary via
subprocess, and parses growth-rate / flux output files.
Reference: Staebler et al., Phys. Plasmas 14 (2007) 055909.
Bases: GKSolverBase
TGLF external solver via GACODE tglf binary.
Parameters¶
binary : str
Path or name of the tglf executable.
work_dir : Path or None
Persistent working directory. If None, uses a tempdir per call.
Bases: RuntimeError
Raised when external TGLF execution or output validation fails.
Render a TGLF namelist string from local plasma parameters.
Parse TGLF output files from run_dir.
Expects out.tglf.run with growth rates and out.tglf.transport
with quasilinear fluxes. If files are missing, returns a zero-flux
unconverged result.
GK TGLF Native¶
Native TGLF-equivalent quasilinear transport model with SAT0/SAT1/SAT2 spectral saturation, E×B shear quench, multi-scale ITG-ETG coupling, and trapped-particle damping. No external binary.
References
- Staebler et al., Phys. Plasmas 14 (2007) 055909 — SAT0/SAT1
- Staebler et al., Phys. Plasmas 24 (2017) 055906 — SAT2
- Maeyama et al., Phys. Rev. Lett. 114 (2015) 255002 — cross-scale
- Waltz et al., Phys. Plasmas 4 (1997) 2482 — E×B shear quench
- Connor et al., Nucl. Fusion 14 (1974) 185 — trapped-particle modes
dataclass
¶
TGLFNativeConfig(sat_model='SAT1', exb_shear_model='linear', multiscale=False, n_ky_ion=16, n_ky_etg=0, n_theta=32, alpha_exb=ALPHA_EXB_DEFAULT, alpha_cs=ALPHA_CS_DEFAULT)
SAT model selection and solver parameters.
dataclass
¶
TGLFNativeResult(chi_i, chi_e, D_e, D_i=0.0, V_e=0.0, k_y=(lambda: empty(0))(), gamma=(lambda: empty(0))(), omega_r=(lambda: empty(0))(), gamma_net=(lambda: empty(0))(), phi_sq=(lambda: empty(0))(), gamma_exb=0.0, dominant_mode='stable', sat_model='SAT0', converged=True, chi_e_etg=0.0)
Transport fluxes and spectral diagnostics.
Bases: GKSolverBase
Native TGLF-equivalent quasilinear transport model.
Pure Python implementation of SAT0/SAT1/SAT2 spectral saturation applied to the native linear GK eigenvalue spectrum.
E×B shearing rate normalised to c_s/a.
Proxy: γ_ExB ≈ |s_hat/q| × ε × R/L_Ti × 0.1 Waltz et al. 1997, Eq. (3).
Growth-rate damping factor from trapped particles.
Connor et al. 1974: scales with f_t × ν*/ω_bounce. Returns multiplicative factor in (0.1, 1.0].
Normalised spectral intensity I_k = (γ_net/k_y) / Σ(γ_net/k_y).
Staebler 2007, Eq. (7).
SAT0: mixing-length saturation with E×B shear quench.
φ²_k = γ_net / (k_y² |ω_r|), where γ_net = max(γ·tp - α|γ_ExB|, 0).
SAT1: spectral saturation with E×B shear quench.
Staebler et al. 2007: amplitude set by peak mode, distributed by spectral weight. φ²_k = I_k × γ_max_net / k_y_max².
SAT2: multi-scale spectral saturation with ITG-ETG cross-scale coupling.
Staebler et al. 2017: ion-scale follows SAT1. ETG modes enhanced by α_cs × (γ_ETG / γ_ITG_max) from Maeyama et al. 2015.
Velocity-space integrated QL weights → physical fluxes.
Returns (chi_i, chi_e, D_e, V_e, chi_e_etg) in [m²/s].
GK Verification Report¶
Per-session verification report for hybrid surrogate+GK transport.
Tracks spot-check statistics, correction factors, OOD trigger rates, and error distributions across a simulation session.
dataclass
¶
VerificationReport(total_steps=0, steps_verified=0, total_spot_checks=0, ood_triggers=0, records=list(), correction_factors=list())
Accumulated verification statistics for one simulation session.
Impurity Transport¶
Impurity transport, tungsten accumulation, and radiative-loss diagnostic utilities.
Parametric cooling rate coefficient L_Z(T_e) [W m³].
Log-normal fits to
W: Putterich et al. 2010, Nucl. Fusion 50, 025012, Fig. 3. C/Ar/Ne: Post et al. 1977, At. Data Nucl. Data Tables 20, 397.
1D transport advance for each impurity species.
Flux: Γ = -(D_anom + D_neo) dn/dr + V n Upwind differencing for convection; Crank-Nicolson implicit for diffusion.
Neoclassical impurity pinch velocity V_neo [m/s].
Full formula (banana regime, trace impurity limit): V_neo,Z = -D_neo [Z ∇n/n + (Z/2 - H_Z) ∇T_i/T_i] where H_Z = 0.5 is the ion screening factor. Hirshman & Sigmar 1981, Nucl. Fusion 21, 1079, Eq. 4.17.
The thermal-gradient term dominates for peaked T_i profiles and produces an inward (negative) pinch when dT_i/dr < 0. Sign convention: V < 0 → inward (toward axis).
D_neo = q² ρ_i² ν_ii / ε^(3/2) (banana regime) Hinton & Hazeltine 1976, Rev. Mod. Phys. 48, 239, Eq. 4.62. Nominal scale D_neo ≈ 0.1 m²/s for ITER-relevant parameters.
Total radiated power P_rad [MW] integrated over plasma volume.
Local power density
p_rad(r) = n_e(r) n_Z(r) L_Z(T_e(r)) [W m^-3]
Post et al. 1977, At. Data Nucl. Data Tables 20, 397.
Volume element for circular cross-section torus
dV = 4π² R0 a² ρ dρ
JAX GK Nonlinear¶
JAX-accelerated nonlinear δf gyrokinetic solver.
Wraps the same physics as gk_nonlinear.py but uses JAX for: - FFT-based E×B bracket (~100× over NumPy on GPU) - jax.checkpoint for memory-efficient RK4
Falls back to the NumPy solver when JAX is unavailable.
JAX-accelerated nonlinear δf solver.
Uses JAX for FFT, lax.scan for the time loop, and jax.checkpoint for RK4 memory efficiency. Falls back to NumPy when JAX is absent.
JAX GK Solver¶
The linear JAX GK solver includes a schema-versioned parity artifact producer
for backend reproducibility. build_jax_gk_parity_artifact() and
write_jax_gk_parity_artifact() bind the native local-dispersion comparison,
backend metadata, dtype/X64 state, solver kwargs, tolerances, and canonical
payload SHA-256 digest while preserving the backend-parity-only claim boundary.
Artifacts also bind case-parameter digests, native/JAX mode spectra, dominant
mode labels, and case acceptance limits for CBC, kinetic-electron TEM, and
low-drive stable-mode parity evidence. validation/validate_jax_gk_parity.py
can require named cases and named backends before admitting an evidence
directory, so archived single-case artifacts cannot be replayed as full parity
coverage.
JAX-accelerated linear GK eigenvalue solver.
Re-implements the native local-dispersion eigenvalue contract from gk_eigenvalue.py using jax.numpy for batched velocity-space quadrature and jax.grad for the transport-stiffness proxy d(chi_i)/d(R_L_Ti). JAX execution is explicit and fails closed when JAX is unavailable.
References
- Dimits et al., Phys. Plasmas 7 (2000) 969 — Cyclone Base Case
- Kotschenreuther et al., Comp. Phys. Comm. 88 (1995) 128
solve_linear_gk_jax(species_list=None, geom=None, vgrid=None, R0=2.78, a=1.0, B0=2.0, q=1.4, s_hat=0.78, Z_eff=1.0, nu_star=0.01, n_ky_ion=16, n_theta=64)
JAX-accelerated linear GK solver. Batched over k_y via vmap.
Parameters match solve_linear_gk from gk_eigenvalue.py. Uses JAX to build the local-dispersion velocity integrals for all k_y simultaneously, then applies the same damped Newton root finder used by the native solver. This keeps backend parity on the physical growth/frequency path while preserving JAX acceleration for the high-volume quadrature assembly.
transport_stiffness_jax(R_L_Ti, R0=2.78, a=1.0, B0=2.0, q=1.4, s_hat=0.78, Z_eff=1.0, nu_star=0.01, n_ky_ion=8, n_theta=32)
d(chi_i_proxy)/d(R_L_Ti) via jax.grad.
Returns the transport stiffness — sensitivity of the ion heat flux proxy to the normalised temperature gradient. Positive above critical gradient, near-zero below.
gk_stiffness_chi_i_profile_jax(R_L_Ti, rho, base_chi_i=0.1, stiffness_scale=1e-06, R0=2.78, a=1.0, B0=2.0, q=1.4, s_hat=0.78, Z_eff=1.0, nu_star=0.01, n_ky_ion=8, n_theta=32)
Map JAX GK stiffness into a bounded ion heat diffusivity profile.
The closure is intentionally conservative: it uses the existing
differentiable transport_stiffness_jax scalar as a local stiffness
amplitude and applies a smooth edge-weighted radial shape. It is a
controller-tuning closure, not a replacement for external GK validation.
build_jax_gk_parity_artifact(*, case='cyclone_base_case', solver_kwargs=None, native_result=None, jax_result=None, gamma_relative_tolerance=0.25, omega_absolute_tolerance=0.25, executed_at=None)
Build a tamper-evident JAX/native gyrokinetic parity artifact.
The artifact is backend-parity evidence for the repository's native local-dispersion formulation. It is intentionally not a cross-code validation artifact and always records that external validation remains required before quantitative facility claims.
write_jax_gk_parity_artifact(output_path, *, case='cyclone_base_case', solver_kwargs=None, gamma_relative_tolerance=0.25, omega_absolute_tolerance=0.25, executed_at=None)
Persist a JAX/native GK parity artifact and return payload plus path.
gk_stiffness_chi_i_profile_jax(R_L_Ti, rho, base_chi_i=0.1, stiffness_scale=1e-06, R0=2.78, a=1.0, B0=2.0, q=1.4, s_hat=0.78, Z_eff=1.0, nu_star=0.01, n_ky_ion=8, n_theta=32)
Map JAX GK stiffness into a bounded ion heat diffusivity profile.
The closure is intentionally conservative: it uses the existing
differentiable transport_stiffness_jax scalar as a local stiffness
amplitude and applies a smooth edge-weighted radial shape. It is a
controller-tuning closure, not a replacement for external GK validation.
build_jax_gk_parity_artifact(*, case='cyclone_base_case', solver_kwargs=None, native_result=None, jax_result=None, gamma_relative_tolerance=0.25, omega_absolute_tolerance=0.25, executed_at=None)
Build a tamper-evident JAX/native gyrokinetic parity artifact.
The artifact is backend-parity evidence for the repository's native local-dispersion formulation. It is intentionally not a cross-code validation artifact and always records that external validation remains required before quantitative facility claims.
write_jax_gk_parity_artifact(output_path, *, case='cyclone_base_case', solver_kwargs=None, gamma_relative_tolerance=0.25, omega_absolute_tolerance=0.25, executed_at=None)
Persist a JAX/native GK parity artifact and return payload plus path.
JAX GS Solver¶
JAX-differentiable fixed-boundary Grad-Shafranov equilibrium solver.
Picard iteration with damped Jacobi inner sweeps, all in pure JAX for
JIT compilation and automatic differentiation. Enables jax.grad
through the full equilibrium solve — closing the autodiff depth gap
with TORAX (JAX) and FUSE (Julia AD).
The GS* operator in cylindrical (R, Z) coordinates:
d²ψ/dR² - (1/R) dψ/dR + d²ψ/dZ² = -μ₀ R J_φ
where J_φ = R p'(ψ_norm) + FF'(ψ_norm) / (μ₀ R).
Key function
jax_gs_solve — Fixed-boundary GS equilibrium, differentiable via jax.grad
NumPy fallback provided when JAX is unavailable.
gs_solve_np(R_min, R_max, Z_min, Z_max, NR, NZ, Ip_target, mu0=4e-07 * pi, n_picard=80, n_jacobi=200, alpha=0.1, omega_j=0.667, beta_mix=0.5)
Fixed-boundary GS solve via Picard iteration (NumPy fallback).
Returns psi on the (NZ, NR) grid with zero Dirichlet boundary.
jax_gs_solve(R_min=0.1, R_max=2.0, Z_min=-1.5, Z_max=1.5, NR=33, NZ=33, Ip_target=1000000.0, mu0=4e-07 * pi, n_picard=80, n_jacobi=200, alpha=0.1, omega_j=0.667, beta_mix=0.5, *, use_jax=True, allow_numpy_fallback=False, allow_legacy_numpy_fallback=False)
Fixed-boundary Grad-Shafranov equilibrium solve.
Picard iteration with damped Jacobi inner sweeps. When JAX is
available, the solve is JIT-compiled and differentiable via
jax.grad. Falls back to NumPy otherwise.
Parameters¶
R_min, R_max : radial domain [m] Z_min, Z_max : vertical domain [m] NR, NZ : grid resolution Ip_target : target plasma current [A] mu0 : vacuum permeability [H/m] n_picard : outer Picard iterations (fixed count for differentiability) n_jacobi : inner Jacobi sweeps per Picard step alpha : Picard under-relaxation factor omega_j : Jacobi damping (2/3 standard, Hackbusch 1985 §4.3) beta_mix : pressure vs current profile weight use_jax : attempt JAX backend
Returns¶
psi : (NZ, NR) poloidal flux, zero on boundary
jax_gs_solve_from_grid(R_grid, psi_init, dR, dZ, Ip_target=1000000.0, mu0=4e-07 * pi, n_picard=80, n_jacobi=200, alpha=0.1, omega_j=0.667, beta_mix=0.5)
GS solve on a pre-existing grid. JAX-only, differentiable.
Parameters¶
R_grid : (NZ, NR) meshgrid of major radius values psi_init : (NZ, NR) initial flux with boundary conditions applied dR, dZ : grid spacing Ip_target : target plasma current [A]
jax_gs_grad_Ip(Ip_target, R_min=0.1, R_max=2.0, Z_min=-1.5, Z_max=1.5, NR=33, NZ=33, mu0=4e-07 * pi, n_picard=40, n_jacobi=100, alpha=0.1, omega_j=0.667, beta_mix=0.5)
Compute d(sum(psi))/d(Ip_target) via JAX autodiff.
Demonstrates full-stack differentiability through the GS solve.
Kinetic EFIT¶
Kinetic-EFIT constraints and reconstruction wrapper for pressure and q-profile coupling.
Anisotropic fast ion pressure model.
dataclass
¶
KineticEFITClaimEvidence(schema_version, source, source_id, diagnostic_source, profile_source, fast_ion_source, mse_calibration_source, model_id, interpolation_geometry, n_te_points, n_ne_points, n_ti_points, n_mse_points, fast_ion_energy_keV, fast_ion_density_fraction, anisotropy_sigma, pressure_consistency, beta_fast, q_axis, q_edge, pressure_relative_error, q_profile_relative_error, anisotropy_abs_error, pressure_relative_tolerance, q_profile_relative_tolerance, anisotropy_abs_tolerance, facility_claim_allowed, claim_status)
Serialisable provenance and reference-comparison evidence for kinetic-EFIT claims.
kinetic_efit_claim_evidence(result, kinetic, fast_ions, *, source, source_id, diagnostic_source, profile_source, fast_ion_source, mse_calibration_source, model_id='bounded_kinetic_efit', reference_pressure=None, reference_q_profile=None, reference_anisotropy_sigma=None, pressure_relative_tolerance=0.05, q_profile_relative_tolerance=0.05, anisotropy_abs_tolerance=0.05)
Build fail-closed kinetic-EFIT evidence for bounded or facility claims.
Facility claims require a recognised facility reference source plus matched pressure, q-profile, and anisotropy references inside declared tolerances. Synthetic regression evidence is retained as bounded controller evidence only and is never promoted to a facility-grade P-EFIT claim.
Raise when kinetic-EFIT evidence is insufficient for a facility claim.
Persist kinetic-EFIT claim evidence as deterministic JSON.
Forward motional-Stark-effect pitch-angle estimate for a tangential beam.
dataclass
¶
KineticEFITClaimEvidence(schema_version, source, source_id, diagnostic_source, profile_source, fast_ion_source, mse_calibration_source, model_id, interpolation_geometry, n_te_points, n_ne_points, n_ti_points, n_mse_points, fast_ion_energy_keV, fast_ion_density_fraction, anisotropy_sigma, pressure_consistency, beta_fast, q_axis, q_edge, pressure_relative_error, q_profile_relative_error, anisotropy_abs_error, pressure_relative_tolerance, q_profile_relative_tolerance, anisotropy_abs_tolerance, facility_claim_allowed, claim_status)
Serialisable provenance and reference-comparison evidence for kinetic-EFIT claims.
kinetic_efit_claim_evidence(result, kinetic, fast_ions, *, source, source_id, diagnostic_source, profile_source, fast_ion_source, mse_calibration_source, model_id='bounded_kinetic_efit', reference_pressure=None, reference_q_profile=None, reference_anisotropy_sigma=None, pressure_relative_tolerance=0.05, q_profile_relative_tolerance=0.05, anisotropy_abs_tolerance=0.05)
Build fail-closed kinetic-EFIT evidence for bounded or facility claims.
Facility claims require a recognised facility reference source plus matched pressure, q-profile, and anisotropy references inside declared tolerances. Synthetic regression evidence is retained as bounded controller evidence only and is never promoted to a facility-grade P-EFIT claim.
Raise when kinetic-EFIT evidence is insufficient for a facility claim.
Persist kinetic-EFIT claim evidence as deterministic JSON.
L-H Transition¶
L-H transition trigger and predator-prey edge-transition controller utilities.
PredatorPreyModel(gamma_L=_KD_GAMMA_L, alpha1=_KD_ALPHA1, alpha2=_KD_ALPHA2, alpha3=_KD_ALPHA3, gamma_damp=_KD_GAMMA_DAMP)
Zonal-flow turbulence model — Kim & Diamond (2003).
Characteristic ZF growth rate γ_ZF ≡ α₃ / (α₁ γ_damp) * γ_L² — Kim & Diamond (2003), Eq. 6.
Used by IPhaseFrequency to estimate the I-phase limit-cycle frequency.
L-H power threshold scaling — Martin et al. (2008).
staticmethod
¶
P_LH [MW] via Martin scaling.
Martin et al., J. Phys.: Conf. Ser. 123, 012033 (2008), Eq. (1).
Parameters¶
ne_19: Line-averaged electron density in units of 10¹⁹ m⁻³. B_T: Toroidal field in T. S_m2: Plasma surface area in m².
Notes¶
The published fit uses n₂₀ = n_e / 10²⁰ m⁻³. Since ne_19 is in 10¹⁹ m⁻³, n₂₀ = ne_19 / 10.
staticmethod
¶
P_LH including the low-density upswing branch.
Below n_min ≈ 0.4 n_GW the threshold increases as (n_min/n_e)² (Ryter et al., Nucl. Fusion 54, 083003 (2014), Fig. 3).
Parameters¶
ne_19: Line-averaged electron density in units of 10¹⁹ m⁻³. B_T: Toroidal field in T. S_m2: Plasma surface area in m². I_p_MA: Plasma current in MA (for Greenwald density). a_m: Minor radius in m (for Greenwald density).
Notes¶
Greenwald density: n_GW [10²⁰ m⁻³] = I_p [MA] / (π a²) [MA m⁻²]. Converted to 10¹⁹ m⁻³: n_GW_19 = 10 * I_p / (π a²).
Detects I-phase limit-cycle oscillations in turbulence traces.
I-phase limit-cycle frequency estimate.
Tynan et al., Nucl. Fusion 53, 073053 (2013), Sec. 4. f_I ≈ γ_ZF / (2π), where γ_ZF is the zonal-flow linear growth rate evaluated at the turbulent saturation level.
Locked Mode¶
Error-field amplification, mode locking, and post-lock island growth chain.
Key references¶
La Haye 2006 : R.J. La Haye, Phys. Plasmas 13, 055501 (2006). [RFA formula Eq. 6; locking condition Eq. 12; MRE in locked-mode regime Eq. 8; C_lock coefficient] Fitzpatrick 1993: R. Fitzpatrick, Nucl. Fusion 33, 1049 (1993). [electromagnetic torque Eq. 28; locking criterion] Rutherford 1973 : P.H. Rutherford, Phys. Fluids 16, 1903 (1973). [classical tearing stability dw/dt ∝ Delta']
Ideal-MHD resonant-field amplification (RFA) model.
B_res = B_err / (1 - beta_N / beta_N_nowall)
La Haye 2006, Phys. Plasmas 13, 055501, Eq. 6. Denominator → 0 as beta_N → beta_N_nowall (no-wall stability boundary).
Toroidal rotation under electromagnetic braking from a resonant error field.
Electromagnetic torque (per unit toroidal angle): T_em = -n² m ψ_ext² sin(Δφ) / (μ₀ R r_s²) Fitzpatrick 1993, Nucl. Fusion 33, 1049, Eq. 28.
Mode-locking condition
|T_em| > T_viscous where T_viscous = ρ χ_φ ω V_island
La Haye 2006, Phys. Plasmas 13, 055501, Eq. 12.
In steady-state the mode locks when the electromagnetic drag exceeds the viscous restoring torque that maintains plasma rotation.
Electromagnetic braking torque on the (m, n) mode.
Derived from Fitzpatrick 1993, Nucl. Fusion 33, 1049, Eq. 28: T_em = -n² m ψ_ext² sin(Δφ) / (μ₀ R r_s²) where ψ_ext ~ B_res × r_s² / n (flux loop estimate) and sin(Δφ) → 1 at the locking threshold. The sign convention here gives a positive magnitude; the caller interprets the torque as opposing rotation.
Viscous restoring torque.
T_viscous = ρ χ_φ ω V_island La Haye 2006, Phys. Plasmas 13, 055501, Eq. 12.
rho_chi = ρ × χ_φ [kg m⁻¹ s⁻¹], V_island ~ 4 π² R₀ r_s (shell volume at the rational surface, thin-shell limit).
Integrate dω/dt = -(ω - ω₀)/τ_visc - T_em / I_eff.
Mode locks when ω ≤ ω_eq (equilibrium with sustained braking). Locking condition: T_em > T_viscous (La Haye 2006, Eq. 12).
Post-locking island growth via the Modified Rutherford Equation.
dw/dt = (η / μ₀ r_s) [r_s Δ' + C_lock r_s / w]
La Haye 2006, Phys. Plasmas 13, 055501, Eq. 8. C_lock accounts for the locked-mode bootstrap drive (same physical origin as the NTM bootstrap term but evaluated at ω = 0). Classical tearing drive: Rutherford 1973, Phys. Fluids 16, 1903.
Integrate MRE for a locked island; stochastic when w > delta_r_mn.
dw/dt = (η / μ₀ r_s) [r_s Δ' + C_lock r_s / w] La Haye 2006, Phys. Plasmas 13, 055501, Eq. 8.
MARFE¶
MARFE radiation-condensation and density-limit stability utilities.
Radiation condensation growth rate.
γ = -(κ_∥ k_∥² + n_e n_Z dL/dT) / (n c_v)
Positive γ ↔ unstable: requires dL/dT < 0 and n_e n_Z |dL/dT| > κ_∥ k_∥². Drake 1987, Phys. Fluids 30, 2429, Eq. 5.
True when dL/dT < 0 and radiation term exceeds parallel conduction damping.
Estimate T_MARFE as the highest T where dL_Z/dT first turns negative.
Lipschultz 1987, J. Nucl. Mater. 145-147, 15: MARFE forms near the temperature where the radiative cooling function has a negative slope. Returns the onset T in eV, or nan if dL/dT is positive everywhere.
n_crit where γ = 0.
From Drake 1987, Eq. 5 with γ = 0: n_e² f_imp |dL/dT| = κ_∥ k_∥²
MARFE criterion: localised cold spot with T_min < 20 eV and T_max > 50 eV. Lipschultz 1987, J. Nucl. Mater. 145-147, 15.
staticmethod
¶
n_GW = I_p / (π a²) [10^20 m^-3]
Greenwald 2002, Plasma Phys. Control. Fusion 44, R27, Eq. 1. I_p in MA, a in m.
staticmethod
¶
Heuristic MARFE-onset density tied to the Greenwald limit.
n_marfe ~ n_GW * sqrt(P_SOL) / (10 * sqrt(f_imp)) Scaling motivated by Drake 1987 radiation condensation criterion with parallel conduction set by P_SOL and impurity fraction f_imp.
MDSplus Acquisition¶
Optional MDSplus shot acquisition with manifest provenance output.
dataclass
¶
Single MDSplus signal acquisition request.
dataclass
¶
Versioned MDSplus acquisition request loaded from JSON.
dataclass
¶
Paths and manifest identity produced by an acquisition run.
Load a versioned MDSplus acquisition request JSON file.
acquire_mdsplus_shot(*, tree, shot, signals, output_npz, manifest_json, source_uri, access_policy, licence, retrieved_at=None, mdsplus_module=None)
Acquire selected signals from MDSplus, write NPZ, and emit a manifest.
Momentum Transport¶
Toroidal momentum-transport and rotation-profile evolution utilities.
Advance rotation profile one step [rad/s].
Momentum equation (cylindrical approximation): ∂(ρ_m R² ω)/∂t + (1/r) ∂/∂r(r Π_φ) = T_tot Π_φ = −χ_φ ∂(ρ_m R² ω)/∂r (pinch neglected) Wesson 2011, "Tokamaks", §4.10.
NBI torque density [N·m/m³].
T_NBI = P_NBI · R_tang / v_beam where R_tang = R₀ sin(θ). Stacey & Sigmar 1985, Phys. Fluids 28, 2800.
Residual-stress torque driving intrinsic counter-current rotation [N·m/m³].
Residual stress Π_res ∝ −∇T_i → counter-current rotation. Empirical Rice scaling: v_φ,intr ∝ W_p / I_p. Rice et al. 2007, Nucl. Fusion 47, 1618, Eq. 3.
E×B shearing rate [rad/s].
ω_E×B = (R B_θ / B) · d/dr (E_r / R B_θ) Burrell 1997, Phys. Plasmas 4, 1499, Eq. 1.
In the rotation-dominated limit E_r ≈ v_φ B_θ = R₀ ω_φ B_θ, so E_r / (R₀ B_θ) ≈ ω_φ and ω_E×B ≈ (R₀ B_θ / B) · dω_φ/dr.
Suppression factor on anomalous transport.
F = 1 / (1 + (ω_E×B / γ_max)²) Biglari, Diamond & Terry 1990, Phys. Fluids B 2, 1.
E_r [V/m] from radial force balance (neoclassical), Z_i = 1.
E_r = (1 / Z_i e n_i) dp_i/dr + v_φ B_θ (v_θ ≈ 0) Hinton & Hazeltine 1976, Rev. Mod. Phys. 48, 239, Eq. 2.3.
Intrinsic toroidal velocity [km/s] from Rice scaling.
v_φ,intr ≈ c_Rice · (W_p / I_p) Rice et al. 2007, Nucl. Fusion 47, 1618, Eq. 3.
c_Rice = 3.5 km s⁻¹ MA MJ⁻¹ (empirical, Ohmic + NBI plasmas, Fig. 4).
Neoclassical¶
Neoclassical transport and bootstrap current models.
Implements Chang-Hinton ion thermal diffusivity, Pfirsch-Schlüter diffusivity, collisionality regime auto-detection, and the full three-coefficient Sauter bootstrap current model for general axisymmetry.
Key references¶
Chang & Hinton 1982 : C.S. Chang, F.L. Hinton, Phys. Fluids 25, 1493 (1982). Hinton & Hazeltine 1976 : F.L. Hinton, R.D. Hazeltine, Rev. Mod. Phys. 48, 239 (1976). Sauter 1999 : O. Sauter et al., Phys. Plasmas 6, 2834 (1999). Wesson 2011 : J. Wesson, "Tokamaks", 4th ed., Oxford (2011). Galeev & Sagdeev 1979 : A.A. Galeev, R.Z. Sagdeev, in Reviews of Plasma Physics, vol. 7 (1979).
Dimensionless neoclassical collisionality nu_star.
nu_star = nu_ii * q * R / (epsilon^1.5 * v_th)
Hinton & Hazeltine 1976, Eq. 4.53.
Parameters¶
n_e_19 : float Electron density [10^19 m^-3]. T_kev : float Species temperature [keV]. q : float Safety factor. R : float Major radius [m]. epsilon : float Inverse aspect ratio r/R. mass_amu : float Species mass in AMU (default 2.0 for D). z_eff : float Effective charge.
Returns¶
float Dimensionless collisionality.
Chang-Hinton banana-regime ion thermal diffusivity.
Chang & Hinton 1982, Phys. Fluids 25, 1493, Eq. 10.
Parameters¶
q : float Safety factor. epsilon : float Inverse aspect ratio r/R. nu_star : float Ion collisionality. rho_i : float Ion Larmor radius [m]. nu_ii : float Ion-ion collision frequency [s^-1].
Returns¶
float Ion thermal diffusivity [m^2/s].
Plateau-regime ion thermal diffusivity.
Valid for 1 < nu_star < q^2 / epsilon^1.5. Hinton & Hazeltine 1976, Rev. Mod. Phys. 48, 239, Eq. 4.66.
Parameters¶
q : float Safety factor. rho_i : float Ion Larmor radius [m]. v_thi : float Ion thermal velocity [m/s]. R : float Major radius [m].
Returns¶
float Ion thermal diffusivity [m^2/s].
Pfirsch-Schlüter regime ion thermal diffusivity.
Valid for nu_star > q^2 / epsilon^1.5 (collisional regime). Wesson 2011, "Tokamaks" 4th ed., Eq. 14.5.7. Galeev & Sagdeev 1979, Reviews of Plasma Physics, vol. 7.
Parameters¶
q : float Safety factor. rho_i : float Ion Larmor radius [m]. nu_ii : float Ion-ion collision frequency [s^-1].
Returns¶
float Ion thermal diffusivity [m^2/s].
Neoclassical ion thermal diffusivity with automatic regime selection.
Regime boundaries follow Hinton & Hazeltine 1976, Rev. Mod. Phys. 48, 239, Section IV and Wesson 2011, Ch. 14.5: - nu_star < epsilon^1.5 : banana (Chang & Hinton 1982, Eq. 10) - epsilon^1.5 <= nu_star < q^2/epsilon^1.5 : plateau (Hinton & Hazeltine 1976, Eq. 4.66) - nu_star >= q^2/epsilon^1.5 : Pfirsch-Schlüter (Wesson 2011, Eq. 14.5.7)
Parameters¶
nu_star : float Dimensionless collisionality. q : float Safety factor. epsilon : float Inverse aspect ratio r/R. rho_i : float Ion Larmor radius [m]. nu_ii : float Ion-ion collision frequency [s^-1]. v_thi : float Ion thermal velocity [m/s]. R : float Major radius [m].
Returns¶
float Ion thermal diffusivity [m^2/s].
Banana-plateau neoclassical diffusivity with Z_eff dependence.
Hinton & Hazeltine 1976, Rev. Mod. Phys. 48, 239, Eq. 4.61.
Parameters¶
q : float Safety factor. epsilon : float Inverse aspect ratio. nu_star : float Collisionality. z_eff : float Effective charge.
Returns¶
float Ion thermal diffusivity (dimensionless units, normalised form).
Sauter bootstrap current density profile with full L31+L32+L34 coefficients.
Sauter et al. 1999, Phys. Plasmas 6, 2834, Eqs. 14–16.
j_bs = -(p_e / B_pol) * [L31 * d(ln p_e)/dr + L32 * d(ln T_e)/dr + L34 * (T_i/T_e) * d(ln T_i)/dr]
Parameters¶
rho : np.ndarray Normalised radial grid. Te, Ti : np.ndarray Electron and ion temperatures [keV]. ne : np.ndarray Electron density [10^19 m^-3]. q : np.ndarray Safety factor profile. R0, a : float Major and minor radii [m]. B0 : float Toroidal magnetic field [T]. z_eff : float Effective charge.
Returns¶
np.ndarray Bootstrap current density [A/m^2].
Neural Turbulence¶
QLKNN-class neural turbulence surrogate, training, and flux prediction utilities.
dataclass
¶
NeuralTurbulenceClaimEvidence(schema_version, model_id, source, source_id, weights_path, weights_sha256, feature_schema, local_sample_count, local_q_i_rmse_gB, local_q_e_rmse_gB, local_gamma_e_rmse_gB, local_flux_relative_mae, local_critical_gradient_accuracy, reference_source, reference_dataset_id, reference_artifact_sha256, reference_sample_count, q_i_rmse_gB, q_e_rmse_gB, gamma_e_rmse_gB, flux_relative_mae, critical_gradient_accuracy, q_i_rmse_tolerance_gB, q_e_rmse_tolerance_gB, gamma_e_rmse_tolerance_gB, flux_relative_mae_tolerance, critical_gradient_accuracy_min, quantitative_claim_allowed, claim_status)
Serialisable admission evidence for neural-turbulence surrogate claims.
Pure NumPy inference for a QLKNN-like neural network. Predicts turbulent fluxes [Q_i, Q_e, Gamma_e] from 10 parameters. van de Plassche et al., Phys. Plasmas 27, 022310 (2020).
Default construction auto-trains on a Jenko et al. (2001) critical gradient model so predictions are physically meaningful out of the box.
staticmethod
¶
Convert physical profiles into the 10 dimensionless QLKNN inputs.
Compare the surrogate against bounded analytic quasilinear targets.
neural_turbulence_claim_evidence(validation_result, *, source, source_id, model_id='neural_turbulence_qlknn_facade', weights_path=None, reference_artifact_path=None)
Build fail-closed evidence for neural-turbulence quantitative claims.
Return evidence only when strict matched-reference admission passed.
Persist neural-turbulence claim evidence as deterministic JSON.
dataclass
¶
NeuralTurbulenceClaimEvidence(schema_version, model_id, source, source_id, weights_path, weights_sha256, feature_schema, local_sample_count, local_q_i_rmse_gB, local_q_e_rmse_gB, local_gamma_e_rmse_gB, local_flux_relative_mae, local_critical_gradient_accuracy, reference_source, reference_dataset_id, reference_artifact_sha256, reference_sample_count, q_i_rmse_gB, q_e_rmse_gB, gamma_e_rmse_gB, flux_relative_mae, critical_gradient_accuracy, q_i_rmse_tolerance_gB, q_e_rmse_tolerance_gB, gamma_e_rmse_tolerance_gB, flux_relative_mae_tolerance, critical_gradient_accuracy_min, quantitative_claim_allowed, claim_status)
Serialisable admission evidence for neural-turbulence surrogate claims.
Compare the surrogate against bounded analytic quasilinear targets.
neural_turbulence_claim_evidence(validation_result, *, source, source_id, model_id='neural_turbulence_qlknn_facade', weights_path=None, reference_artifact_path=None)
Build fail-closed evidence for neural-turbulence quantitative claims.
Return evidence only when strict matched-reference admission passed.
Persist neural-turbulence claim evidence as deterministic JSON.
Orbit Following¶
Guiding-centre orbit-following, classification, and Monte Carlo ensemble utilities.
dataclass
¶
OrbitFollowingClaimEvidence(schema_version, source, source_id, geometry_source, particle_source, collision_model, loss_boundary_source, model_id, q, rho_L_m, epsilon, banana_width_m, first_orbit_loss_fraction, ensemble_particles, ensemble_loss_fraction, ensemble_passing, ensemble_trapped, ensemble_lost, banana_width_relative_error, loss_fraction_abs_error, banana_width_relative_tolerance, loss_fraction_abs_tolerance, external_orbit_claim_allowed, claim_status)
Serialisable provenance and reference-comparison evidence for orbit-following claims.
Guiding-center equations of motion in (R, Z, phi, v_par) for an axisymmetric tokamak.
Drift decomposition follows Boozer 2004, Rev. Mod. Phys. 76, 1071, Eqs. 16–18: dR/dt = v_par * b_R + v_drift_R dZ/dt = v_par * b_Z + v_drift_Z dphi/dt = v_par * b_phi/R + v_drift_phi/R
where v_drift combines grad-B and curvature drifts
v_drift = (mu/q)(b × ∇B)/B + (m v_par²/q)(b × κ)/B = [(m v_par² + mu B) / (q B²)] (B × ∇B)
The combined form uses the identity κ = (b·∇)b = ∇B/B for large-aspect- ratio equilibria; Boozer 2004, Eq. 18.
Stix slowing-down model for fast ions (e.g. fusion alphas).
All formulae from Stix 1972, Plasma Physics 14, 367.
staticmethod
¶
Critical energy E_crit where electron and ion drag are equal.
Stix 1972, Eq. 12 (thermally-averaged form): E_crit = (A_fast / A_e) * T_e * (3√π / 4 * Σ n_j Z_j² / (n_e A_j))^(2/3)
Simplified to leading term for a single-ion background
E_crit = (A_fast m_p / m_e)^(1/3) * (A_fast / (2 A_bg))^(1/3) * (3√π/4)^(2/3) * T_e
We use the standard Cordey 1981 (Nucl. Fusion 21, 1293) simplification: v_c³ = (3√π / 4) * (m_e / m_fast) * Σ_j (n_j Z_j² / n_e) * v_te³
with v_te = √(2 T_e / m_e), giving E_crit = ½ m_fast v_c².
staticmethod
¶
v_c for a DT-plasma background (Z_bg=1, A_bg=2.5). Stix 1972, Eq. 12.
staticmethod
¶
Spitzer electron slowing-down time τ_s for a fast ion.
Stix 1972, Eq. 7 (see also Wesson 2011, Eq. 7.4.7): τ_s = (6 π^(3/2) ε₀² m_e^(1/2) m_fast T_e^(3/2)) / (n_e Z_fast² e⁴ ln Λ * 2^(1/2))
Pre-evaluated for m_fast = 4 m_p (alpha), Z_fast = 2: τ_s [s] ≈ 0.198 A_fast T_e[keV]^(3/2) / (n_e[10²⁰ m⁻³] Z_fast² ln Λ) with A_fast=4, Z_fast=2, ln Λ = 17. Reference: Stix 1972, Plasma Physics 14, 367, Eq. 7.
staticmethod
¶
Stix slowing-down power loss
dE/dt = -(E / τ_s) * (1 + (E_crit / E)^(3/2))
Reference: Stix 1972, Plasma Physics 14, 367, Eq. 12.
Below E_crit the (E_crit/E)^(3/2) term > 1, so ion drag dominates. Above E_crit the term < 1, so electron drag dominates.
staticmethod
¶
Fraction of power deposited into ions vs electrons.
Stix 1972, Eq. 16: f_ion = v_c³ / (v³ + v_c³) f_electron = v³ / (v³ + v_c³)
At v = v_c: equal partition (f_ion = f_electron = 0.5). At v >> v_c: all to electrons. At v << v_c: all to ions.
orbit_following_claim_evidence(*, source, source_id, geometry_source, particle_source, collision_model, loss_boundary_source, q, rho_L_m, epsilon, first_orbit_loss_fraction, ensemble_result=None, model_id='bounded_guiding_centre_orbit_following', reference_banana_width_m=None, reference_loss_fraction=None, banana_width_relative_tolerance=0.05, loss_fraction_abs_tolerance=0.02)
Build fail-closed evidence for bounded or externally referenced orbit claims.
Raise when orbit-following evidence is insufficient for an external-code claim.
Persist orbit-following claim evidence as deterministic JSON.
Prompt loss fraction for birth alphas on their first orbit.
Scaling: f_lost ≈ (ρ_orbit / a)² / I_p Reference: Goldston 1981, J. Comput. Phys. 43, 61, Eq. 15 (orbit width ∝ ρ_L / √ε, loss ∝ (ρ_orbit/a)²/I_p)
ρ_orbit = ρ_L q / √ε — banana orbit width (see banana_orbit_width) Here we use the simpler large-aspect-ratio estimate ρ_orbit ≈ 2 ρ_L valid for passing particles at the loss boundary (Goldston 1981, §3).
dataclass
¶
OrbitFollowingClaimEvidence(schema_version, source, source_id, geometry_source, particle_source, collision_model, loss_boundary_source, model_id, q, rho_L_m, epsilon, banana_width_m, first_orbit_loss_fraction, ensemble_particles, ensemble_loss_fraction, ensemble_passing, ensemble_trapped, ensemble_lost, banana_width_relative_error, loss_fraction_abs_error, banana_width_relative_tolerance, loss_fraction_abs_tolerance, external_orbit_claim_allowed, claim_status)
Serialisable provenance and reference-comparison evidence for orbit-following claims.
orbit_following_claim_evidence(*, source, source_id, geometry_source, particle_source, collision_model, loss_boundary_source, q, rho_L_m, epsilon, first_orbit_loss_fraction, ensemble_result=None, model_id='bounded_guiding_centre_orbit_following', reference_banana_width_m=None, reference_loss_fraction=None, banana_width_relative_tolerance=0.05, loss_fraction_abs_tolerance=0.02)
Build fail-closed evidence for bounded or externally referenced orbit claims.
Raise when orbit-following evidence is insufficient for an external-code claim.
Persist orbit-following claim evidence as deterministic JSON.
Pedestal¶
Pedestal profile model using the modified tanh (mtanh) parameterisation.
Provides H-mode pedestal structures for temperature and density profiles, including width scaling based on the EPED1 model.
dataclass
¶
Parameters for the mtanh pedestal profile.
Parameters¶
f_ped : float Value at the pedestal top (e.g. T_ped in keV). f_sep : float Value at the separatrix (edge). x_ped : float Normalized radial position of the pedestal symmetry point (typically 0.94). delta : float Pedestal width in normalized radius (typically 0.04). a : float Slope parameter for the pedestal top (typically 0.1-0.3).
Modified tanh (mtanh) pedestal profile generator.
Reference: Groebner, R.J. et al., Nucl. Fusion 41, 1789 (2001).
Compute pedestal width using EPED1-like scaling.
Scaling: Delta_psi = 0.076 * sqrt(beta_p_ped) Reference: Snyder, P.B. et al., Phys. Plasmas 16, 056118 (2009).
Parameters¶
beta_p_ped : float Poloidal beta at the pedestal top. delta_psi : float Optional offset or alternative scaling parameter.
Returns¶
float — Predicted pedestal width in poloidal flux (or rho approx).
Pellet Injection¶
Pellet ablation, trajectory, deposition, and fuelling-control utilities.
Plasma Startup¶
Plasma-startup sequence, avalanche, burn-through, and command-planning utilities.
Paschen curve breakdown model for D₂ gas.
Paschen (1889), Wied. Ann. 37, 69. Quantitative coefficients from Lieberman & Lichtenberg (2005), Ch. 14.3. Breakdown voltage: V_bd(pd) = B·pd / (ln(A·pd) - ln(ln(1/γ+1))) Minimum voltage: V_min = e·B/A · ln(1/γ+1) at (pd)_min = e/A · ln(1/γ+1) Reference: Lieberman 2005, Eq. 14.3.2.
Paschen breakdown voltage [V].
Lieberman 2005 Eq. 14.3.1: V_bd = B·(pd) / (ln(A·pd) - C₂) where C₂ = ln(ln(1/γ+1)).
Townsend avalanche and burn-through pre-ionisation model.
Energy balance follows Lieberman & Lichtenberg (2005), Ch. 14.4: d(n_e)/dt = n_e · k_iz · n_0 Ionisation rate coefficient (Maxwellian electrons): k_iz = σ_0 · √(8kT_e/(πm_e)) · exp(−E_iz/T_e) Reference: Janev et al. (1987), Ch. 2, rate coefficient parameterisation.
Rate coefficient k_iz [m³ s⁻¹] × n_neutral [m⁻³] → [s⁻¹].
k_iz = σ_0 · √(8kT_e/(πm_e)) · exp(−E_iz/T_e) Janev et al. (1987) Ch. 2; σ_0 = 1e-20 m², E_iz = 15.4 eV for D₂.
Burn-through energy balance model.
Criterion: P_ohmic > P_radiation_barrier. Energy balance: 3/2 · n_e · dT_e/dt = P_ohmic − P_rad Reference: Wesson (2011), "Tokamaks", 4th ed., §6.4, Eq. 6.4.1–6.4.3. Spitzer resistivity: η = 1.65×10⁻⁹ · Z_eff · ln_Λ / T_e[keV]^{3/2} Reference: Wesson 2011 Eq. 2.5.4.
Plasma Wall Interaction¶
Plasma-wall thermal loading and divertor lifetime-assessment utilities.
Physical sputtering yield Y [atoms ion^{-1}] using Eckstein analytical formula.
Y = Q · S_n(ε) / (1 + Γ · ε^{0.3}) · (1 − (E_th/E)^{2/3}) · (1 − E_th/E)²
Eckstein & Preuss 2003, J. Nucl. Mater. 320, 209. D→W parameters: M1=2, M2=183.8, Z1=1, Z2=74.
Gross and net erosion rates from sputtering.
Γ_erode = Γ_ion · Y(E, θ). Eckstein & Preuss 2003, J. Nucl. Mater. 320, 209.
Γ_erode = Γ_ion · Y(E, θ) [atoms m^{-2} s^{-1}].
Net rate accounting for redeposition [atoms m^{-2} s^{-1}].
Retained T fraction (relative to ion fluence) for W at room temperature.
~1% per 10^24 D m^{-2}. Roth et al. 2009, J. Nucl. Mater. 390-391, 1.
1D explicit heat diffusion into the wall.
Real Data Manifest¶
Validation contracts for real-shot and synthetic-shot data manifests.
The manifest is intentionally strict: synthetic fixtures are useful for CI, but they must never be counted as experimental validation evidence by accident.
Bases: ValueError
Raised when a data manifest cannot support its claimed validation role.
dataclass
¶
Single measured or generated signal entry in a manifest.
dataclass
¶
Acquisition source and provenance for a manifest.
dataclass
¶
Local artefact covered by a manifest checksum.
dataclass
¶
RealDataManifest(schema_version, dataset_id, machine, shot, synthetic, source, signals, retrieved_at=None, checksum_sha256=None, licence=None, synthetic_generator=None, synthetic_seed=None, artifacts=())
Validated manifest separating real-shot evidence from synthetic fixtures.
Load and validate a JSON real-data manifest.
Verify local artefact checksums for manifests that reference local files.
Remote acquisition sources such as MDSplus URIs are provenance records rather
than local files, so they return None here.
Validate a real-shot or synthetic-shot manifest.
Real manifests require stable provenance, licence, units, checksum, and a non-synthetic acquisition source. Synthetic manifests require generator metadata so CI fixtures cannot masquerade as experimental evidence.
Runaway Electrons¶
Runaway-electron generation, avalanche, evolution, and mitigation-assessment utilities.
staticmethod
¶
Electron density [10^20 m^-3] required to raise E_c above E_par.
Derived by inverting critical_field: n_e = E_par (4π ε₀² m_e c²) / (e³ ln Λ). Uses the temperature-dependent Coulomb log (Wesson 2011, Eq. 2.12.4).
staticmethod
¶
Maximum RE energy [MeV] from synchrotron balance.
Martin-Solis et al., Phys. Plasmas 13, 062509 (2006), Eq. 12. Defaults represent a disruption scenario: E_par=10 V/m, post-TQ ne_20=10, Te_keV=1.
staticmethod
¶
Deposited energy density [MJ/m^2] assuming instantaneous loss.
Assumes mean RE energy ≈ E_max / 2 (flat distribution upper bound).
Coulomb logarithm for a thermal plasma.
For T_e > 10 eV: ln Λ = 14.9 - 0.5 ln(n_e / 10^20) + ln(T_e / 10^3) where n_e is in m^-3 and T_e in eV.
Wesson 2011, Tokamaks 4th ed., Eq. 2.12.4.
Dreicer field E_D = n_e e^3 ln Λ / (4π ε₀² T_e).
Connor & Hastie, Nucl. Fusion 15, 415 (1975).
Critical electric field E_c = n_e e^3 ln Λ / (4π ε₀² m_e c²).
Rosenbluth & Putvinski, Nucl. Fusion 37, 1355 (1997). Note: E_D / E_c = m_e c² / T_e ≈ 51 at T_e = 10 keV.
Primary (Dreicer) seed runaway generation rate [m^-3 s^-1].
Connor & Hastie, Nucl. Fusion 15, 415 (1975). Smith et al., Phys. Plasmas 15, 072502 (2008) — hot-tail context uses same Dreicer base rate.
Avalanche multiplication rate [m^-3 s^-1].
Full pitch-angle-scattering formula
dn_RE/dt = n_RE (E/E_c - 1) / (τ_s ln Λ) × F(E/E_c, Z_eff)^{-1/2}
where F = 1 - E_c/E + 4π(Z+1)² / [3(Z+1)(E/E_c)²]
Rosenbluth & Putvinski, Nucl. Fusion 37, 1355 (1997), Eq. 15.
Seed RE density from thermal-quench hot-tail mechanism [m^-3].
Smith et al., Phys. Plasmas 15, 072502 (2008). Parametric fit to Smith (2008) Fig. 3: v_c/v_te scales as τ_q^0.2 at the V_C_V_TE_REF = 4.0 reference (τ_q = 1 ms).
Maximum RE energy [MeV] limited by synchrotron radiation.
E_max = m_e c² (E/E_c)^{1/3} (3 / (4 α_f))^{1/3}
Martin-Solis et al., Phys. Plasmas 13, 062509 (2006), Eq. 12. α_f = 7.297e-3 is the fine-structure constant (CODATA 2018).
Stellarator Geometry¶
Stellarator flux-surface geometry in Boozer coordinates, neoclassical transport, and ISS04 confinement scaling.
References¶
Boozer, A. H., Phys. Fluids 24 (1981) 1999. Magnetic coordinate system with straight field lines; ∇ψ × ∇θ_B = B/B². Grieger, G. et al., Phys. Fluids B 4 (1992) 2081. Wendelstein 7-X design parameters and modular coil geometry. Yamada, H. et al., Nucl. Fusion 45 (2005) 1684. ISS04 empirical scaling: τ_E ∝ a^2.28 R^0.64 P^-0.61 n_e19^0.54 B^0.84 ι^0.41. Nemov, V. V. et al., Phys. Plasmas 6 (1999) 4622. Effective ripple ε_eff via field-line tracing (DKES/NEO-2 basis, Eq. 30). Beidler, C. D. et al., Nucl. Fusion 51 (2011) 076001. Quasi-isodynamic neoclassical transport optimization for W7-X.
Bases: BaseModel
Pydantic v2 schema for stellarator device and magnetic parameters.
Ensure the helical flux surface cannot cross non-positive major radius.
dataclass
¶
StellaratorConfig(N_fp=5, R0=5.5, a=0.53, B0=2.5, iota_0=0.87, iota_a=1.0, mirror_ratio=0.07, helical_excursion=0.05, name='custom')
Stellarator device and magnetic configuration parameters.
Parameters¶
N_fp : int Positive number of toroidal field periods. W7-X: 5 — Grieger et al. 1992, Phys. Fluids B 4, 2081. R0 : float Positive major radius [m]. W7-X: 5.5 m. a : float Positive average minor radius [m]. W7-X: 0.53 m. B0 : float Positive on-axis magnetic field [T]. W7-X standard: 2.5 T. iota_0 : float Positive rotational transform at magnetic axis (s=0). Boozer 1981, Phys. Fluids 24, 1999 — straight field-line coordinate. iota_a : float Positive rotational transform at plasma edge (s=1). mirror_ratio : float Non-negative helical mirror ratio ε_h = δB / B₀. W7-X: ~0.07. helical_excursion : float Non-negative helical axis excursion amplitude [m]; R0 must exceed a plus this excursion.
Wendelstein 7-X standard configuration.
Grieger et al., Phys. Fluids B 4 (1992) 2081 — coil geometry and ι profile. Klinger et al., Nucl. Fusion 59 (2019) 112004 — operational parameters.
Rotational transform ι(s) via linear interpolation.
Boozer 1981, Phys. Fluids 24, 1999: ι = dψ_pol/dψ_tor in Boozer coordinates.
Parameters¶
config : StellaratorConfig s : float or array Normalised toroidal flux label, s ∈ [0, 1].
Raises¶
ValueError If the configuration is non-physical or s leaves [0, 1].
Compute a single flux surface in Boozer coordinates.
Boozer 1981, Phys. Fluids 24, 1999: (θ_B, ζ_B) are 2π-periodic straight field-line angles such that B · ∇θ_B / B · ∇ζ_B = ι.
|B| is represented as B₀(1 − ε_t cos θ − ε_h cos(N_fp ζ − ι θ)). This is the standard two-harmonic Boozer spectrum used for stellarator neoclassical calculations; Grieger et al. 1992, Phys. Fluids B 4, 2081.
Parameters¶
config : StellaratorConfig s : float Normalised toroidal flux, s ∈ (0, 1]. n_theta, n_phi : int Poloidal and toroidal grid resolution.
Returns¶
R : ndarray, shape (n_theta, n_phi) Major radius [m]. Z : ndarray, shape (n_theta, n_phi) Vertical position [m]. B : ndarray, shape (n_theta, n_phi) Magnetic field magnitude [T].
Raises¶
ValueError If the configuration is non-physical, s leaves (0, 1], or either grid resolution is below two points.
Effective helical ripple ε_eff for neoclassical transport.
Nemov et al., Phys. Plasmas 6 (1999) 4622, Eq. 30. Analytic proxy: ε_eff ~ ε_h^(3/2) / √N_fp. Full evaluation requires field-line tracing (DKES/NEO-2).
Parameters¶
config : StellaratorConfig s : float Normalised toroidal flux, s ∈ (0, 1].
Returns¶
float Effective ripple ε_eff ∈ (0, 1).
Raises¶
ValueError If the configuration is non-physical or s leaves (0, 1].
ISS04 energy confinement scaling law for stellarators.
Yamada et al., Nucl. Fusion 45 (2005) 1684, Eq. 4: τ_E = 0.134 · a^2.28 · R^0.64 · P^-0.61 · n_e19^0.54 · B^0.84 · ι_{2/3}^0.41
Parameters¶
config : StellaratorConfig n_e : float Line-averaged electron density [10^19 m^-3]. P_heat : float Absorbed heating power [MW].
Returns¶
float Predicted energy confinement time [s].
Raises¶
ValueError If the configuration, density, or heating power is outside the positive finite ISS04 domain.
Neoclassical thermal diffusivity in the 1/ν regime.
Beidler et al., Nucl. Fusion 51 (2011) 076001 — quasi-isodynamic transport. 1/ν regime: χ_neo ~ ε_eff^(3/2) · v_th² / (ν · R · N_fp)
Collision frequency from Braginskii
ν_ii ~ n_e e⁴ ln Λ / (4π ε₀² m_i² v_th³)
Parameters¶
config : StellaratorConfig s : float Normalised toroidal flux, s ∈ (0, 1]. T_keV : float Ion temperature [keV]. n_e19 : float Electron density [10^19 m^-3].
Returns¶
float Neoclassical χ [m²/s].
Raises¶
ValueError If the configuration is non-physical, s leaves (0, 1], or plasma inputs are not positive finite values.
Tearing Mode Coupling¶
Coupled tearing mode dynamics including Chirikov overlap and inter-mode coupling.
Key references¶
Rutherford 1973 : P.H. Rutherford, Phys. Fluids 16, 1903 (1973). [classical tearing: dw/dt ∝ Δ'] La Haye 2006 : R.J. La Haye, Phys. Plasmas 13, 055501 (2006). [MRE coefficient fits a1 = 6.35 for bootstrap drive] La Haye & Buttery : R.J. La Haye & R.J. Buttery, Phys. Plasmas 16, 022107 2009 (2009). [cross-mode coupling coefficient Eq. 8] Chirikov 1979 : B.V. Chirikov, Phys. Rep. 52, 263 (1979). [overlap criterion σ = (w1 + w2)/(2 Δr) > 1 → stochastic, Eq. 3.1]
Overlap criterion for magnetic island stochasticity.
σ = (w₁ + w₂) / (2 Δr)
σ > 1 → overlapping separatrices → field-line stochasticity. Chirikov 1979, Phys. Rep. 52, 263, Eq. 3.1.
Modified Rutherford Equation for two toroidally coupled tearing modes.
MRE for mode k (k = 1, 2): dw_k/dt = (r_sk / τ_Rk) [r_sk Δ'_k + a1 (j_bs/j_φ)(r_sk/w_k) × w_k²/(w_k² + w_d²) + C_kj w_j² / (r_s1 r_s2)]
Bootstrap drive coefficient a1 = 6.35: La Haye 2006, Phys. Plasmas 13, 055501, Table I. Cross-mode coupling coefficient C_kj: La Haye & Buttery 2009, Phys. Plasmas 16, 022107, Eq. 8. Classical tearing stability Δ': Rutherford 1973, Phys. Fluids 16, 1903.
Cross-mode coupling coefficient C₁₂.
La Haye & Buttery 2009, Phys. Plasmas 16, 022107, Eq. 8: C₁₂ ≈ 0.5 × (a / R₀) This captures the leading-order inverse-aspect-ratio toroidal coupling between the 3/2 and 2/1 modes via the n=1 sideband. Modes with different n do not couple directly at this order.
Effective Δ' for mode 1, modified by presence of mode 2.
Base linear Δ'₀ = −2 m₁ / r_s1 (classically stable). Rutherford 1973, Phys. Fluids 16, 1903, §III. Nonlinear cross-modification: La Haye & Buttery 2009, §II.
Integrate coupled MRE for both modes.
a1 (j_bs/j_φ) (r_s/w) w²/(w²+w_d²)
— La Haye 2006, Phys. Plasmas 13, 055501, Eq. 5; a1 = 6.35 Table I.
Coupling term: C₁₂ w₂²/(r_s1 r_s2) — La Haye & Buttery 2009, Phys. Plasmas 16, 022107, Eq. 8. Stochastic criterion: σ = (w1+w2)/(2Δr) > 1 — Chirikov 1979, Phys. Rep. 52, 263, Eq. 3.1.
Vessel Model¶
Vacuum vessel eddy current model using a lumped-circuit approach.
Models the induction and decay of currents in the conducting vessel walls, providing passive stability effects and flux perturbations.
dataclass
¶
VesselElement(R, Z, resistance, cross_section, inductance, wall_thickness=0.04, conductivity=1350000.0)
Discrete conducting element of the vacuum vessel.
Parameters¶
R : float Major radius of the element center [m]. Z : float Vertical position of the element center [m]. resistance : float Electrical resistance of the toroidal loop [Ohm]. cross_section : float Cross-sectional area of the element [m^2]. inductance : float Self-inductance of the toroidal loop [H]. wall_thickness : float Effective wall thickness for τ_vessel calculation [m]. conductivity : float Electrical conductivity of wall material [S/m]. Stainless steel 316L: σ ≈ 1.35e6 S/m (Wesson 2011, App. C).
Lumped-circuit model for vessel eddy currents.
Solves M dI/dt + R I = −V_ext, where M is the mutual inductance matrix and V_ext is the induced loop voltage.
Return peak halo current I_halo = f_halo × TPF × I_p.
ITER Physics Basis 1999, Ch. 3, §3.8.3: f_halo is the fraction of plasma current flowing in the halo (ITER default 0.3), TPF ≈ 2 is the toroidal peaking factor.
Parameters¶
plasma_current : float Total plasma current I_p [A]. f_halo : float Halo current fraction (default 0.3, ITER Physics Basis).
Returns¶
float — Peak halo current [A].
Return electromagnetic force on the vessel from halo currents.
F = I_halo × B_pol × L
Noll et al. 1993, Fusion Eng. Des. 22, 315 — electromagnetic loads on the first wall and blanket due to disruption halo currents.
Parameters¶
halo_current_a : float Peak halo current [A]. b_poloidal : float Poloidal field at the wall [T]. path_length : float Length of the current path in the wall [m].
Returns¶
float — Electromagnetic force [N].
VMEC Lite¶
Simplified spectral 3D MHD equilibrium solver (VMEC-lite, fixed-boundary).
References¶
Hirshman, S. P. & Whitson, J. C., Phys. Fluids 26 (1983) 3553. Variational moment method for 3D MHD equilibria; spectral representation R(s,θ,ζ) = Σ R_mn cos(mθ − nζ), Z = Σ Z_mn sin(mθ − nζ) [Eqs. 1–2]. Freidberg, J. P., "Ideal MHD" (Cambridge, 2014), Ch. 3. Force balance ∇p = J × B in 3D geometry. Wesson, J., "Tokamaks" 4th ed. (Oxford, 2011), Ch. 3. Toroidal field B_φ ∝ 1/R; poloidal field from rotational transform.
dataclass
¶
VMECLiteClaimEvidence(schema_version, source, source_id, geometry_source, profile_source, current_assumption, model_id, n_s, m_pol, n_tor, n_fp, n_modes, force_residual, iterations, converged, min_major_radius, max_abs_z, pressure_axis, pressure_edge, iota_min, iota_max, q_min, q_max, r_mn_relative_error, z_mn_relative_error, iota_relative_error, residual_tolerance, spectral_relative_tolerance, iota_relative_tolerance, full_vmec_claim_allowed, claim_status)
Serialisable provenance and reference-comparison evidence for VMEC-lite claims.
Fourier basis for VMEC spectral representation.
Hirshman & Whitson 1983, Phys. Fluids 26, 3553, Eqs. 1–2: R(s,θ,ζ) = Σ_{m,n} R_mn(s) cos(mθ − n N_fp ζ) Z(s,θ,ζ) = Σ_{m,n} Z_mn(s) sin(mθ − n N_fp ζ)
The basis requires non-negative poloidal and toroidal truncation orders and a positive integer number of field periods.
Evaluate spectral expansion at (θ, ζ) grid points.
Hirshman & Whitson 1983, Eq. 1–2: Σ C_mn cos(mθ − n N_fp ζ) or Σ C_mn sin(mθ − n N_fp ζ)
Spectral 3D equilibrium solver using steepest descent on the MHD energy.
Hirshman & Whitson 1983, Phys. Fluids 26, 3553 — variational moment method. Force balance ∇p = J × B (Freidberg 2014, Ch. 3) drives the Shafranov shift on the R₀₀ mode; radial tension (finite-difference Laplacian) regularises the spectral coefficients.
n_s must provide at least axis, interior, and boundary surfaces. Profile
inputs are finite one-dimensional arrays; pressure is non-negative and
rotational transform is positive for the reduced q-profile calculation.
Interpolate finite pressure and rotational-transform profiles onto the solver grid.
Steepest-descent equilibrium iteration.
Hirshman & Whitson 1983, Phys. Fluids 26, 3553.
Force balance ∇p = J × B (Freidberg 2014, Ch. 3) drives the pressure-gradient Shafranov shift on R₀₀ via q²(dp/ds)/R₀. Radial tension (d²/ds² finite difference) regularises spectral modes.
max_iter must be a positive integer and tol must be positive and
finite; invalid controls are rejected before iteration.
staticmethod
¶
Low-order Fourier approximation of shaped tokamak boundary.
R = R₀ + a cos(θ + δ sin θ), Z = a κ sin θ Hirshman & Whitson 1983, Phys. Fluids 26, 3553, Eqs. 1–2.
staticmethod
¶
W7-X standard configuration boundary.
Grieger et al., Phys. Fluids B 4 (1992) 2081 — modular coil geometry. N_fp = 5.
vmec_lite_claim_evidence(solver, result, *, source, source_id, geometry_source, profile_source, current_assumption, model_id='bounded_vmec_lite', reference_R_mn=None, reference_Z_mn=None, reference_iota=None, residual_tolerance=0.001, spectral_relative_tolerance=0.05, iota_relative_tolerance=0.02)
Build fail-closed evidence for bounded VMEC-lite or full VMEC claims.
Raise when VMEC-lite evidence is insufficient for a full VMEC claim.
Persist VMEC-lite claim evidence as deterministic JSON.
dataclass
¶
VMECLiteClaimEvidence(schema_version, source, source_id, geometry_source, profile_source, current_assumption, model_id, n_s, m_pol, n_tor, n_fp, n_modes, force_residual, iterations, converged, min_major_radius, max_abs_z, pressure_axis, pressure_edge, iota_min, iota_max, q_min, q_max, r_mn_relative_error, z_mn_relative_error, iota_relative_error, residual_tolerance, spectral_relative_tolerance, iota_relative_tolerance, full_vmec_claim_allowed, claim_status)
Serialisable provenance and reference-comparison evidence for VMEC-lite claims.
vmec_lite_claim_evidence(solver, result, *, source, source_id, geometry_source, profile_source, current_assumption, model_id='bounded_vmec_lite', reference_R_mn=None, reference_Z_mn=None, reference_iota=None, residual_tolerance=0.001, spectral_relative_tolerance=0.05, iota_relative_tolerance=0.02)
Build fail-closed evidence for bounded VMEC-lite or full VMEC claims.
Raise when VMEC-lite evidence is insufficient for a full VMEC claim.
Persist VMEC-lite claim evidence as deterministic JSON.
Phase Modules¶
GK UPDE Bridge¶
Bridge between gyrokinetic transport fluxes and the 8-layer UPDE Kuramoto phase dynamics system.
Maps GK-computed growth rates and diffusivities into adaptive K_nm coupling modulation for layers P0 (microturbulence), P1 (zonal flows), P4 (transport barrier), and P5 (current profile).
Reference layer mappings
P0 ← max(gamma_ITG, gamma_TEM): turbulence drive P1 ← chi_e suppression ratio: zonal flow damping of transport P4 ← chi_i pedestal / chi_i core: transport barrier strength P5 ← bootstrap current contribution (via pressure gradient)
Modulate K_nm based on GK fluxes.
Parameters¶
K_base : array, shape (L, L) Baseline coupling matrix from build_knm_plasma(). gk_output : GKOutput GK solver output (growth rates, fluxes). chi_i_profile : array or None Full chi_i(rho) profile for pedestal ratio calculation. gamma_ref : float Reference growth rate for tanh scaling [c_s/a]. chi_ref : float Reference chi_e for transport modulation [m^2/s].
Adjust layer-0 natural frequency based on GK growth rate.
The turbulence layer's effective frequency increases with the dominant instability growth rate.
SCPN Compiler and Replay Modules¶
FPGA Export¶
Export a CompiledNet to synthesisable HDL (Verilog / VHDL) for FPGA targets.
Generates a complete project directory with top-level module, weight ROM initialisation, timing constraints, and a Makefile for Vivado or Quartus.
dataclass
¶
FPGA synthesis target configuration.
dataclass
¶
HDLExportEvidence(schema_version, generated_utc, controller_artifact_sha256, target, target_part, clock_mhz, lif_bit_width, n_neurons, n_places, bitstream_length, hdl_format, hdl_sha256, weights_mem_sha256, constraints_sha256, makefile_sha256, project_sha256, resource_estimate_sha256, estimated_lut_count, estimated_ff_count, estimated_bram_blocks, weight_rom_depth, estimated_fmax_mhz, synthesis_toolchain, synthesis_report_sha256, synthesis_report_uri, timing_slack_ns, facility_claim_allowed, claim_status, payload_sha256)
Tamper-evident FPGA HDL export and synthesis-admission evidence.
Generate synthesisable Verilog HDL for the LIF neuron array.
Generate synthesisable VHDL for the LIF neuron array.
Estimate FPGA resource usage for the compiled net.
Returns dict with keys: lut_count, ff_count, bram_blocks, estimated_fmax_mhz.
Write a complete FPGA project directory.
Creates
top.v / top.vhd — HDL source weights.mem — weight ROM init constraints.xdc / constraints.sdc — timing Makefile — synthesis flow
hdl_export_evidence(compiled_net, config, project_dir, *, controller_artifact_sha256, target_part, generated_utc=None, synthesis_toolchain=None, synthesis_report_sha256=None, synthesis_report_uri=None, timing_slack_ns=None, facility_claim_allowed=False)
Build tamper-evident FPGA HDL export evidence for a generated project.
Fail closed unless HDL export evidence supports a facility hardware claim.
Persist HDL export evidence as sorted JSON.
Load HDL export evidence with duplicate-key and digest admission.
dataclass
¶
HDLExportEvidence(schema_version, generated_utc, controller_artifact_sha256, target, target_part, clock_mhz, lif_bit_width, n_neurons, n_places, bitstream_length, hdl_format, hdl_sha256, weights_mem_sha256, constraints_sha256, makefile_sha256, project_sha256, resource_estimate_sha256, estimated_lut_count, estimated_ff_count, estimated_bram_blocks, weight_rom_depth, estimated_fmax_mhz, synthesis_toolchain, synthesis_report_sha256, synthesis_report_uri, timing_slack_ns, facility_claim_allowed, claim_status, payload_sha256)
Tamper-evident FPGA HDL export and synthesis-admission evidence.
hdl_export_evidence(compiled_net, config, project_dir, *, controller_artifact_sha256, target_part, generated_utc=None, synthesis_toolchain=None, synthesis_report_sha256=None, synthesis_report_uri=None, timing_slack_ns=None, facility_claim_allowed=False)
Build tamper-evident FPGA HDL export evidence for a generated project.
Fail closed unless HDL export evidence supports a facility hardware claim.
Persist HDL export evidence as sorted JSON.
Load HDL export evidence with duplicate-key and digest admission.
Geometry Neutral Contracts¶
Geometry-neutral control contracts for non-tokamak replay adapters.
dataclass
¶
Magnetic-device metadata independent of tokamak axis observations.
dataclass
¶
ActuatorChannel(name, unit, min_value, max_value, slew_rate_per_s, latency_steps=0, failure_mode='none')
Bounded actuator channel with slew-rate and latency semantics.
dataclass
¶
One measured or replayed diagnostic value.
dataclass
¶
Replay diagnostic frame at one control tick.
dataclass
¶
Named targets, weights, and hard constraints for a replay scenario.
dataclass
¶
ReplayScenario(name, seed, steps, dt_s, magnetic_configuration, actuator_set, objective, initial_frame, fault_schedule)
Deterministic replay scenario contract for compact control adapters.
Geometry Neutral Replay¶
Compact geometry-neutral stellarator replay through the SCPN controller.
dataclass
¶
GeometryNeutralReplayEvidence(schema_version, generated_utc, replay_schema_version, replay_report_sha256, scenario_digest, trace_digest, metrics_digest, thresholds_digest, magnetic_configuration_reference, actuator_calibration, latency_model, fault_model, final_fieldline_spread, improvement_fraction, max_abs_current_A, p95_latency_us, deterministic, passes_thresholds, measured_or_benchmark_artefact_sha256, device_claim_allowed, claim_status, payload_sha256)
Tamper-evident geometry-neutral replay evidence admission object.
Load a bundled geometry-neutral replay JSON Schema document.
Return the bundled v1.1 replay schema after self-identification checks.
Validate a v1 replay report while the v1.1 schema bundle is installed.
build_aer_admission_metadata(*, admission_report, decode_strategy, decode_window_ns, n_features, feature_normalisation='unit', require_monotonic=False, feature_vector=None)
Build replay-safe AER admission metadata from a decoded observation.
Attach digest-bound AER admission metadata to a replay v1.1 report.
Generate a deterministic compact SCPN-control replay report.
Public package alias for :func:generate_report.
Validate the compact replay report contract without external packages.
Public package alias for :func:validate_report.
Persist a validated geometry-neutral replay report as stable JSON.
Load and validate a geometry-neutral replay report with duplicate-key checks.
geometry_neutral_replay_evidence(report, *, generated_utc=None, measured_or_benchmark_artefact_sha256=None, device_claim_allowed=False)
Build tamper-evident replay evidence over an admitted replay report.
Fail closed unless replay evidence supports a device-control claim.
Persist geometry-neutral replay evidence as sorted JSON.
Load geometry-neutral replay evidence with duplicate-key and digest admission.
Public package alias for :func:render_markdown.
dataclass
¶
GeometryNeutralReplayEvidence(schema_version, generated_utc, replay_schema_version, replay_report_sha256, scenario_digest, trace_digest, metrics_digest, thresholds_digest, magnetic_configuration_reference, actuator_calibration, latency_model, fault_model, final_fieldline_spread, improvement_fraction, max_abs_current_A, p95_latency_us, deterministic, passes_thresholds, measured_or_benchmark_artefact_sha256, device_claim_allowed, claim_status, payload_sha256)
Tamper-evident geometry-neutral replay evidence admission object.
geometry_neutral_replay_evidence(report, *, generated_utc=None, measured_or_benchmark_artefact_sha256=None, device_claim_allowed=False)
Build tamper-evident replay evidence over an admitted replay report.
Fail closed unless replay evidence supports a device-control claim.
Persist geometry-neutral replay evidence as sorted JSON.
Load geometry-neutral replay evidence with duplicate-key and digest admission.
CLI¶
scpn-control demo --scenario combined --steps 1000
scpn-control benchmark --n-bench 5000 --json-out
scpn-control validate --json-out
scpn-control validate-release-evidence artifacts/release_evidence_report.json --json-out
scpn-control info --json-out
scpn-control live --port 8765 --zeta 0.5 --layers 16
scpn-control hil-test --shots-dir path/to/shots
| Command | Description |
|---|---|
demo |
Closed-loop control demonstration (PID, SNN, combined) |
benchmark |
PID vs SNN timing benchmark with JSON output option |
validate |
Import hygiene plus data-manifest, JAX GK parity, and physics-traceability gates |
validate-release-evidence |
Admission check for JSON reports emitted by scpn-control validate --json-out, including data manifests, JAX GK parity, physics traceability, multi-shot campaign evidence, and native formal certificate evidence |
info |
Version, Rust backend status, weight provenance, Python/NumPy versions |
live |
Real-time WebSocket phase sync server |
hil-test |
Hardware-in-the-loop test campaign against shot data |
Rust Acceleration¶
When scpn-control-rs is built via maturin, all core solvers use Rust backends automatically:
from scpn_control import RUST_BACKEND
print(RUST_BACKEND) # True if Rust available
# Transparent acceleration — same Python API, Rust execution
kernel = FusionKernel(R0=6.2, a=2.0, B0=5.3)
Build Rust bindings:
PyO3 Bindings¶
| Python Class | Rust Binding | Crate |
|---|---|---|
FusionKernel |
PyFusionKernel |
control-core |
RealtimeMonitor |
PyRealtimeMonitor |
control-math |
SnnPool |
PySnnPool |
control-control |
MpcController |
PyMpcController |
control-control |
Plasma2D |
PyPlasma2D |
control-core |
TransportSolver |
PyTransportSolver |
control-core |
Core — Native Rust Engine Wrapper¶
scpn_control.core.rust_engine is the Python control-plane wrapper for the
optional PyO3 native execution bridge. It configures native campaign execution,
formal-verification mode, runtime admission, transport backend selection, and
emergency telemetry handoff while keeping timing-critical execution inside the
compiled Rust data plane when the extension is available.
This wrapper is an execution boundary, not a physics solver. It does not turn a local workstation run into target-hardware PCS evidence unless the matching runtime-admission and benchmark-context reports pass.
Rust-accelerated control-plane wrapper.
This module exposes a single orchestrator class for launch-time control of the compiled Rust transport primitives.
Responsibility is split by design: - Python owns policy, campaign configuration, and lifecycle boundaries. - Rust owns the per-tick heavy control/transport work and back-pressured UDP path.
Coordinate Rust-native control execution from Python.
Python owns campaign configuration, safety policy, transport policy and telemetry. Execution primitives (SNN, PID, transport publisher) are executed through compiled Rust bindings.
property
¶
Whether transport-native Rust bindings are importable.
Configure control targets and execution-safe scalar bounds.
Supported keys: - target_r / rho_tor_target / R_target - target_z / z_tor_target / Z_target - beta_n_limit - c_gB / c_gB_nominal - optional u_min/u_max safety clamps
Load Kuramoto weights that remain under Python policy control.
Accepts either:
- a mapping of weight names to float coefficients,
- a filesystem path to a JSON object containing weight names.
Record preferred execution affinity for campaign telemetry.
The compiled transport bridge does not currently expose explicit core pinning. These fields are retained for reproducibility, campaign replay, and future native controller migration.
configure_native_formal_verification(*, enabled=True, mode='async_drop', max_marking=100, max_depth=4, dispatch_interval_steps=30, channel_capacity=2)
Configure native formal checking for the fused hardware loop.
Load and apply ITPA gyro-Bohm profile constraints.
Both legacy ({"c_gB": ...}) and upgraded schema
(scaling_parameters.{c_gB_nominal,...}) are supported.
execute_campaign(*, steps=None, runtime_s=None, tick_interval_s=0.001, max_publish_failures=None, initial_state=None, plant_gain=None, core_snn=1, core_z3=2, core_net=3, core_hb=4, execution_backend='auto', pacing_mode='sleep', runtime_admission_policy='warn')
Campaign-facing wrapper used by command-line and external orchestrators.
execute_hardware_loop(*, steps=None, runtime_s=None, tick_interval_s=0.001, initial_state=None, plant_gain=None, max_publish_failures=None, core_snn=1, core_z3=2, core_net=3, core_hb=4, execution_backend='auto', pacing_mode='sleep', runtime_admission_policy='warn')
Run the Rust-accelerated control loop.
Configuration is owned by Python; transport execution is delegated to
compiled Rust primitives (RustSnnController, RustPIDController,
RustUdpTransportBridge).
admit_runtime(*, execution_backend, pacing_mode, tick_interval_s, native_backend_available, policy='warn')
Collect runtime admission evidence and optionally fail closed.
Compatibility alias used by external orchestration scripts.
Best-effort emergency shutdown hook exposed to callers.