Executing the trade studies on a SysML v2 plant model

1.1Exploring the design space: parameters, components, and architecture

Each kind of design change needs different support from the model and its executable program.

Figure 1 sets the three kinds side by side. Each card names a kind of change and describes it; the band above the cards says how we iterate on it.

These choices lead to two ways of iterating. A numerical study runs the existing program with different inputs. Changing component definitions or plant architecture generally means editing the SysML model and generating a new program. Both let us evaluate the consequences of a design decision.

1.2Why build an execution pipeline?

Two reasons. First, when we started, we did not have an out-of-the-box SysML v2 execution tool ready for this workflow. Second, we needed machinery that agents could easily configure and run within a harness. That meant programmatic control over the model's inputs and repeated evaluations, along with the ability to incorporate custom calculations beyond what our SysML execution tools supported.

To support those studies, we built a code-generation pipeline that turns an assembled SysML plant model into an executable Python program. The program accepts design parameters, calculates performance and cost, and checks engineering limits. TEAx runs it repeatedly so a study can compare different design points.

Python also gives us access to more sophisticated computation and flexibility than we can express and execute through our supported SysML subset alone. Later, we will show how a custom calculation, such as a surrogate for a detailed simulation, can become part of the plant model.

1.3From component models to a study

Our demo models and explores stellarator designs. To show how the codegen pipeline works, we will focus on a study of the magnet system: how changing plant radius and coil ampere-turns affects cost and the design's ability to meet engineering limits.

We will follow the component equations and plant connections through code generation, then show how the resulting program evaluates and compares design points.

The worked example follows our stellarator demo through seven steps:

• 2.1 Model a component: separate the magnet's reusable equations from the values and connections chosen for the plant.
• 2.2 Define acceptable designs: express engineering limits as constraints that each design point must pass.
• 2.3 Construct the calculation network: determine where every calculation gets its inputs and which other calculations it depends on.
• 2.4 Implement the mathematics: translate equations into Python and incorporate custom methods, including a simulation surrogate.
• 2.5 Change the component model: give the blanket a more detailed calculation while keeping its place in the plant.
• 2.6 Assemble the executable program: connect the Python calculations with TEAx and expose the inputs a study can vary.
• 2.7 Run a study: evaluate many design points and compare their feasibility and cost.

All code examples come from the stellarator demo as read on 19 September 2026. Excerpts omit surrounding imports and documentation where indicated.

2.1Modeling a component: separate its equations from the plant's choices

We begin by defining a component’s calculations and connecting their inputs and outputs to the plant model.

The stellarator's magnets produce the magnetic field that confines the plasma. For a proposed design, we want to calculate that field from the chosen coil ampere-turns and geometry. Other calculations then use the result to evaluate the plasma and check the magnet's engineering limits.

The magnet model uses a reusable field equation. The plant supplies the design values and connects the magnet to the plasma's radius. The figure below shows how those values reach the equation and how its result becomes available to the rest of the plant.

calc def 'Coil Set Axis Field' {
    in attribute n_coils : Real;
    in attribute I_coil : Real;
    in attribute k_link : Real;
    in attribute R0 : Real;
    in attribute mu0 : Real default 1.25663706212e-6;
    in attribute two_pi : Real default 6.283185307179586;

    out attribute B_axis : Real = mu0 * k_link * n_coils * I_coil / (two_pi * R0);
}

part stellaris : 'MFE Power Plant' {
    part :>> magnet {
        part :>> coil {
            :>> n_coils = 48.0 {}
            :>> I_coil = 15400000.0 {}
            :>> k_link = 0.7731331164622419 {}
        }
    }
}

:>> R0 = plasma.R {}

calc field_calc : 'Coil Set Axis Field' {
    in n_coils = coil.n_coils;
    in I_coil = coil.I_coil;
    in k_link = coil.k_link;
    in R0 = coil.R0;
}

attribute B : Real = field_calc.B_axis;

