Company Website: https://intrafere.com/
Software GitHub that produced this paper: https://github.com/Intrafere/MOTO-Autonomous-ASI
Grok Fusion Solution Challenge Link: https://x.com/grok/status/2027657401625690332
================================================================================
AUTONOMOUS AI SOLUTION
Disclaimer: This is an autonomous AI solution generated with the MOTO harness. This paper was not peer reviewed and was autonomously generated without user oversight or interaction beyond the original user prompt, therefore, this text may contain errors. These papers often contain ambitious content and/or extraordinary claims, all content should be viewed with extreme scrutiny.
(EDITOR NOTE: This single paper does not attempt to solve the user’s prompt entirely, it is meant to be one piece toward the complex solution required for the users prompt – total solutions typically are achieved in later papers) User’s Research Prompt: Deliver a complete, engineering-ready blueprint for a compact stellarator fusion reactor achieving sustained Q>15 net gain by 2030—using only near-term materials, full MHD/plasma stability models, tritium breeding cycle, and <$5B build cost. Include all equations, sim code, and falsifiable tests. AI Model Authors: x-ai/grok-4.1-fast, openai/gpt-5.2, moonshotai/kimi-k2.5 Generated: 2026-02-28 ================================================================================ Paper Title: Contract-Based Digital-Twin Integration for a Compact Q>15 Stellarator: UQ, RAMI, and Cost-Constrained Optimization
Abstract
A compact stellarator digital twin intended to support decisions under stringent performance (e.g., Q ≥ 15), engineering, RAMI, cost, and schedule requirements is only credible if it prevents integration artifacts and issues auditable, solver-independent claims under uncertainty. This paper develops a certificate-oriented mathematical framework for assembling such a system-of-systems twin without asserting the correctness of any particular device scaling law or design. Coupling is formalized by a directed port graph with conserved-quantity incidence structure, unit/dimension tagging, normalized residuals, and explicit tolerance budgets; when contracts are infeasible, dual/Farkas-type witnesses provide checkable conflict explanations. Execution well-posedness is addressed through dissipativity/passivity declarations, scattering-style interconnections for delay-robust co-simulation, and small-gain/contraction conditions that yield convergence and explicit error-amplification bounds, including asynchronous updates with bounded delays.
To eliminate solver-dependent artifacts, a verified-numerics wrapper based on interval arithmetic and Krawczyk-type inclusion tests produces PASS/FAIL/UNKNOWN enclosures for coupled steady states and rigorous interval bounds for constraints, with conservative propagation of UNKNOWN outcomes into downstream claims. Scalable uncertainty propagation is organized around surrogate layers (notably polynomial chaos with explicit surrogate-error envelopes), compressed-sensing coefficient recovery under stated assumptions, and multilevel/multi-fidelity Monte Carlo templates with explicit bias accounting. Feasibility probability is certified using complementary, assumption-explicit tools (moment bounds, scenario methods, Wasserstein DRO, sum-of-squares robust feasibility on polynomial surrogates, and conformal prediction). RAMI and performability are integrated via dependence-aware cut-set bounds, common-cause shock models, and multistate capacity/UGF methods. The framework culminates in certificate-linked cost-constrained continuous–discrete optimization (GP/MICP/GBD) and decision-directed V&V and Bayesian updating with change-control, yielding versioned go/no-go dossiers whose validity is traceable to logged assumptions and artifacts.
I. Introduction
A compact stellarator aimed at high fusion gain (here framed as the requirement \(Q\ge 15\)) must be designed, integrated, and assessed as a tightly coupled system-of-systems: core and edge plasma models interact with divertor heat exhaust, blanket/tritium systems, coils and cryogenics, plant power conversion, reliability/availability/maintainability/inspectability (RAMI), and cost and schedule constraints. A “digital twin” that couples these heterogeneous elements is potentially decision-relevant only if its outputs can be traced to explicit assumptions and if integration artifacts (unit mistakes, sign errors, double counting, solver tolerances, hidden algebraic loops, and untracked surrogate error) are prevented from masquerading as physics or engineering truth.
This paper develops a certificate-oriented mathematical framework for assembling and operating such an integrated twin under uncertainty. The goal is not to claim that any particular stellarator design is feasible or that any specific device-level scaling law is correct. Rather, the goal is to specify an auditable protocol: a collection of mathematical structures and checkable artifacts that (i) enforce physically meaningful interface semantics (conservation, units, sign conventions), (ii) provide stability and well-posedness conditions for modular co-simulation and steady-state coupling, (iii) attach solver-independent verified-numerics enclosures to coupled solves, (iv) propagate uncertainty with explicit accounting for sampling error, surrogate/model-form error, and integration/coupling error, and (v) connect these certificates to RAMI/performability and to cost- and schedule-constrained continuous–discrete optimization.
The central design object throughout is a joint feasibility event \(\mathcal{F}(d,y)\) for continuous design/operation variables \(d\) and discrete architecture variables \(y\): meeting \(Q\ge 15\) together with engineering constraints, RAMI/performability requirements, and cost/schedule bounds. The paper emphasizes that any statement of the form \(\mathbb{P}(\mathcal{F}(d,y))\ge 1-\alpha\) is only meaningful when the semantics of “\(\mathbb{P}\)” and the handling of numerical and surrogate uncertainties are explicit and replayable.
The main contributions are organized as a stack of certificate layers.
1. Contract-based integration spine (Sections II–III). We formalize coupling as a directed port graph with conserved-quantity families and incidence matrices, enforcing unit consistency via dimension vectors and producing dimensionless, normalized residuals with explicit tolerance budgets. When assembled contracts are infeasible, the framework requires auditable explanations via dual/Farkas-type infeasibility witnesses (or conic dual rays), rather than opaque solver failures.
2. Stability and well-posedness of coupled execution (Section IV). We provide assumption-explicit stability conditions for modular execution using dissipativity/passivity declarations, power-preserving interconnection logic, scattering-style couplers for delay-robust passivity, and a small-gain/contraction viewpoint for steady-state fixed-point coupling. These results yield PASS/FAIL/UNKNOWN readiness triage and a quantitative bound on how per-module approximation and coupling errors can amplify into plant-level shared-state error.
3. Verified numerics for coupled steady states and constraint enclosures (Section V). To avoid solver-dependent artifacts, we attach interval-based inclusion tests (Krawczyk-type operators) to certify existence/uniqueness of a coupled steady state in a box \(X\), and then certify constraint satisfaction by interval-enclosing each constraint over \(X\). The output is a replayable dossier with a tri-valued PASS/FAIL/UNKNOWN classification that propagates conservatively into downstream uncertainty quantification.
4. Surrogates and scalable uncertainty propagation with explicit error envelopes (Section VI). We present polynomial chaos expansions (PCE) as an algebraic surrogate layer; compressed-sensing recovery for high-dimensional coefficient estimation, stated with assumption-conditional bounds; multilevel and multi-fidelity Monte Carlo templates with explicit bias/variance/cost assumptions; and active-subspace diagnostics with auditable eigen-spectrum/tail-energy reporting. A recurring principle is that surrogate error must appear as an explicit envelope \(\Delta\) (deterministic or probabilistic) rather than being hidden in validation plots.
5. Probability-of-feasibility and safety certificates without unverifiable Monte Carlo claims (Section VII). We compile complementary certificate mechanisms: distribution-free moment bounds (Cantelli/union allocations), scenario methods with finite-sample guarantees (conditional on convexity and i.i.d. sampling), Wasserstein DRO templates (conditional on Lipschitz bounds and metric/radius choices), SOS robust-feasibility certificates for polynomial/semialgebraic representations, and conformal prediction tubes for black-box surrogate stacks under exchangeability. The treatment of verified-numerics UNKNOWN outcomes is made explicit, since it directly affects conservatism and validity of probability claims.
6. RAMI and performability under dependence and partial degradation (Section VIII). We incorporate cut-set bounds that remain valid without independence assumptions, introduce common-cause-aware modeling via explicit shock variables (with conservative bounding when parameters are not validated), and model partial degradation via multistate capacity and universal generating functions (UGF) with quantization-induced conservatism tracked in the certificate.
7. Optimization, V\&V prioritization, and updating under contracts (Sections IX–XIII). We standardize how component-level physics/risk modules must declare validity domains and error envelopes to be admissible in the coupled twin. We then frame cost- and schedule-constrained co-design as producing not only candidate designs but also optimization evidence (e.g., MICP gaps, Benders cut provenance) that is explicitly linked to the feasibility-certificate layers. Finally, we describe decision-directed V\&V planning (BOED, sensitivity diagnostics, conformal tube tightening) and contract-conditioned Bayesian calibration with change-control rules, so that certificates remain meaningful under versioned updates.
A guiding methodological stance is explicit skepticism: many quantities required by the above certificates (Lipschitz/gain bounds, passivity indices, surrogate envelopes, dependence parameters, monotonicity regimes, tail models) are not automatic consequences of running a code. Accordingly, the paper attaches “gap flags” whenever a certificate pathway depends on assumptions that must be empirically justified for a specific module and operating regime.
Roadmap. Section II fixes notation and formalizes the plant-level constraints and uncertainty objects. Section III builds the contract-based coupling spine and infeasibility explanations. Section IV develops stability and well-posedness conditions for modular execution (including asynchronous delays). Section V introduces verified-numerics enclosures for coupled steady states and constraints. Section VI develops surrogate and scalable UQ machinery with explicit error tracking. Sections VII and VIII turn those ingredients into probability-of-feasibility and RAMI/performability certificates. Sections IX–XIII describe how component modules, optimization, V\&V planning, and Bayesian updating fit into a single auditable workflow.
II. Preliminaries, Notation, and Problem Formalization
This section defines the plant-level quantities of interest, constraints, design and architecture variables, and a contract vocabulary for coupling heterogeneous modules in a digital-twin “system-of-systems”. The emphasis is on statements that are mathematically well-posed and auditable. Any specific physics limits or numerical values mentioned below are treated as placeholders until externally validated for the target device.
II.A. Plant-level quantities of interest and constraints
We consider a plant digital twin that, for a given design/architecture choice and uncertain inputs, returns a vector of quantities of interest (QoIs) and constraint functions.
1. Design, uncertainty, and outputs.
Let
\[
d \in \mathcal{D} \subset \mathbb{R}^{n_d}
\]
be continuous design/operating variables (geometry scalars, field levels, set-points, capacity margins), and let
\[
y \in \mathcal{Y} \subset \mathbb{Z}^{n_y}
\]
be discrete architecture variables (redundancy counts, number of modules, spare counts, maintenance resources).
Let
\[
\xi \in \Xi
\]
denote uncertain inputs (aleatoric and epistemic) defined precisely in II.B. The coupled digital twin returns a vector
\[
z(d,y,\xi) \in \mathbb{R}^{n_z}
\]
that includes (at minimum) fusion gain, engineering margins, RAMI measures, and cost/schedule outputs.
2. Fusion-gain requirement.
Let \(Q(d,y,\xi)\) denote fusion gain. The core performance requirement is
\[
Q(d,y,\xi) \ge 15
\]
interpreted as a constraint that must hold with a specified probability level (Section VII) or robustly with respect to an uncertainty set, depending on the certification layer used.
3. Engineering constraints (illustrative).
Engineering and safety constraints are encoded as functions
\[
g_j(d,y,\xi) \le 0, \qquad j=1,\dots,J_{eng}.
\]
Examples (purely illustrative; the framework does not assume their correctness for any device):
– Divertor/SOL heat-load limit: \(g_{q}(d,y,\xi) := q_{\perp}^{peak}(d,y,\xi) – q_{\perp}^{lim}\le 0\), where \(q_{\perp}^{peak}\) has units of \([\mathrm{W}/\mathrm{m}^2]\) and \(q_{\perp}^{lim}\) is an allowable limit.
– Coil current/field margin constraint: \(g_{I}(d,y,\xi) := I_{op}(d,y,\xi) – I_{allow}(d,y,\xi)\le 0\) (units \([\mathrm{A}]\)), or an equivalent margin form.
– Tritium breeding constraint: \(g_{TBR}(d,y,\xi) := TBR_{min} – TBR(d,y,\xi)\le 0\), where \(TBR\) is dimensionless.
– Inventory limits: \(g_{inv}(d,y,\xi) := I_T(d,y,\xi) – I_T^{lim} \le 0\) with tritium inventory \(I_T\) in \([\mathrm{g}]\) or \([\mathrm{mol}]\), depending on declared units.
All constraints are written in the convention “\(\le 0\) is safe”.
4. RAMI constraints.
Let \(U(d,y,\xi)\in[0,1]\) denote unavailability (steady-state or over a specified mission time) and \(A(d,y,\xi):=1-U(d,y,\xi)\) availability. RAMI requirements can be encoded as
\[
g_{RAMI}(d,y,\xi) := U(d,y,\xi) – U_{max} \le 0,
\]
or as multi-criteria constraints on planned maintenance unavailability, repair resources, and spares policies (Section VIII). Since dependence assumptions strongly affect RAMI estimates, any certificate must explicitly state whether it assumes independence, bounded dependence, or a common-cause model.
5. Cost and schedule constraints.
Let \(C(d,y,\xi)\) denote cost (e.g., in dollars) and \(S(d,y,\xi)\) a schedule metric (e.g., completion year or duration). Hard constraints are
\[
C(d,y,\xi) \le C_{max}, \qquad S(d,y,\xi) \le S_{max}.
\]
The outline motivates \(C_{max}=5\,\mathrm{B\$}\) and \(S_{max}=2030\), but the mathematics below only requires constants \(C_{max},S_{max}\) with declared units and reference frames.
II.B. Decision variables and uncertainty models
We require a formal separation of (i) design variables under control, (ii) stochastic/uncertain inputs, and (iii) numerical/integration error sources.
1. Continuous decision variables.
The continuous decision vector \(d\) may include:
– geometry scalars (major radius proxy, aspect ratio proxy, blanket thickness proxies),
– field/operation set-points (on-axis field proxy, rotational transform control parameters, heating power set-points),
– capacity margins (installed auxiliary heating capacity, cryogenic capacity margin, pumping capacity margin),
– control set-points (e.g., density/temperature targets) when modeled algebraically.
We keep \(\mathcal{D}\) abstract; later sections impose convex/GP-compatible parameterizations when needed.
2. Discrete architecture variables.
The integer vector \(y\) may include:
– redundancy counts (e.g., number of parallel compressor trains, number of heating sources),
– spare inventories and logistics choices,
– modular segmentation counts (e.g., number of blanket sectors),
– maintenance-resource levels (remote-handling tool sets, crew shifts),
again leaving the set \(\mathcal{Y}\) abstract until the optimization formulation in Section X.
3. Uncertain inputs: aleatoric vs. epistemic.
We model uncertainty via \(\xi\), optionally decomposed as
\[
\xi = (\xi_a, \xi_e)
\]
where \(\xi_a\) represents aleatoric variability (intrinsic fluctuations) and \(\xi_e\) epistemic uncertainty (limited knowledge, model-form uncertainty). The certification layer must declare whether it uses:
– a probability model \(\mathbb{P}\) on \(\Xi\) (scenario methods, conformal prediction, Monte Carlo/MLMC),
– an ambiguity set \(\mathcal{P}\) of distributions (distributionally robust optimization), or
– a deterministic uncertainty set \(\Xi\) (robust optimization).
Examples of uncertain inputs (illustrative only): transport multipliers, confinement efficiency factors, blanket breeding multipliers, coil fabrication tolerance descriptors, cross-section multipliers, failure and repair-time parameters.
II.C. Interface and contract vocabulary
The document uses “contracts” to make module coupling explicit, checkable, and auditable.
1. Modules, ports, and shared state.
A module is a map that consumes inputs and produces outputs, typically written as
\[
M_i: (u_i,\xi_i) \mapsto (y_i, s_i),
\]
where \(u_i\) are inputs received from other modules (through ports), \(y_i\) are exported outputs, and \(s_i\) is an internal state (possibly empty) that is not exported except through declared ports.
A port is a declared interface variable (or vector) with:
– a physical unit (or dimension vector),
– a sign convention (especially for conserved flows),
– an admissible range and/or tolerance,
– and a semantic label (e.g., “net electric power”, “tritium flow”, “heat to divertor”).
The shared state \(x\) of the coupled twin is the collection of all port variables that participate in inter-module coupling (Section III), together with any globally defined algebraic states used by multiple modules.
2. Contract residuals, tolerances, and margins.
A contract is a set of declared conditions that must hold at the interface of a module and/or for an interconnection, typically expressed as inequalities and equalities
\[
h_{i,k}(x,d,y,\xi) \le \tau_{i,k}, \qquad r_{i,\ell}(x,d,y,\xi)=0.
\]
Here \(h_{i,k}\) is a residual function (e.g., conservation mismatch, unit consistency check, or envelope violation) and \(\tau_{i,k}\ge 0\) is an auditable tolerance budget. A margin is a signed slack variable that converts a constraint into a “distance to violation”, e.g.
\[
m_j(d,y,\xi) := -g_j(d,y,\xi) \quad (\text{so } m_j\ge 0 \text{ is safe}).
\]
When certificates are computed, it is often more meaningful to certify lower bounds on margins than to certify raw values.
3. Certificates and falsifiers.
A certificate is a verifiable artifact (e.g., a dual witness, an interval enclosure, a probabilistic lower bound with specified assumptions) whose validity can be checked independently of the particular solver run that produced it.
A falsifier is an input scenario, test case, or counterexample that violates a contract or invalidates an assumption (e.g., a parameter regime where Lipschitz bounds fail, or a unit/sign mismatch).
4. Distinguishing error sources.
We explicitly distinguish:
(i) Model-form uncertainty: discrepancy between physics and the model class (e.g., missing physics).
(ii) Numerical/discretization error: error due to finite resolution or solver tolerance.
(iii) Integration/coupling error: error due to module interconnection (lagging, inconsistent shared variables, algebraic loops, unit/sign mismatches).
These must not be conflated in UQ claims; later sections treat them as separate contributions to conservative constraint tightening.
II.D. Standardization and gap flag
The internal research context that motivated this outline contains many proposed theorems and claims with inconsistent numbering and unverified assumptions. In this paper, theorem statements are re-derived from standard results (e.g., contraction mapping, small-gain, scenario bounds, interval inclusion tests) with all assumptions stated explicitly.
Accordingly, throughout Sections III–XII we will attach “gap flags” to any step requiring empirical validation for a specific integrated module, including but not limited to:
– numerical values of Lipschitz/passivity bounds,
– surrogate approximation error envelopes,
– dependence/independence structure in RAMI,
– tail-behavior assumptions for probabilistic lifetime/damage models,
– validity of any low-degree/sparsity assumptions used for polynomial surrogates.
This completes the formal problem setup used by the contract-based integration, stability analysis, verified numerics, UQ/certification, RAMI, and optimization layers that follow.
III. Contract-Based Digital-Twin Integration Spine (Conservation + Units + Infeasibility Explanations)
This section defines a contract-based coupling spine for a heterogeneous “system-of-systems” digital twin. The purpose is (i) to prevent nonphysical exchanges caused by sign/unit mismatches and double counting, and (ii) to produce auditable explanations when a coupled set of contracts is infeasible. The mathematics below is independent of any particular stellarator physics model; it only requires that each module declares ports (typed by conserved-quantity family and units) and, when needed, tolerance budgets.
III.A. Port-graph representation of the system-of-systems
III.A.1. Directed port graph and conserved-quantity families
Let \(G=(V,E)\) be a directed graph representing module couplings. Each node \(i\in V\) is a module (core plasma, SOL/divertor, blanket, tritium plant, auxiliary heating, coils/cryo, power conversion, RAMI, cost/schedule, etc.). Each directed edge \(e=(i\to j)\in E\) represents an exchanged flow variable from module \(i\) to module \(j\).
We distinguish conserved-quantity families \(q\in\mathcal{Q}\) (e.g., energy/power, particle species, tritium mass). In practice the same physical connection may carry multiple conserved families; this can be represented by either (a) separate edges per family, or (b) one edge with a vector-valued flow whose components are separately unit-tagged. For mathematical clarity, we use per-family edge sets \(E_q\) and write \(G_q=(V,E_q)\).
For each family \(q\), define the signed node-edge incidence matrix \(B_q\in\mathbb{R}^{|V|\times |E_q|}\) by
\[
(B_q)_{i,e} :=
\begin{cases}
+1, & \text{if edge } e \text{ enters node } i,\\
-1, & \text{if edge } e \text{ leaves node } i,\\
0, & \text{otherwise.}
\end{cases}
\]
A flow vector \(f_q\in\mathbb{R}^{|E_q|}\) assigns a signed flow to each edge, with the convention that \(f_{q,e}>0\) means the flow travels in the declared direction of \(e\).
III.A.2. Node inventories, sources/sinks, and algebraic vs dynamic closure
At node \(i\) for family \(q\), we allow (a) an inventory/state \(n_{q,i}(t)\) (e.g., tritium holdup, stored thermal energy), and (b) an internal source/sink term \(s_{q,i}(t)\) (e.g., fusion power injected to the plasma node; external grid export; tritium decay/processing), both declared by the module.
A generic conservation law on \(G_q\) can be written in dynamic form as
\[
\dot n_q(t) + B_q f_q(t) – s_q(t) = 0,
\]
where \(n_q(t),s_q(t)\in\mathbb{R}^{|V|}\) stack the node inventories and sources. In steady-state couplings (typical for plant-level design loops), inventories are constant (or absent), so \(\dot n_q=0\) and conservation becomes an algebraic closure
\[
r_q(f_q,s_q) := B_q f_q – s_q = 0.
\]
The coupling spine treats \(r_q=0\) as a global interconnection constraint rather than a “per-module check”. The incidence structure is the mechanism that prevents (by construction) double counting of a shared flow: each edge appears exactly once with opposite signs at its two endpoint nodes.
III.B. Dimensional analysis embedded in contracts
III.B.1. Dimension vectors and unit consistency checks
To make unit checking machine-auditable, we represent physical dimensions by exponent vectors in \(\mathbb{Z}^7\) using the SI base dimensions (Mass, Length, Time, Electric current, Temperature, Amount, Luminous intensity). Let \(\dim(\cdot): \text{(quantity)} \to \mathbb{Z}^7\) assign a dimension vector. Examples:
– Power \([\mathrm{W}]=[\mathrm{kg\,m^2\,s^{-3}}]\) corresponds to \((1,2,-3,0,0,0,0)\).
– Energy \([\mathrm{J}]\) corresponds to \((1,2,-2,0,0,0,0)\).
– Mass flow \([\mathrm{kg/s}]\) corresponds to \((1,0,-1,0,0,0,0)\).
A port declaration for an edge flow \(f_{q,e}\) must include its dimension vector \(\dim(f_{q,e})\). A scalar algebraic expression \(\varphi\) is unit-consistent if all terms in sums have identical dimension vectors and if products/quotients satisfy dimension-vector addition/subtraction. Concretely, if a contract contains a residual equation of the form
\[
\sum_{k=1}^K a_k = 0,
\]
then the contract is rejected as ill-posed unless \(\dim(a_1)=\cdots=\dim(a_K)\). (In implementation, this is a compile-time or model-build-time check.)
III.B.2. Dimensionless normalized residuals for auditability
For auditing, it is often preferable to store and reason about dimensionless residuals. For each equality-type conservation residual component \((r_q)_i\) with declared physical units, choose a positive scaling \(\sigma_{q,i}>0\) with matching units and define
\[
\hat r_{q,i} := \frac{(r_q)_i}{\sigma_{q,i}},
\]
so \(\hat r_{q,i}\) is dimensionless.
Similarly, for an inequality constraint \(g_j(d,y,\xi)\le 0\), select a positive scale \(\sigma_j\) in the same units as \(g_j\) and define \(\hat g_j := g_j/\sigma_j\). The choice of \(\sigma\) is part of the audit record; typical scalings include allowable limits, nominal operating values, or explicitly chosen “engineering significance” scales.
III.C. Global conservation residual requirements
III.C.1. Residual bounds and tolerance budgeting
Even when the physical law is an equality, numerical time stepping, surrogate error, and asynchronous co-simulation introduce mismatch. We therefore impose contract bounds on normalized residuals:
\[
\|\hat r_q(t)\|_\infty \le \tau_q
\]
for each conserved family \(q\), where \(\tau_q\ge 0\) is a declared tolerance budget. The budget must be traceable to identified contributions, e.g.
\[
\tau_q = \tau^{\mathrm{disc}}_q + \tau^{\mathrm{surr}}_q + \tau^{\mathrm{cpl}}_q,
\]
corresponding to discretization error, surrogate/model-form envelopes at ports, and coupling/execution error (e.g., lagging). The contract spine does not assume any particular formula for these terms; it requires only that the decomposition be explicitly recorded so that tightening or loosening \(\tau_q\) has an auditable justification.
III.C.2. Prevention of sign and double-counting errors
The following property is immediate from the incidence formulation.
Proposition III.1 (Edge-wise antisymmetry of accounting).
For any \(q\) and any edge flow vector \(f_q\), the total graph sum of conservation residuals satisfies
\[
\mathbf{1}^\top B_q f_q = 0.
\]
In particular, if \(s_q\equiv 0\), then \(\sum_i (r_q)_i = 0\) identically.
Proof.
Each column of \(B_q\) has exactly one \(+1\) and one \(-1\); hence its column sum is zero, so \(\mathbf{1}^\top B_q=0\) and the claim follows. \(\square\)
Operationally, Proposition III.1 means that if the physical source/sink vector \(s_q\) is separately unit-checked, then any mismatch in global closure cannot be “hidden” by inconsistent sign conventions across modules: the spine forces a single sign convention per edge and guarantees that every inter-module exchange is counted exactly once.
III.D. Infeasibility certificates for coupled constraints
When the coupled set of module contracts is infeasible, engineering use requires more than “solver failed”: it requires an explanation identifying a minimal conflicting subset or a quantitative witness of conflict. Producing exact minimal infeasible subsets is generally combinatorial; instead we focus on auditable certificates derived from standard convex duality / Farkas-type alternatives in cases where the assembled constraints are (or are conservatively approximated by) linear or convex.
III.D.1. Linear-inequality infeasibility and Farkas witnesses
Consider a linearized or convex-relaxed constraint system in variables \(x\in\mathbb{R}^n\) (representing shared ports, slack variables, and selected module internal interface variables):
\[
A x \le b,
\]
where each row corresponds to a normalized contract residual bound, a port admissible range, or a linearized conservation closure. (This system may arise directly, or as a relaxation used solely for diagnosis.)
Theorem III.2 (Farkas infeasibility certificate).
Exactly one of the following statements holds:
(i) There exists \(x\) such that \(A x \le b\).
(ii) There exists a vector \(\lambda\ge 0\) such that
\[
\lambda^\top A = 0, \qquad \lambda^\top b < 0.
\]
Any \(\lambda\) satisfying (ii) certifies infeasibility of \(A x\le b\).
Proof sketch.
This is a standard theorem of alternatives (Farkas' lemma) for linear inequalities. If (ii) holds and \(x\) satisfied \(A x\le b\), then \(0=\lambda^\top A x \le \lambda^\top b<0\), a contradiction. Conversely, if the system is infeasible, such a \(\lambda\) exists by the theorem of alternatives. \(\square\)
Auditable explanation format.
Given a witness \(\lambda\), the set of active/conflicting contracts can be reported as the index set
\[
\mathcal{I}(\lambda) := \{ i : \lambda_i > 0\},
\]
together with the scalar contradiction margin \(-\lambda^\top b>0\). This produces a reproducible “because” explanation: the positive-weighted combination of those inequalities implies \(0< -\lambda^\top b\) while also implying \(0\le 0\), hence inconsistency.
Gap flag (diagnostic vs ground truth).
If \(A x\le b\) is a relaxation/linearization of the true (possibly nonlinear) contract set, then an infeasibility witness explains infeasibility of the relaxed system, not necessarily the original. The workflow should therefore label certificates as either:
- exact (the assembled constraints are truly linear/convex as written), or
- diagnostic (witness for a conservative relaxation/linearization used to localize the conflict).
III.D.2. Convex infeasibility certificates (separating hyperplanes)
A similar explanation mechanism exists for convex constraints. Suppose the coupled contract set is expressed as membership in a closed convex set \(\mathcal{K}\subset\mathbb{R}^n\), i.e., find \(x\in\mathcal{K}\). If \(\mathcal{K}=\emptyset\), then for many structured representations (e.g., conic form), dual infeasibility certificates exist and can be logged as part of the solver output. Because details depend on the chosen conic modeling layer, we state only the audit requirement:
Audit requirement III.3 (Dual logging for convex contract solves).
If a convex solver is used to test coupled contract feasibility, the integration spine must log a checkable dual infeasibility certificate whenever the solver reports infeasible, including:
- the explicit dual variables,
- the dual feasibility residuals,
- the dual objective value demonstrating inconsistency,
- and a mapping from dual variables to contract rows/blocks.
III.D.3. Minimal conflicting subsets: practical approach
Exact minimal infeasible subsets are NP-hard in general. A practical, auditable compromise is:
1. Compute a certificate \(\lambda\) (linear) or a dual ray (conic) from a feasibility solve.
2. Rank contracts by weight (e.g., \(\lambda_i\) magnitude).
3. Extract a small subset of high-weight constraints and re-solve feasibility restricted to that subset; iterate until the subset is irreducible under single-constraint removal.
This procedure does not guarantee global minimality, but it produces a small, reproducible, solver-checkable explanation that can be fed into redesign, model debugging, or V\&V prioritization.
III.D.4. Contract spine outputs
For each coupled run (deterministic, UQ sample, or optimization iteration), the spine is expected to emit a machine-checkable dossier containing at least:
- the declared port list (graph edges, direction, physical dimension vectors),
- the incidence matrices \(B_q\) (or hashes thereof) for each conserved family,
- the scaled residuals \(\hat r_q\) and their tolerance budgets \(\tau_q\),
- unit-consistency check results (pass/fail with the first offending expression), and
- if infeasible: an infeasibility certificate artifact (Farkas witness or convex dual ray) and an extracted short explanation set.
The remaining sections build on this spine: stability/well-posedness conditions (Section IV) and verified numerics (Section V) assume the coupling manifold is well-defined by the port graph and its contracts, while UQ and optimization (Sections VI–XII) assume every constraint evaluation can be traced back to unit-checked ports and, when necessary, conservative tolerance budgets.
IV. Stability and Well-Posedness of Coupled Execution
This section provides mathematically checkable conditions under which a coupled digital twin (assembled from heterogeneous modules connected through the contract spine of Section III) is well-posed and stable under modular execution. We separate three concerns:
(i) physics-consistent interconnection: the coupling should not create or destroy conserved quantities except through declared sources/sinks;
(ii) numerical coupling stability: co-simulation/time stepping and communication delays should not inject nonphysical energy or destabilize the closed loop;
(iii) steady-state coupling solvability: fixed-point/Gauss\u2013Seidel style modular solves should converge to a unique coupled solution, and per-module approximation errors should admit an auditable amplification bound.
Throughout, we emphasize that stability certificates are only as credible as their assumptions (passivity indices, Lipschitz bounds, delay bounds). These parameters are therefore treated as contract-declared quantities whose estimation must be supported by separate V\&V evidence (see gap flags in IV.E).
IV.A. Port-Hamiltonian network modeling layer
IV.A.1. Effort/flow ports and power variables
For each conserved-quantity family q \u2208 \u211a\u211a (notationally distinct from the fusion gain Q), we represent inter-module exchange using paired port variables (effort, flow). For an energy-like family, effort\u00d7flow has units of power. The exact physical interpretation (temperature\u00d7entropy flow, voltage\u00d7current, etc.) depends on the module; what matters for the interconnection is that the port product has units and sign conventions consistent with the Section III incidence graph.
For module i, let u_i denote the vector of incoming port variables and y_i denote outgoing port variables. In a passivity framework, it is convenient to collect these into an input\u2013output pair (u_i, y_i) with a declared supply rate
\[
w_i(u_i,y_i) := y_i^\top u_i,
\]
interpreted as instantaneous power supplied to module i through its ports (after the unit/sign conventions have been fixed).
IV.A.2. Dissipativity/passivity as a contract-declared property
Definition IV.1 (Incremental dissipativity).
A dynamical module \(\Sigma_i\) with state \(x_i\), input \(u_i\), and output \(y_i\) is incrementally dissipative with respect to the supply rate \(w_i(u_i,y_i)=y_i^\top u_i\) if there exists a nonnegative storage function \(V_i(x_i,\bar x_i)\) such that for any two trajectories \((x_i,u_i,y_i)\) and \((\bar x_i,\bar u_i,\bar y_i)\),
\[
\dot V_i(x_i,\bar x_i) \le (y_i-\bar y_i)^\top (u_i-\bar u_i) - \psi_i(x_i,\bar x_i),
\]
where \(\psi_i\ge 0\) quantifies dissipation (possibly identically zero).
If \(\psi_i\) is positive definite in \(x_i-\bar x_i\) one obtains an incremental stability property. In practice, modules may only satisfy a \"shortage of passivity\" inequality of the form
\[
\dot V_i \le (\tilde y_i)^\top (\tilde u_i) + \nu_i \|\tilde u_i\|^2 - \psi_i,
\]
for some declared \(\nu_i\ge 0\); such indices are common in robust passivity analysis. We do not assume a particular index form; instead, the contract for module i must explicitly declare the inequality it claims to satisfy and the units/norms in which it is stated.
Audit requirement IV.2 (Passivity/dissipativity declaration).
If a stability certificate in this paper is invoked, each participating module i must provide in its contract dossier:
(a) the precise supply rate and norms used;
(b) the storage function class (or a verified numeric bound sufficient for the analysis);
(c) the numerical value(s) of passivity/dissipativity indices (e.g., \(\nu_i\), lower bounds on \(\psi_i\));
(d) an empirical or analytic justification trace (test cases, linearizations, or a posteriori estimates) supporting those values.
IV.A.3. Power-preserving interconnection and closed-loop dissipativity
Let the interconnection be power-preserving in the sense that, after wiring all port variables according to the Section III port graph (including sign conventions), the net supplied power sums to the declared external sources/sinks:
\[
\sum_{i\in V} y_i^\top u_i = P_{ext},
\]
where \(P_{ext}\) aggregates only explicitly declared injections/extractions (e.g., auxiliary heating, grid export, radiative losses modeled as sinks). A sufficient algebraic representation is that the interconnection constraints define a linear subspace \(\mathcal{D}\) of admissible (u,y) pairs such that \(\sum_i y_i^\top u_i = 0\) for internal ports; this is the classical \"Dirac structure\" condition in port-Hamiltonian interconnections.
Theorem IV.3 (Network dissipativity under power-preserving coupling).
Assume each module i is incrementally dissipative with storage \(V_i\) and dissipation term \(\psi_i\ge 0\) as in Definition IV.1. If the interconnection is power-preserving for internal ports (so that the internal contribution to \(\sum_i (\tilde y_i)^\top(\tilde u_i)\) is zero), then the interconnected network is incrementally dissipative with storage \(V := \sum_i V_i\) and dissipation \(\psi := \sum_i \psi_i\), up to the contribution of declared external ports.
Proof.
Sum the incremental dissipativity inequalities over i. The right-hand side contains \(\sum_i (\tilde y_i)^\top(\tilde u_i)\), which cancels for internal ports by the power-preserving interconnection property, leaving only declared external supply terms and \(-\sum_i \psi_i\). \(\square\)
Interpretation.
Theorem IV.3 does not assert that the physical plant is stable; it asserts that if each module contract provides a correct dissipativity inequality and if the wiring is power-preserving, then the integration itself does not introduce spurious net power injection. This addresses a common co-simulation failure mode: unstable numerical coupling caused by inconsistent port conventions.
IV.B. Structure-preserving co-simulation via scattering interconnection
IV.B.1. Motivation: decoupled time stepping and delay tolerance
In practice, modules run on different time bases (multi-rate stepping), may update asynchronously, and may communicate through sampled/held values. Even if the continuous-time interconnection is power-preserving, naive discrete-time coupling can destroy passivity. A standard remedy is the scattering (wave-variable) transformation, which converts a power-conjugate port (effort, flow) into wave variables so that the interconnection with communication delays remains passive.
We state this at an abstract level suitable for contracts: each energy-like port must declare an associated scattering parameter \(\gamma>0\) (with units chosen so that the transformation is dimensionally consistent).
IV.B.2. Scattering transformation and passivity with delays (audit-level statement)
Consider a single power port with effort e and flow f such that the instantaneous power is \(p = e^\top f\). For a chosen \(\gamma>0\), define wave variables
\[
a := \frac{1}{\sqrt{2\gamma}}(e + \gamma f), \qquad b := \frac{1}{\sqrt{2\gamma}}(e – \gamma f).
\]
Then
\[
e^\top f = \|a\|^2 – \|b\|^2.
\]
In particular, if a communication channel transmits a and returns b with arbitrary (finite) delay, and if the remote side is passive in wave variables, then the channel does not generate net energy. This identity is the basis for passivity-preserving co-simulation; details depend on channel structure and discrete-time implementation.
Audit requirement IV.4 (Scattering parameter logging and adaptation).
If scattering is used, the coupled run must log:
(a) \(\gamma\) values per port and any adaptation law;
(b) the discrete-time sampling/hold schedule per port;
(c) the measured discrete-time energy balance \(\sum (\|a_k\|^2-\|b_k\|^2)\) to verify that numerical coupling does not inject energy beyond the declared tolerance budgets of Section III.
Gap flag.
Although scattering methods are classical in networked control/co-simulation, the choice of \(\gamma\) and the numerical implementation can materially affect accuracy and convergence speed. For each integrated port, \(\gamma\) must be justified (e.g., impedance matching arguments, conditioning studies) and treated as an auditable integration parameter, not a hidden solver tweak.
IV.C. Small-gain fixed-point certificate for modular steady-state coupling
IV.C.1. Modular steady state as a fixed point
Many plant-level digital twins evaluate a steady operating point by iterating modules to consistency (e.g., block Gauss\u2013Seidel): each module computes its outputs given the latest available inputs from other modules, until shared ports stop changing.
Let x \u2208 \u211d^n stack all shared port variables (and any algebraic shared variables). Write the coupled steady-state condition as a fixed point
\[
x = F(x; d,y,\xi),
\]
where F is the \”one-sweep\” operator defined by the chosen module update schedule (ordering, relaxation, etc.). This is distinct from the Section III conservation residual r(x)=0 viewpoint; both are compatible because a fixed point of F is intended to satisfy the coupled contract manifold.
IV.C.2. Contraction/small-gain condition
Theorem IV.5 (Contraction certificate and error amplification bound).
Fix (d,y,\xi). Suppose there exists a norm \(\|\cdot\|\) on \(\mathbb{R}^n\) and a constant \(\kappa\in[0,1)\) such that
\[
\|F(x)-F(x’)\| \le \kappa \|x-x’\| \quad \text{for all } x,x’ \text{ in a closed set } \Omega\subset\mathbb{R}^n.
\]
Then:
(a) (existence/uniqueness) F has a unique fixed point \(x^\star\in\Omega\).
(b) (iterative convergence) the synchronous iteration \(x^{k+1}=F(x^k)\) converges to \(x^\star\) for any \(x^0\in\Omega\), with \(\|x^k-x^\star\|\le \kappa^k\|x^0-x^\star\|\).
(c) (implementation error amplification) if an implemented sweep \(\widetilde F\) satisfies a uniform error bound
\[
\|\widetilde F(x)-F(x)\| \le \varepsilon \quad \text{for all } x\in\Omega
\]
(for example, due to surrogate error at module interfaces, incomplete convergence of module internal solves, or contract tolerance budgets), then any fixed point \(\tilde x\) of \(\widetilde F\) in \(\Omega\) satisfies
\[
\|\tilde x – x^\star\| \le \frac{\varepsilon}{1-\kappa}.
\]
Proof.
(a)\u2013(b) are the Banach fixed-point theorem on \(\Omega\). For (c), write \(\tilde x – x^\star = \widetilde F(\tilde x)-F(x^\star)\) and add/subtract F(\tilde x):
\(\|\tilde x-x^\star\|\le \|\widetilde F(\tilde x)-F(\tilde x)\| + \|F(\tilde x)-F(x^\star)\|\le \varepsilon + \kappa\|\tilde x-x^\star\|\), yielding the bound. \(\square\)
IV.C.3. Jacobian bounding matrices as an auditable gain model
To apply Theorem IV.5 in a modular way, we can bound the sensitivity of each module’s outputs with respect to its inputs. Let x be partitioned into blocks corresponding to module ports, x=(x_1,\dots,x_M). Suppose the i-th block update map satisfies a Lipschitz bound
\[
\|F_i(x)-F_i(x’)\| \le \sum_{j=1}^M G_{ij} \|x_j-x’_j\|
\]
for some nonnegative matrix G (a \”gain matrix\”). Standard arguments imply that if the induced matrix norm \(\|G\|<1\) (in a compatible block norm) then F is a contraction. One auditable sufficient condition is
\[
\rho(G) < 1,
\]
where \(\rho(\cdot)\) is the spectral radius, together with a declared norm under which the bound is enforced.
Audit requirement IV.6 (Coupling-gain contract).
If the small-gain certificate is used, the dossier must include:
(a) the chosen block structure and norm(s);
(b) the numeric G and the method used to bound each entry (analytic derivative bounds, interval slopes, finite-difference envelopes, or a posteriori sensitivity estimates);
(c) the computed contraction metric (e.g., \(\|G\|\) or \(\rho(G)\)) and PASS/FAIL/UNKNOWN status.
PASS/FAIL/UNKNOWN triage.
If \(\rho(G)\) can be proven <1, report PASS (contraction certified). If a lower bound proves \(\rho(G)\ge 1\), report FAIL (no contraction; redesign or different execution scheme needed). If only inconclusive bounds are available, report UNKNOWN and route to V\&V prioritization (Section XI) aimed at tightening the dominant gain estimates.
IV.D. Certified asynchronous contraction mapping with communication delays
IV.D.1. Asynchronous updates as a realistic execution model
Even for steady-state solves, modules may update on different schedules (e.g., expensive blanket neutronics surrogate runs less frequently than a 0D power balance update), and they may use stale values of other ports. This motivates an asynchronous fixed-point model.
Let F: \(\Omega\to\Omega\) be as in IV.C, with a unique fixed point \(x^\star\). Consider iterations where, at each integer time k, only a subset of components is updated and updates may use delayed values:
\[
x_i^{k+1} =
\begin{cases}
F_i\big(x_1^{\tau_1^i(k)},\dots,x_M^{\tau_M^i(k)}\big), & i\in\mathcal{U}(k),\\
x_i^k, & i\notin\mathcal{U}(k),
\end{cases}
\]
where \(\mathcal{U}(k)\) is the set of updated blocks and the delay indices satisfy \(k-D\le \tau_j^i(k)\le k\) for some finite delay bound D.
IV.D.2. Convergence under a contraction assumption
Theorem IV.7 (Asynchronous convergence for contractive maps; bounded delays).
Assume:
(A1) F is a contraction on \(\Omega\) in a norm \(\|\cdot\|\) with modulus \(\kappa<1\) (Theorem IV.5).
(A2) (bounded delays) there exists an integer D\ge 0 such that all delays satisfy \(k-D\le \tau_j^i(k)\le k\).
(A3) (persistent updates) each block i is updated infinitely often: for every i, \(\{k: i\in\mathcal{U}(k)\}\) is infinite.
Then the asynchronous iterates converge to the unique fixed point \(x^\star\).
Proof sketch.
This is a classical result in asynchronous iterations for contraction mappings: bounded delays and persistent updates preserve convergence to the unique fixed point. A complete proof can be constructed by bounding the error sequence \(e^k:=\|x^k-x^\star\|\) using the contraction property and the bounded delay to show that the error decreases geometrically at a rate depending on \(\kappa\) and D. \(\square\)
Remark (Why this is non-redundant with synchronous Banach arguments).
Theorem IV.7 justifies convergence even when modules use stale information and update at different rates, which is the typical execution mode of a multi-fidelity digital twin. Synchronous contraction alone does not address this regime.
Audit requirement IV.8 (Delay-bound declarations).
If asynchronous certificates are invoked, the run must log:
(a) the maximum observed communication delay D (in iteration counts or in physical time, with unit declarations);
(b) the update schedule evidence that each block updated persistently;
(c) a justification that the contraction modulus \(\kappa\) used in IV.C is valid for the asynchronous execution (same norm and domain \(\Omega\)).
IV.E. Gap flags and empirical validation requirements
The stability results above rely on parameters that must be validated for a specific integrated digital twin. The framework therefore treats the following as gap-flag items requiring explicit evidence before any certificate is interpreted as engineering reality:
1. Passivity/dissipativity indices and storage inequalities (IV.A): must be justified per module, including units/norms.
2. Power-preserving wiring (IV.A\u2013IV.B): relies on correct sign conventions, effort\u2013flow pairing, and consistent unit declarations (Section III). Unit correctness is mechanically checkable; correct physical pairing is not and requires domain review.
3. Scattering implementation parameters (IV.B): choice of \(\gamma\), sampling/hold schemes, and the discrete-time energy accounting must be demonstrated on representative integration tests.
4. Lipschitz/gain bounds and contraction moduli (IV.C): bounding G entries or \(\kappa\) is often the hardest step; when only empirical bounds exist, the certificate must be labeled as conditional-on-tests.
5. Delay bounds and persistent updates (IV.D): must be enforced by runtime scheduling and logged; otherwise asynchronous convergence claims are not auditable.
With these caveats, Section IV provides a contract-compatible systems-theory layer: it ties together physically meaningful interface constraints (Section III) with execution-level stability checks and explicit error-amplification bounds that are later consumed by verified numerics (Section V), UQ/certification layers (Sections VI\u2013VII), and optimization (Section X).
V. Verified Numerics for the Coupled Steady State and Constraint Enclosures
This section addresses solver-dependence and hidden algebraic-loop artifacts by attaching a solver-independent numerical certificate to a coupled steady-state solution. The goal is not to produce the most efficient nonlinear solver, but to produce a checkable artifact stating (when possible): (i) a true solution exists in a declared neighborhood, (ii) it is unique there, and (iii) all constraints are satisfied within rigorous numerical bounds.
The core tools are interval arithmetic and inclusion tests (Krawczyk/interval-Newton type). These results are classical; we present them in a form compatible with the contract spine (Section III) and the convergence/error-amplification layer (Section IV).
V.A. Motivation and problem statement
V.A.1. Coupled steady state as a nonlinear equation
For fixed (d,y,\xi), let the coupled steady-state closure be represented as a system of equations
\[
r(x;d,y,\xi)=0, \qquad r: \mathbb{R}^n \to \mathbb{R}^n,
\]
where x stacks the shared port variables and any global algebraic states (Section II.C). The residual r may include (a) contract-spine conservation closures \(r_q\) (Section III), and (b) additional algebraic equalities encoding steady-state consistency between modules (e.g., power-balance and auxiliary system setpoints).
The twin also defines inequality constraints
\[
g_j(x;d,y,\xi) \le 0, \qquad j=1,\dots,J,
\]
including the Q constraint (written as \(g_Q(x)=15-Q(x)\le 0\)) and engineering/RAMI/cost/schedule constraints, all with declared units and sign conventions.
V.A.2. Why a verified layer is non-redundant
Even when a modular fixed-point iteration converges (Section IV), the limit may correspond to:
- a numerical fixed point that is not close to any true solution of \(r(x)=0\) (e.g., due to internal solver tolerances or inconsistent coupling updates),
- one of multiple distinct solutions (hysteresis/multi-root regimes), or
- a solution whose constraint status depends on floating-point roundoff or solver tolerances.
Interval inclusion tests are designed to produce an auditable PASS/FAIL/UNKNOWN classification that is independent of the initial solver used to propose a candidate solution.
V.B. Interval preliminaries and contract compatV.B.1. Interval vectors and interval extensions
An interval (box) in \(\mathbb{R}^n\) is written
\[
X = [\\underline x,\overline x] := \{x\in\mathbb{R}^n : \\underline x_i \le x_i \le \overline x_i,\ i=1,\dots,n\}.
\]
An interval extension of a function \(\varphi: \mathbb{R}^n\to\mathbb{R}^m\) is a map \([\varphi]: \mathbb{IR}^n \to \mathbb{IR}^m\) such that
\[
\varphi(X) := \{\varphi(x): x\in X\} \subseteq [\varphi](X).
\]
In practice, \([\varphi](X)\) is computed by interval arithmetic (with outward rounding) and yields guaranteed enclosures.
Audit requirement V.1 (Typed intervals).
All interval computations must respect the unit/scale declarations of Section III.B: interval bounds are carried in the same units as the corresponding variables/residuals, and any normalization to dimensionless form must be logged (scales \(\sigma\)).
V.B.2. Interval Jacobians and slope bounds
Let \(J(x):=Dr(x)\) denote the Jacobian of r. A verified method needs an enclosure \([J](X)\) such that \(J(x)\in [J](X)\) for all \(x\in X\). This may be obtained by:
- analytic derivative bounds with interval evaluation,
- algorithmic differentiation with interval arithmetic,
- or, in limited cases, finite-difference envelopes (gap-flagged unless numerically justified).
Since Jacobian bounding is often the bottleneck, the dossier must record how \([J](X)\) was obtained and what assumptions were used.
V.C. Krawczyk-type inclusion tests (existence/uniqueness)
V.C.1. Krawczyk operator
Fix a point \(\hat x\in\mathbb{R}^n\) (typically a floating-point solution guess) and an interval box X containing \(\hat x\). Choose a nonsingular preconditioning matrix \(C\in\mathbb{R}^{n\times n}\), intended as an approximate inverse of \(J(\hat x)\). Define the Krawczyk operator
\[
K(X;\hat x,C) := \hat x - C\,r(\hat x) + \bigl(I - C\,[J](X)\bigr)(X-\hat x),
\]
where matrix\(\times\)interval products are interpreted as interval evaluations producing an enclosure.
V.C.2. Existence/uniqueness theorem
Theorem V.2 (Krawczyk inclusion test).
Assume r is continuously differentiable on X and \([J](X)\) encloses \(J(x)\) for all \(x\in X\). Let \(\hat x\in X\) and C be nonsingular.
(a) If
\[
K(X;\hat x,C) \subseteq \operatorname{int}(X),
\]
then there exists a unique \(x^\star\in X\) such that \(r(x^\star)=0\).
(b) If
\[
K(X;\hat x,C) \cap X = \emptyset,
\]
then \(r(x)=0\) has no solution in X.
Otherwise the test is inconclusive (UNKNOWN).
Proof sketch.
The operator is constructed so that, under the inclusion condition, it defines a self-map on X whose fixed point corresponds to a zero of r; the inclusion into the interior yields a contraction-like argument ensuring uniqueness. The exclusion condition implies that no fixed point (hence no root) can lie in X. Full proofs are standard in verified numerics references on Krawczyk/interval Newton methods. \(\square\)
V.C.3. Practical PASS/FAIL/UNKNOWN logic
Given (d,y,\xi), a coupled run produces a candidate \(\hat x\). The verified wrapper then:
1. Constructs an initial box X (e.g., \(\hat x\pm\delta\) where \(\delta\) is based on solver tolerances, contract residual tolerances, and/or the amplification bound \(\varepsilon/(1-\kappa)\) from Theorem IV.5).
2. Builds \([J](X)\) and computes K(X).
3. Returns:
- PASS (unique solution certified) if K(X) \(\subseteq\) int(X),
- FAIL (no solution in box) if K(X) \(\cap X=\emptyset\),
- UNKNOWN otherwise, triggering refinement (V.C.4).
Audit requirement V.3 (Inclusion-log artifact).
A PASS/FAIL outcome must be accompanied by the data needed to re-check it: \(\hat x\), the full interval box X, the chosen preconditioner C, the enclosure \([J](X)\), and the computed set bounds for K(X).
V.C.4. Refinement loops and branching
UNKNOWN outcomes may be resolved by:
- shrinking X (when the candidate \(\hat x\) is believed accurate),
- improving C (better approximate inverse, better scaling),
- tightening \([J](X)\) (improve derivative bounds), or
- splitting X into sub-boxes \(X=\cup_k X_k\) and applying Theorem V.2 to each (branch-and-bound).
Because splitting is exponential in dimension, this layer is primarily intended to certify solutions near a candidate operating point, not to globally solve arbitrary nonlinear systems.
V.D. Rigorous interval evaluation of constraints
V.D.1. Constraint enclosures
Assume a unique root \(x^\star\in X\) has been certified (PASS). For each constraint \(g_j(x;d,y,\xi) \le 0\), compute a rigorous enclosure \([g_j](X)\) via interval arithmetic.
- If \(\overline{g_j}(X) \le 0\), then \(g_j(x^\star)\le 0\) is certified.
- If \(\\underline{g_j}(X) > 0\), then \(g_j(x^\star)>0\) is certified (violation).
– Otherwise, constraint status is UNKNOWN; one may refine X or tighten the constraint evaluation.
This produces a solver-independent feasibility statement for the coupled plant constraints at (d,y,\xi).
V.D.2. Conservative treatment when uniqueness is not certified
If uniqueness/existence is UNKNOWN, the safe posture for decision-making depends on context. For safety-critical or go/no-go decisions, a conservative rule is:
– treat UNKNOWN as not certified (i.e., do not count the sample as feasible in a probability-of-feasibility lower bound unless additional evidence is provided),
and route the case to targeted model refinement or V\&V experiments.
V.E. Multiple roots, hysteresis, and robust optimization implications
V.E.1. Detecting multiple roots
Multiple operating points can arise from nonlinear couplings (e.g., different steady detachment states, different thermal equilibria, or different feasible auxiliary heating allocations). Verified numerics provides an auditable way to detect this possibility:
– If two disjoint boxes \(X_1\) and \(X_2\) both satisfy the PASS condition of Theorem V.2, then the system has at least two distinct solutions.
Conversely, inability to certify uniqueness globally does not prove multiplicity; it only flags that the local certificate is insufficient.
V.E.2. Continuation and hysteresis audits
To study parameter-dependent hysteresis, consider a scalar continuation parameter \(\lambda\) (e.g., density setpoint or edge-transport multiplier) and a family \(r(x;\lambda)=0\). One can perform continuation with verified steps:
– at each \(\lambda_k\), use a previous certified box as a predictor to initialize a new box,
– attempt a PASS certificate; if FAIL, bracket and split to locate branch termination.
The resulting audit artifact is a set of certified boxes \(\{(\lambda_k,X_k)\}\) that document where the coupled twin has a unique steady state and where it does not.
V.E.3. Implications for robust design
If multiple certified solutions exist for the same (d,y,\xi), then a design is not well-defined unless an operational selection rule is specified (e.g., a control policy selecting a particular branch). In a conservative optimization loop, one may instead enforce constraints for the worst-case over all certified solutions, which can be phrased as:
\[
\max_{x\in\mathcal{S}(d,y,\xi)} g_j(x;d,y,\xi) \le 0,
\]
where \(\mathcal{S}(d,y,\xi)\) is the (possibly set-valued) steady-state solution set. Section X discusses how such conservatism interacts with convexification and discrete decisions.
V.F. Integration with UQ sampling and feasibility probability estimation
For sampling-based UQ layers (Sections VI–VII), the verified solver wrapper yields, for each sample \(\xi^{(k)}\), a tri-valued outcome for each constraint: certified safe / certified violated / unknown. A conservative lower bound on the feasible-event indicator uses
\[
\mathbf{1}_{\mathrm{cert\_feas}}^{(k)} := \prod_{j=1}^J \mathbf{1}\{\overline{g_j}(X^{(k)})\le 0\},
\]
counting a sample as feasible only if all constraints are certified safe. This can be combined with scenario/DRO/conformal methods, provided the sampling/certification logic is recorded (in particular, the handling of UNKNOWN samples must be explicit because it affects any claimed probability lower bounds).
Audit requirement V.4 (Per-sample certificate attachment).
For any probability-of-feasibility statement derived from samples, the dossier must include, per sample or per hashable batch:
– PASS/FAIL/UNKNOWN status for the root inclusion test,
– the certified box X (if PASS) or the attempted boxes (if UNKNOWN/FAIL),
– interval bounds for each constraint \([g_j](X)\), and
– the rule used to map UNKNOWN to the final probability bound.
V.G. Output artifact format (design dossier expectations)
A verified steady-state evaluation at a design point (d,y,\xi) should export a compact, replayable artifact containing:
1. Candidate solution \(\hat x\) (floating-point) and solver metadata.
2. Verified box X and inclusion-test data (Theorem V.2): C, \([J](X)\), K(X) bounds, and the PASS/FAIL/UNKNOWN decision.
3. Interval bounds for all plant-level constraints \([g_j](X)\) and margins \([m_j](X)=-[g_j](X)\).
4. A linkage to the contract spine: the port graph hash, units/dimension-vector checks, and tolerance budgets used to construct X.
This verified-numerics layer is intentionally orthogonal to any particular physics module: it only requires (i) a residual function r for the coupled closure and (ii) interval evaluators for r, its Jacobian enclosure, and the constraints g. The next sections use these artifacts as building blocks for scalable surrogate modeling and certified uncertainty propagation.
VI. Surrogates and Certified Uncertainty Propagation at Scale
This section describes a surrogate-and-UQ layer intended to (i) reduce inner-loop cost for design and certification, while (ii) preserving auditable, conservative statements by explicitly tracking surrogate error and numerical/integration uncertainty. The emphasis is on constructions that can produce checkable error bounds or conservative envelopes. Any claim that a given plant map is “low-degree” or “sparse” is treated as a gap-flag item requiring device-specific validation.
VI.A. Polynomial chaos expansions (PCE) as a common surrogate layer
VI.A.1. Setup and notation
Let \(\xi\in\Xi\subset\mathbb{R}^{n_\xi}\) denote uncertain inputs equipped with a reference probability measure \(\mu\) (e.g., a product measure when independence is assumed). Let \(f(\xi)\) denote a scalar quantity returned by the coupled twin (e.g., \(Q\), a margin, or a constraint function \(g_j\)). A generalized polynomial chaos expansion approximates \(f\) by
\[
f(\xi) \approx f_P(\xi) := \sum_{\alpha\in\mathcal{A}_P} c_\alpha\,\psi_\alpha(\xi),
\]
where \(\{\psi_\alpha\}\) are multivariate polynomials orthonormal in \(L^2(\mu)\) and \(\mathcal{A}_P\) is an index set determined by a truncation rule (total degree \(\le P\), hyperbolic truncation, etc.). The coefficient vector \(c\) is what must be learned from data.
For vector outputs \(f\in\mathbb{R}^m\) we use componentwise expansions. When \(f\) is a constraint, recall the sign convention from Section II: \(g_j\le 0\) is safe.
VI.A.2. Moments and constraints from PCE coefficients
Because \(\psi_0\equiv 1\) and the basis is orthonormal,
\[
\mathbb{E}_\mu[f_P] = c_0, \qquad \mathrm{Var}_\mu(f_P)= \sum_{\alpha\ne 0} c_\alpha^2.
\]
These identities are algebraic and auditable once the basis \(\{\psi_\alpha\}\), truncation \(\mathcal{A}_P\), and coefficients \(c\) are logged.
Gap flag (basis/measure mismatch).
If the true uncertainty distribution differs from the reference \(\mu\) (e.g., dependence, heavier tails, or distribution shift), PCE-derived moments are moments under \(\mu\), not necessarily under reality. Later probability-of-feasibility layers (Section VII) must therefore declare whether they are (a) conditional on \(\mu\), (b) distributionally robust (e.g., Wasserstein DRO), or (c) purely finite-sample (scenario/conformal).
VI.A.3. Surrogate/model-form error as a deterministic robustifier
Let \(g_j(\xi)\) denote the true (expensive) constraint function evaluation produced by the coupled twin (possibly including the verified solve wrapper of Section V). Let \(\widehat g_j(\xi)\) be a surrogate (e.g., PCE).
Define a deterministic error envelope \(\Delta_j \ge 0\) such that for all \(\xi \in \Xi_0\) (the domain of validity),
\[
|g_j(\xi) – \widehat g_j(\xi)| \le \Delta_j.
\]
The surrogate-safe constraint is then
\[
\widehat g_j(\xi) + \Delta_j \le 0 \;\Rightarrow\; g_j(\xi)\le 0.
\]
**Audit requirement VI.1 (surrogate envelope provenance).** Any use of \(\Delta_j\) must log (i) the definition of \(\Xi_0\), (ii) whether the bound is deterministic or probabilistic, and (iii) the evidence used to compute it (cross-validation residual bounds, interval arithmetic bounds, reduced-basis a posteriori bounds, etc.). If only empirical evidence exists, label \(\Delta_j\) as conditional-on-tests.
We introduce an explicit deterministic envelope \(\Delta_j\ge 0\) such that
\[
|g_j(\xi)-\widehat g_j(\xi)| \le \Delta_j \quad \text{for all } \xi\in\Xi_0,
\]
for a declared uncertainty domain \(\Xi_0\) (or with high probability, if only probabilistic surrogate error control is available).
A conservative surrogate-safe constraint is then
\[
\widehat g_j(\xi) + \Delta_j \le 0 \;\Rightarrow\; g_j(\xi)\le 0.
\]
Audit requirement VI.1 (surrogate envelope provenance).
Any use of \(\Delta_j\) must log (i) the definition of \(\Xi_0\), (ii) whether the bound is deterministic or probabilistic, and (iii) the evidence used to compute it (cross-validation residual bounds, interval arithmetic bounds, reduced-basis a posteriori bounds, etc.). If only empirical evidence exists, label \(\Delta_j\) as conditional-on-tests.
VI.B. Compressed sensing recovery of sparse PCE coefficients
A truncated PCE can be high-dimensional: \(|\mathcal{A}_P|\) grows combinatorially with \(n_\xi\) and \(P\). A common mitigation is to assume that the coefficient vector \(c\) is sparse or compressible in the chosen basis and to recover it from \(N\) samples with \(N \ll |\mathcal{A}_P|\) using \(\ell_1\)-regularized regression.
VI.B.1. Linear observation model
Given sample points \(\{\xi^{(k)}\}_{k=1}^N\), define the design matrix \(\Psi\in\mathbb{R}^{N\times M}\) with \(M:=|\mathcal{A}_P|\) by
\[
\Psi_{k,\alpha} := \psi_\alpha(\xi^{(k)}).
\]
Let \(y\in\mathbb{R}^N\) be observations
\[
y_k = f(\xi^{(k)}) + \eta_k,
\]
where \(\eta\) aggregates (a) numerical error, (b) integration/coupling error, and (c) any measurement noise. In vector form,
\[
y = \Psi c + \eta.
\]
VI.B.2. Basis pursuit denoising and a standard recovery bound
Define the \(\ell_1\)-regularized recovery
\[
\widehat c \in \arg\min_{c\in\mathbb{R}^M} \|c\|_1 \quad \text{s.t.}\quad \|\Psi c – y\|_2 \le \varepsilon,
\]
where \(\varepsilon\) is a declared residual tolerance (in the same units as \(f\); if \(f\) is normalized to be dimensionless, then \(\varepsilon\) is dimensionless).
Theorem VI.2 (A canonical compressed-sensing error bound; conditional).
Assume the sampling and basis are such that \(\Psi\) satisfies a restricted isometry property (RIP) of order \(2s\) with constant \(\delta_{2s}\) below a standard sufficient threshold, and that \(\|\eta\|_2\le \varepsilon\). Let \(c_s\) be the best \(s\)-sparse approximation to \(c\) in \(\ell_1\) (i.e., keep the \(s\) largest magnitudes). Then there exist constants \(C_0,C_1\) (depending only on \(\delta_{2s}\)) such that
\[
\|\widehat c – c\|_2 \le C_0\,\frac{\|c-c_s\|_1}{\sqrt{s}} + C_1\,\varepsilon.
\]
Consequently, for any linear functional \(L\) on coefficients (e.g., mean, variance contributions), the induced error can be bounded by operator norms times \(\|\widehat c-c\|_2\).
Gap flag (RIP/coherence verification).
For polynomial chaos matrices, RIP-type conditions are nontrivial to verify and can fail under poor sampling or high degree. Therefore, any use of Theorem VI.2 must be labeled as conditional on a validated sampling design (e.g., randomized sampling matched to \(\mu\), preconditioning, and empirical coherence diagnostics).
VI.B.3. Propagation of coefficient error to moments and tail surrogates
From orthonormality, the mean error satisfies
\[
|\widehat{\mathbb{E}}[f]-\mathbb{E}[f]| = |\widehat c_0 – c_0| \le \|\widehat c-c\|_2.
\]
For the variance, write \(v(c):=\sum_{\alpha\ne 0} c_\alpha^2\). Then
\[
|v(\widehat c)-v(c)| = |\langle \widehat c_{\neq 0}+c_{\neq 0},\;\widehat c_{\neq 0}-c_{\neq 0}\rangle|
\le (\|\widehat c_{\neq 0}\|_2+\|c_{\neq 0}\|_2)\,\|\widehat c-c\|_2,
\]
which is auditable provided a bound or estimate on \(\|c_{\neq 0}\|_2\) is recorded (e.g., from \(\|\widehat c_{\neq 0}\|_2\) plus the recovery error).
These moment bounds can be used by Section VII.A moment-based probability bounds, but any resulting feasibility probability statement must log that it is conditional on (i) the PCE truncation and (ii) the compressed-sensing recovery assumptions.
VI.C. Certified multilevel Monte Carlo (MLMC) for feasibility metrics
PCE is not always appropriate (non-smooth responses, discontinuous feasibility boundaries, heavy tails). Sampling-based estimators remain essential, but naive Monte Carlo can be too expensive at high fidelity. MLMC reduces cost by combining many cheap low-fidelity samples with a few expensive high-fidelity samples.
VI.C.1. MLMC estimator for a scalar QoI
Let \(P\) be a scalar QoI (e.g., a smoothed feasibility indicator, or net power). Let \(P_\ell\) denote a hierarchy of approximations (\(\ell=0,1,\dots,L\)) of increasing fidelity/cost, with \(P_L\) the highest fidelity available in the workflow (possibly including the verified solve wrapper of Section V).
Using the telescoping identity
\[
\mathbb{E}[P_L] = \mathbb{E}[P_0] + \sum_{\ell=1}^L \mathbb{E}[P_\ell – P_{\ell-1}],
\]
an MLMC estimator is
\[
\widehat{\mathbb{E}}[P_L] := \frac{1}{N_0}\sum_{k=1}^{N_0} P_0^{(k)} + \sum_{\ell=1}^L \frac{1}{N_\ell}\sum_{k=1}^{N_\ell} \bigl(P_\ell^{(k)}-P_{\ell-1}^{(k)}\bigr),
\]
where coupled pairs \((P_\ell^{(k)},P_{\ell-1}^{(k)})\) are evaluated using the same underlying sample \(\xi^{(k)}\) to reduce variance.
VI.C.2. A standard complexity statement (assumption-explicit)
Proposition VI.3 (Canonical MLMC cost scaling; conditional).
Assume there exist positive constants \(\alpha,\beta,\gamma\) such that:
(i) Bias decay: \(|\mathbb{E}[P_L-P]| \le c_1 2^{-\alpha L}\).
(ii) Variance decay: \(\mathrm{Var}(P_\ell-P_{\ell-1}) \le c_2 2^{-\beta \ell}\).
(iii) Cost growth: the expected cost of one coupled sample at level \(\ell\) is \(\le c_3 2^{\gamma \ell}\).
Then, by choosing \(L\) to meet a target bias and selecting \(N_\ell\) proportional to \(\sqrt{\mathrm{Var}(P_\ell-P_{\ell-1})/\mathrm{Cost}_\ell}\), one can achieve mean-square error \(\mathbb{E}[|\widehat{\mathbb{E}}[P_L]-\mathbb{E}[P]|^2]\le \epsilon^2\) with total expected cost scaling as:
– \(\mathcal{O}(\epsilon^{-2})\) if \(\beta>\gamma\),
– \(\mathcal{O}(\epsilon^{-2}\log^2\epsilon^{-1})\) if \(\beta=\gamma\),
– \(\mathcal{O}(\epsilon^{-2-(\gamma-\beta)/\alpha})\) if \(\beta<\gamma\).
This statement is conditional on the bias/variance/cost models and must be empirically validated for the chosen hierarchy.
Gap flag (hierarchy design).
For a coupled stellarator plant twin, level definitions are architecture choices (0D/+1D power balance, add-on surrogates, reduced-basis blankets, limited high-fidelity calls). The exponents \(\alpha,\beta,\gamma\) are not universal constants; they must be estimated from pilot runs and logged. If they cannot be reliably estimated, MLMC must be treated as a heuristic acceleration rather than a certificate-bearing estimator.
VI.D. Multi-fidelity MLMC (MF-MLMC) with bias control for feasibility events
VI.D.1. The feasibility event and smooth indicators
The plant feasibility event for fixed \((d,y)\) is
\[
\mathcal{F} := \{\xi : Q(d,y,\xi)\ge 15 \ \wedge\ g_j(d,y,\xi)\le 0 \ \forall j\}.
\]
The raw indicator \(\mathbf{1}_{\mathcal{F}}\) is discontinuous, which can degrade variance decay. A common device is to replace it by a smooth surrogate \(\phi_\delta\) applied to the maximum normalized violation:
\[
v(\xi) := \max\Bigl\{\frac{15-Q(d,y,\xi)}{s_Q},\ \frac{g_1(d,y,\xi)}{s_1},\dots,\frac{g_J(d,y,\xi)}{s_J}\Bigr\},
\]
with declared positive scales \(s_\cdot\) (Section III.B), and a smooth monotone function \(\phi_\delta\) such that
\[
\mathbf{1}\{v\le 0\} \le \phi_\delta(v) \le \mathbf{1}\{v\le \delta\}.
\]
Then \(\mathbb{E}[\phi_\delta(v)]\) upper-bounds \(\mathbb{P}(\mathcal{F})\) only in a one-sided sense; it also introduces a smoothing bias controlled by \(\delta\). For conservative lower bounds on \(\mathbb{P}(\mathcal{F})\), one instead uses a smoothed *inner* approximation, e.g., \(\psi_\delta(v)=\mathbf{1}\{v\le -\delta\}\) or a smooth function satisfying \(\mathbf{1}\{v\le -\delta\}\le \psi_\delta(v)\le \mathbf{1}\{v\le 0\}\).
Audit requirement VI.4 (indicator approximation choice).
Any probability-of-feasibility calculation using smooth indicators must record (i) the definition of \(v\) (including scales), (ii) the smoothing function and parameter \(\delta\), and (iii) whether the resulting estimate is a lower bound, upper bound, or neither.
VI.D.2. Multi-fidelity control variates and explicit bias accounting
Multi-fidelity estimators generalize MLMC by using correlated low-fidelity models as control variates for a high-fidelity target. A conservative, auditable approach is to separate two errors:
- Sampling error: controlled by concentration/finite-sample bounds (Sections VII.B and VII.E).
- Model bias: controlled by explicit bias bounds \(|\mathbb{E}[P_L-P]|\le b_L\) or by verified comparisons on a held-out validation set. If no bound exists, the bias must be treated as an additional epistemic uncertainty, not ignored.
In the present framework, bias control is contractually tied to model hierarchy: each level declares either (a) an a posteriori error estimator (e.g., reduced-basis bounds), or (b) a validation-derived envelope \(\Delta\) that can be propagated into conservative constraint tightening.
VI.E. Certified active-subspace dimension reduction
High-dimensional uncertainty and design spaces motivate dimension reduction. Active-subspace methods seek a low-dimensional linear subspace in which a scalar risk functional varies most.
VI.E.1. Decision-relevant risk functional
Define a scalar, dimensionless risk functional for fixed \((d,y)\), for example
\[
R(d,y,\xi) := \operatorname{lse}_\tau\Bigl(\frac{15-Q(d,y,\xi)}{s_Q},\ \frac{g_1(d,y,\xi)}{s_1},\dots,\frac{g_J(d,y,\xi)}{s_J}\Bigr),
\]
where \(\operatorname{lse}_\tau(a_1,\dots,a_m):=\tau\log\big(\sum_{i=1}^m e^{a_i/\tau}\big)\) is a smooth approximation of \(\max\{a_i\}\) and \(\tau>0\) is chosen and logged. Small \(\tau\) makes \(R\) close to max-violation; larger \(\tau\) improves smoothness.
VI.E.2. Active-subspace matrix and a conservative surrogate notion
Assume \(R\) is differentiable in \(\xi\) and that gradients can be computed (e.g., by automatic differentiation through surrogates, or by finite differences with logged step sizes and error checks). Define
\[
C := \mathbb{E}_\mu\big[\nabla_\xi R\, (\nabla_\xi R)^\top\big].
\]
Let \(C=W\Lambda W^\top\) be an eigen-decomposition with eigenvalues \(\lambda_1\ge\cdots\ge\lambda_{n_\xi}\ge 0\). For a chosen dimension \(r\), set \(W_1\in\mathbb{R}^{n_\xi\times r}\) to the leading eigenvectors and project \(\xi\mapsto y:=W_1^\top\xi\). A ridge surrogate \(\widetilde R(d,y,\cdot)\) is then fit so that
\[
R(d,y,\xi) \approx \widetilde R(d,y,W_1^\top\xi).
\]
VI.E.3. What can be certified
A fully rigorous, distribution-free certification of active-subspace error for general black-box \(R\) is not available without strong assumptions. However, two auditable partial guarantees are often feasible:
(a) Matrix estimation error: use concentration bounds for empirical covariance matrices (e.g., matrix Bernstein) to bound \(\|\widehat C-C\|\) given bounded-gradient assumptions. This yields a certified uncertainty on estimated eigenvalues/eigenspaces.
(b) Tail-energy reporting: report the residual gradient energy \(\sum_{i>r} \widehat\lambda_i\) as an auditable diagnostic. If this tail is not small relative to \(\sum_i \widehat\lambda_i\), then no low-dimensional reduction should be claimed.
Audit requirement VI.5 (active-subspace gap flag).
Any use of active subspaces must log (i) the risk functional \(R\) and smoothing parameters, (ii) the gradient computation method, (iii) the sample size used to estimate \(C\), (iv) the estimated spectrum \(\widehat\lambda_i\), and (v) a declared decision rule for selecting r (e.g., eigenvalue gap threshold or tail-energy tolerance). Claims that reduced-dimensional UQ is “certified” must be explicitly conditioned on the assumptions used to bound \(\|\widehat C-C\|\) and on evidence that gradients are reliable.
VI.F. Interface to the verified numerics layer
The surrogate and UQ tools above are most defensible when they sit on top of the verified-numerics wrapper (Section V) rather than replacing it. Concretely:
– Training/evaluation data should record whether each sample evaluation was PASS/FAIL/UNKNOWN in the inclusion test.
– If a surrogate is fit to constraint values, the training set must specify how UNKNOWN samples were handled (excluded, treated as censored intervals, or conservatively labeled unsafe).
– Any surrogate-based feasibility claim must incorporate explicit envelopes \(\Delta_j\) that include (at least) surrogate fit error and any unresolved verified-numerics UNKNOWN regions.
With these interfaces, Section VI supplies the computational machinery used later in Section VII (probability-of-feasibility statements), Section X (optimization), and Section XI (V\&V prioritization). The next section turns these surrogate/UQ outputs into explicit, auditable probability-of-feasibility and safety certificates.
VII. Probability-of-Feasibility and Safety Certificates (Without Unverifiable Monte Carlo Claims)
This section defines several complementary ways to issue auditable statements of the form
\[
\mathbb{P}(\mathcal{F}(d,y)) \ge 1-\alpha,
\]
where \(\mathcal{F}(d,y)\) is the joint feasibility event (Q constraint, engineering constraints, RAMI, cost, schedule). Each method has a different set of assumptions and therefore a different role in a decision dossier. When the assumptions required for a method cannot be defended for a specific module, the method must be marked as conditional and replaced by a more conservative alternative (e.g., scenario or conformal).
VII.A. Distribution-free probability lower bounds from moments
Moment-based bounds convert (estimated) means and variances of margins into conservative probability statements without requiring a parametric distribution (Gaussianity, tail models). They can be computed cheaply from PCE coefficients (Section VI.A) when the reference measure \(\mu\) is trusted, or from samples.
VII.A.1. One-sided Cantelli bound for a single constraint
Let \(g(\xi)\) be a scalar constraint with safe convention \(g(\xi)\le 0\). Define the margin random variable
\[
M(\xi) := -g(\xi),
\]
so feasibility of the constraint is the event \(\{M\ge 0\}\). Assume \(\mu_M:=\mathbb{E}[M]\) exists and \(\sigma_M^2:=\mathrm{Var}(M)<\infty\).
Proposition VII.1 (Cantelli one-sided bound, distribution-free).
If \(\mu_M>0\), then
\[
\mathbb{P}(M\ge 0) \;\ge\; 1-\frac{\sigma_M^2}{\sigma_M^2+\mu_M^2} \;=\; \frac{\mu_M^2}{\sigma_M^2+\mu_M^2}.
\]
If \(\mu_M\le 0\), this bound gives no nontrivial guarantee.
Proof.
Cantelli’s inequality states that for any random variable \(X\) with finite mean \(\mu\) and variance \(\sigma^2\), and any \(a>0\),
\[
\mathbb{P}(X-\mu\le -a) \le \frac{\sigma^2}{\sigma^2+a^2}.
\]
Set \(X=M\), \(a=\mu_M\), so \(\mathbb{P}(M\le 0)=\mathbb{P}(M-\mu_M\le -\mu_M)\le \sigma_M^2/(\sigma_M^2+\mu_M^2)\). Rearranging yields the claim. \(\square\)
Audit requirement VII.2 (moment provenance).
Any use of Proposition VII.1 must log:
– how \(\mu_M\) and \(\sigma_M^2\) were computed (PCE moments under \(\mu\), sample moments, or conservative enclosures),
– whether the result is conditional on a reference distribution \(\mu\) (PCE) or is purely empirical,
– how numerical/integration/surrogate error envelopes (Sections III and VI) were incorporated (e.g., by tightening \(M\) to \(M-\Delta\)).
VII.A.2. Joint constraints by union bound (risk budgeting)
Let constraints be \(g_j(\xi)\le 0\) with margins \(M_j=-g_j\). Define the joint feasibility event
\[
\mathcal{F} := \bigcap_{j=1}^J \{M_j\ge 0\}.
\]
Then, without any independence assumptions,
\[
\mathbb{P}(\mathcal{F}) = 1-\mathbb{P}\Big(\bigcup_{j=1}^J \{M_j<0\}\Big) \ge 1-\sum_{j=1}^J \mathbb{P}(M_j<0).
\]
Combining with Proposition VII.1 yields an auditable but potentially loose bound.
A practical way to issue a target \(1-\alpha\) joint certificate is to allocate per-constraint risk budgets \(\alpha_j\ge 0\) with \(\sum_j \alpha_j\le \alpha\), and require (using any single-constraint method)
\[
\mathbb{P}(M_j\ge 0) \ge 1-\alpha_j \quad \forall j.
\]
The dossier must record the allocation \(\alpha_j\) and the rationale (engineering criticality, tightness of available bounds, or dual sensitivities from optimization in Section X).
Gap flag (tightness).
Union bounds can be very conservative when constraints are strongly correlated. Any design decision that is sensitive to this conservatism should be cross-checked with a dependence-aware method (scenario or conformal) or with explicit dependence modeling (Section VIII for RAMI; otherwise empirical validation).
VII.B. Scenario optimization with finite-sample guarantees
Scenario methods provide finite-sample, solver-checkable guarantees for chance constraints under i.i.d. sampling from a specified uncertainty distribution. Unlike moment bounds, they do not require low-order moment sufficiency; unlike pure Monte Carlo, they provide explicit out-of-sample violation probability bounds.
VII.B.1. Scenario program for joint feasibility
Fix a design \((d,y)\). Let \(g(d,y,\xi)\in\mathbb{R}^J\) stack constraints with \(g\le 0\) safe. Draw i.i.d. samples \(\xi^1,\dots,\xi^N\) from the declared distribution model for \(\xi\). Consider the scenario-feasibility condition
\[
g(d,y,\xi^k) \le 0, \qquad k=1,\dots,N.
\]
In optimization mode (Section X), one solves an objective subject to these constraints; in certification mode, one checks feasibility at a fixed \((d,y)\).
Crucially, in this paper's workflow, each sampled evaluation should be performed through the verified wrapper (Section V) when feasible, so that the scenario constraints correspond to certified-safe checks rather than solver-dependent values.
VII.B.2. A standard sample-complexity bound (assumption explicit)
We state a canonical form of the scenario guarantee for convex programs; it is included as a blueprint for auditability, but its applicability must be checked case-by-case.
Theorem VII.3 (Scenario bound, convex case; conditional).
Consider an optimization problem with decision variable \(x\in\mathbb{R}^{n_x}\) (here \(x\) may represent a convexified parameterization of \((d,y)\) or continuous sub-decisions) and convex constraint functions \(g(x,\xi)\le 0\). Assume:
- \(\xi^k\) are i.i.d. samples from the declared distribution;
- the scenario program is convex and attains a (measurable) optimizer \(\hat x_N\);
- the number of support constraints \(s\) is almost surely bounded by \(n_x\) (a standard assumption for nondegenerate convex programs).
Then for any target violation probability \(\varepsilon\in(0,1)\),
\[
\mathbb{P}\big( V(\hat x_N) > \varepsilon \big) \le \sum_{i=0}^{n_x-1} \binom{N}{i} \varepsilon^i (1-\varepsilon)^{N-i},
\]
where \(V(\hat x_N):=\mathbb{P}(g(\hat x_N,\xi)\nleq 0)\) is the true violation probability under the sampling distribution.
Discussion (audit use).
Given a confidence parameter \(\beta\), one can choose \(N\) such that the right-hand side is \(\le \beta\), yielding an explicit statement: with confidence \(1-\beta\) over the random draw of the scenario set, the returned design has violation probability \(\le \varepsilon\).
Gap flag (nonconvexity and discrete variables).
For nonconvex models and mixed-integer decisions, classical scenario bounds do not apply directly. In such cases, scenario constraints are still useful as a conservative test set, but any probability guarantee must be derived from a method whose assumptions match the solver/model class (e.g., conformal prediction for black-box predictors, or a conservative relaxation with logged conditions).
Audit requirement VII.4 (replayable scenario sets).
Any scenario-based certificate must log: random seed, sampling distribution version/hash, the full sample set (or a cryptographic hash with a retrieval path), and the verified-numerics PASS/FAIL/UNKNOWN outcomes for each sampled constraint evaluation.
VII.C. Wasserstein distributionally robust chance constraints (DRO)
DRO aims to certify feasibility under distribution shift: the true distribution \(\mathbb{P}\) may differ from the nominal/sampled distribution. A practical approach is to define an ambiguity set \(\mathcal{P}\) of distributions close to the empirical distribution \(\widehat{\mathbb{P}}_N\) in Wasserstein distance, and enforce safety for all \(\mathbb{Q}\in\mathcal{P}\).
VII.C.1. Wasserstein balls and Lipschitz requirements
Fix a metric \(d(\cdot,\cdot)\) on \(\Xi\), and define a 1-Wasserstein ball
\[
\mathcal{P}_{\rho} := \{\mathbb{Q} : W_1(\mathbb{Q},\widehat{\mathbb{P}}_N) \le \rho\}.
\]
Let \(g(\xi)\) be a scalar constraint (safe if \(g\le 0\)). A key quantity is a Lipschitz constant \(L_g\) such that
\[
|g(\xi)-g(\xi’)| \le L_g\, d(\xi,\xi’) \quad \forall \xi,\xi’\in\Xi_0.
\]
Because \(L_g\) is rarely available for complex coupled twins, this is a gap-flag item: any DRO certificate must log how \(L_g\) was bounded (analytic bounds, interval slopes, or conservative envelopes).
VII.C.2. A conservative DRO chance-certificate via worst-case CVaR
A sufficient condition for a chance constraint is \(\mathrm{VaR}_{\alpha}(g)\le 0\), which is in turn implied by \(\mathrm{CVaR}_{\alpha}(g)\le 0\). We therefore propose a conservative, auditable DRO certificate by upper bounding the worst-case CVaR over \(\mathcal{P}_{\rho}\).
Proposition VII.5 (Wasserstein-robust CVaR upper bound; conditional).
Assume \(g\) is \(L_g\)-Lipschitz on \(\Xi_0\) under metric \(d\). For any \(t\in\mathbb{R}\), define \(\ell_t(\xi):=(g(\xi)-t)_+\). Then \(\ell_t\) is also \(L_g\)-Lipschitz, and
\[
\sup_{\mathbb{Q}\in\mathcal{P}_{\rho}} \mathbb{E}_{\mathbb{Q}}[\ell_t] \le \frac{1}{N}\sum_{k=1}^N (g(\xi^k)-t)_+ + \rho L_g.
\]
Consequently,
\[
\sup_{\mathbb{Q}\in\mathcal{P}_{\rho}} \mathrm{CVaR}^{\mathbb{Q}}_{\alpha}(g)
\le \inf_{t\in\mathbb{R}} \Big( t + \frac{1}{\alpha}\Big[\frac{1}{N}\sum_{k=1}^N (g(\xi^k)-t)_+ + \rho L_g\Big]\Big).
\]
Therefore, if the right-hand side is \(\le 0\), then \(\mathrm{CVaR}^{\mathbb{Q}}_{\alpha}(g)\le 0\) for all \(\mathbb{Q}\in\mathcal{P}_{\rho}\), hence \(\mathrm{VaR}^{\mathbb{Q}}_{\alpha}(g)\le 0\) for all \(\mathbb{Q}\in\mathcal{P}_{\rho}\), and thus
\[
\inf_{\mathbb{Q}\in\mathcal{P}_{\rho}} \mathbb{Q}(g\le 0) \ge 1-\alpha.
\]
Proof sketch.
The expectation bound is a standard Kantorovich\u2013Rubinstein argument: the worst-case shift in expectation over a 1-Wasserstein ball is bounded by \(\rho\) times the Lipschitz constant of the integrand. Apply this to \(\ell_t\), then use the variational form of CVaR: \(\mathrm{CVaR}_{\alpha}(g)=\inf_t\{t + \frac{1}{\alpha}\mathbb{E}[(g-t)_+]\}\). \(\square\)
Audit requirement VII.6 (DRO parameters).
Any DRO-based certificate must log: metric \(d\), radius \(\rho\), Lipschitz bound \(L_g\), the optimization over \(t\) (including attained \(t^\star\)), and how surrogate/verified-numerics UNKNOWN outcomes were mapped into the sample values \(g(\xi^k)\).
VII.D. Sum-of-squares (SOS) certificates for polynomial feasible regions
SOS methods apply when constraints can be represented (or conservatively approximated) as polynomials in a bounded uncertainty set. In this paper, the most plausible pathway is: a PCE surrogate yields polynomial approximations \(\widehat g_j(\xi)\) and the uncertainty support is a semialgebraic set \(K\) (e.g., a box for standardized germs). Then SOS can certify *robust* feasibility on \(K\), which implies probabilistic feasibility for any distribution supported on \(K\).
VII.D.1. Robust feasibility over a semialgebraic support
Let
\[
K := \{\xi\in\mathbb{R}^{n_\xi}: \kappa_m(\xi)\ge 0,\ m=1,\dots,M\}
\]
be a compact basic closed semialgebraic set (e.g., \([-1,1]^{n_\xi}\) can be written via \(1-\xi_i^2\ge 0\)). Suppose \(\widehat g_j(\xi)\) are polynomials.
Proposition VII.7 (SOS robust-feasibility certificate; conditional).
If for each j there exist SOS polynomials \(\sigma_{j,0},\sigma_{j,1},\dots,\sigma_{j,M}\) such that
\[
-\widehat g_j(\xi) = \sigma_{j,0}(\xi) + \sum_{m=1}^M \sigma_{j,m}(\xi)\, \kappa_m(\xi),
\]
then \(\widehat g_j(\xi)\le 0\) for all \(\xi\in K\). If additionally the true constraint satisfies \(g_j(\xi)\le \widehat g_j(\xi)+\Delta_j\) on K (Section VI.A.3), then
\[
\widehat g_j(\xi)+\Delta_j\le 0\ \forall \xi\in K \;\Rightarrow\; g_j(\xi)\le 0\ \forall \xi\in K.
\]
Audit requirement VII.8 (SOS scaling and conservatism).
Any SOS certificate must log: the defining polynomials \(\kappa_m\), degree limits used in the SOS search, the resulting SOS polynomials (or a hash with reproducible solver settings), and the surrogate envelopes \(\Delta_j\). The dossier must also explicitly state that the certificate is robust-over-K and therefore may be conservative relative to the true uncertainty distribution.
Gap flag (scaling).
SOS/SDP computations can become intractable as \(n_\xi\) or degrees grow. Therefore, SOS is best viewed here as a local-to-moderate dimension tool (e.g., after active-subspace reduction, or for a small set of dominant uncertainties), not a universal plant-level method.
VII.E. Conformal prediction tubes for joint feasibility
Conformal prediction provides finite-sample, distribution-free coverage guarantees under an exchangeability assumption, and is applicable to black-box predictors (including heterogeneous surrogate stacks). It therefore serves as a certificate layer when distributional modeling is uncertain or when convexity assumptions needed for scenario/DRO do not hold.
VII.E.1. Prediction sets for a scalar score encoding joint feasibility
Let the coupled surrogate stack produce predictions \(\widehat z(\xi)\) for a vector of outputs sufficient to compute constraints. Define a scalar, dimensionless score that encodes joint violation; one convenient choice is the maximum normalized violation computed from the *true* coupled twin output \(z(\xi)\):
\[
s(\xi) := \max\Bigl\{\frac{15-Q(\xi)}{s_Q},\ \frac{g_1(\xi)}{s_1},\dots,\frac{g_J(\xi)}{s_J}\Bigr\},
\]
with declared positive scales \(s_\cdot\) (Section III.B). Then the joint feasibility event is \(\mathcal{F}=\{s\le 0\}\).
Given a predictor \(\widehat s(\xi)\) (computed from \(\widehat z\)), define nonconformity scores on calibration data \(\{\xi^i\}_{i=1}^{n_{cal}}\) by
\[
r_i := s(\xi^i) – \widehat s(\xi^i).
\]
Let \(\widehat q\) be the \((1-\alpha)\)-quantile of \(\{r_i\}\) using the standard split-conformal rank rule (with the \(\lceil(1-\alpha)(n_{cal}+1)\rceil\)-th order statistic).
Theorem VII.9 (Split conformal upper tube for the joint-violation score).
If the calibration pair \((\xi^i,s(\xi^i))\) and a future pair \((\xi^{new},s(\xi^{new}))\) are exchangeable, then
\[
\mathbb{P}\big(s(\xi^{new}) \le \widehat s(\xi^{new}) + \widehat q\big) \ge 1-\alpha.
\]
Corollary VII.10 (Auditable lower bound on joint feasibility, conditional on exchangeability).
If for a given \((d,y)\) and a prospective deployment distribution of \(\xi\) one has
\[
\widehat s(\xi) + \widehat q \le 0,
\]
then, with probability at least \(1-\alpha\) over the randomness in \(\xi\) and the exchangeability model,
\[
s(\xi)\le 0 \quad \Rightarrow\quad Q\ge 15 \ \wedge\ g_j\le 0\ \forall j.
\]
Interpretation.
The conformal tube does not require a parametric distribution for \(\xi\); it requires that the calibration data and future data be drawn in an exchangeable manner from the same mechanism. When distribution shift is a concern, the dossier must explicitly state that exchangeability may fail and either (a) widen the tube conservatively using domain-shift heuristics, or (b) adopt a DRO certificate (VII.C) with logged assumptions.
Audit requirement VII.11 (conformal provenance).
Any conformal certificate must log: training/calibration split, the score definition (including scales), the calibration residuals (or their hash), the quantile index used, and performance checks on a held-out validation set (coverage estimated with binomial confidence intervals).
VII.F. How certificates interact with verified numerics and UNKNOWN outcomes
All probability statements above implicitly assume that constraint evaluations are meaningful. In this paper’s workflow, the Section V verified-numerics layer produces tri-valued outputs (safe / violated / unknown). The probability-of-feasibility layer must therefore specify a conservative mapping from tri-valued outcomes to the quantities used in bounds.
A baseline conservative rule is:
– For moment bounds: treat UNKNOWN as violated by setting margins to worst-case values consistent with interval enclosures.
– For scenario programs: count UNKNOWN samples as constraint violations (or exclude them and explicitly weaken the claim to a conditional-on-PASS region).
– For conformal calibration: treat UNKNOWN as large positive violation scores.
These rules ensure that probability lower bounds remain conservative at the expense of potentially large conservatism; Section XI’s V\&V prioritization is then responsible for reducing UNKNOWN frequency in the regions that dominate design decisions.
This completes the probability-of-feasibility certificate layer. The next section develops RAMI and performability certificates, including dependence and partial degradation, and then couples those plant-availability constraints back into the joint feasibility event.
VIII. RAMI and Performability Under Dependence and Partial Degradation
This section formalizes RAMI (reliability, availability, maintainability, inspectability) constraints in a way that is compatible with the contract spine (Section III) and the probability-of-feasibility layer (Section VII). The key outputs are (i) auditable bounds on unavailability/availability under explicit dependence assumptions, and (ii) performability (probability of meeting performance thresholds while degraded), avoiding the need to enumerate large continuous-time Markov chains.
Throughout, the system “down” event and degraded-performance events are defined with respect to declared requirements. For example, availability constraints may appear as
\[
U(d,y,\xi) \le U_{max}\quad\text{or}\quad A(d,y,\xi)=1-U(d,y,\xi)\ge A_{min},
\]
with all time references (mission duration, planned maintenance windows) declared and unit-checked.
VIII.A. Minimal cut sets and availability/unavailability bounding
VIII.A.1. Binary component states and the top event
Let components (or subsystems) be indexed by \(i\in\{1,\dots,n\}\). For a fixed design/architecture \((d,y)\) and uncertainty sample \(\xi\), define a binary down indicator
\[
X_i = \mathbf{1}\{\text{component } i \text{ is unavailable over the declared mission definition}\} \in \{0,1\}.
\]
Let \(E_i:=\{X_i=1\}\) be the corresponding event. The system-down (top) event \(D\) is a monotone Boolean function of \(X\), typically represented via a fault tree or reliability block diagram.
A (minimal) cut set is a subset \(C\subseteq\{1,\dots,n\}\) such that simultaneous unavailability of all components in \(C\) implies system down:
\[
\bigcap_{i\in C} E_i \subseteq D,
\]
and minimality means no proper subset has this property. Let \(\mathcal{C}=\{C_1,\dots,C_m\}\) be a set of minimal cut sets.
Then
\[
D = \bigcup_{k=1}^m \bigcap_{i\in C_k} E_i
\]
for a complete minimal-cut representation. In practice, \(\mathcal{C}\) may be truncated; any truncation must be logged because it changes bounds (see VIII.A.3).
VIII.A.2. Distribution-free bounds on unavailability from cut sets
Define cut-set events \(D_k := \bigcap_{i\in C_k} E_i\). Then
\[
U := \mathbb{P}(D) = \mathbb{P}\Big(\bigcup_{k=1}^m D_k\Big).
\]
Without any independence assumptions, the following bounds are immediate:
– Lower bound (max bound):
\[
U \ge \max_{k} \mathbb{P}(D_k).
\]
– Upper bound (Boole/union bound):
\[
U \le \sum_{k=1}^m \mathbb{P}(D_k).
\]
These bounds are auditable once each \(\mathbb{P}(D_k)\) is computed or bounded (Sections VIII.B\u2013VIII.C address dependence models for that task).
A second-order (Bonferroni) lower bound can also be stated:
\[
U \ge \sum_{k=1}^m \mathbb{P}(D_k) – \sum_{1\le k<\ell\le m} \mathbb{P}(D_k\cap D_\ell),
\]
but it is only useful when pairwise intersection probabilities can be bounded conservatively.
VIII.A.3. Truncation and “bounding direction” audit
If only a subset \(\widetilde{\mathcal{C}}\subseteq \mathcal{C}\) is used (e.g., only cut sets up to size \(r\), or only the top-K most probable cut sets), then
\[
\widetilde D := \bigcup_{C\in\widetilde{\mathcal{C}}} \bigcap_{i\in C} E_i \subseteq D.
\]
Hence \(\mathbb{P}(\widetilde D)\) is a lower bound on \(U\) (it under-approximates the down event). Conversely, using the union bound on a truncated set yields an upper bound on \(\mathbb{P}(\widetilde D)\), not necessarily on \(U\).
Audit requirement VIII.1 (RAMI cut-set completeness/truncation).
Any unavailability number must declare whether it is:
- exact for the declared fault-tree model,
- a lower bound on \(U\) due to cut-set truncation (and how truncation was chosen), or
- an upper bound on \(U\) due to Boole/Bonferroni-type bounding.
VIII.A.4. Repairable components: rates and units
When components are repairable, one often parameterizes with failure rates \(\lambda_i\) and repair rates \(\mu_i\), with units \([\lambda_i]=[\mu_i]=\mathrm{time}^{-1}\). Under a two-state exponential (CTMC) model, the steady-state unavailability of component \(i\) is
\[
q_i := \mathbb{P}(X_i=1) = \frac{\lambda_i}{\lambda_i+\mu_i}.
\]
This formula is only valid under explicit modeling assumptions (memoryless up/down times). In this paper, it is acceptable as a module-internal RAMI surrogate only when the assumption is declared and validated for the component class; otherwise \(q_i\) should be treated as an uncertain input inferred from data.
VIII.B. Common-cause-aware RAMI certificates
Independence between component down events is rarely defensible for shared infrastructures (controls, power supplies, cryogenic distribution, common vendor batches). We therefore encode dependence by augmenting the fault model with explicit common-cause “shock” events, which are treated as additional components in the cut-set logic.
VIII.B.1. Virtual shock components (beta-factor style) as an auditable dependence model
Let \(g\) index a common-cause group (e.g., all gyrotrons sharing a modulator bus, or all cryoplant compressor trains sharing a contamination vulnerability). Introduce a group-level shock indicator
\[
Y_g \in \{0,1\}, \quad C_g := \{Y_g=1\}.
\]
For a component \(i\) in group \(g(i)\), define an independent-cause indicator \(Z_i\in\{0,1\}\) and set
\[
X_i = \max\{Y_{g(i)}, Z_i\}.
\]
Then any cut-set probability can be bounded in terms of \(\mathbb{P}(C_g)\) and the independent-cause probabilities \(\mathbb{P}(Z_i=1)\), without assuming full independence of the \(X_i\).
Mapping to beta-factor parameterization (conditional statement).
If one adopts the classical beta-factor idea that a fraction \(\beta_g\) of the marginal unavailability of each component in group \(g\) is attributable to common cause, then one may relate
\[
q_i = \mathbb{P}(X_i=1),\qquad q_g^C := \mathbb{P}(Y_g=1),\qquad q_i^I := \mathbb{P}(Z_i=1)
\]
by
\[
q_i = q_g^C + (1-q_g^C)q_i^I.
\]
If additionally \(q_g^C = \beta_g q_i\) is imposed for components with identical marginals inside the group, then \(q_i^I\) is determined by \(q_i\) and \(\beta_g\). This is a modeling choice and must be logged as such.
Audit requirement VIII.2 (common-cause declaration).
Any use of a shock model must record:
- group definitions (which components are in each \(g\)),
- the meaning of \(Y_g\) (what it represents physically),
- how \(q_g^C\) and/or \(\beta_g\) were chosen or inferred,
- whether shocks are assumed independent across groups (a further dependence assumption).
VIII.B.2. Bounding without validated dependence parameters
If only marginal unavailabilities \(q_i\) are known (or bounded) but dependence is not validated, then any computed plant unavailability is necessarily conditional. A conservative approach is to provide Fr\u00e9chet-type bounds for key intersections.
For any two events \(E,F\),
\[
\max\{0,\mathbb{P}(E)+\mathbb{P}(F)-1\} \le \mathbb{P}(E\cap F) \le \min\{\mathbb{P}(E),\mathbb{P}(F)\}.
\]
Applied to cut-set events \(D_k\), these bounds can be used to bound \(\mathbb{P}(D_k)\) even when component dependencies are not modeled, albeit often very conservatively.
VIII.C. Dependence-aware unavailability via Stein–Chen Poisson approximation (rare-event regime)
This subsection provides an optional certificate mechanism for regimes where down events are rare but dependent, and a Poisson approximation can be justified with an auditable error bound.
VIII.C.1. Dependent Bernoulli sums and Poisson approximation
Let \(\{I_a\}_{a\in\mathcal{A}}\) be indicator random variables for a family of “bad” events (e.g., specific basic failures, or specific minimal cut-set occurrences after common-cause expansion). Define
\[
W := \sum_{a\in\mathcal{A}} I_a, \qquad \lambda := \mathbb{E}[W] = \sum_{a\in\mathcal{A}} p_a,\quad p_a:=\mathbb{P}(I_a=1).
\]
Then \(\mathbb{P}(W\ge 1)\) upper bounds the system-down probability when \(D\subseteq\{W\ge 1\}\), and equals it when the events \(\{I_a=1\}\) exactly represent the down mechanism.
Assume each event a has a declared dependency neighborhood \(\Gamma(a)\subseteq\mathcal{A}\) such that \(I_a\) is independent of \(\{I_b: b\notin\Gamma(a)\}\). Define
\[
b_1 := \sum_{a}\sum_{b\in\Gamma(a)} p_a p_b,\qquad
b_2 := \sum_{a}\sum_{b\in\Gamma(a)\setminus\{a\}} \mathbb{P}(I_a=1, I_b=1).
\]
A classical Stein–Chen bound (stated here as a template; constants depend on the specific theorem variant used) yields a total-variation error bound
\[
d_{TV}(\mathcal{L}(W),\mathrm{Poisson}(\lambda)) \le \min\{1,\lambda^{-1}\}(b_1+b_2).
\]
Consequently,
\[
\big|\mathbb{P}(W=0) - e^{-\lambda}\big| \le d_{TV}(\mathcal{L}(W),\mathrm{Poisson}(\lambda)),
\]
and therefore
\[
\mathbb{P}(W\ge 1) = 1-\mathbb{P}(W=0) \in \Big[1-e^{-\lambda}-\epsilon_{SC},\ 1-e^{-\lambda}+\epsilon_{SC}\Big],
\]
with \(\epsilon_{SC}\) given by the right-hand side bound.
Audit requirement VIII.3 (Stein–Chen eligibility).
If this approximation is used, the dossier must include:
- the explicit index set \(\mathcal{A}\) (what events are being counted),
- the dependency neighborhoods \(\Gamma(a)\),
- computed values (or conservative bounds) for \(p_a\) and \(\mathbb{P}(I_a I_b=1)\), and
- the resulting \(\epsilon_{SC}\) and whether it is small enough to be decision-relevant.
Gap flag.
The Stein–Chen method is only compelling when events are rare and dependencies are weak/local; otherwise \(\epsilon_{SC}\) will be too large to certify anything useful.
VIII.D. Multistate RAMI and performability via universal generating functions (UGF)
Binary availability (up/down) can be inadequate for fusion plants where partial degradation materially changes physics feasibility (e.g., losing 1 of \(N\) heating sources may still allow operation but with reduced Q margin). We therefore model performability by discrete capacity states.
VIII.D.1. Multistate subsystem model
For a subsystem \(i\) (heating bank, cryogenic trains, pumping, tritium processing), define a nonnegative random capacity
\[
S_i \in \{c_{i,0},c_{i,1},\dots,c_{i,K_i}\}, \qquad 0\le c_{i,0}<\cdots
\[
0\le S_i – \widetilde S_i < \Delta \quad \text{almost surely}.
\]
If \(H\) is monotone nondecreasing in each argument and 1-Lipschitz with respect to the \(\ell_1\) norm (true for sum and min), then
\[
0\le S_{sys} - \widetilde S_{sys} \le m\Delta,
\]
where \(\widetilde S_{sys}:=H(\widetilde S_1,\dots,\widetilde S_m)\).
Therefore, for any threshold \(t\),
\[
\mathbb{P}(\widetilde S_{sys}\ge t) \le \mathbb{P}(S_{sys}\ge t) \le \mathbb{P}(\widetilde S_{sys}\ge t-m\Delta).
\]
The left inequality shows that using rounded-down capacities yields a conservative lower bound on performability.
Audit requirement VIII.4 (UGF quantization).
Any UGF-based performability statement must log \(\Delta\), whether rounding was down (conservative) or not, and the implied threshold slack \(m\Delta\).
VIII.E. Coupling performability to physics feasibility thresholds
Performability must connect to physics feasibility (Q and engineering margins). A key pattern is: compute a required auxiliary capacity threshold from a (monotone) physics surrogate, then propagate the probability that the plant can supply at least that capacity.
VIII.E.1. Auxiliary heating threshold from a monotone Q surrogate
Fix \((d,y,\xi)\) and consider the dependence of fusion gain on auxiliary heating set-point \(P_{aux}\ge 0\):
\[
Q = Q(d,y,\xi; P_{aux}).
\]
Assumption VIII.5 (monotonicity; gap flag).
For the operating regime of interest, \(Q\) is nondecreasing in \(P_{aux}\). This is a modeling assumption and must be validated for the surrogate used.
Define the required auxiliary power as the minimal value achieving the Q target:
\[
P_{aux}^{min}(d,y,\xi) := \inf\{P\ge 0 : Q(d,y,\xi;P)\ge 15\}.
\]
If monotonicity holds and \(Q\) is continuous in \(P\), then \(P_{aux}^{min}\) can be approximated by bisection.
Incorporating surrogate error envelopes.
Suppose a surrogate \(\widehat Q(P)\) satisfies \(|Q(P)-\widehat Q(P)|\le \Delta_Q\) uniformly in \(P\in[P_{lo},P_{hi}]\). Then a conservative sufficient condition for meeting Q>15 is
\[
\widehat Q(P) – \Delta_Q \ge 15.
\]
Hence one can define a conservative threshold
\[
\widehat P_{aux}^{min} := \inf\{P\in[P_{lo},P_{hi}] : \widehat Q(P)-\Delta_Q \ge 15\},
\]
which guarantees that any delivered auxiliary power \(P_{aux}\ge \widehat P_{aux}^{min}\) implies \(Q\ge 15\) within the declared envelope.
Audit requirement VIII.6 (threshold inversion dossier).
Any derived \(P_{aux}^{min}\) (or conservative \(\widehat P_{aux}^{min}\)) must log: the monotonicity assumption, the inversion method and bracket \([P_{lo},P_{hi}]\), and any surrogate/verified-numerics envelopes used.
VIII.E.2. Performability constraints from capacity states
Let \(S_{heat}\) denote the delivered heating capacity state (random due to failures/degradation). A physics-conditioned performability constraint can be written as
\[
\mathbb{P}(S_{heat} \ge P_{aux}^{min}(d,y,\xi)) \ge 1-\alpha_{heat},
\]
or in a conservative two-stage form:
1. Compute \(\widehat P_{aux}^{min}(d,y,\xi)\) using a surrogate with envelopes.
2. Require \(\mathbb{P}(S_{heat} \ge \widehat P_{aux}^{min})\ge 1-\alpha_{heat}\).
Analogous constructions apply to other support systems by defining thresholds for required cryogenic capacity, pumping capacity, tritium processing throughput, or thermal margins.
VIII.E.3. Coupling to joint feasibility (Sections VII and X)
The joint feasibility event \(\mathcal{F}\) (Section VII) can be extended to include performability by adding constraints of the form
\[
g_{perf}(d,y,\xi) := \widehat P_{aux}^{min}(d,y,\xi) – S_{heat}(y,\xi) \le 0,
\]
with \(S_{heat}\) treated as an additional uncertain input derived from the RAMI model. This makes “physics feasibility under degradation” part of the same certificate pipeline as Q and engineering margins.
Gap flag (closing the loop).
When thresholds depend on \(\xi\) and when \(S_{heat}\) depends on \(\xi\) (e.g., repair-time uncertainty), the dependence structure must be explicit. If not validated, the result should be stated conditionally (e.g., conditional independence given group shocks).
VIII.F. Summary of RAMI/performability outputs
A RAMI/performability evaluation compatible with this paper’s audit requirements should emit:
– a cut-set list (or hash) and whether it is complete or truncated;
– declared component/group unavailability models with unit-checked rate parameters (if used);
– explicit dependence assumptions (common-cause groups and shock variables), and conservative bounds when parameters are not validated;
– for performability: UGF definitions, quantization step \(\Delta\), and threshold probabilities with conservative rounding choices;
– a mapping showing how availability/performability constraints enter the plant-level joint feasibility event \(\mathcal{F}(d,y)\).
The next sections build on this RAMI layer: Section IX discusses component-level physics and risk modules to be unified by contracts; Section X integrates the joint feasibility and RAMI/performability constraints into cost- and schedule-constrained continuous–discrete optimization.
IX. Component-Level Physics/Risk Modules Claimed in the Database (to be Unified by Contracts)
This section describes examples of component-level physics and risk modules that may be coupled into the plant digital twin. The role of this section is not to assert that any particular stellarator-specific scaling law, surrogate, or tail model is correct. Rather, it standardizes how such modules must be presented so that they (i) are compatible with the contract spine (Sections II–III), and (ii) produce outputs with explicit uncertainty/error envelopes that can be consumed by verified numerics (Section V), UQ (Section VI), feasibility certificates (Section VII), RAMI/performability (Section VIII), and optimization (Section X).
IX.A. General contract template for a component-level module
A component-level module \(\mathcal{M}\) is accepted into the integrated twin only if it declares, in addition to ordinary input/output ports, an **error and validity contract**:
1. Typed inputs/outputs.
A module consumes a subset of global ports and design variables \((x,d,y)\), and a subset of uncertainties \(\xi_{\mathcal{M}}\), and produces outputs \(o_{\mathcal{M}}\) used elsewhere:
\[
o_{\mathcal{M}} = \mathcal{M}(x,d,y,\xi_{\mathcal{M}}).
\]
Every component of \(o_{\mathcal{M}}\) must carry declared units/dimension vectors and sign conventions (Section III.B).
2. Validity domain.
A declared domain \(\Omega_{\mathcal{M}}\subset \mathcal{D}\times\mathcal{Y}\times\Xi_{\mathcal{M}}\) in which the module is intended to be used. Outside \(\Omega_{\mathcal{M}}\), the module must return UNKNOWN (or refuse evaluation) rather than silently extrapolate.
3. Error envelope.
An explicit bound (deterministic or probabilistic) on the discrepancy between the module output and a reference definition.
– Deterministic envelope: \(\|o_{\mathcal{M}}-\widehat o_{\mathcal{M}}\|\le \Delta_{\mathcal{M}}\) for all \((d,y,\xi)\in\Omega_{\mathcal{M}}\).
– Probabilistic envelope: \(\mathbb{P}(\|o_{\mathcal{M}}-\widehat o_{\mathcal{M}}\|\le \Delta_{\mathcal{M}})\ge 1-\alpha_{\mathcal{M}}\).
The dossier must state which interpretation is used, because it changes how Section VII aggregates risks.
4. Monotonicity/Lipschitz declarations (optional but powerful).
If the module claims a monotonicity property (used to derive conservative thresholds as in VIII.E) or a Lipschitz constant (used in small-gain or DRO layers), it must declare:
\[
\|o_{\mathcal{M}}(\theta)-o_{\mathcal{M}}(\theta’)\| \le L_{\mathcal{M}}\,\|\theta-\theta’\| \quad \text{on } \Omega_{\mathcal{M}},
\]
with the chosen norm and units.
These requirements are the mechanism by which “physics modules” become auditable: their claims are reduced to checkable interface objects (units, envelopes, and optionally monotonicity/Lipschitz indices) rather than informal narrative.
IX.B. Coil manufacturing tolerances \(\to\) field-spectrum perturbations \(\to\) confinement degradation
IX.B.1. Input/output ports
A tolerance module consumes manufacturing descriptors (e.g., geometric deviations, alignment errors) and produces a representation of the induced magnetic-field perturbation. One generic representation is via a truncated Fourier/harmonic coefficient vector \(b\in\mathbb{R}^m\) (e.g., coefficients of a field-strength perturbation in a chosen angular basis). The choice of basis is a modeling decision; what matters for integration is:
– the coefficient vector \(b\) is dimensionless or explicitly unit-tagged,
– truncation order \(m\) is recorded,
– and a truncation remainder bound \(\|b-b^{(m)}\|\le \Delta_{trunc}\) is provided (gap flag unless independently validated).
IX.B.2. Propagating to a performance-relevant scalar via a Lipschitz contract
Let \(\varepsilon\) denote a scalar “ripple/perturbation metric” derived from \(b\), for instance \(\varepsilon := \|W b\|\) for some declared weight matrix \(W\). Let \(Q\) be the plant fusion gain. A defensible contract for coupling is not “\(Q\) degrades by X%”, but rather:
Assumption IX.1 (local Lipschitz degradation; gap flag).
There exists \(L_Q\ge 0\) such that for all \((d,y,\xi)\) in a declared regime,
\[
|Q(d,y,\xi;\varepsilon)-Q(d,y,\xi;\varepsilon’)| \le L_Q\,|\varepsilon-\varepsilon’|.
\]
If \(\varepsilon\) is uncertain because of tolerance uncertainty and truncation error, we can form a conservative Q bound:
\[
Q(d,y,\xi;\varepsilon) \ge Q(d,y,\xi;\widehat\varepsilon) – L_Q\,|\varepsilon-\widehat\varepsilon|.
\]
This inequality becomes actionable only when \(L_Q\) and \(|\varepsilon-\widehat\varepsilon|\) have auditable bounds. Estimating them is a V\&V task (Section XI).
IX.C. High-beta stability certification as a margin-under-uncertainty problem
IX.C.1. Stability margin definition and contract form
Let \(\beta\) denote a dimensionless plasma pressure parameter and let \(\mathcal{S}(d,y,\xi)\) denote a scalar stability margin returned by an equilibrium/MHD stability module, with convention:
\[
\mathcal{S}(d,y,\xi) \ge 0 \;\text{stable/acceptable}, \qquad \mathcal{S}(d,y,\xi) < 0 \;\text{unstable/unacceptable}.
\]
The plant constraint is \(g_{\beta}(d,y,\xi):=-\mathcal{S}(d,y,\xi)\le 0\).
IX.C.2. Perturbation-to-margin bounding via sensitivity envelopes
Many stability surrogates are ultimately functions of an underlying eigenvalue or generalized eigenvalue problem. Without committing to a particular model class, we can standardize the integration requirement as:
Audit requirement IX.2 (stability-margin envelope).
The module must provide an interval enclosure or deterministic error envelope
\[
\mathcal{S}(d,y,\xi) \in [\\underline{\mathcal{S}},\overline{\mathcal{S}}]
\]
or \(|\mathcal{S}-\widehat{\mathcal{S}}|\le \Delta_{\mathcal{S}}\), with provenance (solver tolerances alone are not an acceptable provenance).
If an interval enclosure is available, then the stability constraint is certified safe whenever \(\\underline{\mathcal{S}}\ge 0\). This directly matches the Section V logic for constraint certification.
Gap flag.
If the stability module cannot provide a posteriori error bounds (common for complex physics codes), then any stability constraint should be treated as scenario/conformal-certified rather than “verified” (Section VII), unless an independent verification route is established.
IX.D. Divertor/edge heat-exhaust safety under epistemic uncertainty
IX.D.1. Contracted interface variables and conservation closure
Let \(q_{\perp}^{peak}(d,y,\xi)\) be peak divertor heat flux (units \([\mathrm{W}/\mathrm{m}^2]\)). The associated plant constraint is \(g_q:=q_{\perp}^{peak}-q_{\perp}^{lim}\le 0\). A heat-exhaust module must couple consistently to the Section III energy-family graph: the net power crossing the edge/divertor boundary in the port graph must match the integral power loads used by the heat-flux computation, up to the declared conservation residual tolerance.
Audit requirement IX.3 (edge/divertor energy consistency).
The module dossier must identify which port-graph edge flows contribute to the divertor power balance and provide a normalized residual check of the form \(\|\hat r_E\|_\infty\le \tau_E\) (Section III.C) for the mapping from edge power to divertor loads.
IX.D.2. Monotone-enclosure pattern for safety
Often a dominant uncertainty is an edge-transport multiplier \(\chi\) (dimensionless) influencing \(q_{\perp}^{peak}\). If the module can justify monotonicity in \(\chi\) on a regime (gap flag), i.e.
\[
\chi_1\le \chi_2 \implies q_{\perp}^{peak}(\chi_1)\le q_{\perp}^{peak}(\chi_2),
\]
then safety under \(\chi\in[\\underline\chi,\overline\chi]\) reduces to a single evaluation at \(\overline\chi\) (plus envelopes):
\[
q_{\perp}^{peak}(\chi) \le q_{\perp}^{peak}(\overline\chi) + \Delta_q.
\]
When monotonicity cannot be justified, interval evaluation or sample-based certification should be used. In particular, if the heat-flux calculation is embedded in a nonlinear steady-state solve, the Section V verified wrapper can be applied to certify bounds on \(q_{\perp}^{peak}\) for a given \(\xi\).
IX.E. Blanket neutronics–thermal–tritium surrogate modeling and dynamic tritium-cycle safety
IX.E.1. Reduced-basis surrogate with a posteriori bounds (static quantities)
Let \(TBR(d,y,\xi)\) be tritium breeding ratio (dimensionless) and \(\sigma_{th}(d,y,\xi)\) a thermal stress metric (units declared). If a reduced-basis (RB) surrogate is used, the key integration feature is the existence of a rigorous a posteriori error bound:
\[
|TBR-\widehat{TBR}| \le \Delta_{TBR}^{RB}, \qquad |\sigma_{th}-\widehat\sigma_{th}|\le \Delta_{\sigma}^{RB}.
\]
These yield conservative constraints
\[
\widehat{TBR} - \Delta_{TBR}^{RB} \ge TBR_{min}, \qquad \widehat\sigma_{th}+\Delta_{\sigma}^{RB} \le \sigma_{lim}.
\]
Such bounds are “certificate-compatible” because they are checkable without rerunning the high-fidelity solver, provided the RB theory’s assumptions (coercivity, residual evaluation accuracy) are satisfied.
IX.E.2. Tritium inventory transient safety as a contract-compatible inequality
Let \(I_T(t)\) be tritium inventory (units \([\mathrm{g}]\) or \([\mathrm{mol}]\), declared). A minimal abstract model for coupling is a balance law
\[
\dot I_T(t) = b(t) - e(t),
\]
where \(b(t)\) is breeding/production inflow and \(e(t)\) is extraction/outflow, both derived from other modules and connected to the tritium-mass conserved family in the port graph (Section III).
To make a safety claim \(I_T(t)\le I_T^{lim}\) auditable under uncertainty and delays, one can require the tritium-cycle module to provide a **comparison-system bound**: upper envelopes \(\overline b(t)\), lower envelopes \(\\underline e(t)\) (including delay effects) such that \(b\le \overline b\), \(e\ge \\underline e\), implying
\[
I_T(t) \le I_T(0) + \int_0^t (\overline b(s)-\\underline e(s))\,ds.
\]
This is conservative but contract-compatible: it turns a dynamic safety claim into an explicit inequality whose inputs are unit-checked flows and whose conservatism is auditable via the choice of envelopes.
Gap flag.
Any more detailed stochastic-delay or queueing model of extraction must be validated for the actual plant architecture; absent that, the envelope approach should be treated as the default certificate layer.
IX.F. Neutron wall loading \(\to\) damage accumulation \(\to\) lifetime/schedule coupling
IX.F.1. Unit-checked mapping from fusion power to wall loading
Let \(P_{fus}(t)\) be fusion power (units \([\mathrm{W}]\)). A wall-loading module defines a wall heat/particle/neutron load metric \(w(t)\) (units declared, e.g. \([\mathrm{W}/\mathrm{m}^2]\) for power flux) via a contractually declared mapping
\[
w(t) = \mathcal{W}(P_{fus}(t),d,y,\xi_{wall}),
\]
with explicit geometric and conversion factors and unit checks. This mapping must be consistent with the Section III energy-family accounting (the same \(P_{fus}\) cannot be injected twice into different subsystems without an explicit split).
IX.F.2. Damage accumulation as an integral functional with conservative enclosures
Let damage be an accumulated functional
\[
D(T) := \int_0^T \varphi(w(t),\xi)\,dt,
\]
where \(\varphi\ge 0\) converts load to damage rate. If \(w(t)\) and \(\varphi\) are available only through surrogates with envelopes, one can compute an interval enclosure for \(D(T)\) by bounding \(\varphi\) above and integrating:
\[
D(T) \le \int_0^T \overline\varphi(t)\,dt.
\]
This yields an auditable sufficient condition for staying below a damage limit \(D_{lim}\): if \(\int_0^T \overline\varphi(t)dt \le D_{lim}\), then \(D(T)\le D_{lim}\).
IX.F.3. Lifetime as a probabilistic statement (explicitly conditional)
If a probabilistic lifetime model is introduced (e.g., a parametric survival distribution for the random threshold \(D_{fail}\)), then any tail probability claim \(\mathbb{P}(D(T)\le D_{fail})\ge 1-\alpha\) must be explicitly labeled as conditional on that distributional model and on the calibration dataset used. When such a distributional model is not defensible, Section VII’s distribution-free tools (scenario/conformal) should be used directly on observed or simulated failure indicators.
IX.G. HTS coil margin and quench reliability under uncertainty
IX.G.1. Margin definition
Let \(I_{op}(d,y,\xi)\) be operating current and \(I_c(d,y,\xi)\) an effective critical current (both in \([\mathrm{A}]\)). Define a current margin
\[
M_I(d,y,\xi) := I_c(d,y,\xi) - I_{op}(d,y,\xi).
\]
The safe constraint is \(M_I\ge 0\) (equivalently \(g_I:= -M_I\le 0\)).
IX.G.2. Concentration-bound pattern for aggregate defect counts (conditional)
Many reliability arguments for complex coils reduce to bounding the probability that an aggregate defect count exceeds a threshold. If one can model a defect count \(N\) as a sum of bounded or Bernoulli-like indicators, then concentration inequalities (Hoeffding/Chernoff variants) can yield auditable upper bounds on \(\mathbb{P}(N\ge n_{crit})\) from mean-rate assumptions.
Gap flag.
These bounds are only as credible as the underlying independence/boundedness assumptions. Therefore, any such quench reliability certificate must:
- state the probabilistic model (including dependence structure),
- supply parameter provenance (inspection statistics, vendor data, test coupons),
- and, if dependence cannot be justified, revert to conservative dependence-aware RAMI treatments (Section VIII.B) or finite-sample empirical guarantees (Section VII).
IX.H. How these modules connect back to plant-level certificates
Each component-level module in this section is intended to produce one or more of the following **contract-consumable objects**:
1. A scalar/vector output that appears directly in plant constraints \(g_j\) with unit-checked definitions.
2. A deterministic envelope \(\Delta_j\) suitable for conservative constraint tightening (Section VI.A.3).
3. A verified interval enclosure \([g_j](X)\) when wrapped inside the Section V verified solve.
4. A monotonicity or Lipschitz declaration that enables conservative threshold inversion (VIII.E) or distributionally robust bounds (VII.C).
The next section (Section X) shows how these contract-consumable outputs are assembled into a unified cost- and schedule-constrained continuous–discrete co-design problem, and how optimization algorithms can return not only candidate designs but also auditable feasibility/optimality certificates and dual-based explanations.
X. Unified Cost-Constrained Optimization: Continuous–Discrete Co-Design with Certificates
This section assembles the preceding layers into optimization problems whose *outputs are not only candidate designs* \((d,y)\) but also *auditable certificates* (or explicitly conditional claims) supporting (i) feasibility under uncertainty, (ii) cost/schedule compliance, and (iii) explainable redesign actions when infeasible.
We emphasize that the integrated digital twin is generally nonconvex and mixed-integer. Therefore, any claim of global optimality or probabilistic feasibility must be tied to a specific *restricted optimization class* (e.g., GP/MICP, scenario convex programs, or conservative relaxations) and must explicitly log the assumptions under which the certificate is valid.
X.A. Plant-level optimization statement
X.A.1. Joint feasibility event and objective components
Fix a continuous design vector \(d\in\mathcal{D}\) and discrete architecture vector \(y\in\mathcal{Y}\). Define the joint feasibility event (Section VII)
\[
\mathcal{F}(d,y) := \{\xi:\ Q(d,y,\xi)\ge 15\ \wedge\ g_j(d,y,\xi)\le 0\ \forall j=1,\dots,J\ \wedge\ C(d,y,\xi)\le C_{max}\ \wedge\ S(d,y,\xi)\le S_{max}\},
\]
where \(g_j\) may include engineering margins and RAMI/performability constraints encoded as inequality constraints (Section VIII.E).
We consider multi-objective plant performance of the schematic form
\[
\max_{d\in\mathcal{D},\ y\in\mathcal{Y}} \; \Big(\mathbb{E}[P_{net}(d,y,\xi)],\; \mathbb{E}[A(d,y,\xi)],\; -\mathbb{E}[C(d,y,\xi)]\Big)
\]
subject to risk and hard constraints. For algorithmic development, we typically scalarize via weights \(w\ge 0\):
\[
\max_{d,y} \; J(d,y) := w_P\,\mathbb{E}[P_{net}] + w_A\,\mathbb{E}[A] - w_C\,\mathbb{E}[C] - w_S\,\mathbb{E}[\max\{0,S-S_{max}\}],
\]
with all expectations and penalties interpreted under an explicitly declared uncertainty model (probabilistic, DRO ambiguity set, or finite-sample statement).
X.A.2. Certificate-compatible constraint formats
We target constraints in one (or more) of the following certificate-compatible forms.
1. Deterministic robust constraints (SOS / robust GP / verified numerics):
\[
g_j(d,y,\xi) \le 0 \quad \forall \xi\in\Xi_0.
\]
2. Chance constraints with explicit certificate type:
\[
\mathbb{P}(g_j(d,y,\xi)\le 0) \ge 1-\alpha_j \quad \text{(with stated assumptions and estimator)}.
\]
3. Joint chance constraints, either by allocation \(\sum_j \alpha_j\le \alpha\) or by a joint score (Section VII.E):
\[
\mathbb{P}(\mathcal{F}(d,y)) \ge 1-\alpha.
\]
4. Hard budget constraints (enforced pointwise in \(\xi\) or robustly depending on interpretation):
\[
C(d,y,\xi)\le C_{max},\qquad S(d,y,\xi)\le S_{max}.
\]
Audit requirement X.1 (optimization-to-certificate linkage).
Any optimization run must log, for each constraint class, which certificate mechanism is used (verified/robust, scenario, DRO, conformal, moment-based), and must link to the corresponding evidence artifact (Sections V–VII) rather than only reporting a single numeric “probability of success”.
X.B. Convexifiable chance-constrained geometric programming (GP) approach
X.B.1. GP-compatible surrogate structure and log-domain convexity
A geometric program (GP) in standard form minimizes a posynomial objective subject to posynomial inequality constraints and monomial equalities. In this paper’s context, a GP-like subproblem is appropriate only after choosing surrogate parameterizations for physics/cost/RAMI relationships that are posynomial (or can be conservatively upper-bounded by posynomials).
Assumption X.2 (GP representability; gap flag).
For a selected subset of constraints and objectives, there exist positive decision variables \(x\in\mathbb{R}_{++}^n\) (a reparameterization of parts of \(d\)) such that each constraint can be written as
\[
p_k(x) \le 1,\qquad k=1,\dots,K,
\]
where \(p_k\) are posynomials, possibly after conservative upper-bounding and after embedding deterministic error envelopes \(\Delta\) (Section VI.A.3).
Under Assumption X.2, the standard log change of variables \(u=\log x\) transforms the GP into a convex program in \(u\). Any such transformation must be logged because it changes the interpretation of Lipschitz constants and tolerance budgets.
X.B.2. Conservative incorporation of uncertainty: three auditable patterns
We list three patterns that preserve auditability. Each can be used to produce an optimizer-ready constraint without requiring unverifiable Monte Carlo claims.
Pattern 1 (Robustification by deterministic envelopes).
If a surrogate constraint \(\widehat g_j\) is used with a deterministic envelope \(\Delta_j\) such that \(g_j\le \widehat g_j+\Delta_j\), then the conservative constraint
\[
\widehat g_j(d,y,\xi) + \Delta_j \le 0
\]
implies the true constraint. If \(\widehat g_j\) is posynomial in GP variables, this yields a GP-compatible conservative constraint.
Pattern 2 (Moment-based chance constraints as deterministic sufficient conditions).
Let \(M_j(\xi):=-g_j(d,y,\xi)\) denote a margin random variable. If \(\mu_j:=\mathbb{E}[M_j]\) and \(\sigma_j^2:=\mathrm{Var}(M_j)\) are available (from PCE moments, Section VI.A, or samples), then Proposition VII.1 provides
\[
\mathbb{P}(g_j\le 0)=\mathbb{P}(M_j\ge 0) \ge \frac{\mu_j^2}{\mu_j^2+\sigma_j^2}.
\]
Thus, the sufficient condition
\[
\frac{\mu_j^2}{\mu_j^2+\sigma_j^2} \ge 1-\alpha_j
\]
is a deterministic inequality in \((\mu_j,\sigma_j)\). This is not automatically GP-compatible; however, it can be used as a screening constraint or embedded via conservative convex inequalities once \(\mu_j\) and \(\sigma_j\) are themselves expressed by convex surrogates.
Pattern 3 (Scenario constraints on a convex GP subproblem).
If the GP subproblem is convex in log-variables and sampled constraints remain convex for each scenario \(\xi^k\), then the scenario program (Section VII.B) can be solved and equipped with a finite-sample bound conditional on the convexity/support-constraint assumptions.
Audit requirement X.3 (no silent mixing of probability semantics).
If moment bounds, scenario constraints, DRO, or conformal tubes are used inside an optimizer, the design dossier must record whether the resulting design is certified under (i) the nominal distribution, (ii) an ambiguity set, or (iii) exchangeability-only assumptions, and how surrogate envelopes \(\Delta\) and verified-numerics UNKNOWN outcomes were handled.
X.B.3. Output artifacts from GP-based screening/optimization
A GP-based optimization layer is most defensible as an *inner-loop screening* tool that outputs:
- a candidate \((d,y)\) and operating point \(x\),
- the explicit set of GP constraints used (posynomials, scalings, and envelopes),
- a conservative feasibility classification (PASS/FAIL/UNKNOWN) at representative samples using the verified wrapper (Section V), and
- if chance constraints are claimed: the exact certificate form (Cantelli allocation, scenario bound parameters, or conformal coverage statement).
X.C. Mixed-integer convex programming (MICP) for discrete architecture decisions
X.C.1. MICP template and auditability expectations
Discrete architecture choices (redundancy counts, number of parallel trains, spare inventories) motivate a mixed-integer model. In this section we restrict attention to a *convex* continuous relaxation so that branch-and-bound can produce a global optimality certificate in the standard optimization sense.
Definition X.4 (MICP subproblem).
A design problem is an MICP if it can be written as
\[
\begin{aligned}
\min_{x\in\mathbb{R}^n,\ y\in\mathbb{Z}^m}\quad & f(x,y)\\
\text{s.t.}\quad & (x,y)\in\mathcal{K},
\end{aligned}
\]
where \(\mathcal{K}\) is a closed convex set described by linear, second-order cone, exponential-cone, or semidefinite constraints (or intersections thereof) for each fixed integer assignment y.
In our context, \(x\) can include continuous design and margin variables, while \(y\) encodes discrete counts. Cost and schedule are naturally handled as linear or convex constraints once suitable surrogates are declared.
Audit requirement X.5 (global optimality evidence).
If an MICP solver is used, the dossier must log:
- the best feasible objective value (upper bound) \(UB\),
- the best bound on the global optimum (lower bound) \(LB\) reported by the solver,
- the final optimality gap \((UB-LB)/(1+|UB|)\),
- and the exact solver tolerances and termination criteria.
This produces a mathematically standard optimality certificate for the *MICP model as written* (not for the underlying physical plant), which must be clearly distinguished.
X.C.2. Linking MICP feasibility to plant feasibility certificates
MICP feasibility is feasibility of a surrogate/relaxed optimization model. To prevent “certificate laundering”, any MICP-feasible design \((x,y)\) must be post-checked by:
1. unit and conservation contracts (Section III),
2. coupling stability/well-posedness checks where applicable (Section IV),
3. verified numerics at representative uncertainty points (Section V),
4. and a chosen probability-of-feasibility certificate (Section VII).
This produces a two-level statement:
- Optimization certificate: global optimality (or bounded gap) for the convexified mixed-integer surrogate.
- Feasibility certificate: probabilistic/robust statements for the underlying digital twin, conditional on the declared UQ and verification layers.
X.D. Generalized Benders Decomposition (GBD) for architecture + margins
X.D.1. Decomposition structure
GBD is appropriate when, for each fixed integer architecture \(y\), the continuous subproblem in \(x\) is convex (or can be made convex by conservative approximation), but the full \((x,y)\) problem is hard.
We consider the template
\[
\begin{aligned}
\min_{y\in\mathcal{Y}}\quad & \Phi(y) := \min_{x\in\mathcal{X}(y)} f(x,y)\\
\text{s.t.}\quad & \mathcal{X}(y) := \{x: \; h(x,y)\le 0\},
\end{aligned}
\]
where, for each y, the set \(\mathcal{X}(y)\) is convex and the constraint mapping \(h\) is convex in x.
In this paper’s context, x can include explicit *margins* (Section II.C.2) and slack variables used to absorb conservative envelopes. The aim is to make feasibility of \(\mathcal{X}(y)\) correspond to an auditable “architecture + margin” realizability statement.
X.D.2. Feasibility and optimality cuts from dual witnesses
Assume the y-fixed subproblem admits strong duality under standard constraint qualifications (gap flag: Slater-like conditions must be checked for the chosen convexification). Then dual variables provide two valuable artifacts:
1. Feasibility cuts: when \(\mathcal{X}(y)\) is empty, a dual ray or Farkas-type witness (analogous to Section III.D) yields a cut excluding the current architecture y and (often) a region of similar architectures. This produces an auditable explanation of *which contract blocks are incompatible* under that architecture.
2. Optimality cuts: when \(\mathcal{X}(y)\) is nonempty, dual multipliers produce a lower bound on \(\Phi(y)\), yielding a cut in the master problem that drives architecture search.
Audit requirement X.6 (GBD cut provenance).
Every Benders cut added to the master problem must be logged with:
- the associated y at which it was generated,
- primal/dual solver status for the subproblem,
- the dual variables (or dual ray) used,
- and a mapping back to contract rows/blocks (engineering interpretability).
X.D.3. “What to change” explanations as part of the certificate
A core benefit of dual-derived feasibility cuts is explainability. For example, if the continuous subproblem contains constraints of the form
\[
a_k(y)^\top x \le b_k(y),
\]
then a dual infeasibility certificate \(\lambda\ge 0\) implies that the weighted sum \(\sum_k \lambda_k a_k(y)\) cannot be satisfied with the corresponding \(\sum_k \lambda_k b_k(y)\). Reporting the dominant indices \(k\) with largest \(\lambda_k\) yields an auditable shortlist of “which constraints conflict”, and, via sensitivity of \(a_k(y),b_k(y)\) to discrete decisions, suggests “which architecture variables to increase/decrease” (e.g., add redundancy, add spares, increase capacity margin) to restore feasibility.
X.E. Shadow-price / dual-multiplier interpretations for engineering action ranking
X.E.1. Shadow prices for margin and cost trade-offs
When a convexified continuous subproblem is feasible and dual variables \(\pi\) exist for key constraints (cost, schedule, peak heat flux, coil margin, RAMI bounds), the values of \(\pi\) provide local sensitivity information:
- A large multiplier on a constraint indicates that tightening that constraint increases objective value rapidly (i.e., the constraint is “active and costly”).
- In a budget allocation setting (risk allocation \(\alpha_j\), tolerance budgets \(\tau_q\), margin allocations), dual variables can be used to rank which constraints are worth reducing uncertainty in (distinct from variance-based BOED ranking).
Audit requirement X.7 (shadow-price interpretability limits).
Dual multipliers are interpretable only for the convexified model used in the solve and only locally. Therefore, any “engineering action ranking” derived from multipliers must be explicitly labeled as:
- model-conditional (depends on the convexification/surrogate),
- local (near the optimum found), and
- subject to verification against the coupled twin via Section V (and, where relevant, via stability assumptions in Section IV).
X.E.2. Distinguishing engineering-action ranking from BOED
We distinguish two ranking outputs:
1. Engineering-action ranking (this section): derived from dual multipliers and discrete-choice trade-offs; yields “what to buy/build/change” actions (spares, redundancy, capacities, margins) under cost/schedule constraints.
2. V\&V/experiment ranking (Section XI): derived from reducing epistemic uncertainty in parameters that dominate feasibility probability or constraint tightness.
Both must appear in an auditable go/no-go dossier, but they answer different decision questions.
X.F. Summary: optimization outputs as certificate-bearing dossiers
For each candidate design \((d,y)\) produced by optimization, the expected artifact bundle includes:
- the optimization model specification (GP/MICP/GBD; convexifications and envelopes),
- solver evidence (optimality gap for MICP; cut logs for GBD),
- post-checks on the coupled twin: unit/conservation contracts (III), stability/well-posedness checks (IV), verified steady-state enclosure where feasible (V),
- and a probability-of-feasibility certificate of an explicitly declared type (VII), including treatment of UNKNOWN outcomes.
The next section (Section XI) uses these artifacts to plan 12–18 month V\&V programs: selecting falsifiers and experiments that most reduce uncertified assumptions and shrink conservative envelopes in the constraints that dominate feasibility and cost.
XI. V\&V Prioritization and Go/No-Go Decision Matrix (12–18 Month Test Planning)
This section turns the conditional assumptions and UNKNOWN outcomes produced by Sections III–X into an explicit verification-and-validation (V\&V) plan. The goal is not generic model improvement; it is **decision-directed** uncertainty reduction: shrinking exactly those envelopes and conditional assumptions that dominate (i) the certified lower bound on \(\mathbb{P}(\mathcal{F}(d,y))\) and (ii) cost/schedule feasibility.
XI.A. Bayesian optimal experimental design (BOED) via mutual information
XI.A.1. Parameterization of epistemic uncertainty and experiments
Let \(\Theta\in\mathbb{R}^{n_\Theta}\) denote epistemic parameters appearing in one or more modules (transport multipliers, confinement efficiency scalars, blanket/TBR multipliers, manufacturing tolerance statistics, failure/repair-rate parameters, etc.). Let \(D\) denote the current V\&V dataset (experiments, mockups, and/or high-fidelity simulations).
An experiment design choice is denoted \(e\in\mathcal{E}\). Given \(e\) and \(\Theta\), the experiment produces an outcome \(Y\) with likelihood
\[
p(y\mid \Theta,e)\,.
\]
The posterior after running experiment \(e\) and observing \(Y=y\) is
\[
p(\Theta\mid D, y, e) \propto p(y\mid \Theta,e)\,p(\Theta\mid D).
\]
XI.A.2. Mutual-information objective and “top falsifiers”
A standard BOED utility is the mutual information between \(\Theta\) and the prospective outcome \(Y\) conditional on current data \(D\):
\[
U(e) := I(\Theta;Y\mid D,e)
= \mathbb{E}_{\Theta\sim p(\cdot\mid D)}\,\mathbb{E}_{Y\sim p(\cdot\mid \Theta,e)}\Big[\log\frac{p(Y\mid\Theta,e)}{p(Y\mid D,e)}\Big],
\]
where \(p(y\mid D,e)=\int p(y\mid\theta,e)p(\theta\mid D)\,d\theta\).
An experiment \(e\) is called a *top falsifier* for a certificate claim if it is expected to maximally reduce (in a well-defined sense) uncertainty along failure modes that currently drive a conservative feasibility certificate. Concretely, define a (dimensionless) joint violation score for a design \((d,y)\) (as in Sections VI.E and VII.E)
\[
s(d,y,\xi,\Theta) := \max\Big\{\frac{15-Q(d,y,\xi,\Theta)}{s_Q},\ \frac{g_1(d,y,\xi,\Theta)}{s_1},\dots,\frac{g_J(d,y,\xi,\Theta)}{s_J}\Big\}.
\]
Then BOED may be targeted at reducing uncertainty in the tail of \(s\) that controls the joint feasibility claim. One auditable way to formalize this is to choose a risk functional \(R(\Theta)\) (e.g., a CVaR of \(s\) under \(\xi\), or a conservative lower bound on \(\mathbb{P}(s\le 0)\) computed by one of Section VII’s methods) and to replace \(U(e)\) with an information gain about \(R(\Theta)\). The exact choice must be logged because it changes the meaning of “top”.
Audit requirement XI.1 (BOED reproducibility).
Any BOED ranking must log: the prior/posterior version \(p(\Theta\mid D)\), the experiment design space \(\mathcal{E}\), the likelihood model \(p(y\mid \Theta,e)\) including unit/sign conventions for residuals (Section III), the Monte Carlo estimator used for \(U(e)\) (sample sizes, random seeds), and the final ranked list with uncertainty bars (e.g., standard errors of \(\widehat U(e)\)).
XI.A.3. Decision rules linking test outcomes to feasibility-certificate improvement
Let \(\mathfrak{C}(D)\) denote a certificate object computed from data \(D\), for example:
- a conformal tube parameter \(\widehat q(D)\) for the joint violation score (Section VII.E),
- a scenario-bound design guarantee (Section VII.B), or
- a surrogate error envelope \(\Delta_j(D)\) used for constraint tightening (Section VI.A.3).
A test plan should specify an auditable go/no-go logic of the form:
1. If outcome \(y\) yields \(\mathfrak{C}(D\cup\{(e,y)\})\) with improved bound \(\mathbb{P}(\mathcal{F})\ge 1-\alpha\) at the target design(s), then proceed to the next integration stage;
2. If outcome \(y\) invalidates a core assumption (e.g., violates a declared Lipschitz bound or monotonicity region used in Sections IV, VII.C, or VIII.E), then mark the corresponding certificate pathways as FAIL and route redesign or model-class change.
XI.B. Certified Sobol variance decomposition for coupled outputs
Sobol indices provide a decomposition of output variance across uncertain inputs and their interactions. In this paper, Sobol analysis is not used as a standalone truth claim; it is used as an auditable ranking signal to decide which uncertainties most affect decision-relevant quantities (margins, joint violation score, availability proxies).
XI.B.1. Definitions
Assume a reference probability model for \(\xi=(\xi_1,\dots,\xi_{n_\xi})\) under which Sobol indices are defined (typically requiring independence of components; if independence is not defensible, the dossier must explicitly treat Sobol indices as *model-conditional diagnostics*).
Let \(f(\xi)\in L^2\) be a scalar QoI, e.g. \(f=s(d,y,\xi)\) or a single margin. Let \(\mathrm{Var}(f)>0\). The total-order Sobol index for variable \(\xi_i\) is
\[
S^{\mathrm{tot}}_i := 1-\frac{\mathrm{Var}(\mathbb{E}[f\mid \xi_{\sim i}])}{\mathrm{Var}(f)} = \frac{\mathbb{E}[\mathrm{Var}(f\mid \xi_{\sim i})]}{\mathrm{Var}(f)},
\]
where \(\xi_{\sim i}\) denotes all components except \(\xi_i\).
XI.B.2. Exact computation for PCE surrogates (when valid)
If \(f\) is represented by an orthonormal PCE under the reference measure \(\mu\) (Section VI.A), then Sobol indices can be computed algebraically from coefficients by grouping polynomial terms by their active variable sets. This is fully auditable **provided**:
(i) the PCE basis and measure \(\mu\) are logged;
(ii) coefficient uncertainty (from compressed sensing, Section VI.B) is propagated into bounds on \(S_i^{\mathrm{tot}}\).
A conservative bound can be produced by using coefficient error bounds \(\|\widehat c-c\|_2\le \delta_c\) to bound the numerator/denominator in the Sobol ratio; the dossier must state exactly how this propagation was done (e.g., interval arithmetic on polynomial coefficient group norms), because ratios can be unstable when \(\mathrm{Var}(f)\) is small.
XI.B.3. Certified Monte Carlo bounds (distribution-free, finite-sample)
When PCE is not trusted, Sobol indices can be estimated by Monte Carlo estimators (Saltelli/Jansen-type). For auditability, any reported \(\widehat S_i^{\mathrm{tot}}\) must include an uncertainty statement. One conservative option is to treat the estimator as a bounded random variable (after normalization/clipping of \(f\) to an auditable range) and attach a Hoeffding-type confidence interval; another is to use bootstrap with explicit caveats (bootstrap is not distribution-free and must be labeled as such).
Audit requirement XI.2 (Sobol audit).
Any Sobol-based ranking must log: the uncertainty model for \(\xi\) (including any independence assumptions), the estimator used, sample size and seeds, and an uncertainty quantification for \(\widehat S_i^{\mathrm{tot}}\). If the estimator relies on surrogate outputs, it must also log how surrogate envelopes \(\Delta\) and verified-numerics UNKNOWN outcomes were handled.
XI.C. Conformal coverage improvement as a V\&V metric
Section VII.E yields a split-conformal tube
\[
s(\xi) \le \widehat s(\xi)+\widehat q
\]
with marginal coverage \(\ge 1-\alpha\) under exchangeability. A test campaign should measure progress in terms of **coverage and tightness**, not only pointwise error.
XI.C.1. Tube-width and feasibility-margin metrics
Define a tube width functional (dimensionless)
\[
W(D) := \widehat q(D)
\]
and, for a given design \((d,y)\), define a conservative feasibility gap under the surrogate stack
\[
G(D;d,y) := \sup_{\xi\in\Xi_{eval}} \big(\widehat s(d,y,\xi)+\widehat q(D)\big),
\]
where \(\Xi_{eval}\) is an auditable evaluation set (e.g., a scenario set or a certified uncertainty box, depending on context). If \(G\le 0\), then the conformal tube certifies joint feasibility on \(\Xi_{eval}\) at level \(1-\alpha\) (conditional on exchangeability and the evaluation-set semantics). Reduction in \(W\) (or in \(G\)) is therefore a directly interpretable V\&V success metric.
XI.C.2. Integration with scenario/DRO/PCE
Conformal tubes can wrap any underlying predictor, including PCE and multi-fidelity stacks; scenario and DRO certificates can be used in parallel as robustness cross-checks. The test plan must specify which certificate is primary for each decision gate, because different certificates fail under different assumption violations (exchangeability vs Lipschitz vs convexity).
XI.D. Example test concepts (formalized as certificate-affecting experiments)
The following test concepts are framed only as examples of how an experiment connects to certificate parameters. They are not claims that these particular tests are sufficient.
XI.D.1. Divertor tile puff-testbed
Goal: tighten an envelope \(\Delta_q\) on peak heat-flux predictions and/or validate monotonicity in an edge-transport parameter \(\chi\) used for conservative enclosure (Section IX.D).
Contract-facing outputs:
– updated error envelope \(\Delta_q(D)\) or an interval enclosure routine for \(q_\perp^{peak}\);
– a falsification report if a claimed monotonic regime fails.
XI.D.2. Neutronics/tritium coupling validation
Goal: validate the mapping from blanket state to TBR and tritium-flow ports, reducing \(\Delta_{TBR}\) and tightening transient inventory envelopes (Section IX.E).
Contract-facing outputs:
– improved bounds \(\Delta_{TBR}^{RB}\) (if RB is used) or calibrated conformal tubes for \(TBR\);
– verification of unit-consistent port mapping in the tritium conserved-family graph.
XI.D.3. Full-twin “dry run” for certificate consistency
Goal: demonstrate that (i) unit/conservation checks pass, (ii) stability/well-posedness checks return PASS/UNKNOWN with logged assumptions, and (iii) verified-numerics PASS/UNKNOWN rates are measured and tied to runtime monitoring thresholds.
Contract-facing outputs:
– a reproducible integration dossier template with PASS/FAIL/UNKNOWN definitions executed end-to-end.
XI.E. Go/No-Go decision matrix (certificate-centered)
A practical decision matrix should be organized by certificate layers and their failure modes. A minimal structure is:
1. **Integration correctness gate (Sections III–V):** unit/sign checks PASS; conservation residual budgets met; verified steady-state solve PASS for a declared operating envelope (or else UNKNOWN count below a declared threshold with a mitigation plan).
2. **Stability/execution gate (Section IV):** passivity/scattering energy accounting within declared budgets; contraction/small-gain checks PASS or UNKNOWN with explicit V\&V plan to reduce gain uncertainty.
3. **Feasibility probability gate (Section VII):** at least one primary certificate method reaches \(\mathbb{P}(\mathcal{F})\ge 1-\alpha\) (with logged assumptions); secondary methods do not contradict it beyond logged conservatism.
4. **RAMI/performability gate (Section VIII):** availability/performability constraints certified under explicitly declared dependence assumptions; common-cause parameters either calibrated or bounded conservatively.
5. **Cost/schedule gate (Section X):** optimization outputs include solver certificates (gap/duals) and post-checks against the coupled twin.
XII. Implementation Workflow and Auditability Requirements
This section specifies the minimum metadata, logging, and pipeline structure needed so that the mathematical certificates of Sections III–XI are checkable artifacts rather than narrative claims.
XII.A. Data/metadata requirements for contracts
XII.A.1. Port schema
Every port variable must carry:
– a unique identifier (module, port name, graph edge),
– a conserved-family tag (energy/particle/tritium/other),
– a dimension vector \(\dim(\cdot)\in\mathbb{Z}^7\) and explicit unit string,
– a sign convention (direction on the port graph),
– an admissible range and/or tolerance budget \(\tau\),
– a scaling \(\sigma\) if dimensionless normalization is used.
XII.A.2. Contract residuals and tolerance budgets
For every residual inequality \(h\le \tau\) or equality \(r=0\) (implemented as \(|\hat r|\le \tau\)), the dossier must include:
– the symbolic definition of \(h\) or \(r\),
– normalization scalings \(\sigma\),
– the decomposition of \(\tau\) into declared sources (discretization/surrogate/coupling),
– the runtime values of the residuals and whether they triggered escalation.
XII.A.3. Provenance identifiers
All certificate-relevant objects must be versioned:
– module versions and configuration hashes,
– parameter priors/posteriors hash (Section XIII),
– solver versions and tolerance settings,
– random seeds and sample-set hashes for scenario/conformal/MLMC computations.
XII.B. End-to-end certification pipeline
XII.B.1. Pipeline stages and required artifacts
A minimal end-to-end run at \((d,y)\) should proceed:
1. **Contract compilation:** unit/dimension checks; port-graph construction; incidence matrices \(B_q\) generated and hashed.
2. **Execution readiness:** stability prerequisites logged (passivity indices, scattering parameters, delay bounds) and PASS/UNKNOWN computed under explicit assumptions.
3. **Coupled solve:** compute a candidate \(\hat x\) and run verified enclosure (Section V) to yield PASS/FAIL/UNKNOWN with inclusion logs.
4. **Constraint enclosure:** compute \([g_j](X)\) or envelope-tightened surrogate constraints, producing certified margins.
5. **UQ/certification:** compute one or more certificate statements (moment/scenario/DRO/SOS/conformal), logging treatment of UNKNOWN.
6. **RAMI/performability:** compute availability/performance-tail bounds with dependence declarations.
7. **Optimization linkage (if applicable):** attach solver optimality/feasibility certificates, dual explanations, and post-checks.
XII.B.2. PASS/FAIL/UNKNOWN triage and escalation rules
The workflow must define machine-enforceable escalation thresholds, e.g.:
– if unit/sign mismatch: immediate FAIL (integration bug);
– if conservation residual exceeds budget: FAIL or re-solve with tightened coupling;
– if verified solve UNKNOWN rate exceeds a declared fraction in decision-relevant regions: trigger V\&V/mesh refinement/surrogate upgrade;
– if certificate assumptions invalidated (e.g., Lipschitz bound falsified): mark certificate pathway FAIL and switch to an assumption-light method (scenario/conformal).
XII.C. Computational budgeting aligned with optimization needs
The computational plan must state where cost is spent and how acceleration methods are justified:
– compressed sensing for PCE coefficients (Section VI.B) only when coherence/sparsity diagnostics are logged;
– MLMC/MF-MLMC (Section VI.C–VI.D) only when bias/variance/cost pilot estimates are logged;
– active-subspace reduction (Section VI.E) only when tail-energy diagnostics are logged;
– verified numerics (Section V) applied to representative operating points and to samples that dominate certificate tightness.
Gap flag.
Any claimed speedup or wall-clock feasibility must be supported by reproducible benchmarks: hardware, code versions, and the exact hierarchy definition.
XIII. Bayesian Calibration and Digital-Twin Updating Under Contracts (Posterior UQ, Error Budgets, and Decision Provenance)
This section defines how the digital twin is updated as new data arrive, while preserving contract auditability and avoiding untracked shifts in certificate meaning.
XIII.A. Contract-conditioned Bayesian updating of epistemic parameters
XIII.A.1. Likelihood construction with explicit error separation
Let \(r(D;\Theta)\) denote a vector of unit-consistent residuals between observed data and model predictions, constructed from contract-declared port variables and normalized by declared scales \(\sigma\) (Section III.B). We explicitly separate three error sources (Section II.C.4) by a hierarchical noise model:
\[
r = r_{model}(\Theta) + \varepsilon_{num} + \varepsilon_{cpl} + \varepsilon_{meas},
\]
where each term has a declared interpretation and (if probabilistic) a declared distribution family/parameters.
A simple auditable Gaussian likelihood (used only as a working model, not a truth claim) is
\[
p(r\mid\Theta) \propto \exp\Big(-\tfrac12 r^\top \Sigma^{-1} r\Big),
\]
with \(\Sigma\) decomposed (and logged) as
\[
\Sigma = \Sigma_{num} + \Sigma_{cpl} + \Sigma_{meas} + \Sigma_{disc},
\]
where \(\Sigma_{disc}\) represents explicit discrepancy/model-form uncertainty if such a term is included. If \(\Sigma_{disc}\) is omitted, the dossier must mark the posterior as conditional on “no discrepancy model”.
XIII.A.2. Conditioning on hard constraints vs soft constraints
Contracts include hard constraints (unit consistency, sign conventions, admissible ranges) and soft constraints (tolerance residuals). Bayesian conditioning must distinguish:
– **Hard support restrictions:** \(\Theta\in\Theta_{adm}\) where \(\Theta_{adm}\) encodes physically or contractually impossible regions. This is implemented as prior support restriction.
– **Soft penalties/likelihood terms:** violations of approximate equalities or conservative envelopes are incorporated through likelihood terms with explicit noise levels.
Audit requirement XIII.1 (conditioning semantics).
For every constraint-like datum, the dossier must specify whether it enters as a hard support restriction or as a likelihood penalty, and must report how this choice affects certificates (e.g., narrower posterior may increase certificate tightness but risks overconfidence if a hard restriction is unjustified).
XIII.A.3. Updating dependence structure assumptions
RAMI common-cause parameters (Section VIII.B) and other dependence structures may be updated from evidence. This must be done explicitly: dependence parameters must be included in \(\Theta\) and posterior updates must log whether dependence assumptions changed. Any such change invalidates previous availability certificates unless re-issued.
XIII.B. Constraint-aware uncertainty representations for downstream certificates
XIII.B.1. Posterior objects usable by scenario/DRO/conformal layers
Different certificate layers require different uncertainty objects:
– scenario methods require a sampling procedure for \(\xi\) (which may include posterior draws of \(\Theta\));
– DRO requires an ambiguity set radius and a metric, plus (often) Lipschitz bounds;
– conformal requires exchangeability of calibration and deployment data.
Therefore, the calibrated twin must export a **posterior uncertainty bundle** containing:
– posterior samples \(\{\Theta^{(k)}\}\) (or a parametric posterior),
– a declared sampling method for \(\xi\) conditional on \(\Theta\),
– and the mapping to the uncertainty objects required by each certificate layer.
XIII.B.2. Handling verified-numerics UNKNOWN as censoring/outer approximation
If Section V returns UNKNOWN for a sample or operating point, then that evaluation cannot be treated as a precise scalar constraint value. Two conservative, auditable treatments are:
1. **Censoring:** treat the observation as interval-valued \(g_j\in[\\underline g_j,\overline g_j]\) and use an interval likelihood (or a conservative likelihood bound) when updating \(\Theta\).
2. **Outer approximation:** treat UNKNOWN as a constraint-violation possibility and propagate it as an increased envelope \(\Delta_j\) for downstream certificates.
The dossier must specify which treatment is used for each constraint and dataset.
XIII.C. Goal-oriented a posteriori error estimation as calibration and certification input
XIII.C.1. Abstract goal-oriented estimator statement
Suppose a high-fidelity module solves a PDE or large-scale discretized system \(\mathcal{R}(u;\Theta)=0\) for state \(u\), and a plant-level quantity of interest is \(J(u)\) (e.g., a stability margin). Let \(u_h\) be a numerical approximation and let \(z\) solve the adjoint (dual) problem associated with \(J\).
A standard goal-oriented a posteriori error estimate has the form
\[
J(u)-J(u_h) \approx \langle \mathcal{R}(u_h;\Theta), z\rangle,
\]
with computable upper bounds
\[
|J(u)-J(u_h)| \le \eta_J(u_h,z),
\]
where \(\eta_J\) depends on residual norms and stability constants.
Integration into contracts.
If \(J\) enters a plant constraint \(g_j\), then \(\eta_J\) becomes a contract-level deterministic envelope \(\Delta_j\) (Section VI.A.3) that can be consumed by Sections V and VII. The dossier must log the estimator definition and the assumptions under which it is valid (coercivity/stability bounds, residual evaluation accuracy).
XIII.D. Versioned decision dossiers and change-control logic
XIII.D.1. Immutable records
Each certificate issuance must be accompanied by an immutable record containing:
– priors/posteriors (or their hashes),
– contract versions and tolerance budgets,
– solver configurations and verified-numerics logs,
– sampling sets and certificate computation artifacts.
XIII.D.2. Re-certification triggers
The pipeline must declare re-certification triggers, including:
– any change in module code or surrogate definition,
– any change in port schema, unit definitions, sign conventions, or tolerance budgets,
– any posterior update that materially changes \(p(\Theta\mid D)\) in decision-relevant directions,
– any discovered falsifier that violates an assumption used by a certificate pathway.
These triggers are essential to prevent stale certificates from being treated as current evidence.
With Sections XI–XIII, the body of the paper provides a complete certificate-bearing workflow: (i) contract-based integration (III), (ii) stability and execution readiness (IV), (iii) verified numerics (V), (iv) scalable UQ and surrogate layers (VI), (v) probability-of-feasibility certificates (VII), (vi) RAMI/performability under dependence (VIII), (vii) standardized component module interfaces (IX), (viii) certificate-linked optimization (X), and (ix) V\&V and calibration/update mechanisms (XI–XIII).
Conclusion
This paper developed a certificate-oriented mathematical framework for integrating a compact stellarator “digital twin” as a heterogeneous system-of-systems under joint physics, engineering, RAMI, cost, and schedule constraints. The central theme was that an executable integrated model is decision-relevant only to the extent that it can emit solver-independent and assumption-explicit artifacts: unit- and sign-consistent interfaces; conservation residual budgets; stability and well-posedness checks; verified numerical enclosures; and probability-of-feasibility statements whose semantics are unambiguous and auditable.
At the integration level (Sections II–III), we formalized the coupled twin as a directed port graph whose edges carry typed, unit-tagged exchanged quantities grouped by conserved families. The incidence-matrix representation provides a structural guard against sign and double-counting errors, while contract residuals with explicit tolerance budgets make coupling error an auditable quantity rather than an implicit solver byproduct. When the assembled contracts are inconsistent, we emphasized infeasibility explanations via checkable dual/Farkas-type witnesses (or dual rays in convex settings), together with a practical, reproducible procedure for extracting small conflicting subsets.
At the execution level (Section IV), we separated physically consistent interconnection from numerically stable co-simulation. A port-Hamiltonian/dissipativity viewpoint links “do not inject nonphysical power” to explicit module-declared inequalities, while scattering-style interconnections provide an auditable mechanism for preserving passivity under multi-rate stepping and delays (provided parameters and energy balances are logged). For steady-state modular solves, we recast coupling as a fixed point and stated contraction/small-gain conditions that (i) imply existence/uniqueness in a declared domain and (ii) yield a quantitative amplification bound mapping per-module approximation/integration errors to plant-level shared-state error. An asynchronous iteration model with bounded delays was included to match realistic digital-twin scheduling, with explicit logging requirements for delay bounds and update persistence.
To prevent solver-dependent artifacts from being mistaken for plant truths (Section V), we attached a verified-numerics wrapper to the coupled steady-state solve. Interval extensions and a Krawczyk-type inclusion test provide a tri-valued PASS/FAIL/UNKNOWN classification with replayable evidence (candidate point, interval box, Jacobian enclosure, and operator bounds). When PASS holds, all constraints are evaluated with rigorous interval enclosures, producing certified margins; when UNKNOWN persists, the workflow forces an explicit decision rule (conservative treatment or targeted refinement) rather than silent overconfidence. This verified layer also provides an auditable pathway to detect (or at least flag) multi-solution regimes that complicate both robust design and certificate interpretation.
For scalability (Section VI), we presented surrogate and uncertainty-propagation machinery designed to preserve auditability by explicitly tracking surrogate/model-form envelopes. Polynomial chaos expansions serve as a common algebraic surrogate layer, but all low-degree/sparsity assumptions were treated as gap-flag items requiring validation. Compressed-sensing recovery was stated with conditional recovery bounds and explicit warnings about unverifiable RIP/coherence assumptions. MLMC and multi-fidelity variants were presented via assumption-explicit bias/variance/cost models and with logging requirements to prevent implicit bias from being misreported as statistical uncertainty reduction. Active-subspace reduction was framed as a diagnostic/certification-adjacent tool: what can be audited are matrix-estimation uncertainties and tail-energy diagnostics, not blanket claims of low-dimensional truth.
Building on these ingredients, Section VII assembled multiple, complementary probability-of-feasibility certificate layers whose assumptions are distinct and therefore must not be conflated: distribution-free moment bounds (useful but potentially loose for joint events), scenario-based finite-sample guarantees (conditional on convexity and i.i.d. sampling), Wasserstein DRO certificates (conditional on auditable Lipschitz bounds and metric/radius choices), SOS-based robust feasibility over semialgebraic supports (powerful but scaling-limited and conservative), and conformal prediction tubes (distribution-free under exchangeability and applicable to black-box surrogate stacks). A key methodological point was the explicit handling of verified-numerics UNKNOWN outcomes: probability statements remain meaningful only if the mapping from tri-valued evaluations to the bound’s inputs is conservative and logged.
Section VIII integrated RAMI/performability into the same certificate pipeline. Minimal cut sets yield immediate distribution-free bounds that remain valid without independence assumptions, while common-cause “shock” variables provide an auditable way to represent dependence when supported by evidence. We also included dependence-aware rare-event approximations (Stein–Chen style) as an optional, explicitly conditioned tool. For partial degradation, multistate capacity modeling via universal generating functions gives a tractable performability layer with quantization-induced conservatism that can be bounded and logged. A concrete coupling pattern was provided: invert monotone physics surrogates to obtain required auxiliary-capacity thresholds (with envelopes), then certify the probability that degraded plant capacity exceeds those thresholds.
Sections IX–XIII completed the decision-focused synthesis. Component-level physics/risk modules were not assumed correct; instead they were required to present contract-consumable outputs: unit-checked quantities, declared validity domains, and explicit error envelopes (deterministic or probabilistic). Optimization under cost and schedule constraints was then framed as producing dossiers rather than merely candidate designs: GP/MICP/GBD formulations can yield solver-side optimality/feasibility evidence for the surrogate model, but must be paired with post-checks against the coupled twin via contracts, stability prerequisites, verified numerics, and an explicitly chosen feasibility-probability certificate. V\&V planning was formulated as decision-directed uncertainty reduction using BOED, sensitivity diagnostics, and conformal coverage/tube tightening, with a go/no-go matrix organized by certificate layer failure modes. Finally, Bayesian calibration was cast as contract-conditioned updating with explicit separation of model-form, numerical, coupling, and measurement errors, together with change-control rules that prevent stale certificates from persisting across model revisions.
The main limitation of the framework is not mathematical completeness but assumption credibility: passivity indices, Lipschitz/gain bounds, surrogate envelopes, dependence structures in RAMI, and any monotonicity/tail-model assertions must be empirically justified for each concrete module and operating regime. Accordingly, the paper’s outputs should be read as a rigorous *protocol* for producing and auditing digital-twin evidence, not as a guarantee that any particular compact Q>15 stellarator design is feasible. The immediate open tasks are therefore operational: (i) instantiate the contract schema for a specific module set; (ii) establish credible, module-level evidence for the declared bounds/indices; (iii) benchmark the verified and probabilistic certificate layers for tightness and computational cost; and (iv) demonstrate end-to-end decision dossiers whose PASS/FAIL/UNKNOWN logic remains stable under versioned updates. Subject to these validations, the framework provides a principled pathway to integrate heterogeneous models into a single, auditable decision engine for performance, safety, availability, and cost/schedule feasibility under uncertainty.
[HARD CODED END-OF-PAPER MARK — ALL CONTENT SHOULD BE ABOVE THIS LINE]
================================================================================
MODEL CREDITS
This autonomous solution attempt was generated with the Intrafere LLC AI Harness,
MOTO, and the following model(s):
– x-ai/grok-4.1-fast (67 API calls)
– openai/gpt-5.2 (30 API calls)
– moonshotai/kimi-k2.5 (8 API calls)
Total AI Model API Calls: 105
================================================================================