Private API — kernel backend

KernelCodegen

Compile one ModelingToolkit System with declared inputs and outputs into plain callable functions. generate_io_function is the only entry point; everything else is upstream NetworkDynamics support code kept close to its original form so it stays diffable against a fresh checkout.

generate_io_function(sys, inputs, outputs; verbose=false, cse=true)

Compile sys into the callables one component type needs. inputs and outputs name variables of sys (array-valued ones are scalarized). Returns a named tuple:

• f — state derivatives, nothing when the component has no state
• g — the declared outputs
• obs — every remaining observed variable, nothing when there are none
• mass_matrix, and the symbolic states, inputs, outputs, obsstates, params, callable_params in the order the buffers use
• reads, which inputs and states each output and each state derivative depends on, as dependency_indices lists

All three callables take the same argument list, (target, u, input, numeric, callables, instances, batch, t), and run over every instance named in batch, as compile_batched wraps them.

compile_batched(func_expr, target_field) -> RuntimeGeneratedFunction

Wrap a build_function expression in a loop over component instances, giving

(target, u, input, numeric, callables, instances, batch, t)

which writes target[inst.target_field] for every instance named in batch. target_field is :states for a state-derivative map, :outputs for an output map and :observables for an observed map.

The body appears once, so what is compiled grows with the number of component types, not with the number of instances; calling the scalar body per instance instead leaves a real call boundary, since kernel bodies are too large to inline.

The loop's own names are fixed rather than gensymed, so two compiles of the same component share a type: a RuntimeGeneratedFunction's type is a hash of its argument names and body. The # in each name keeps it from colliding with the body's own variables.

simplify_with_mtkcompile(sys, allinputs, alloutputs; verbose)

Run mtkcompile with the declared inputs left unbound and the declared outputs kept, then scalarize what it left as array symbolics and repair the metadata it stripped. Returns (sys, eqs, obseqs_sorted, states, params) with states ordered to match eqs.

drop_unused_inputs(inputs, used) -> Set

The inputs to hide from mtkcompile because no equation reads them. An array is all-in or all-out here: an array some of whose components are read stays, and mtkcompile_inputs then passes it whole.

io_array_base(sym)

The array a scalarized element belongs to (pos[2] → pos), or the symbol itself when it is not an array element.

mtkcompile_inputs(inputs, used) -> Vector

The input list mtkcompile accepts, given the inputs left after drop_unused_inputs and the set of variables the equations actually read. mtkcompile refuses part of a declared array — "the entire array must be an input" — and equally refuses an input that appears in no equation, so an array only some of whose components are read is passed as the whole array symbolic and everything else element by element. Our components read parts of vectors routinely: a wrench reads pos[3] for the air density at its height, and a flap hinge about [0, 1, 0] leaves six of its nine frame entries unread.

pick_best_alias_names(eqs, obseqs, states, outputs, inputs; verbose)

Consolidate the alias chains simplification leaves behind. Pure alias equations a ~ b in obseqs are grouped into transitively connected clusters; one main representative is picked per cluster (differential states first, then inputs, then outputs, then the least deeply nested name); the substitution nonmain → main is applied throughout eqs, obseqs and states, and canonical nonmain ~ main observations are reinserted.

_alias_connected_components(pairs)

Given alias pairs [(a,b), (c,d), ...], returns groups of transitively connected variables as a Vector{Vector}. E.g. pairs (a,b), (b,c) → [[a, b, c]].

_scalarize_eqs(eqs) -> Vector{Equation}

Expand array-valued equations (lhs::Array ~ rhs::Array) into their scalar component equations; scalar equations pass through unchanged. Used so array I/O variables can be scalarized without vector equations leaking into the scalar-only codegen pipeline.

_scalarize_system(sys) -> System