The walk below follows the numbered boxes in Figure 3: the plant's values (1), the calculation usage that selects them (2), the equation they feed (3), and the attribute that exposes its result (4).

Start with the plant's choices. At the reference point, it specifies 48 coils and 15.4 million ampere-turns per reference coil. The coil's radius attribute, coil.R0, is connected to the plasma's major radius, plasma.R. Changing the plant's plasma radius therefore also changes the radius used by the magnet.

The calculation definition declares the field equation and its inputs. Inside the magnet, the calculation usage, field_calc, connects the coil's attributes to those inputs. For example, it binds the equation's R0 input to coil.R0, completing the connection from the plant's plasma radius to the field equation. The equation's constants, mu0 (vacuum permeability) and two_pi (twice π), come from their declared defaults.

The equation returns the calculated field as B_axis. The magnet exposes that result through an attribute named B, which other components and calculations can read. Downstream calculations use it to determine peak conductor field and plasma beta, the ratio of plasma pressure to magnetic pressure. A new radius or coil ampere-turn value therefore affects both the field calculation and the calculations that depend on its result.

2.2Constraints test whether a design point meets engineering limits

Next, we add constraints that test whether a design meets the modeled engineering limits.

Later, we will show how a study sweeps across a wide range of input combinations, such as plant radius and coil ampere-turns. Each complete set of input values defines a design point. The calculations tell us what performance and cost each point produces, but some combinations will violate engineering limits. We need a way to identify those failures as the study runs.

We express engineering limits as constraints: expressions that evaluate to true or false for a design point. A constraint can compare a calculated quantity with a limit, or relate several quantities that must be consistent with one another. In SysML, an assert constraint applies the condition to the model so the generated program can check it during evaluation.

For example, the magnet's conductor is arranged into a winding pack, which must fit inside its casing. A larger pack may provide the needed current capacity, but it also needs enough physical space.

The magnet's fit calculation compares the pack dimensions with the available space, accounting for insulation and assembly allowances. It calculates clearance in two directions and exposes the smaller margin as fit_minimum_margin.

Figure 4 shows what that calculation compares, in a cross-section of one coil. Read it from the outside in: the casing wall, the cavity it leaves, then the required envelope, which is the pack plus its insulation and assembly clearance on every face. The shaded band between the envelope and the cavity is the margin, one value in the radial direction and one in the transverse direction; the check reads the smaller of the two. The proportions are schematic, and the layers are drawn thicker than they are.

calc wp_fit : 'Winding Pack Casing Fit' {
    in wp_side = winding_pack.wp_side;
    in aspect_ratio = winding_pack.fit_aspect_ratio;
    in internal_fraction_x = winding_pack.internal_build_x;
    in internal_fraction_y = winding_pack.internal_build_y;
    in ground_insulation = winding_pack.ground_insulation;
    in radial_allocation = coil.coil_t;
    in interior_y = casing.interior_y;
    in wall_thickness = casing.wall_thickness;
    in assembly_clearance = casing.assembly_clearance;
}

attribute fit_minimum_margin : Real = wp_fit.minimum_margin;

The reusable constraint definition models/library/analyses/mfe_winding_pack_fit.sysml states the condition: the minimum clearance must be nonnegative.

constraint def 'Winding Pack Fits Casing' {
    in attribute minimum_margin_in : Real;
    minimum_margin_in >= 0.0
}

The stellarator plant applies that condition to its magnet's calculated clearance:

assert constraint wp_fit_ok : 'Winding Pack Fits Casing' {
    in minimum_margin_in = magnet.fit_minimum_margin;
}

A negative clearance means the winding pack does not fit. The program still reports the design's calculated performance and cost, alongside the failed check. That lets a study show what a proposed change improves and which engineering limits it violates.

A design point passes the modeled feasibility checks only when all applicable assertions have been evaluated and satisfied. Results identify which assertions were assessed; an unassessed check cannot count as a pass. Even a point that passes every check still needs evidence for engineering limits the model does not yet represent.

2.3Reconstructing the calculation network

Codegen now identifies the calculations used by the assembled plant and determines where each one gets its inputs.

Those connections determine what must be calculated first. The program cannot check winding-pack clearance until it knows the pack dimensions, and in the version shown here it could not calculate those dimensions until it knew how much conductor the current needed. TEAx will use these dependencies to order the calculations when it runs the program.

SysML organizes the model around physical components. Codegen constructs a second view organized around calculations: each calculation is a node, and a connection from one calculation's output to another's input is an edge. This is the computational graph. It tells the execution system which results must be available before a calculation can run. The following figure shows the magnet branch of the generated stellarator graph.

Read the graph from the top. The calculated field and magnet geometry determine the peak field at the conductor. That field affects how much current the conductor can carry, so it feeds the calculation that sizes the conductor supply. The sizing result then feeds the calculation of winding-pack dimensions.

The pack dimensions feed two branches: one checks whether the pack fits, while the other calculates material quantities and their procurement cost. A separate check determines whether the resulting conductor supply can carry the required current. Following the arrows shows how a change in field can affect both the physical fit and the cost of the magnet.

This graph is the model as it stood on 19 September. The next day we removed the two sizing calculations, because they chose the winding's size from the current instead of evaluating a size the designer supplies (section 2.7.1). The pack side is now a supplied input: the fit and volume calculations read it directly, the cost calculations follow from the volume, and the current check reports whether that winding carries its current. Codegen resolves and orders those connections the same way.

This execution order depends on the model author declaring each calculation's inputs and outputs. If two calculations each need the other's result before they can run, neither can start. The pipeline therefore requires a network without circular dependencies (a directed acyclic graph, or DAG). It does not rearrange equations or build a solver automatically; a coupled calculation needs a numerical method supplied by the author. The direction of the connections also fixes which values a study can set: a value one calculation produces cannot also be supplied as a design choice.

Codegen constructs this network in three steps.

2.4Turning modeled mathematics into executable Python

The graph tells us which calculations are needed and how they depend on one another. Each calculation also needs executable code that turns its input numbers into output numbers. For equations written in the supported SysML expression language, codegen generates that code directly, so the author does not have to maintain a second copy of the equation in Python.

Consider the magnetic-field equation from section 2.1. SysIDE represents its arithmetic as nested operations: multiply the numerator's terms, multiply the denominator's terms, then divide. Codegen visits those operations and writes the equivalent Python expression. The equation's I_coil input, for example, becomes inputs.I_coil in Python. The generated function exploration/stellarator_e2e/generated/handwritten/mfe_magnet_field/coil_set_axis_field_impl.py returns:

return ((((inputs.mu0 * inputs.k_link) * inputs.n_coils) * inputs.I_coil)
        / (inputs.two_pi * inputs.R0))

Only the line break has changed for readability. A generated wrapper checks the inputs and packages the result so TEAx can call the calculation.

def run(
    self, two_pi: float, I_coil: float, R0: float, n_coils: float, k_link: float, mu0: float    ) -> ModuleResult[Float]:
    # Validate inputs
    validated_inputs = self.validate_and_fill_default(two_pi, I_coil, R0, n_coils, k_link, mu0)

    # Import handwritten implementation
    from stellarator_tea.handwritten.mfe_magnet_field.coil_set_axis_field_impl import (
        run_coil_set_axis_field,
    )

    # Execute implementation - returns single value
    B_axis = run_coil_set_axis_field(validated_inputs)

    # Single output - return Float directly (RootModel[float])
    # TEAx assigns entire return value to the one channel declared in YAML
    return ModuleResult(data=Float(B_axis))

The compiler supports a subset of SysML mathematics. The generated code does not automatically convert units, so equations and supplied values must use consistent units.

2.4.1Incorporating calculations beyond the expression translator