Rebuild sys with every array-valued unknown/parameter/equation expanded to its scalar components, so the scalar-only codegen sees one consistent set of scalar symbolics (avoids array/scalar identity mismatches between our scalarization and mtkcompile's). A no-op when the system already has no array symbolics.

_scalarize_io_syms(syms) -> Vector

Expand any array-valued symbolics in syms into their scalar components, leaving scalars untouched, so array I/O variables (pos(t)[1:3]) become pos[1], pos[2], pos[3] for the otherwise scalar-only codegen pipeline.

_match_scalar_element(invalids, valid)

When invalids is a scalar array element (pos[1]) whose metadata-carrying resolution valid came back as the whole array, return the matching scalarized element of valid (Symbolics.scalarize(valid)[i]); otherwise return valid.

Build the build_function argument list for eqs: one Let block holding every observed assignment eqs needs plus the equation assignments, returning the first output, then the remaining outputs.

is_nonnumeric_param(p) -> Bool

Whether a parameter symbolic is nonnumeric, i.e. a callable operator (an interpolant/polar declared as (name::T)(..)). Its symtype is a SymbolicUtils.FnType. Numeric scalars and scalarized array elements have a Number symtype and stay in the flat parameter buffer.

scalar_name(sym) -> Symbol

Unique label for a symbolic. A scalarized array element pos(t)[2], whose getname collapses to :pos, becomes :pos_2, so every slot of a buffer has its own name.

param_default_value(defaults, p)

The build-time default of parameter p in defaults (a system's initial conditions). A scalarized array element p[k] falls back to element k of its whole-array default, which is where _scalarize_system keeps it. Returns nothing when the parameter has no default.

eq_type(eq::Equation)

Checks the type of the equation. Returns:

• (:explicit_diffeq, lhs_variable) for explicit differential equations
• (:implicit_diffeq, nothing) for implicit differential equations
• (:explicit_algebraic, lhs_variable) for explicit algebraic equations
• (:implicit_algebraic, nothing) for implicit algebraic equations

match_diff_states(eqs, states)

Match a set of unknowns to a vector of equations. Returns a vector of states with diff states at the correct places and non-diff states filled in arbitrarily.

getproperty_symbolic(sys, var; might_contain_toplevel_ns=true)

Like getproperty but works on a greater variety of "var"

• var can be Num or Symbolic (resolved using getname)
• strip namespace of sys if present (don't strip if might_contain_top_level_ns=false)
• for nested variables (foo₊bar₊baz) resolve them one by one

The concrete BasicSymbolic type every symbolic in these pipelines has.

The set of D(x) terms appearing on the right-hand side of eqs.

Diagonal mass matrix: 1 for each explicit differential equation, 0 for each implicit algebraic one. Returns I when every equation is differential.

The variables in exprs that lost their symbolic metadata during simplification.

Substitute metadata-carrying symbolics back into invalid_eqs in place, resolving each stripped variable by name against sys. Scalarized array elements never resolve, so the residual report is behind verbose rather than a warning.

Like get_variables but replaces any D(x) term with the inner variable x, which get_variables intentionally treats as one atomic symbol.

Insert newobs (already in topological order) into obseqs, each just after the last equation it depends on.

The (a, b) pair of an equation that is a pure alias a ~ b, else nothing.

dependency_indices(eqs, obs_subs, symbols) -> Vector{Vector{Int}}

For each equation, the positions in symbols its right-hand side reads once every observed substitution is resolved. obs_subs is in definition-before-use order, so each observed variable is resolved once and reused, where _all_dependencies re-descends the whole chain per call.

The positions term reads, taking an already-resolved observed variable's own.

The symbols term depends on once every substitution in dict is resolved recursively.

SlotMap(symbols)

Name → position map for one of a kernel's buffers. Array-valued symbolics are scalarized before codegen, so pos(t)[2] is stored under :pos_2; the group :pos additionally resolves to all of pos_1 … pos_n in order, which is how a whole vector is wired in one call. A name that carries no trailing index is its own one-element group. A group is ordered by the element index its name carries, never by buffer position, because simplification reorders variables freely.

group_name(name) -> (base, component)

Split a scalarized slot name into the vector it belongs to and its element index: :pos_2 → (:pos, 2). One trailing _<digits> is stripped, which is exactly what scalarizing a vector variable appends. A name without such a suffix is returned unchanged with element index 0.

slots(map, name) -> Vector{Int}

The buffer positions name occupies: the whole vector when name is a group (:pos → the slots of pos_1 … pos_n), otherwise the single named slot.

Whether map has any slot under name, as a slot or as a group.

KernelReads(output_input, output_state, dstate_input, dstate_state)

Which of a kernel's inputs and states each of its outputs, and each of its state derivatives, actually reads, as slot lists. This is the component's own Jacobian sparsity; walked out through the wiring it gives the model's (state_dependencies), and collapsed to one bool per input it gives the ComponentKernel's input_feeds_output.

Whether each of count inputs is read by any output of a kernel with reads.

ComponentKernel

One compiled component type. All three maps take (target, u, input, numeric, callables, instances, batch, t) and run over every instance named in batch: rhs! integrates the component's own states and is nothing for a stateless component, out! writes its declared outputs and obs! its remaining observed variables. The SlotMaps name each buffer's positions, input_feeds_output marks the inputs out! reads — all the schedule needs to know about the component's internal dependencies — and reads resolves that per slot.

compile_kernel(system, inputs, outputs; name, verbose=false)

Compile one component System into a ComponentKernel. inputs and outputs name variables of system; array-valued ones are scalarized, so an output named :pos becomes the slots pos_1 … pos_3. Any parameter without a build-time default is an error — the kernel's parameter buffer is seeded from those defaults and only then overwritten per instance.

required_default(kernel, param, value)

value as a SimFloat, or an error naming the parameter that has none. Every parameter needs a build-time default: the kernel's buffer is seeded from them and only then overwritten per instance, so a missing one would leave a silent zero.

ComponentInstance

One occurrence of a kernel. position counts the instance among its own kernel's, which is how it finds its callable parameters; every other field is this instance's contiguous slice of the corresponding global buffer — states of u/du, params of the numeric parameter vector, and so on.

Wiring(sources, weights, offsets)

Flattened compressed-row map from output slots to input slots: input slot k is the weighted sum of the output-buffer slots sources[offsets[k]:offsets[k + 1] - 1]. Most weights are 1; the exception is a weighted reference point — a wing frame fitted from a blend of structural points — which therefore needs no component of its own. An input with no sources reads zero.

GatherPlan(copy_targets, copy_sources, general)

One gather pass, split by how much work each slot needs. A slot fed by a single unit-weighted source is a plain copy, which nearly all of them are — on a replicated four-kite model the wiring averages 1.45 sources per slot, so the loop setup guarding the sum costs more than the sum. The rest keep the weighted form in general. A slot with no source appears in neither list: only a gather writes an input, so it keeps the zero its buffer was built with and never needs touching.

gather_plan(wiring, slots) -> GatherPlan

Sort slots into the copies, the weighted sums and the untouched, once, so that gather! does not rediscover the shape of the wiring on every call.

gather!(input, output, wiring, plan)

Sum every wired output slot into the input slots plan covers. This is the whole coupling model: a segment writing force into a point, a ride point forwarding a moment into its body and a winch driving a tether length all go through it, and because it indexes by slot there is no uniform I/O width to pad to. A plan holds what the next pass is about to read rather than the whole buffer, and each slot belongs to exactly one plan — the schedule orders every producer before the earliest consumer, so a slot that two layers read is already final when the first of them runs.

KernelSystem

An assembled model. kernels is a vector, not a tuple: a concrete tuple names every kernel type in the model, and carrying that type on a field here makes inference of the code reaching a KernelSystem grow explosively in the number of kernel types. The evaluation loop still needs a tuple to unroll over, and gets it from KernelRHS.

A batch is the instance list of one kernel: layers[l][k] are the instances of kernel k whose outputs run in layer l, and state_batches[k] the stateful instances of kernel k. layer_inputs[l] is the GatherPlan that layer needs and final_inputs what the derivative and observable passes need on top, so each gather! touches only the slots about to be used and no slot is gathered twice in a call. The _active lists hold the positions of the non-empty batches of each pass, which is what run_batches! walks, skipping the kernels idle in that layer.

SystemBuilder()

Accumulates kernels, instances and connections, then build_systems them into a KernelSystem. Buffer slices are handed out as instances are added, so a connection can be recorded as soon as both endpoints exist.

add_kernel!(builder, kernel, seconds) -> Int

Register a compiled component type and return its index in the kernel table. seconds is what compiling it cost, kept for the build-time breakdown.

add_instance!(builder, kernel_idx) -> Int

Add one occurrence of a registered kernel, reserving its slice of every buffer. Returns the instance index.

Reserve count slots at the end of the buffer counted by field.

connect!(builder, source, out_name, target, in_name; weight=1)

Add weight times instance source's output out_name to instance target's input in_name, slot by slot. Both names may be whole vectors (:pos covering pos_1 … pos_3) and must then have the same width. Several sources may feed one input; they sum.

build_system(builder) -> KernelSystem

Turn the recorded connections into a Wiring and a layered schedule.

A slot a layer gathers is left out of final_inputs: the schedule puts every producer of a feeding input before the layer that reads it, so such a slot is already final when its layer gathered it and the derivative pass can read it as is. What remains for the final gather is the inputs no output map reads — a point's aggregated force and its like — which is where the coupling comes back down.

build_wiring(targets, sources, weights, n_inputs) -> Wiring

Invert the recorded (target input slot, source output slot, weight) triples into the compressed-row form gather! walks.

build_schedule(builder, wiring) -> Vector{NTuple{N, Vector{Int}}}

Layer the instances so that every output map runs after the maps it reads. An instance depends on another when one of the inputs its own outputs read is fed by that instance's outputs; inputs the outputs ignore (a point's aggregated force, say) impose no order, which is what lets the force flow back down the same chain. Errors on a cycle, naming the kernels involved.

ScheduleCycleError(kernels)

Thrown when the output-dependency graph has a cycle — a genuine algebraic loop between components. kernels names the component types still unscheduled.

state_sparsity(builder, wiring, layers) -> SparseMatrixCSC{Bool, Int}

The Jacobian's nonzero pattern: entry (i, j) is set when du[i] reads u[j]. Handing it to the integrator turns a dense finite-difference Jacobian — O(n) right- hand side evaluations and a dense factorization — into a coloured sparse one. The diagonal is always stored, read or not, because the solver's W = M/γ - J writes there and can only write where the pattern has room.

state_dependencies(builder, wiring, layers) -> Vector{Vector{Int}}

For each state, the states its derivative reads: what its own instance reads directly, plus the output_reach of everything wired into the inputs it reads.

output_reach(model, wiring, layers) -> Vector{BitSet}

For each output slot, the states it depends on. An output reaches its instance's own states plus everything reaching the outputs wired into the inputs it reads, so one sweep in schedule order resolves the model; the sweep repeats until nothing grows, which costs one confirming pass and converges an instance that feeds itself. model is a SystemBuilder while assembling and a KernelSystem after, and only its instances and kernels are read.

Grow one instance's outputs' reachable states; true if any of them grew.

Union into into the states inst reads directly, plus those its inputs reach.

global_mass_matrix(builder)

Assemble the state mass matrix from the kernels' own. Returns I when every component is a plain ODE, which is the case for all of them today.

KernelParams(numeric, callables)

The parameter object the assembled ODEProblem carries. numeric is the flat buffer the generated kernels index — its layout is ours, (instance, slot) → flat index, so the same field read on two instances can never collide. callables holds the polars and rigidity laws, which cannot live in a numeric buffer: one entry per kernel, each a vector holding one instance's callables per element as a tuple, so every slot keeps its own type and calling a polar dispatches statically.

CallableSlot(store, position, slot)

One instance's callable at slot, as an address that can be written. store is its kernel's per-instance callable tuples and position the instance's index among them. The write goes through Base.setindex because the tuple is immutable — and it is a tuple rather than a vector so that each slot keeps its own concrete type.

Put value in the slot target addresses.

CallableGroup(targets, readers)

The callable readers sharing one concrete type with the CallableSlots they fill, the callable counterpart of ReaderGroup and there for the same reason: one dispatch per group rather than one per slot.

group_callables(targets, readers) -> Vector{CallableGroup}

Split the parallel targets/readers vectors into one CallableGroup per (target, reader) type pair, keeping the original order within each group.

Write every callable of group into the slot it addresses.

KernelParamSync(groups, callable_groups)

Copies live SystemStructure fields into a KernelParams: the ReaderGroups supply numeric, the CallableGroups the callables. This is the KernelBackend's parameter sync, applied every step through the same ProbWithAttributes machinery as the monolith's, so it is on the hot path twice per step and is grouped by type to stay allocation-free.

KernelInitialSync(model)

Pushes the struct's initial conditions onto a problem's u0, the KernelBackend's counterpart of the monolith's InitialSync. It exists for the same reason: a serialized problem carries the state its build-time struct had, and the model hash is structural, so the same bin is reused for another set of positions and rest lengths.

KernelBuffers(input, output, observable)

The scratch buffers one evaluation needs, at one element type. Cached per eltype(u) so a ForwardDiff dual pass gets its own set.

KernelRHS(system)

The callable (du, u, p, t) an ODEProblem integrates. Each schedule layer gathers the inputs and runs its batches' output maps; once every layer has run the inputs are complete and the stateful kernels write their derivatives.

kernels is the system's kernels as a tuple, so the loop unrolls over them and every kernel call is statically dispatched. It lives here rather than on the KernelSystem because only this type is on the evaluation path; the system is what assembly, parameter binding and the state getter all carry, and the tuple's type on it is what made their inference explode.

Scratch buffers for element type T, allocated on first use. The solver's own element type is a field rather than a lookup, so the hot path allocates nothing.

settled_slots(builder, wiring, layers, layer_slots) -> Set{Int}

The slots a layer gathered whose value cannot change afterwards, so the final gather can leave them alone. Almost every feeding input qualifies — the schedule orders a producer before the layer that reads its output. The exception is an instance whose own output feeds its own input, which build_schedule deliberately does not order against itself; that slot is only complete once its layer has run, so it stays in final_inputs.

layer_gather_slots(builder, layers) -> Vector{Vector{Int}}

The input slots each layer has to gather before its output maps run: those the maps actually read, minus the ones an earlier layer already gathered. An input its instance's outputs ignore — a point's aggregated force, say — is left to final_inputs; and since every producer is ordered before the earliest consumer, a slot two layers read is already final when the first of them runs.

output_read_slots(builder, instances) -> Vector{Int}

The input slots instances read in their output maps, ascending. read_slots is the same question for the instance as a whole.

read_slots(builder, instances) -> Vector{Int}

The input slots instances read, in ascending order so the gather walks the buffer and the wiring forwards.

Gather and run every schedule layer, leaving the inputs the derivative and observable passes read complete.

run_batches!(select, kernels, active, batches, callables, instances, input,
             target, u, numeric, t)

Run the map select picks out of every kernel that has work, one batch at a time. active holds the positions of the non-empty batches, so a kernel idle in this pass costs nothing; walking the whole kernel tuple instead costs a call frame per kernel per pass.

The buffers arrive as arrays rather than as the KernelSystem, KernelBuffers and KernelParams holding them: reaching through those structs inside the loop makes every instance reload the field (a store into target may alias the struct) and puts a huge type on a call the compiler then declines to inline.

dispatch_batch!(select, kernels, callables, batches, k, instances, input, target,
                u, numeric, t)

Run batch k with the types its entries have. The kernels differ in type, so indexing the tuple at a position only known at run time would erase that type and copy the entry into a box; a ladder of comparisons against literal positions keeps every call concrete instead, at a handful of integer compares per batch.

produces_output(builder, i) -> Bool

Whether instance i's kernel declares an output. One that declares none — a wing's aerodynamic sum, whose results are all observables — has nothing for the output pass to write, so it is left out of the schedule and only its observed map runs.

Split instances into one list per kernel, in kernel-table order.

active_batches(batches) -> Vector{Int}

The positions of the batches that hold an instance, in order, so a pass skips the kernels with no work in it instead of calling into each of them to find none.

The map writing a kernel's declared outputs.

The map writing a kernel's state derivatives.

The map writing a kernel's remaining observed variables.

refresh_outputs!(rhs, u, p, t) -> KernelBuffers

Re-run every output map and every observed map for the state u, and return the buffers. The state getter reads component results out of these after a step, where the integrator's last RHS evaluation need not correspond to u.

ScaledReader(reader, factor)

Reads reader and multiplies by factor. Used where a parameter is a fixed share of a struct field — an unwinched tether segment's rest length is its tether's length over the segment count.

kernel_param_slots(kernel) -> Dict{Symbol, Vector{Int}}

Map each of kernel's numeric parameter names to every slot it occupies, so a registry entry can be matched to the slots it survived as. Names, not symbolics: a parameter a nested component declared comes back namespaced (aero₊p_wings_1_aero_v_ind_1_2) while the registry holds the bare symbol it was created under, so the two are different objects for the same field.

name_slots(names) -> Dict{Symbol, Vector{Int}}

name => slots for one buffer's names, where a namespaced name is registered under its leaf (the part after the last ₊) as well as in full. One field can occupy several slots: a component that both reads a field itself and passes it to a nested subsystem declares it twice, once bare and once namespaced, and only the namespaced copy is the one the subsystem's equations read. Writing every match is right rather than lucky — the names agree because the field is the same, so the value is too.

The slots a parameter symbolic occupies in a name_slots map, by the name it scalarizes to; empty when the kernel does not carry it.

No slots — the shared empty result of matched_slots.

The part of a namespaced name after the last ₊, or the name itself.

kernel_callable_slots(kernel) -> Dict{Symbol, Vector{Int}}

As kernel_param_slots, for the kernel's callable parameters.

instance_readers(kernel, registry, index_map)

The live readers for one instance of kernel. Returns (numeric, callable), each a vector of (slot, reader) pairs: reader(sys_struct) is the current value of the parameter at that slot of the instance's slice. index_map gives this instance's index in each container the kernel reads (Dict(:points => 7)); parameters the kernel minted itself (a pulley's damping, say) are not in the registry and keep their build-time default.

entry_reader(entry, index_map, element)

The reader for one registry entry on the instance described by index_map: the recorded path with its container index swapped, plus element when the entry is an array whose components occupy separate slots. A computed entry (not a plain field read) is instance-independent and is reused as it stands.

remap_path(path, index_map) -> Tuple

path with the index following each container named in index_map swapped for that instance's. The container may sit at any depth, so a panel addressed as wings[1].aero.panels[3].v_ind is remapped on panels while wings stays put.

record_build_index!(built_from, entry, index_map)

Record which component of each remapped container the kernel was built from, and error on a second one. Two components of one container in a single kernel would both be swapped to the same instance index, silently reading the wrong one; the components are written so this cannot happen, and this says so if that changes.

KernelJacobian

jac(J, u, p, t) for a KernelSystem, composed from per-kernel local Jacobians and the constant wiring instead of differentiated globally.

Per instance, rows holds its outputs' Y block over the state columns cols, and nzindex maps its derivative block over dstate_cols onto the nonzeros of matrix. Both column sets are the instance's own reach, so every block is small and dense and the substitution is a sweep of little matrix products rather than sparse algebra over the whole model.

blocks and shape restate what the KernelDuals hold, at concrete types: a workspace's element type carries its dual width, so a vector of them is abstract and reaching through one inside the sweep dispatches dynamically. The workspaces are touched once per kernel and the blocks once per instance, so only the latter has to be concrete.

Each component is differentiated numerically, so a registered leaf without a symbolic derivative — the wind profile, an aerodynamic polar — is differentiated by running it on duals, where the monolith's symbolic Jacobian fails.

build_jacobian(rhs; prn=true) -> KernelJacobian or nothing

Plan the analytical Jacobian of rhs: the per-kernel dual workspaces, each instance's column support and input gather, and the map from its derivative block onto the nonzeros of the assembled matrix.

Returns nothing, with a warning, when an instance's own outputs feed an input its outputs read. The schedule deliberately does not order such an instance against itself, so its block is an algebraic loop that one sweep would not settle; those models keep the solver's own Jacobian.

KernelDuals{W, T}

One kernel's forward-mode workspace, at that kernel's own width W = states + inputs. Every instance of a kernel owns a disjoint slice of u, of input and of output, so seeding partial j on local slot j of every instance at once returns all of their blocks from a single pass without the seeds mixing — a compressed differentiation whose colouring is the local slot index.

The buffers are packed, one instance after another, and instances is the kernel's instance list rebased onto them, so no kernel holds a set of duals over the scattered global buffers. blocks receives the local Jacobians as (row, slot, position), the rows being the kernel's outputs followed by its state derivatives.

The slice packed buffers give the position-th instance of a width-slot buffer.

fill_blocks!(duals, kernel, system, u, input, numeric, callables, t)

Run one kernel's output and derivative maps over duals seeded on its own slots, leaving every instance's local Jacobian in duals.blocks. One pass per kernel type, at that kernel's width, and nothing per model state.

InputGather(slots, offsets, producer, row, weight)

Where one instance's input rows of G·Y come from, in compressed-row form: the kth gathered row is local input slot slots[k], and sums weight times row row of the Y block of instance producer over offsets[k]:offsets[k + 1] - 1. The same wiring gather! walks for values, resolved to producing instances so the substitution never searches for one.

slots holds only the inputs the pass in hand reads — those the output maps read for the sweep, those the derivatives read for the composition — since an input the pass ignores contributes a row of zeros.

self_feeding_instances(system) -> Vector{Int}

The instances whose own outputs feed an input their output map reads. build_schedule skips such an edge, so the schedule does not order the instance against itself and one sweep of Y = A + B·G·Y would leave that block short of its own contribution.

output_owners(system) -> (owner, local_row)

Per output slot, the instance that writes it and the slot's position within that instance's outputs — which is the row of its Y block.

instance_gather(system, owner, local_row, index, which) -> InputGather

Resolve one instance's wiring into the producing instances and rows the substitution reads, one compressed row per input slot its output maps (which = :output) or its state derivatives (which = :dstate) read.

instance_columns(system, reach, index, which) -> Vector{Int}

The state columns one instance's blocks span: its own states plus everything reaching the inputs its outputs read (which = :output) or its state derivatives read (which = :dstate). One column set per instance rather than per row keeps each block dense, and the rows that span less of it hold structural zeros.

nonzero_indices(matrix, rows, columns) -> Matrix{Int}

Where each entry of one instance's derivative block lands in matrix's value array, 0 for an entry the pattern does not hold. Those entries are structurally zero: a derivative reads a column only through an input whose reach the pattern already covers, so skipping them writes nothing that is not zero.

settle_outputs!(jac, system, index)

Fill one instance's Y block: gather G·Y over the rows its outputs read, multiply by the instance's own B, and add A on its own state columns. Called in schedule order, so every producer it reads is already settled.

block_product!(target, blocks, row_offset, column_offset, position, n_row, slots,
               gathered, n_column)

target = B · gathered, where row k of gathered is input slot slots[k] and B is the corresponding window of one instance's local Jacobian at (row_offset, column_offset). Written out rather than left to mul! because the blocks are a handful of rows wide: a BLAS call costs more than the product, and the loop can skip the zeros gathered is full of.

compose_derivatives!(jac, system, index)

Add one instance's rows of J: C on its own state columns plus D times the gathered G·Y over the inputs its derivatives read, scattered onto the assembled matrix's nonzeros.

gather_blocks!(jac, entries, columns) -> Matrix

G·Y for one instance, in the leading rows × columns corner of the shared scratch, one row per entry of entries.slots. Marks columns in column_of so a producer's own column set can be mapped onto this one by lookup; a producer column this instance does not span carries a structural zero in the row being read, so dropping it loses nothing. The caller clears the marks.

Undo gather_blocks!'s marks, so column_of stays all zero between uses.

store_jacobian!(target, assembled)

Copy the assembled sparse Jacobian into whatever matrix the solver handed us: its value array directly when it carries the same pattern, entry by entry otherwise.

Write every stored entry of assembled into target by index.

point_variables()

The variables a point component may declare, as a named tuple: its world pos/vel outputs, the summed force_in/mass_in of its incident segments, their drag_in share, and the observed total_drag, wind_vec and net_force. Each component lists only the ones it uses — a static point has no net_force, carrying no force input to build one from.

point_wind_eqs(s, params, idx, io)

Bind wind_vec to the point's point_wind_source at its own height — the monolith's wind_at_point, scattered into point.wind_vec — and total_drag to the point's own aerodynamic drag against it plus the share its segments deliver, which the state getter scatters into point.drag_force.

Particle(s, params, idx; name)

A free point mass: integrates pos/vel under the force gathered from its incident segments, its own drag, gravity and world-frame damping, through the shared point_acceleration. Its translational mass is extra_mass plus the incident segments' half-masses (mass_in).

Anchor(s, params, idx; name)

A ground-anchored point: pos is pinned to params.points[idx].pos_w and vel is zero, matching the STATIC branch of point_eqs!. Its segments' force and mass still arrive and are still summed into the observed net_force — nothing moves in response, but that is how the load the anchor carries is read off.

pulley_rope_mass(params, pulley_idx, segment_idx)

The rope mass sum_len · ρ · π (d/2)² driving a pulley's rope split, built in equation from the pulley's sum_len and one of its segments' material, so it follows the struct live.

PulleyParticle(s, params, idx, pulley_idx, segment_idx; name)

A Particle that also owns a pulley's rope split: D(pulley_len) = pulley_vel and D(pulley_vel) = tension_in / rope_mass − damping · pulley_vel, with tension_in the imbalance spring[seg1] − spring[seg2] its two segments deliver and line_in their mean, which the friction scales with. pulley_len_out exposes the split so those segments read it as their rest length.

WinchAnchor(s, params, winch, idx; name)

A reeling winch at the STATIC point idx. It owns the motor speed winch_vel and one tether_len_k state per connected tether, each an output its tether's segments read as their rest length, alongside a tether_vel_k output carrying that length's rate for their dampers. The motor law is winch_component reused verbatim; it reads the tension gathered at the winch point (the vector sum of the spring forces there, as winch_eqs! forms it), the mean tether length, the control set_value and the brake, and returns the drum acceleration.

segment_variables()

The variables every segment kernel declares, as a named tuple: the two endpoints' pos/vel inputs, the force on each endpoint (positive force-on-point sign), the half_mass and half_drag each endpoint carries, and the spring_force/len/l0 diagnostics. Mass and drag are one output each because both endpoints receive the same value; the wiring delivers it twice.

segment_eqs(s, params, idx, io, rest_len; with_drag=true, rest_len_rate=0.0)

The equations every segment kernel shares: the shared segment_load_terms evaluated at rest_len, bound to the segment_variables outputs and diagnostics. Returns (eqs, loads, wind_params); loads carries the scalar and vector spring tension a pulley or tether segment emits on top of these and wind_params the parameters segment_wind_params minted, which the kernel must declare. rest_len_rate is d(rest_len)/dt for the kernels whose rest length is driven by another component's state.

SpringSegment(s, params, idx; name, with_drag=true)

A spring-damper segment at a fixed rest length (params.segments[idx].l0). With with_drag = false it is the drag-free wing-structural link — a distinct compiled type rather than a zeroed drag coefficient.

PulleySegment(s, params, idx, pulley_idx; name)

One of a pulley's two segments. Its rest length comes from the pulley point's pulley_len_out, read as rest_len: the first segment takes it directly, the second sum_len − rest_len, selected by the assembly-set pulley_side (±1). The matching pulley_vel_out arrives as rest_vel and gives the damper the same sign, so the split it drives also damps it. It emits pulley_side · spring as tension, so the pulley point aggregates spring[seg1] − spring[seg2], and half its spring force as line, so the same point aggregates the mean of the two leg tensions for the friction to scale with.

TetherSegment(s, params, idx; name)

A segment of a winched tether. Its rest length is the winch's tether_len_k output divided by the assembly-set segment_count, and its rate the matching tether_vel_k over the same count, so reeling relaxes the segment instead of straining it. Both endpoints' spring force are emitted as src_tension/dst_tension; the assembly wires only the one at the winch point, which is how the winch sees the force winch_eqs! sums there.

body_variables()

The variables a rigid-body component declares, as a named tuple: the aggregated force_in/moment_in at and about its COM, the pose outputs a point riding the body reads (pos, vel, frame — R_b_to_w column-major — com, com_velocity, omega_w), and the observed acc/orientation/omega_b/alpha_b the state getter scatters into the struct or feeds to the reported scalars.

vector with its component along the unit axis removed.

Only vector's component along the unit axis.

body_integration(params, idx, com_w, com_vel, omega_p, alpha_p, com_acc,
                 orientation_p; frozen=false, wing_frame=nothing,
                 wing_vel=nothing)

The four integration overrides rigid_body_pose_expressions takes. fix_sphere confines the body to a sphere about the world origin by keeping only the radial part of its COM velocity and acceleration and dropping the radial part of its spin. alpha_p/com_acc are the caller's torn variables, so the overrides can name the accelerations they correct without a cycle.

Damping is folded in here so both backends share one definition: angular_damping on the absolute spin, world_frame_damping on the COM velocity on the world axes, and body_frame_damping through the same body_frame_damp_accel a wing node uses — the velocity relative to the parent wing, resolved on the wing's axes, so it damps deformation and not rigid flight. wing_frame/wing_vel carry that parent's frame and velocity; a body with no parent wing (wing_idx == 0) gets its own frame and a resting parent.

frozen holds all four derivatives at zero, clamping a body that is integrated but must not move; a backend that gives such a body no state leaves it false. A body not frozen at build time reads the fix_static parameter instead, the same clamp under runtime control.

RigidBody(s, params, idx; name)

Free 6-DOF rigid body: integrates the 13-state principal pose (com_w, com_vel, Q_p_to_w, ω_p) under gravity, the wrench gathered at force_in/moment_in, the external wrench (ext_force_w/ext_force_b/ext_moment_b) and its damping, through rigid_body_pose_expressions. fix_sphere stays a parameter. A clamped body is StaticBody instead.

parented declares the wing_frame/wing_velocity inputs a body with a parent wing needs, so its body_frame_damping resists motion relative to that wing rather than absolute motion; the assembly connects them from the parent instance.

StaticBody(s, params, idx; name)

Clamped (STATIC) rigid body: no state at all, just the fixed pose params.bodies[idx] already holds, with zero velocity and spin. The monolith freezes such a body by holding its thirteen states at their initial values; here there is nothing to hold, so the pose is emitted directly. Its force_in/moment_in are declared, so joints and ride points may deliver to it, and ignored.

wing_frame_variables()

The reference-point inputs a fitted (KINEMATIC) wing body reads: the two z-axis and two y-axis structural reference positions its frame is built from, plus the origin's position and velocity. Each is a weighted blend of real points, which the wiring supplies — see Wiring.

KinematicBody(s, params, idx; name)

A KINEMATIC (particle-wing) body: no state at all. Its orientation is fitted from four structural reference points through the shared wing_frame_columns, and its origin pose is read from its origin reference point, exactly as the KINEMATIC branch of wing_eqs! does. The principal frame is aliased to the body frame and the spin is zero, because such a wing has no rigid rotation of its own — the particles carry the motion.

body_frame_damp_accel(vel, body_damp, orientation, wing_vel)

Body-frame damping acceleration R·(coeff ⊙ (Rᵀ·(vel − wing_vel))), the term point_damping_accel adds for a point that damps against its wing's frame rather than the world.

WingNodePoint(s, params, idx; name, with_damping=true)

A particle belonging to a fitted wing: the shared Particle motion plus the body-frame damping point_eqs! adds for a point that damps against its wing's frame rather than the world, which reads the wing_frame/wing_velocity the fitted body supplies. Its aerodynamic force arrives at force_in from ParticleWingAero, like any other load. with_damping = false drops the damping term whose coefficient the struct leaves unset, so no zero-valued parameter is generated, matching the monolith.

ParticleWingAero(s, params, idx; name)

The live aerodynamic force on a fitted wing's structural points: every point's world pos/vel and the wing's own pose in, the world force on each point out. It builds each point's apparent wind and height from the pose, feeds the wing's aero component through particle_wing_aero_wiring and rotates the per-point force it returns back into the world. A frozen mode compiles down to one parameter read per point; a live mode solves in place. The wing's lumped aero_force_b/aero_moment_b are observed.

WagnerLag(s, params, wing_idx; name)

The wing's two-state Wagner lift lag as one component, holding the states the whole wing shares and handing every panel the angle-of-attack deficiency to subtract. Its va_in is the wing's mean body-frame apparent wind, gathered over the wing's nodes by the Wiring exactly as AeroInflow gathers a panel group's. The physics is the shared wagner_lag_eqs, so this emits what a whole-wing system emits.

AeroInflowPoint(s, params, idx; name)

What one structural point contributes to its wing's aerodynamics: its world pos and vel and the wing's pose in; its body-frame position, apparent wind and air density out. These are exactly the per-point quantities the PARTICLE_DYNAMICS branch of aero_eqs! builds, and the only thing that distinguishes point idx from any other is the parameter its point_wind_source reads — remapped per instance — so one compiled kernel serves every aerodynamic point of every wing.

AeroInflow(; name)

The mean apparent wind and density of one group of aerodynamic points, materialized so the panels that share it can read it as an output. The averaging itself is the weighted gather in the Wiring — this only turns the gathered input into an output, which is what lets a wing-wide mean be summed once rather than once per panel.

AeroPanel(s, params, wing_idx, panel_idx, orient; name, with_flap)

One refined VSM panel's aerodynamic load: its two sections' leading and trailing edges (each already the gathered strut interpolation), their apparent wind and density, and its flap deflection in; its body-frame force and the couple its mode's scatter places (scatter_couple) out. The physics is the shared panel_force_eqs on a single column, so the expressions are those a whole-wing system emits for this panel. orient is the panel's ±1 span sign, baked in because it costs a second kernel and saves a parameter on every instance; the chord blend weight cannot be, because it differs per panel and would cost a kernel each. with_flap selects the (α, δ) polars.

A wing with flow_curvature_enabled takes two more inputs, its sections' trailing minus leading edge apparent wind, gathered at strut_pitch_weights. A wing with wagner_enabled takes one more, the lag deficiency its WagnerLag hands to every panel.

AeroPointForce(s, params, wing_idx, point_idx; name)

The aerodynamic force on one wing point: its body-frame position, the wing's frame and the panels' gathered body-frame load in; the world force it delivers to the point out, plus the body-frame force and its moment about the wing origin for the wing's lumped wrench. Whatever constant the mode adds on top of the panels — a frozen surface traction — enters here as aero_point_offset.

WingAeroSum(; name)

A wing's lumped aerodynamic wrench: the body-frame force and moment its points deliver, gathered and observed as aero_force_b and aero_moment_b. It carries nothing else — a PARTICLE_DYNAMICS wing's loads act on its points, so the wrench is a readout, and summing it over points rather than panels is what keeps it comparable with ParticleWingAero.

WingAero(s, params, idx; name)

The aerodynamic wrench of a RIGID_DYNAMICS wing: its body's pose in, the world force and the moment about its COM out, plus one twist moment per station it carries. It builds the wing's apparent wind from that pose, feeds the wing's aero component through rigid_wing_aero_wiring, and transports the returned body-frame wrench to the COM as create_sys does. Like ParticleWingAero it is a component of its own: a body's pose does not depend on the force it receives, so aero after pose and before the body's derivative is just another schedule layer.

StationDOF(s, params, idx; name)

The added twist degree of freedom of a DYNAMIC station: a thin plate hinged at its leading edge, driven by the aerodynamic moment its wing's aero returns and the bridle couple its points deliver, restrained by the surface's own stiffness and damping. Its inertia ⅓·m·L² takes the mass from those same points as an input, so the component reads only its own surface's parameters. The monolith's fix_wing freeze is not carried over: it is a parameter nothing ever sets.

station_diagnostics()

The three quantities get_all_state copies out of a station beyond its twist: the bridle couple's tether_force and tether_moment about the hinge, and the aero_moment its wing's aero returns. Observed, never read by any equation.

PrescribedTwist(s, params, idx; name)

A STATIC station's prescribed section twist: no state and no inputs, just the twist its parameters hold. It exists so a node reading a twist angle reads one whether its surface twists dynamically (StationDOF) or not.

TwistNodePoint(s, params, idx; name, surface_idx=0)

The kinematics of a structural node on a RIGID_DYNAMICS wing: its body's pose and its surface's twist in, its world pos, its moment arm about the body COM and its height out. Such a node is placed from the body's COM by a twist-deformed body-frame offset (twist_deformed_offset) rather than by a fixed anchor, and point_eqs! gives it no velocity of its own, which is why it is not a RidePoint.

TwistNodeWrench(s, params, idx; name, surface_idx=0, gated=false)

The statics of a structural node on a RIGID_DYNAMICS wing: the load its segments, its own drag and its external force deliver — no gravity, because the wing body already carries the node's mass — the moment that load makes about the body COM and, when the node belongs to a station, the bridle couple it exerts on that surface's hinge and the mass it lends to the surface's inertia. gated is the wing's group_points_moment = false, which drops an in-surface node's moment on the body.

twist_deformed_offset(params, idx, surface_idx, angle)

The body-frame offset of a wing's structural node, rotated about its surface's leading edge by the section twist angle — the placement point_eqs! gives a wing node whose surface twists. Without a surface (surface_idx == 0) it is the node's own pos_undeformed_b.

A full rotation about the (unit) spanwise axis, so a node offset along the span keeps that offset: a station's nodes need not share one chordwise line, and dropping the axial term would shrink their spanwise spread by cos(angle).

twist_bridle_couple(surface, pos_b, twist, orientation)

The geometry of the couple a structural node at the body-frame offset pos_b exerts on its twist surface's hinge, at section twist twist and body orientation orientation (R_b_to_w). Returns (; axis, offset, arm, direction): axis along the surface chord, offset from the moment reference le_pos + moment_frac·chord to the node, its chordwise arm = offset ⋅ axis, and the world direction the twisted section normal points against. The node's load projected on direction, times arm, is the hinge moment it delivers. Used by TwistNodeWrench and station_eqs!.

rigid_body_point_velocity(pose, lever)

World velocity of a point rigidly attached to a body: com_velocity + ω × lever, where lever is the point's offset from that body's COM.

A scalar output variable named name.

FlapDelta(s, params, idx; name)

The live flap deflection of a flapped KINEMATIC station: its two flap bodies' orientations in, δ out, through the shared flap_delta_expression. Deflection is all such a surface contributes — it carries no twist DOF of its own — and the wing's aero component reads it per panel.

PointFlapDelta(s, params, idx; name)

The live deflection of a point-flap KINEMATIC station: its three flap points' positions and the owning wing's frame in, δ out, through the shared point_flap_delta_expression. The point counterpart of FlapDelta.

A flap's two frames, nine scalars each rather than two vectors: a hinge axis leaves some entries unread, and mtkcompile will not take part of a declared array.

A point flap's three chord positions and the wing frame its hinge axis is in.

flap_delta_expression(station, R_main, R_flap)

The signed live deflection δ of a flapped station: the angle between its two flap bodies' reference chords about the world hinge axis, referenced to rest. The axis, the reference chords and the rest angle are frozen rest geometry, so δ is a function of the two bodies' orientations alone.

point_flap_delta_expression(station, fore, hinge, aft, R_wing)

Deflection δ [rad] of a point flap from three structural positions: the signed angle the aft segment hinge→aft makes with the fore segment fore→hinge about the hinge axis, less the rest angle the CAD pose holds. Positive δ is a trailing edge deflected down, the sign the polars are tabulated on.

R_wing takes the axis from the wing's frame to the frame the positions are in, so the same expression reads world positions at run time and CAD positions at build. Nothing here is a body: a chord bending over several beam elements still reads a deflection, and it is read off the very points the aerodynamics is built on rather than off a pair of orientations standing in for a hinge.

hinge_angle(vec_a, vec_b, axis)

Signed angle [rad] from vec_a to vec_b about the unit axis, measured between their components in the plane normal to it. An atan, so it holds over the full polar δ range rather than only near zero. Shared by every form of the flap deflection — body pair or point triple, Julia or symbolic — so they cannot drift apart.

The projections are never built: the axial parts drop out of the cross product's own axial component, and out of the dot product once their product is subtracted. The two forms agree exactly, and this one is a third of the expression the RHS would otherwise carry per station.

flap_delta_eqs(wing, subsys, delta_of) -> Vector{Equation}

Bind every per-panel delta connector of wing's aero component to the deflection delta_of returns for the station that panel deflects with. A panel mapped to no surface is bound to zero. Empty when the aero exposes no delta connector.

flap_delta_inputs(wing, subsys) -> (; vars, eqs)

One flap_delta_g input per station some panel of wing maps to, and the flap_delta_eqs binding the aero component's per-panel delta connector to them. Empty when the wing's aero exposes no delta connector.

wing_flap_surfaces(wing) -> Vector{Int}

The stations wing's aero panels deflect with, sorted and unique, or empty when its aero mode has no per-panel flap coupling. One flap_delta input and one FlapDelta instance exist per entry. The 0 of an unmapped panel is not a surface and is left out.

station_aero_driven(station) -> Bool

Whether a wing's aero drives this station's hinge moment. A STATIC surface with no aero sections has nothing to drive it with, and station_eqs! binds its moment to zero instead.

particle_wing_aero_wiring(s, subsys; orientation, origin, positions,
                          velocities, apparent_winds, heights)
    -> (; eqs, force_b, moment_b)

Feed a PARTICLE_DYNAMICS wing's aero component: each point's world position and velocity rotated into the wing frame orientation about origin, its body-frame apparent wind, and the air_density at its height. force_b and moment_b are the body-frame wrench the component's per-point forces sum to, the moment taken about origin. apparent_winds and heights are supplied by the caller so each backend passes the variables it already carries rather than the expressions behind them; heights are raw world z, since air_density applies the ground clamp.

rigid_wing_aero_wiring(s, subsys, wing; apparent_wind_b, height, frame,
                       omega_b, twist_angles, twist_rates)
    -> (; eqs, force_b, moment_b, twist_moments)

Feed a RIGID_DYNAMICS wing's aero component: its body-frame apparent wind, the air_density at height, its body-to-world frame as a column-major nine-vector, its body-frame angular velocity, and one twist angle and rate per station it carries. Returns the body-frame wrench and one hinge moment per station, in the order of wing.station_idxs; which of those a backend binds is its own choice (station_aero_driven).

indexed_vector_variables(base, count; input=false) -> Vector

count three-component variables named base_1 … base_count, declared as inputs or as outputs. A component whose I/O width follows the model — a wing's aero reads one position per structural point — needs them named individually, because the wiring addresses base_k as one vector.

indexed_scalar_variables(base, count; input=false)

count scalar variables named base_1 … base_count, declared as inputs or as outputs. A slot map groups them under base exactly as it groups a scalarized vector, so the wiring still addresses them as one vector — but each is a variable in its own right, so mtkcompile may drop any one nothing reads. A declared array is all-or-nothing to it, and a component that reads only part of a matrix it is handed needs this: a flap hinge about [0, 1, 0] reads six of its nine frame entries.

A scalar input variable named name.

A length-n input variable named name.

body_pose_variables()

The pose inputs a component riding a rigid body reads off that body's outputs: the body origin pose_pos, its orientation pose_frame (R_b_to_w, column-major), its pose_com/pose_com_velocity and its world spin pose_omega.

RidePoint(s, params, idx; name)

The kinematics of a BODY_STATIC point: its body's pose in, and its world pos/vel, its moment arm (pos − com) and its height out, exactly as body_ride_eqs places it. It is deliberately only half of a ride point — the load it feeds back to the body is RideWrench — because the two halves sit on opposite sides of the same dependency chain: the position must be known before the incident segments can compute their forces, and the force is only known after. One component holding both would be a cycle; two are a four-layer schedule.

RideWrench(s, params, idx; name, with_gravity=true)

The statics of a BODY_STATIC point: everything that flows back to its body. It takes the point's own height/vel/arm from RidePoint and the load ride_load builds, and emits force_out and moment_out = arm × force_out into the body's force_in/moment_in. total_drag is observed, as for any other point.

HermiteRidePoint(s, params, idx; name)

The kinematics of a point riding a Timoshenko beam's deflected centerline: both end bodies' poses in, and the point's world pos/vel, its two moment arms (arm_a/arm_b, the offsets from each end body's COM) and its height out. The beam-anchored counterpart of RidePoint, and split for the same reason; the placement itself is beam_hermite_ride_expressions. The velocity is evaluated at the already-bound pos output, so the heavy element-frame subtree is built once.

HermiteRideWrench(s, params, idx; name)

The statics of a point riding a Timoshenko beam: the load ride_load builds, split along the element by the point's axial fraction beam_frac — (1 − s) onto the first end body and s onto the second — with each half's moment about that body's COM. The beam-anchored counterpart of RideWrench.

ride_wrench_variables()

The variables shared by every anchored point's statics half, as a named tuple: its own height and vel from the kinematics half, the force_in/mass_in/ drag_in its segments deliver, and the observed total_drag/wind_vec/net_force — the same three point_variables declares, so the struct reads the same fields back whether a point rides something or flies free. The moment arms differ per anchor kind and are declared by the component itself.

ride_load(s, params, idx, io; with_gravity) -> (; load, drag, wind)

The world load an anchored point delivers to whatever carries it: the force its segments deliver, its own aerodynamic drag at its height, its gravity and its external force — the monolith's point_force. with_gravity = false is the point that rides its own wing body, whose mass is already counted at that body's COM (rides_own_wing in point_eqs!). The drag and the wind at the point's own height come back separately because they are the other two quantities ride_wrench_eqs reports.

ride_wrench_eqs(io, ride) -> Vector{Equation}

The three observables every anchored point's statics half reports, whatever it rides: its total_drag (its own plus its segments' share), the wind_vec at its own height and its net_force. point_eqs! writes all three for every point, so binding them here is what lets the struct read the same fields back on either backend.

joint_variables()

Everything a body-to-body joint component declares: both end bodies' poses in (joint_pose_variables) and the restoring wrench on each out (joint_wrench_variables).

joint_pose_variables()

The two end bodies' poses as inputs (a_pos/a_frame/a_com/a_com_velocity/ a_omega and the b_… set), as a named tuple. Declared by everything that spans an element: the joints themselves and a point riding one.

joint_wrench_variables()

The wrench an element delivers to each of its two end bodies, as outputs.

joint_poses(io)

The two end bodies' poses in the argument form the shared joint wrench builders take, reading them off a joint_variables tuple. The world spin arrives directly as a_omega/b_omega, so unlike the monolith — which carries ω_b and rotates it — there is nothing to rotate here.

joint_wrench_eqs(io, ex)

Bind a joint component's four wrench outputs from a shared wrench builder's result ex, together with the torn equations it needs.

ElasticJointComponent(s, params, idx; name)

Lumped 6-DOF elastic joint between two bodies: reads both poses and emits the restoring wrench on each, through elastic_joint_wrench.

TimoshenkoJointComponent(s, params, idx; name)

Corotational Timoshenko beam element between two bodies: reads both poses and emits the restoring wrench on each, through timoshenko_element_wrench. The element frame, the two nodes' chord-relative rotations and the element forces are torn, so the shared frame subtree is built once.

PointRole(kind, pulley_idx, segment_idx, winch_idx, body_idx, joint_idx)

How one point is realised: kind is :particle, :anchor, :pulley, :winch, :wing_node, :ride or :hermite. A :pulley point carries the pulley it splits and one of that pulley's segments (whose material gives the rope mass); a :winch point carries its winch; a :ride point carries the body it is anchored to and a :hermite point the Timoshenko joint whose beam it rides.

classify_points(sys_struct) -> Vector{PointRole}

Decide each point's component type, following point_eqs!'s anchor rule: a point anchored to a beam rides its Timoshenko element, one anchored to a body rides that body, a point that splits a pulley is a pulley particle, a point carrying a winch is a winch anchor, and the rest follow their DynamicsType.

pulley_point_index(sys_struct, pulley) -> Int

The point a pulley splits its rope over: the single point its two segments share.

SegmentRole(kind, pulley_idx, pulley_side, tether_idx, segment_count, winch_point)

How one segment is realised: kind is :spring (fixed rest length), :structural (drag-free wing link), :pulley (rest length from a pulley split) or :tether (rest length from a winch). pulley_side is +1 for a pulley's first segment and −1 for its second; segment_count is its tether's segment count and winch_point the winch point it touches, or 0.

classify_segments(sys_struct) -> Vector{SegmentRole}

Decide each segment's component type. A pulley's two segments take their rest length from the pulley split; a winched tether's segments take theirs from the winch. An unwinched tether's rest length never changes, so its segments stay plain springs whose l0 parameter the assembler reads from the tether instead.

KernelEntry(index, registry, source)

A compiled component type: its position in the runtime's kernel table, the ParamRegistry its build filled, and the component index it was built from — the index every instance's parameter paths are remapped away from.

kernel!(builder, table, sam, key, source, make, inputs, outputs) -> KernelEntry

Return the kernel registered under key, compiling it on first use. make(params) builds the component System from a fresh parameter view, so every field the component reads is recorded against source and can be remapped to each instance.

callable_field_key(key, component, fields) -> Symbol

key narrowed by which of fields the component holds a callable in rather than a Real. A component's equations branch on that when they are built — a Real rigidity becomes a numeric parameter and a callable one a callable slot — so two components that differ there cannot share a kernel even though they are the same kind, and writing one's law into the other's slot fails on type.

The joint fields whose equations differ between a Real and a callable, so a kernel is shared only by joints that agree on them (callable_field_key).

KernelModel

An assembled model: the runtime system, its initial state and parameters, the parameter sync, and the instance index of every point, segment and body, which is all the state getter and the control setter need to find their values. A BODY_STATIC point has two: point_instances holds its kinematics and wrench_instances the load it feeds back (0 for every other point). aero_instances holds each wing's aero instance and twist_instances each twist surface's, 0 where there is none. aero_force_instances holds each point's AeroPointForce and inflow_instances its AeroInflowPoint, whose va_b output is the apparent wind the aero is solved on; 0 for a point with neither.

assemble(sam) -> KernelModel

Translate sam.sys_struct into a KernelModel. Kernels are compiled once per component type, one instance is added per component, the wiring is recorded, and the parameters and initial state are bound from the struct afterwards, when every instance's buffer slice is final.

add_point!(builder, table, bindings, sam, idx, role, bodies, wrenches, twists)
    -> Int

Add the instance for point idx and record which container indices its parameters must be remapped to. A :ride point becomes two instances — the kinematics, whose index is returned so segments wire to it, and the wrench, recorded in wrenches and wired to its body here.

add_segment!(builder, table, bindings, sam, idx, role) -> Int

Add the instance for segment idx and record its parameter remapping.

add_body!(builder, table, bindings, sam, idx) -> Int

Add the instance for body idx. Each of the three kinds is its own compiled type, because which one a body is, is topology: a DYNAMIC body integrates, a STATIC one is clamped, and a KINEMATIC one is fitted from reference points.

add_kinematic_body!(builder, table, bindings, sam, idx) -> Int

Add a KINEMATIC wing body. It has no state: its frame is fitted from four structural reference points and its origin pose read from a fifth, so the whole component is wiring. Each reference is a weighted blend of real points, delivered by the weights in the Wiring.

wire_kinematic_body!(builder, sys_struct, idx, body_instance, point_instances)

Wire a fitted wing's reference points into it: the four frame references and the origin's position and velocity, each a weighted blend of point outputs.

kinematic_wing_of(sys_struct, point) -> Int

The fitted (KINEMATIC) wing whose frame point needs, or 0. point_eqs! damps any DYNAMIC point against its wing's frame, not only aerodynamic surface nodes, so a steering or pulley point with body_frame_damping needs it too.

add_ride_point!(builder, table, bindings, sam, idx, role, bodies, wrenches)

Add the two instances a BODY_STATIC point needs and wire them to its body: the RidePoint reads the body's pose, the RideWrench reads the ride point's pos/vel/arm and delivers force and moment back. Returns the ride point, which is what the incident segments connect to.

add_hermite_ride_point!(builder, table, bindings, sam, idx, role, bodies, wrenches)

Add the two instances a beam-anchored point needs and wire them to the two end bodies of the Timoshenko element it rides: the HermiteRidePoint reads both poses, the HermiteRideWrench reads the ride point's pos/vel and two moment arms and delivers a share of its load to each body. Returns the ride point, which is what the incident segments connect to.

add_wing_aero!(builder, table, bindings, sam, wing, bodies, points, flaps,
               twists, wrenches) -> Int

Add a wing's aero instance and wire it, dispatching on how the wing carries its loads: a PARTICLE_DYNAMICS wing's aero delivers a force to each structural point (ParticleWingAero), a RIGID_DYNAMICS wing's delivers one wrench to its body (add_rigid_wing_aero!).

Either way the aero is a component of its own, not a term inside the points or the body: nothing about it is a cycle, so the schedule runs the structure, then the wing frame, then the aero, then the derivatives.

add_panel_wing_aero!(builder, table, bindings, sam, wing, bodies, points, flaps,
                     wrenches) -> Int

Add a PARTICLE_DYNAMICS wing's aerodynamics as small repeated components: one AeroInflowPoint per structural point, one AeroInflow per inflow group, one AeroPanel per refined panel, one AeroPointForce per point again, and one WingAeroSum for the readouts. Everything between them — the strut interpolation, the inflow average, the load scatter — is a constant weight, so it is wiring rather than equations. Returns the sum's instance, which is where the wing's readouts hang.

The one-component equivalent, ParticleWingAero, is superlinear in the wing's size; this decomposition is not, which is what makes a large wing buildable.

add_aero_inflow_points!(builder, table, bindings, sam, nodes, points, body)

One AeroInflowPoint per structural point of a wing, wired from that point's kinematics and its wing body's pose. The only per-component field the kernel reads is the point's prescribed wind, remapped per instance, so every point of every wing shares it. Returns the instances, indexed as nodes is.

add_aero_inflows!(builder, table, bindings, sam, groups, inflow_points) -> Vector{Int}

One AeroInflow per group of aero_inflow_groups, gathering its points' apparent wind and density at the group's weights.

add_aero_panels!(builder, table, bindings, sam, wing, nodes, inflow_points,
                 inflows, pitches, section_group, flaps) -> Vector{Int}

One AeroPanel per refined panel, reading its two sections' corners from the strut interpolation (aero_geometry_entries), their inflow from the group they belong to, and its flap deflection from the station it deflects with. Panels differ only in their ±1 span sign, so a wing needs at most two kernels.

pitches are the aero_pitch_groups gathers, or nothing on a wing without flow_curvature_enabled, which then has no such inputs to connect. wagner is the wing's add_wagner_lag! instance, or nothing in the same way.

add_wagner_lag!(builder, table, bindings, sam, wing, inflow_points) -> Int

The wing's one WagnerLag instance, its va_in gathered at equal weight over every node of the wing so the lag rides the wing's mean apparent wind. Returns the instance the panels read their deficiency from.

wire_panel_flaps!(builder, mode, panels, flaps)

Connect each panel's flap deflection to the FlapDelta of the station it deflects with. A panel mapped to no surface keeps its input unconnected, which reads zero — the same deflection flap_delta_inputs gives it.

add_aero_point_forces!(builder, table, bindings, sam, wing, nodes, inflow_points,
                       body, points, wrenches) -> Vector{Int}

One AeroPointForce per structural point, gathering the panels' scattered body-frame load and delivering the world force to the point (or to its wrench half, when the point rides a body).

The output slots of point idx's body-frame aerodynamic force, or none when it has no AeroPointForce — a point outside a panel-decomposed wing.

add_wing_aero_sum!(builder, table, bindings, sam, forces) -> Int

The WingAeroSum gathering every point's body-frame force and moment, whose observables are the wing's aero_force_b and aero_moment_b readouts.

add_rigid_wing_aero!(builder, table, bindings, sam, wing, bodies, twists) -> Int

Add a RIGID_DYNAMICS wing's WingAero and wire it: the wing body's pose and each of its stations' angle and rate in, the world wrench about the body COM back into the body, and each surface's aerodynamic hinge moment on to that surface. Returns the instance. aero_eqs! drives every surface's hinge moment except a prescribed one with no aero sections, which has nothing to drive it with.

add_stations!(builder, table, bindings, sam) -> Vector{Int}

Add one twist instance per station that has a section twist to report: a StationDOF for a DYNAMIC surface, whose twist is a state, and a PrescribedTwist for a STATIC one, whose twist is a parameter. A KINEMATIC surface has neither and gets 0; its deflection is a FlapDelta instead.

add_twist_node!(builder, table, bindings, sam, idx, role, bodies, wrenches, twists)

Add the two instances a RIGID_DYNAMICS wing's structural node needs and wire them: the TwistNodePoint places it from the wing body's pose and its surface's twist, the TwistNodeWrench sends its load and moment back to the body and its bridle couple and mass on to the surface. Returns the node, which is what the incident segments connect to. Only a DYNAMIC surface takes the couple; a prescribed one holds its twist whatever the nodes pull.

rigid_wing_node(sys_struct, point) -> Bool

Whether point is a structural node of a RIGID_DYNAMICS wing. point_eqs! takes such a node out of the anchor rule entirely — the wing is the rigid body, so the node is placed by a twist-deformed body-frame offset instead of riding anything — and so does classify_points.

station_of(sys_struct, idx) -> Int

The one station point idx belongs to, or 0. Two would make its section twist ambiguous, which point_eqs! also rejects.

add_flap_deltas!(builder, table, bindings, sam, bodies, points) -> Dict{Int, Int}

Add one flap-deflection kernel per flapped station and wire what it reads: a point flap's three chord positions and its wing's frame (PointFlapDelta), or a body flap's two orientations (FlapDelta). Returns the instance of each such surface; a surface with no flap is absent, and the aero input it would feed stays unconnected and so reads the zero station_delta_eqs! binds it to.

add_joint!(builder, table, bindings, sam, joint, bodies, container, make)

Add one body-to-body joint and wire it: both end bodies' poses in, the restoring wrench on each back out. container names the joint collection its parameters are remapped over and make builds the component.

wire_segment!(builder, sys_struct, idx, role, points, segments)

Connect segment idx to its two endpoints: their pose into its inputs, its endpoint force, mass share and drag share back into theirs — a clamped endpoint included, which does not move in response but does report the load it carries. A pulley segment also exchanges the rope split and its tension with the pulley point; a tether segment takes its rest length from the winch and delivers its tension there.

The index of the winch sitting at point point_idx.

load_target(points, wrenches, idx) -> Int

The instance a point's loads are delivered to: its wrench half when it is anchored to a body or a beam, otherwise the point itself. An anchored point is two components, and only the second one accepts force, mass and drag.

buffer_slots(system, instance, buffer, name) -> Vector{Int}

The global indices of name in instance's slice of buffer (:states, :inputs, :outputs, :observables or :params).

winch_state_slots(system, instance, name, alias) -> Vector{Int}

The state slots of name on a winch anchor, taking alias when name is not the unknown that survived. WinchAnchor binds motor.vel to its own winch_vel and, for a single tether, motor.len to tether_len_1; either member of such a pair can be the one mtkcompile keeps, so the reader may not assume. alias is nothing where no pair exists.

bind_params(system, sys_struct, bindings, segment_roles, segment_instances)

Fill the parameter buffer with every kernel's build-time defaults, then record the live reader for each parameter that is a struct field read, remapped to its own instance. Returns (KernelParams, KernelParamSync).

callable_store(system) -> Tuple

The callable parameters, seeded from each kernel's build-time defaults: one entry per kernel, holding one tuple per instance of that kernel. The tuple is what keeps a polar call inside a generated kernel statically dispatched — a vector would widen to the join of the slots' types, and a wing panel's cl, cd and cm are three different types, so every coefficient lookup would become a dynamic dispatch.

retarget_tether_rest_lengths!(readers, segment_roles)

Point an unwinched tether segment's rest-length reader at its tether's length instead of the segment's own l0, matching segment_eqs!'s l0 = tether_len / n_segs for a tether with no winch to integrate that length. The reader is found by the path it reads, which names the segment and so is unique.

apply_constants!(params, system, segment_roles, segment_instances)

Write the parameters fixed by the topology straight into the buffer: a pulley segment's side and a tether segment's segment count. They never change, so no reader syncs them each step.

initial_state!(u0, system, sys_struct, point_roles, point_instances,
               body_instances, twist_instances)

Fill u0 with the initial state: each DYNAMIC station's twist and rate, each body's principal pose, each particle's pos/vel, each pulley's split pulley_len/pulley_vel and each winch's winch_vel and per-tether lengths, read from the struct. Called again by KernelInitialSync whenever a problem is reused, so a slot the struct does not reach has to be zeroed here.

PointReadout(point, pos, vel, drag, wind, force, mass, va, va_frame, aero)

Where one point's results live: its pos/vel in the output buffer and its total_drag, wind_vec and net_force in the observable buffer — all three from the instance that owns its loads, which for a point riding a body is its wrench half. A name that instance does not observe gets no slots and is left alone.

mass is the slot its segments sum their half-masses into, which with extra_mass is the monolith's point_mass.

va is the point's va_b from its AeroInflowPoint, the apparent wind the aero is solved on, and is empty for a point that has no aero instance. va_frame is then the body whose frame to rebuild it in (va_frame_body).

SegmentReadout(segment, spring_force, len, l0)

Where one segment's diagnostics live in the observable buffer.

PulleyReadout(pulley, len, vel)

Where one pulley's rope split lives in the state vector.

WinchReadout(winch, vel, force, acc, friction, set_value, tether_lengths)

Where one winch's results live: vel and its per-tether lengths in the state vector, the force/acceleration/friction in the observable buffer, and the control setpoint in the parameter buffer.

BodyReadout(body, pos, vel, acc, com, com_velocity, orientation, omega_b,
            principal, spin_p)

Where one body's results live: its pose in the output buffer, its acceleration, body-frame orientation and spin in the observable buffer, and its principal attitude/spin in the state vector.

KinematicWingReadout(body, z1, z2, y1, y2, origin, aero_points)

The reference points a fitted (KINEMATIC) wing's kinematics are rebuilt from after a step. Such a wing has no state of its own — its frame, apparent wind and per-point apparent wind are functions of where its points ended up, which is what wing_kinematics_from_points! computes. Each reference is the wing's own WeightedRefPoints, so a blend of several points is rebuilt at the same weights the Wiring feeds the kinematic body kernel.

RigidWingReadout(wing, frame, alpha_b, base_point)

wing is the body index, as Body.idx is. Where a RIGID_DYNAMICS wing's reported scalars come from: its body's frame output and its observed angular acceleration, plus the transform base point the elevation and azimuth are measured against. Its pose, spin and apparent wind are already scattered by the body and aero readouts, so these are all that is missing to fill the same scalars the monolith reads out of the integrator.

KernelStateGetter(model)

Callable (integrator, sys_struct) that scatters the runtime's results back into the struct, mirroring the monolith's get_all_state: point positions, velocities, apparent wind and drag, segment tension/length/rest length, pulley splits, winch state, tether lengths, and every wing's reported scalars. It re-runs the output and observable maps for the integrator's current state first, since the last RHS evaluation need not correspond to it. A point's total_mass comes from the half-masses its segments gather into it, so it follows their rest lengths as the monolith's point_mass does.

A point's va_b is copied from its AeroInflowPoint after the fitted-wing pass, so the value the aero is actually solved on is the compiled model's rather than a refit of it. A point with no such instance gets the monolith's fallback, its apparent wind in its own wing's frame, or the first body's if it is not a wing node.

One WinchReadout per winch, resolving its per-tether length slots.

One BodyReadout per body. Only a DYNAMIC body integrates, so a clamped or fitted one keeps whatever principal attitude and spin the struct holds — which for a fitted wing is what its own pose output already implies.

kinematic_wing_readouts(sys_struct)

One KinematicWingReadout per KINEMATIC body, holding its weighted reference points and its aerodynamic surface points, resolved once.

rigid_wing_readouts(model, sys_struct)

One RigidWingReadout per RIGID_DYNAMICS wing. A fitted wing is covered by kinematic_wing_readouts instead.

WingAeroReadout(wing, force, moment, apparent, wind)

wing is the body index, as Body.idx is. Where one wing's aero component observes the quantities the struct carries: its lumped body-frame aero_force_b/aero_moment_b, and — for a rigid wing, whose apparent wind the component computes rather than the getter — its va_b and wind_vec. A name the component does not observe gets no slots and is left alone.

StationReadout(surface, angle, rate, tether_force, tether_moment, aero_moment)

Where one station's results live: its twist and rate in the output buffer, and the bridle couple and aerodynamic hinge moment in the observable buffer. A surface whose twist is prescribed reports the same names, with the couple bound to zero.

wing_aero_readouts(model, sys_struct) -> Vector{WingAeroReadout}

One WingAeroReadout per wing that has an aero instance, resolving whichever of the four names that instance observes.

station_readouts(model) -> Vector{StationReadout}

One StationReadout per station that has a twist instance. A KINEMATIC surface has none — its deflection is a FlapDelta, which the aero reads directly and the struct does not carry.

The slots name occupies in instance's observables, or none if it has no such observable.

The output slots of point idx's va_b, or none when it has no AeroInflowPoint. This is the apparent wind the compiled model solves the aero on; reading it keeps the struct on the same value rather than a refit.

The slot instance gathers its incident segments' half-masses into, or 0 if it has no mass_in input. With extra_mass this is the monolith's point_mass, which moves as the segments' rest lengths do.

The body whose frame point idx's va_b is expressed in when the point has no AeroInflowPoint to read it from: its own wing, or the first body if it is not a wing node, which is the fallback point_eqs! uses. 0 for a model with no wing, where the monolith leaves va_b at zero. Resolved once here because the choice is fixed for the life of the model.

The instance that observes point idx's total_drag: its own, or its wrench half when the point rides a body.

write_stretched_lengths!(sys_struct)

Each tether's stretched length is the sum of its segments' current lengths, which is what tether_eqs! computes symbolically.

Copy buffer[slots] into target, component by component.

KernelControlSetter(model)

Callable (target, values) writing each winch's control setpoint into the parameter buffer, mirroring the monolith's set_set_values. values[winch.idx] is the setpoint for that winch.

write_total_mass!(sys_struct)

Write each point's effective translational mass — extra_mass plus the half-mass of every incident segment — into point.total_mass, so validate_sys_struct and the diagnostics see it. The monolith gets the same value from its point_mass observable; the batched backends compute it here because no single component owns it.

segment_load_terms(s, src_pos, src_vel, dst_pos, dst_vel, unit_stiffness,
                   unit_damping, compression_frac, compression_damping_frac,
                   l0, diameter, density, cd_tether, wind_source;
                   with_drag=true, nonlinear=false, rest_len_rate=0.0)

Every load term a segment produces, as a named tuple: the geometry (segment_vec, len, unit_vec, spring_vel), the signed scalar spring-damper tension spring and its vector spring_vec, the half_mass and half_drag each endpoint carries, and the total force_on_src/force_on_dst in the positive force-on-point sign. With nonlinear the unit_stiffness is a callable force law of strain (segment_nonlinear_force) rather than a linear rate; with with_drag = false the tether drag is dropped entirely, so cd_tether/wind_source are unused.

rest_len_rate is d(l0)/dt for a segment whose rest length is a state — a pulley leg or a winched tether member. The damper resists the rate of change of the extension len - l0, so spring_vel is the endpoint closing speed plus rest_len_rate; zero there would leave the rest-length degree of freedom undamped while its segments' dampers still drive it, which lets the damping inject energy.

segment_nonlinear_force(len, l0, spring_vel, force_law, unit_damping)

Scalar force along a segment whose unit_stiffness is a callable force law of the strain (len − l0) / l0 (the callable branch of segment_eqs!). The law owns the whole curve, slack and compression included, so there is no compression_frac; the damping term (unit_damping / len) · spring_vel is the same as in the linear case.

workload_fixture() -> String

Copy the 2-plate fixture into a scratch directory and return its path. The workload builds real models and init! serialises one, so it must not write into the package's own data.

workload_model(fixture, geometry, system_name, backend) -> SymbolicAWEModel

Load the 2-plate structure from geometry and build it on backend, ready for init!. Reloads the structure per call because building consumes it.

run_workload(fixture)

Build, initialise and step the 2-plate models on both backends. Called from this package's workload and again from the Makie extension's, where the same code is compiled against the method table Makie leaves behind — loading Makie invalidates around 19500 method instances, and a cold start is almost entirely inference.