Some engineering calculations need numerical algorithms, external libraries, or a fast approximation to an expensive simulation. For these methods, the model declares the calculation's inputs, outputs, and plant connections. Codegen generates the Python interface, and the author supplies the computational method. TEAx can then run that method alongside the automatically translated equations. An untranslated calculation cannot run until its implementation is supplied.

The blanket surrounding the plasma uses fusion neutrons to breed tritium fuel. We want to explore how its thickness affects breeding, but running a detailed neutron-transport simulation at every design point would be expensive.

Instead, the model uses a response table from prior neutron-transport calculations models/designs/stellarator_09/breeding_response.json. A handwritten Python function exploration/stellarator_e2e/generated/handwritten/mfe_tritium_breeding/blanket_tritium_breeding_impl.py estimates the breeding response between the stored results by interpolation. That function is a surrogate: a cheaper approximation to the detailed simulation, which the plant can call at each design point.

Figure 7 shows the five stored results and the response the function returns between them. Each dot is one transport calculation; the line between two dots is the straight interpolation the function uses. The shaded band spans the stored results.

The table contains results at five blanket thicknesses spanning 0.60–1.00 m, with the remaining geometry and scenario assumptions fixed. The function checks that the inputs stay within those conditions and marks the response undefined outside them.

index = min(bisect.bisect_right(thicknesses, thickness) - 1, len(nodes) - 2)
left, right = nodes[index:index+2]
fraction = (thickness - left['thickness_m']) / (right['thickness_m'] - left['thickness_m'])
li6 = left['tbr_li6'] * (1-fraction) + right['tbr_li6'] * fraction
li7 = left['tbr_li7'] * (1-fraction) + right['tbr_li7'] * fraction
mean = li6 + li7
variance = ((1-fraction)*left['std_error'])**2 + (fraction*right['std_error'])**2

Both translated equations and custom methods expose declared inputs and outputs, so codegen can connect them into the same plant calculation network. The implementation note archive/write-up/sysml-codegen-assets/model-evidence.md, "Implementation preservation", covers how Python bodies are retained when regenerating a package.

2.5Selecting a different component model

Changing numerical inputs lets us explore a component model’s existing behavior. To explore a different component or add more detailed behavior, we need to select a different definition for that component in the plant.

The blanket provides an example of adding detail. Its base model supplies a fixed breeding estimate, so varying blanket thickness does not change the breeding result. We want to replace that assumption with the thickness-dependent calculation introduced in section 2.4.

The stellarator design inherits a common plant assembly that already includes a blanket and its connections. The more detailed definition extends the base blanket, and the stellarator selects it for the existing blanket usage. Other parts of the plant can continue reading the breeding result from that same component.

The three declarations below express those choices. They come from separate source locations; the bodies are omitted:

part def 'Blanket' :> 'CAS22 Power Core' {
    attribute tbr : Real default 1.0;
    // ...
}

// Specialized component definition
part def 'Transport Calculated Blanket' :> 'Blanket' {
    // ...
}

// Generic plant: the inherited component usage
part blanket : 'Blanket' {
    // ...
}

// Stellarator design: select the specialization
part :>> blanket : 'Transport Calculated Blanket' {
    // ...
}

part def 'Transport Calculated Blanket' :> 'Blanket' {
    attribute breeding_R : Real;
    attribute breeding_a : Real;
    // ... five more geometry attributes the plant binds
    calc breeding : 'Blanket Tritium Breeding' {
        in blanket_t_in = blanket_t;
        in R_in = breeding_R;
        in a_in = breeding_a;
        // ... the remaining geometry inputs
    }
    :>> tbr = breeding.tbr_mean;
    attribute tbr_lower : Real = breeding.tbr_lower;
    attribute breeding_defined : Real = breeding.defined_flag;
    // ...
}

After this model change, we regenerate the program. Codegen includes the selected blanket's calculation and connects it to the plant, so the breeding result now comes from the surrogate. A study can vary blanket thickness and evaluate its effect on breeding alongside the other plant calculations.

2.6Putting together the Python computational model

We now have the calculation network and its Python functions, including any custom implementations. To run the plant, we need instructions for passing values between those functions and a way to supply the external inputs. The generated package captures both: a TEAx specification describes the connections, and input files hold the externally supplied parameters.

2.6.1Wiring the functions into an executable calculation network

A Python function for winding-pack size knows how to calculate a size from current density and ampere-turns. It does not decide which plant values to use. Codegen writes those choices into a YAML specification using the connections resolved in section 2.3. TEAx reads that specification to assemble and run the program.

Here is the complete pack-sizing entry from the stellarator pipeline:

stellarator_09__stellaris__magnet__wp_sizing:
  module_type: mfe_magnet_field.Winding_Pack_SizingModule
  inputs:
    j_wp: float stellarator_09__stellaris__magnet__current_sizing__selected_effective_density
    I_coil: float stellarator_plant_params.stellarator_09__stellaris__magnet__coil__I_coil
  outputs:
    root: RootModel[float] stellarator_09__stellaris__magnet__wp_sizing__wp_side

The module_type selects the winding-pack sizing calculation. Its j_wp input receives the density calculated by current sizing, while I_coil receives the supplied coil ampere-turns. The output entry names the calculated pack side so other calculations can read it. The long names identify each value's location in the plant. Generated pipeline: exploration/stellarator_e2e/studies/20260917-pre-reveal-feasible-neighborhood/preparation/package-pipelines/pipeline.yaml.

Figure 8 draws the same entry as connections. Read it left to right: the two inputs and where each comes from, the module in the middle, and the two calculations that read its output, taken from the bindings behind Figure 5.

stellarator_09__stellaris__magnet__field_calc:
  module_type: mfe_magnet_field.Coil_Set_Axis_FieldModule
  inputs:
    two_pi: float mfe_magnet_field_params.stellarator_09__stellaris__magnet__field_calc__two_pi
    I_coil: float stellarator_plant_params.stellarator_09__stellaris__magnet__coil__I_coil
    R0: float stellarator_plant_params.stellarator_09__stellaris__plasma__R
    n_coils: float stellarator_plant_params.stellarator_09__stellaris__magnet__coil__n_coils
    k_link: float stellarator_plant_params.stellarator_09__stellaris__magnet__coil__k_link
    mu0: float mfe_magnet_field_params.stellarator_09__stellaris__magnet__field_calc__mu0
  outputs:
    root: RootModel[float] stellarator_09__stellaris__magnet__field_calc__B_axis

TEAx uses these connections to run each calculation after the results it needs are available. It gathers the inputs, calls the Python function, and makes the outputs available to the next calculations. In the 17 September package, it computed field before current-based sizing, sizing before pack dimensions, and dimensions before the fit check.

2.6.2Exposing the free input parameters

The connected calculations still need starting values. In our example, the program calculates magnetic field, but the 17 September study supplied radius and coil ampere-turns; the current model calculates the ampere-turns from the supplied turns and current per turn. These externally supplied values form the program's input interface. Changing them lets us evaluate a new design point without regenerating the program.

Codegen identifies these externally supplied values while resolving the graph and exposes them as public inputs. The model's authored values provide their initial settings. Values supplied by another calculation remain connected to that calculation's output.

Here is an exact subset of the stellarator's generated plant-input JSON:

{
  "stellarator_09__stellaris__plasma__R": 12.7,
  "stellarator_09__stellaris__magnet__coil__I_coil": 15400000.0,
  "stellarator_09__stellaris__magnet__coil__n_coils": 48.0,
  "stellarator_09__stellaris__magnet__coil__k_link": 0.7731331164622419
}

These settings supply the plasma radius, coil ampere-turns, coil count, and linkage factor. Changing the radius changes the value received by every calculation connected to that plant attribute. The magnetic field remains a calculated result. Full generated input file: exploration/stellarator_e2e/studies/20260917-pre-reveal-feasible-neighborhood/preparation/package-inputs/stellarator_plant_params.json.

The long JSON keys identify each input's location in the plant. File-based runs load these settings from JSON; studies can supply the values directly in memory. The input-generation note archive/write-up/sysml-codegen-assets/model-evidence.md, "Input-generation details", covers file grouping, defaults, and inputs that need values before evaluation.

A study chooses which inputs to vary and which to hold fixed. In the radius–current study, radius and ampere-turns vary while coil count and linkage factor remain fixed.

2.7Running a study: compare feasibility and cost across design points

We can now use the generated program to ask an engineering question: which combinations of geometry and operating parameters meet the modeled limits, and what do they cost?

Figure 9 follows one study from left to right: the inputs it chooses, the program it loads once, the design points it evaluates, and the record it keeps. Below them are the study's two uses on this page, each a section of its own.

The current study runner evaluates grids and prepared lists of candidates. An external optimizer can call the same evaluator; the study-runner note archive/write-up/sysml-codegen-assets/study-evidence.md, "Study runner capabilities", distinguishes that option from the built-in study strategies.

2.7.1See how a local change affects the rest of the design

With the tape performance we assume, the published Stellaris winding cannot carry its current: its tape would run at about 1.7 times its estimated critical current, where we allow 0.8. The winding also does not fit the radial space our model allocates for it. The 15 September joint magnet-sizing study exploration/stellarator_e2e/studies/20260915-joint-magnet-sizing/report.md evaluated two changes intended to fix this: enough tape to carry the current, then more space for the larger winding.

In that study, the model chose the amount of tape from the current requirement. We removed that option on 20 September, because a calculation that sizes equipment takes the choice away from the designer (Part 4b, section 2). The model now evaluates the winding a designer supplies, both its size and the space allocated for it, and reports any shortfall instead of resizing it. The three cases below are the September study's results: the original design, a design with enough tape to carry the current, and a design with more space for that tape.

Read Figure 10 left to right. Each card is one evaluated case: its electricity cost, the conductor it needs, and four of its check results. Between the cards is the one change made, with what it fixed and what it cost or broke. The small drawing on each card uses the vocabulary of Figure 4: the clear cavity with the required envelope inside it, or spilling out of it.

In the study, sizing the conductor supply to carry the required current fixed the current-capacity problem but increased pack size and cost. Enlarging the allocated space fixed the fit problem, but the changed radial geometry raised peak conductor field to 25.30 T, above the retained 24.9 T limit.

The 24.9 T limit is Stellaris's published peak field, held as a ceiling. The model is calibrated so the published design sits exactly on it, so any rise in field fails the check. The field calculation also responds only to coil position and leaves out a term that would lower the field for a larger winding, so the model does not establish which way the real field moves.

Evaluating the connected model shows why fixing one limit does not necessarily produce an acceptable design. Each proposed change must be checked against its effects on the rest of the plant.

2.7.2Map feasibility and cost together

The 17 September neighborhood study exploration/stellarator_e2e/studies/20260917-pre-reveal-feasible-neighborhood/report.md provides the broader radius–current sweep introduced earlier. In the results shown below, major radius and coil ampere-turns vary while the other inputs remain fixed. Each mark represents one evaluated design point.

The passing samples form a narrow band. At fixed ampere-turns, reducing radius increases the modeled field and can violate its limit. Increasing radius reduces field and changes the plasma power balance, eventually running into divertor heat or installed-heating limits. The cost plot explains why minimizing price alone would not answer the design question: low calculated costs also occur at points that fail engineering checks.

This pattern is conditional on the other inputs held fixed for the sweep. The study evidence archive/write-up/sysml-codegen-assets/study-evidence.md, "Held inputs for the radius–current map", summarizes those settings and links to the complete record.

These studies turn design choices into engineering feedback. We can see what a proposed change improves, which limits it violates, and what to investigate next. When that investigation calls for a different component model or additional physics, we revise the model, regenerate the program, and run the next study.