1008 lines
38 KiB
Python
1008 lines
38 KiB
Python
"""Proof-gated tangent columns for supported piston branch networks.
|
|
|
|
This module is deliberately narrower than the generic algebraic solver. It
|
|
only compiles a tangent provider after proving the state layout, component
|
|
types, physical connections, and causal execution plan used by every selected
|
|
piston branch. The legacy three-piston entry point remains available for its
|
|
committed fixture, while the topology-driven entry point discovers any number
|
|
of branches without depending on component names. A failed proof leaves the
|
|
ordinary seed-0 numerical Jacobian in control; a runtime mode boundary requests
|
|
the same one-build fallback through :class:`ExactColumnsUnavailable`.
|
|
"""
|
|
|
|
from __future__ import annotations
|
|
|
|
from collections.abc import Mapping, Sequence
|
|
from dataclasses import dataclass
|
|
from math import isfinite
|
|
from typing import TYPE_CHECKING
|
|
|
|
import numpy as np
|
|
|
|
from app.simulation.solvers.jacobian import ExactColumnsUnavailable
|
|
from app.simulation.solvers.mechanical import MechanicalConstraintGroup
|
|
from app.simulation.systems.network import Endpoint
|
|
|
|
if TYPE_CHECKING:
|
|
from app.simulation.systems.generic import GenericFluidSystem
|
|
|
|
|
|
_TARGET_BRANCH_NAMES = (
|
|
(
|
|
"mass_friction_endstops_10",
|
|
"pn_brp2_8",
|
|
"pn_c1_8",
|
|
"pneumatic_69",
|
|
"elasticendstop_8",
|
|
),
|
|
(
|
|
"mass_friction_endstops_11",
|
|
"pn_brp2_9",
|
|
"pn_c1_9",
|
|
"pneumatic_68",
|
|
"elasticendstop_9",
|
|
),
|
|
(
|
|
"mass_friction_endstops_12",
|
|
"pn_brp2_10",
|
|
"pn_c1_10",
|
|
"pneumatic_66",
|
|
"elasticendstop_10",
|
|
),
|
|
)
|
|
|
|
|
|
@dataclass(frozen=True)
|
|
class ThreePistonBranch:
|
|
mass: object
|
|
piston: object
|
|
chamber: object
|
|
pipe: object
|
|
contact: object
|
|
chamber_connection_port: str
|
|
velocity_index: int
|
|
position_index: int
|
|
|
|
|
|
@dataclass(frozen=True)
|
|
class ThreePistonTangentCompilation:
|
|
eligible: bool
|
|
reason: str | None
|
|
columns: tuple[int, ...] = ()
|
|
provider: "ThreePistonTangentProvider | None" = None
|
|
reached_assignment_count: int = 0
|
|
|
|
def diagnostics(self) -> dict[str, object]:
|
|
return {
|
|
"eligible": self.eligible,
|
|
"fallbackReason": self.reason,
|
|
"columns": list(self.columns),
|
|
"columnCount": len(self.columns),
|
|
"reachedAssignmentCount": self.reached_assignment_count,
|
|
}
|
|
|
|
|
|
@dataclass(frozen=True)
|
|
class _PrimalContext:
|
|
time: float
|
|
state: np.ndarray
|
|
connected_h: dict[str, dict[str, float]]
|
|
|
|
|
|
def _failed(reason: str) -> ThreePistonTangentCompilation:
|
|
return ThreePistonTangentCompilation(False, reason)
|
|
|
|
|
|
class ThreePistonTangentProvider:
|
|
"""Batched selected-branch provider compiled for one system instance."""
|
|
|
|
def __init__(
|
|
self,
|
|
system: "GenericFluidSystem",
|
|
branches: tuple[ThreePistonBranch, ...],
|
|
columns: tuple[int, ...],
|
|
state_offsets: Mapping[str, tuple[int, int]],
|
|
reached_assignment_count: int,
|
|
) -> None:
|
|
self.system = system
|
|
self.branches = branches
|
|
self.columns = columns
|
|
self.state_offsets = dict(state_offsets)
|
|
self.reached_assignment_count = reached_assignment_count
|
|
self._context: _PrimalContext | None = None
|
|
self._capture_requested = False
|
|
|
|
def request_primal_capture(self) -> None:
|
|
"""Capture exactly the next completed RHS primal closure."""
|
|
|
|
self._capture_requested = True
|
|
|
|
def cancel_primal_capture(self) -> None:
|
|
"""Discard a pending capture after an interrupted base RHS."""
|
|
|
|
self._capture_requested = False
|
|
|
|
def record_primal(
|
|
self,
|
|
time: float,
|
|
state: Sequence[float],
|
|
connected_h: Mapping[str, Mapping[str, float]],
|
|
) -> None:
|
|
if not self._capture_requested:
|
|
return
|
|
self._capture_requested = False
|
|
required_components = {
|
|
item.name
|
|
for branch in self.branches
|
|
for item in (branch.chamber, branch.pipe)
|
|
}
|
|
self._context = _PrimalContext(
|
|
time=float(time),
|
|
state=np.asarray(state, dtype=float).copy(),
|
|
connected_h={
|
|
name: {port: float(value) for port, value in values.items()}
|
|
for name, values in connected_h.items()
|
|
if name in required_components
|
|
},
|
|
)
|
|
|
|
def __call__(
|
|
self,
|
|
time: float,
|
|
state: np.ndarray,
|
|
normalized_indexes: tuple[int, ...],
|
|
) -> np.ndarray:
|
|
if normalized_indexes != self.columns:
|
|
raise ValueError("Compiled tangent columns were requested out of order.")
|
|
context = self._context
|
|
values = np.asarray(state, dtype=float)
|
|
if (
|
|
context is None
|
|
or float(time) != context.time
|
|
or values.shape != context.state.shape
|
|
or not np.array_equal(values, context.state)
|
|
):
|
|
raise ExactColumnsUnavailable("stalePrimalContext")
|
|
return self._evaluate(context)
|
|
|
|
@staticmethod
|
|
def _zeros(width: int) -> tuple[float, ...]:
|
|
return (0.0,) * width
|
|
|
|
@staticmethod
|
|
def _has_signal(vector: Sequence[float]) -> bool:
|
|
return any(float(value) != 0.0 for value in vector)
|
|
|
|
@staticmethod
|
|
def _require_valid(value: object, prefix: str) -> None:
|
|
if not bool(getattr(value, "valid", False)):
|
|
reason = getattr(value, "reason", None) or "invalid"
|
|
raise ExactColumnsUnavailable(f"{prefix}:{reason}")
|
|
|
|
def _evaluate(self, context: _PrimalContext) -> np.ndarray:
|
|
# The primal RHS may disable the causal fast path after a residual
|
|
# audit. Recheck after that closure and before replaying its compiled
|
|
# assignments so a dynamic downgrade uses the ordinary full FD build.
|
|
if not self.system.pressure_flow_solver.causal_fast_path_enabled:
|
|
raise ExactColumnsUnavailable("causalFastPathDisabled")
|
|
width = len(self.columns)
|
|
zero = self._zeros(width)
|
|
out = np.zeros((len(context.state), width), dtype=float)
|
|
seeds = {
|
|
column: tuple(float(index == seed_index) for index in range(width))
|
|
for seed_index, column in enumerate(self.columns)
|
|
}
|
|
solver = self.system.pressure_flow_solver
|
|
tangents: dict[str, tuple[float, ...]] = {}
|
|
|
|
def tangent(key: str) -> tuple[float, ...]:
|
|
return tangents.get(key, zero)
|
|
|
|
def negative_sum(
|
|
vectors: Sequence[Sequence[float]],
|
|
) -> tuple[float, ...]:
|
|
return tuple(
|
|
-sum(float(vector[index]) for vector in vectors)
|
|
for index in range(width)
|
|
)
|
|
|
|
for variable in ("x", "v"):
|
|
plan = solver._causal_effort_plan_by_variable.get(variable)
|
|
if plan is None:
|
|
raise ExactColumnsUnavailable("causalEffortPlanUnavailable")
|
|
for assignment in plan:
|
|
matched = [
|
|
branch
|
|
for branch in self.branches
|
|
if any(
|
|
unknown.component == branch.mass.name
|
|
for unknown in assignment.members
|
|
)
|
|
]
|
|
if len(matched) > 1:
|
|
raise ExactColumnsUnavailable("coupledTargetMechanicalSeeds")
|
|
vector = zero
|
|
if matched:
|
|
branch = matched[0]
|
|
column = (
|
|
branch.position_index
|
|
if variable == "x"
|
|
else branch.velocity_index
|
|
)
|
|
vector = seeds[column]
|
|
for member in assignment.members:
|
|
tangents[member.id] = vector
|
|
|
|
geometry: dict[str, object] = {}
|
|
chamber_properties: dict[str, object] = {}
|
|
for branch in self.branches:
|
|
piston = branch.piston
|
|
linearization = piston.linearize_geometry_and_force(
|
|
tangent(f"{piston.name}.port_4.x"),
|
|
tangent(f"{piston.name}.port_5.x"),
|
|
tangent(f"{piston.name}.port_4.v"),
|
|
tangent(f"{piston.name}.port_5.v"),
|
|
zero,
|
|
)
|
|
self._require_valid(linearization, "pistonGeometry")
|
|
chamber = branch.chamber
|
|
volume = float(chamber.total_volume())
|
|
if volume <= float(chamber.cvol0) / 100.0:
|
|
raise ExactColumnsUnavailable("volumeFloor")
|
|
primal = chamber.medium.properties_from_mU(
|
|
chamber.state.m,
|
|
chamber.state.U,
|
|
volume,
|
|
)
|
|
properties = chamber.medium.linearize_properties_from_mU(
|
|
chamber.state.m,
|
|
chamber.state.U,
|
|
volume,
|
|
zero,
|
|
zero,
|
|
linearization.volume_tangent,
|
|
properties=primal,
|
|
)
|
|
self._require_valid(properties, "chamberProperties")
|
|
geometry[piston.name] = linearization
|
|
chamber_properties[chamber.name] = properties
|
|
|
|
pressure_plan = solver._causal_effort_plan_by_variable.get("p")
|
|
if pressure_plan is None:
|
|
raise ExactColumnsUnavailable("causalPressurePlanUnavailable")
|
|
for assignment in pressure_plan:
|
|
matched = [
|
|
branch
|
|
for branch in self.branches
|
|
if any(
|
|
unknown.component == branch.chamber.name
|
|
for unknown in assignment.members
|
|
)
|
|
]
|
|
if len(matched) > 1:
|
|
raise ExactColumnsUnavailable("coupledTargetPressureSeeds")
|
|
vector = (
|
|
chamber_properties[matched[0].chamber.name].tangents.p
|
|
if matched
|
|
else zero
|
|
)
|
|
for member in assignment.members:
|
|
tangents[member.id] = tuple(vector)
|
|
|
|
contacts: dict[str, object] = {}
|
|
for branch in self.branches:
|
|
piston = branch.piston
|
|
linearization = piston.linearize_geometry_and_force(
|
|
tangent(f"{piston.name}.port_4.x"),
|
|
tangent(f"{piston.name}.port_5.x"),
|
|
tangent(f"{piston.name}.port_4.v"),
|
|
tangent(f"{piston.name}.port_5.v"),
|
|
tangent(f"{piston.name}.port_1.p"),
|
|
)
|
|
self._require_valid(linearization, "pistonPressureForce")
|
|
geometry[piston.name] = linearization
|
|
contact = branch.contact
|
|
contact_linearization = contact.linearize_contact_force(
|
|
tangent(f"{contact.name}.port_1.x"),
|
|
tangent(f"{contact.name}.port_2.x"),
|
|
tangent(f"{contact.name}.port_1.v"),
|
|
tangent(f"{contact.name}.port_2.v"),
|
|
)
|
|
self._require_valid(contact_linearization, "contactMode")
|
|
contacts[contact.name] = contact_linearization
|
|
|
|
equations = {
|
|
equation.id: equation for equation in solver.equation_templates
|
|
}
|
|
branch_component = {
|
|
item.name: branch
|
|
for branch in self.branches
|
|
for item in (
|
|
branch.mass,
|
|
branch.piston,
|
|
branch.chamber,
|
|
branch.pipe,
|
|
branch.contact,
|
|
)
|
|
}
|
|
pipe_flows: dict[str, object] = {}
|
|
for stage in solver._explicit_flow_plan:
|
|
pending: list[tuple[str, tuple[float, ...]]] = []
|
|
for assignment in stage.assignments:
|
|
equation = equations.get(assignment.equation_id)
|
|
if equation is None:
|
|
raise ExactColumnsUnavailable("causalEquationMissing")
|
|
dependencies = [
|
|
tangent(variable)
|
|
for variable in equation.variables
|
|
if variable != assignment.unknown.id
|
|
]
|
|
if not any(
|
|
self._has_signal(vector) for vector in dependencies
|
|
):
|
|
pending.append((assignment.unknown.id, zero))
|
|
continue
|
|
if equation.relation == "sumToZero":
|
|
vectors = [
|
|
tangent(variable)
|
|
for variable in equation.variables
|
|
if variable != assignment.unknown.id
|
|
and variable in solver._unknowns_by_id
|
|
and solver._unknowns_by_id[variable].role == "flow"
|
|
]
|
|
pending.append(
|
|
(assignment.unknown.id, negative_sum(vectors))
|
|
)
|
|
continue
|
|
|
|
component = assignment.component
|
|
model = getattr(component, "MODEL_TYPE", None)
|
|
if model == "amesim_pnl0001":
|
|
branch = branch_component.get(component.name)
|
|
if branch is None or component is not branch.pipe:
|
|
raise ExactColumnsUnavailable("unexpectedPipeReach")
|
|
flow_linearization = pipe_flows.get(component.name)
|
|
if flow_linearization is None:
|
|
properties = component.medium.properties_from_mU(
|
|
component.state.m,
|
|
component.state.U,
|
|
component.volume,
|
|
)
|
|
flow_linearization = component.linearize_mass_flow(
|
|
component.port_1.p,
|
|
component.port_2.p,
|
|
properties.T,
|
|
)
|
|
self._require_valid(
|
|
flow_linearization,
|
|
"pipeMassFlow",
|
|
)
|
|
pipe_flows[component.name] = flow_linearization
|
|
vector = tuple(
|
|
flow_linearization.partial_p_1 * first
|
|
+ flow_linearization.partial_p_2 * second
|
|
for first, second in zip(
|
|
tangent(f"{component.name}.port_1.p"),
|
|
tangent(f"{component.name}.port_2.p"),
|
|
strict=True,
|
|
)
|
|
)
|
|
elif model == "amesim_lstp00a":
|
|
contact_linearization = contacts.get(component.name)
|
|
if contact_linearization is None:
|
|
raise ExactColumnsUnavailable(
|
|
f"unexpectedContactReach:{component.name}"
|
|
)
|
|
sign = (
|
|
1.0
|
|
if assignment.unknown.port == "port_1"
|
|
else -1.0
|
|
)
|
|
vector = tuple(
|
|
sign * value
|
|
for value in contact_linearization.force_tangent
|
|
)
|
|
elif model == "amesim_pnrp17":
|
|
piston_geometry = geometry.get(component.name)
|
|
if piston_geometry is None:
|
|
raise ExactColumnsUnavailable(
|
|
f"unexpectedPistonReach:{component.name}"
|
|
)
|
|
other = next(
|
|
(
|
|
tangent(variable)
|
|
for variable in equation.variables
|
|
if variable != assignment.unknown.id
|
|
and variable.endswith(".f")
|
|
),
|
|
zero,
|
|
)
|
|
sign = (
|
|
-1.0
|
|
if equation.id.endswith(
|
|
"piston_side_force_balance"
|
|
)
|
|
else 1.0
|
|
)
|
|
vector = tuple(
|
|
-force + sign * pressure
|
|
for force, pressure in zip(
|
|
other,
|
|
piston_geometry.pressure_force_tangent,
|
|
strict=True,
|
|
)
|
|
)
|
|
else:
|
|
raise ExactColumnsUnavailable(
|
|
f"unsupportedReach:{model or 'connection'}"
|
|
)
|
|
pending.append((assignment.unknown.id, tuple(vector)))
|
|
for key, vector in pending:
|
|
tangents[key] = vector
|
|
|
|
outflow: dict[Endpoint, tuple[float, ...]] = {}
|
|
for branch in self.branches:
|
|
enthalpy_tangent = tuple(
|
|
chamber_properties[branch.chamber.name].tangents.h
|
|
)
|
|
for port_name in branch.chamber.ports:
|
|
outflow[Endpoint(branch.chamber.name, port_name)] = (
|
|
enthalpy_tangent
|
|
)
|
|
# PNRP17 mirrors the connected chamber enthalpy on its sole
|
|
# pneumatic port at the already closed primal point.
|
|
outflow[Endpoint(branch.piston.name, "port_1")] = enthalpy_tangent
|
|
|
|
connected_tangent: dict[
|
|
str,
|
|
dict[str, tuple[float, ...]],
|
|
] = {name: {} for name in self.system.network.components}
|
|
for connection in self.system.network.connections:
|
|
if connection.kind != "physical":
|
|
continue
|
|
first, second = connection.endpoints
|
|
connected_tangent[first.component][first.port] = outflow.get(
|
|
second,
|
|
zero,
|
|
)
|
|
connected_tangent[second.component][second.port] = outflow.get(
|
|
first,
|
|
zero,
|
|
)
|
|
|
|
for component_name, port_values in connected_tangent.items():
|
|
component = self.system.network.components[component_name]
|
|
if not getattr(component, "PRESSURE_FLOW_DEPENDS_ON_STREAM", False):
|
|
continue
|
|
if any(
|
|
self._has_signal(vector) for vector in port_values.values()
|
|
):
|
|
raise ExactColumnsUnavailable(
|
|
f"streamSensitiveReach:{component_name}"
|
|
)
|
|
|
|
for branch in self.branches:
|
|
chamber = branch.chamber
|
|
chamber_linearization = chamber.linearize_state_derivative(
|
|
context.connected_h[chamber.name],
|
|
state_mass_tangent=zero,
|
|
state_energy_tangent=zero,
|
|
external_volume_tangent=(
|
|
geometry[branch.piston.name].volume_tangent
|
|
),
|
|
external_volume_rate_tangent=(
|
|
geometry[branch.piston.name].volume_flow_tangent
|
|
),
|
|
port_mass_flow_tangents={
|
|
name: tangent(f"{chamber.name}.{name}.m_flow")
|
|
for name in chamber.ports
|
|
},
|
|
connected_h_tangents=connected_tangent[chamber.name],
|
|
property_linearization=chamber_properties[chamber.name],
|
|
)
|
|
self._require_valid(chamber_linearization, "chamberDerivative")
|
|
offset, size = self.state_offsets[chamber.name]
|
|
if size != 2:
|
|
raise ValueError("Target chamber state layout changed.")
|
|
out[offset, :], out[offset + 1, :] = (
|
|
chamber_linearization.tangents
|
|
)
|
|
|
|
pipe = branch.pipe
|
|
pipe_properties = pipe.medium.linearize_properties_from_mU(
|
|
pipe.state.m,
|
|
pipe.state.U,
|
|
pipe.volume,
|
|
zero,
|
|
zero,
|
|
zero,
|
|
)
|
|
self._require_valid(pipe_properties, "pipeProperties")
|
|
pipe_linearization = pipe.linearize_state_derivative(
|
|
context.connected_h[pipe.name],
|
|
state_mass_tangent=zero,
|
|
state_energy_tangent=zero,
|
|
port_mass_flow_tangents={
|
|
name: tangent(f"{pipe.name}.{name}.m_flow")
|
|
for name in pipe.ports
|
|
},
|
|
connected_h_tangents=connected_tangent[pipe.name],
|
|
property_linearization=pipe_properties,
|
|
)
|
|
self._require_valid(pipe_linearization, "pipeDerivative")
|
|
offset, size = self.state_offsets[pipe.name]
|
|
if size != 2:
|
|
raise ValueError("Target pipe state layout changed.")
|
|
out[offset, :], out[offset + 1, :] = pipe_linearization.tangents
|
|
|
|
target_pneumatic = {
|
|
item.name
|
|
for branch in self.branches
|
|
for item in (branch.chamber, branch.pipe)
|
|
}
|
|
for entry in self.system.mechanical_state_reducer.state_entries:
|
|
if (
|
|
isinstance(entry, MechanicalConstraintGroup)
|
|
or entry.name in target_pneumatic
|
|
):
|
|
continue
|
|
for port_name, port in entry.ports.items():
|
|
definition = port.definition
|
|
if (
|
|
definition is not None
|
|
and definition.kind == "physical"
|
|
and definition.domain == "pneumatic"
|
|
and self._has_signal(
|
|
tangent(f"{entry.name}.{port_name}.m_flow")
|
|
)
|
|
):
|
|
raise ExactColumnsUnavailable(
|
|
f"unsupportedDynamicReach:{entry.name}"
|
|
)
|
|
|
|
mass_seeds = {
|
|
branch.mass.name: (
|
|
seeds[branch.velocity_index],
|
|
seeds[branch.position_index],
|
|
)
|
|
for branch in self.branches
|
|
}
|
|
for group in self.system.mechanical_state_reducer.groups:
|
|
if len(group.components) != 1:
|
|
raise ExactColumnsUnavailable("reachableRigidMassGroup")
|
|
mass = group.representative
|
|
force_1 = tangent(f"{mass.name}.port_1.f")
|
|
force_2 = tangent(f"{mass.name}.port_2.f")
|
|
velocity, position = mass_seeds.get(mass.name, (zero, zero))
|
|
if not any(
|
|
self._has_signal(vector)
|
|
for vector in (force_1, force_2, velocity, position)
|
|
):
|
|
continue
|
|
fixed = (
|
|
mass._constraint_acceleration == 0.0
|
|
and mass._constraint_velocity == 0.0
|
|
)
|
|
if fixed and (
|
|
abs(group.total_unconstrained_force())
|
|
<= 1.0e-12
|
|
* max(
|
|
abs(mass.port_1.f),
|
|
abs(mass.port_2.f),
|
|
1.0,
|
|
)
|
|
):
|
|
raise ExactColumnsUnavailable("mechanicalReleaseBoundary")
|
|
mass_linearization = mass.linearize_state_derivative(
|
|
force_1,
|
|
force_2,
|
|
velocity,
|
|
position,
|
|
constraint_mode="current" if fixed else "free",
|
|
)
|
|
self._require_valid(
|
|
mass_linearization,
|
|
"mechanicalDerivative",
|
|
)
|
|
offset, size = self.state_offsets[mass.name]
|
|
if size != 2:
|
|
raise ValueError("Mechanical state layout changed.")
|
|
out[offset, :], out[offset + 1, :] = mass_linearization.tangents
|
|
|
|
if not np.all(np.isfinite(out)):
|
|
raise ExactColumnsUnavailable("nonFiniteTangentColumns")
|
|
return out
|
|
|
|
|
|
@dataclass(frozen=True)
|
|
class _PistonBranchSpec:
|
|
names: tuple[str, str, str, str, str]
|
|
chamber_connection_port: str
|
|
|
|
|
|
def _compile_named_piston_tangent_provider(
|
|
system: "GenericFluidSystem",
|
|
branch_specs: Sequence[_PistonBranchSpec],
|
|
) -> ThreePistonTangentCompilation:
|
|
"""Compile a named, topology-proven set of supported piston branches."""
|
|
|
|
solver = system.pressure_flow_solver
|
|
if not solver.causal_fast_path_eligible:
|
|
return _failed("causalFastPathIneligible")
|
|
if not solver.causal_fast_path_enabled:
|
|
return _failed("causalFastPathDisabled")
|
|
closure = system._thermofluid_closure_plan
|
|
if closure.uses_conservative_global_solver:
|
|
return _failed("conservativeThermofluidClosure")
|
|
if system.pneumatic_storage_reducer.groups:
|
|
return _failed("coupledPneumaticStorage")
|
|
|
|
state_offsets: dict[str, tuple[int, int]] = {}
|
|
cursor = 0
|
|
group_by_name: dict[str, MechanicalConstraintGroup] = {}
|
|
for entry in system.mechanical_state_reducer.state_entries:
|
|
if isinstance(entry, MechanicalConstraintGroup):
|
|
if len(entry.components) != 1:
|
|
cursor += 2
|
|
continue
|
|
component = entry.representative
|
|
state_offsets[component.name] = (cursor, 2)
|
|
group_by_name[component.name] = entry
|
|
cursor += 2
|
|
else:
|
|
state_offsets[entry.name] = (cursor, int(entry.state_size))
|
|
cursor += int(entry.state_size)
|
|
|
|
expected_types = (
|
|
"amesim_mecmas21",
|
|
"amesim_pnrp17",
|
|
"amesim_pnch012",
|
|
"amesim_pnl0001",
|
|
"amesim_lstp00a",
|
|
)
|
|
branches: list[ThreePistonBranch] = []
|
|
for branch_spec in branch_specs:
|
|
names = branch_spec.names
|
|
try:
|
|
components = tuple(system.network.components[name] for name in names)
|
|
except KeyError:
|
|
return _failed("targetComponentMissing")
|
|
if tuple(getattr(item, "MODEL_TYPE", None) for item in components) != expected_types:
|
|
return _failed("targetComponentTypeMismatch")
|
|
mass, piston, chamber, pipe, contact = components
|
|
if mass.name not in group_by_name or mass.name not in state_offsets:
|
|
return _failed("targetMechanicalStateLayout")
|
|
offset, size = state_offsets[mass.name]
|
|
if size != 2:
|
|
return _failed("targetMechanicalStateLayout")
|
|
branches.append(
|
|
ThreePistonBranch(
|
|
mass=mass,
|
|
piston=piston,
|
|
chamber=chamber,
|
|
pipe=pipe,
|
|
contact=contact,
|
|
chamber_connection_port=(
|
|
branch_spec.chamber_connection_port
|
|
),
|
|
velocity_index=offset,
|
|
position_index=offset + 1,
|
|
)
|
|
)
|
|
|
|
required_pairs: set[frozenset[Endpoint]] = set()
|
|
for branch in branches:
|
|
required_pairs.update(
|
|
{
|
|
frozenset((Endpoint(branch.mass.name, "port_1"), Endpoint(branch.piston.name, "port_2"))),
|
|
frozenset(
|
|
(
|
|
Endpoint(branch.piston.name, "port_1"),
|
|
Endpoint(
|
|
branch.chamber.name,
|
|
branch.chamber_connection_port,
|
|
),
|
|
)
|
|
),
|
|
frozenset((Endpoint(branch.chamber.name, "port_1"), Endpoint(branch.pipe.name, "port_1"))),
|
|
frozenset((Endpoint(branch.piston.name, "port_5"), Endpoint(branch.contact.name, "port_1"))),
|
|
}
|
|
)
|
|
actual_pairs = {
|
|
frozenset(connection.endpoints)
|
|
for connection in system.network.connections
|
|
if connection.kind == "physical"
|
|
}
|
|
if not required_pairs <= actual_pairs:
|
|
return _failed("targetTopologyMismatch")
|
|
|
|
required_methods = (
|
|
("piston", "linearize_geometry_and_force"),
|
|
("chamber", "linearize_state_derivative"),
|
|
("pipe", "linearize_mass_flow"),
|
|
("pipe", "linearize_state_derivative"),
|
|
("contact", "linearize_contact_force"),
|
|
("mass", "linearize_state_derivative"),
|
|
)
|
|
for branch in branches:
|
|
for owner, method in required_methods:
|
|
if not callable(getattr(getattr(branch, owner), method, None)):
|
|
return _failed(f"missingTangentPrimitive:{owner}.{method}")
|
|
if not callable(getattr(branch.chamber.medium, "linearize_properties_from_mU", None)):
|
|
return _failed("missingTangentPrimitive:medium.linearize_properties_from_mU")
|
|
|
|
if any(
|
|
len(group.components) != 1
|
|
for group in system.mechanical_state_reducer.groups
|
|
):
|
|
return _failed("rigidMassAggregation")
|
|
|
|
target_mass_names = {branch.mass.name for branch in branches}
|
|
target_chamber_names = {branch.chamber.name for branch in branches}
|
|
target_pipe_names = {branch.pipe.name for branch in branches}
|
|
target_piston_names = {branch.piston.name for branch in branches}
|
|
target_contact_names = {branch.contact.name for branch in branches}
|
|
reached_ids: set[str] = set()
|
|
for variable in ("x", "v"):
|
|
for assignment in solver._causal_effort_plan_by_variable.get(
|
|
variable,
|
|
(),
|
|
):
|
|
if any(
|
|
member.component in target_mass_names
|
|
for member in assignment.members
|
|
):
|
|
reached_ids.update(member.id for member in assignment.members)
|
|
for assignment in solver._causal_effort_plan_by_variable.get("p", ()):
|
|
if any(
|
|
member.component in target_chamber_names
|
|
for member in assignment.members
|
|
):
|
|
reached_ids.update(member.id for member in assignment.members)
|
|
|
|
equations = {
|
|
item.id: item for item in solver.equation_templates
|
|
}
|
|
reached_assignments = []
|
|
allowed_reached_models = {
|
|
"amesim_pnl0001",
|
|
"amesim_pnrp17",
|
|
"amesim_lstp00a",
|
|
}
|
|
for stage in solver._explicit_flow_plan:
|
|
stage_reached = []
|
|
for assignment in stage.assignments:
|
|
equation = equations.get(assignment.equation_id)
|
|
if equation is None:
|
|
return _failed("causalEquationMissing")
|
|
dependencies = set(equation.variables) - {
|
|
assignment.unknown.id
|
|
}
|
|
if not dependencies.intersection(reached_ids):
|
|
continue
|
|
component = assignment.component
|
|
if equation.relation != "sumToZero":
|
|
model_type = getattr(component, "MODEL_TYPE", None)
|
|
if model_type not in allowed_reached_models:
|
|
return _failed(f"unsupportedReach:{model_type}")
|
|
expected_names = {
|
|
"amesim_pnl0001": target_pipe_names,
|
|
"amesim_pnrp17": target_piston_names,
|
|
"amesim_lstp00a": target_contact_names,
|
|
}[model_type]
|
|
if component.name not in expected_names:
|
|
return _failed(
|
|
f"unexpectedReach:{model_type}:{component.name}"
|
|
)
|
|
if bool(
|
|
getattr(
|
|
component,
|
|
"PRESSURE_FLOW_DEPENDS_ON_STREAM",
|
|
False,
|
|
)
|
|
):
|
|
return _failed(f"streamSensitiveReach:{component.name}")
|
|
stage_reached.append(assignment)
|
|
reached_assignments.extend(stage_reached)
|
|
reached_ids.update(item.unknown.id for item in stage_reached)
|
|
|
|
# The only non-zero h_outflow seeds are the target PNCH ports. Each
|
|
# direct neighbour must consume it in a supported target balance, mirror
|
|
# it through the one-port piston, or terminate at a fixed PNPL01 cap.
|
|
# Secondary pressure blocks may contain the same causal flow coordinates,
|
|
# so membership alone is not evidence of a stream derivative. The direct
|
|
# enthalpy reach proof below, plus the runtime dynamic-owner gate, is the
|
|
# relevant condition for this proof-gated branch program.
|
|
neighbor_by_endpoint: dict[Endpoint, Endpoint] = {}
|
|
for connection in system.network.connections:
|
|
if connection.kind != "physical":
|
|
continue
|
|
first, second = connection.endpoints
|
|
neighbor_by_endpoint[first] = second
|
|
neighbor_by_endpoint[second] = first
|
|
permitted_h_neighbors = {
|
|
"amesim_pnl0001",
|
|
"amesim_pnrp17",
|
|
"amesim_pnpl01",
|
|
}
|
|
for branch in branches:
|
|
for port_name in branch.chamber.ports:
|
|
neighbor = neighbor_by_endpoint.get(
|
|
Endpoint(branch.chamber.name, port_name)
|
|
)
|
|
if neighbor is None:
|
|
return _failed("targetStreamBindingMissing")
|
|
component = system.network.components[neighbor.component]
|
|
if (
|
|
getattr(component, "MODEL_TYPE", None)
|
|
not in permitted_h_neighbors
|
|
):
|
|
return _failed(
|
|
f"unsupportedStreamReach:{component.name}"
|
|
)
|
|
if bool(
|
|
getattr(
|
|
component,
|
|
"PRESSURE_FLOW_DEPENDS_ON_STREAM",
|
|
False,
|
|
)
|
|
):
|
|
return _failed(f"streamSensitiveReach:{component.name}")
|
|
|
|
columns = tuple(
|
|
sorted(
|
|
index
|
|
for branch in branches
|
|
for index in (branch.velocity_index, branch.position_index)
|
|
)
|
|
)
|
|
provider = ThreePistonTangentProvider(
|
|
system,
|
|
tuple(branches),
|
|
columns,
|
|
state_offsets,
|
|
reached_assignment_count=len(reached_assignments),
|
|
)
|
|
return ThreePistonTangentCompilation(
|
|
True,
|
|
None,
|
|
columns,
|
|
provider,
|
|
reached_assignment_count=len(reached_assignments),
|
|
)
|
|
|
|
|
|
def _physical_neighbor_map(
|
|
system: "GenericFluidSystem",
|
|
) -> dict[Endpoint, Endpoint]:
|
|
"""Return the one-to-one physical connector map proved by the network."""
|
|
|
|
neighbors: dict[Endpoint, Endpoint] = {}
|
|
for connection in system.network.connections:
|
|
if connection.kind != "physical":
|
|
continue
|
|
first, second = connection.endpoints
|
|
# SimulationNetwork already rejects multiply connected physical ports.
|
|
# Retain a defensive gate because this compiler may also be called by
|
|
# custom network builders in tests or downstream applications.
|
|
if first in neighbors or second in neighbors:
|
|
raise ValueError("Physical endpoint has more than one connection.")
|
|
neighbors[first] = second
|
|
neighbors[second] = first
|
|
return neighbors
|
|
|
|
|
|
def _discover_supported_piston_branch_specs(
|
|
system: "GenericFluidSystem",
|
|
) -> tuple[_PistonBranchSpec, ...] | ThreePistonTangentCompilation:
|
|
"""Discover every complete catalog piston branch by type and port topology.
|
|
|
|
A PNRP17 is the unambiguous root: its mechanical piston-side port must be
|
|
driven by a singleton MECMAS21 coordinate, its pneumatic port must feed a
|
|
PNCH012 whose first port feeds PNL0001, and its rod-side port must meet an
|
|
LSTP00A contact. If even one PNRP17 is only partially supported, reject the
|
|
batch with a stable reason instead of silently omitting derivative columns.
|
|
"""
|
|
|
|
try:
|
|
neighbors = _physical_neighbor_map(system)
|
|
except ValueError:
|
|
return _failed("unsupportedPistonBranchTopology:multipleConnection")
|
|
|
|
components = system.network.components
|
|
|
|
def model_type(endpoint: Endpoint | None) -> str | None:
|
|
if endpoint is None:
|
|
return None
|
|
return getattr(components[endpoint.component], "MODEL_TYPE", None)
|
|
|
|
pistons = tuple(
|
|
component
|
|
for component in components.values()
|
|
if getattr(component, "MODEL_TYPE", None) == "amesim_pnrp17"
|
|
)
|
|
if not pistons:
|
|
return _failed("supportedPistonBranchMissing")
|
|
|
|
specs: list[_PistonBranchSpec] = []
|
|
for piston in pistons:
|
|
mass_endpoint = neighbors.get(Endpoint(piston.name, "port_2"))
|
|
if (
|
|
model_type(mass_endpoint) != "amesim_mecmas21"
|
|
or mass_endpoint is None
|
|
or mass_endpoint.port != "port_1"
|
|
):
|
|
return _failed("unsupportedPistonBranchTopology:mass")
|
|
|
|
chamber_endpoint = neighbors.get(Endpoint(piston.name, "port_1"))
|
|
if model_type(chamber_endpoint) != "amesim_pnch012":
|
|
return _failed("unsupportedPistonBranchTopology:chamber")
|
|
assert chamber_endpoint is not None
|
|
|
|
pipe_endpoint = neighbors.get(
|
|
Endpoint(chamber_endpoint.component, "port_1")
|
|
)
|
|
if (
|
|
model_type(pipe_endpoint) != "amesim_pnl0001"
|
|
or pipe_endpoint is None
|
|
or pipe_endpoint.port != "port_1"
|
|
):
|
|
return _failed("unsupportedPistonBranchTopology:pipe")
|
|
|
|
contact_endpoint = neighbors.get(Endpoint(piston.name, "port_5"))
|
|
if (
|
|
model_type(contact_endpoint) != "amesim_lstp00a"
|
|
or contact_endpoint is None
|
|
or contact_endpoint.port != "port_1"
|
|
):
|
|
return _failed("unsupportedPistonBranchTopology:contact")
|
|
|
|
specs.append(
|
|
_PistonBranchSpec(
|
|
names=(
|
|
mass_endpoint.component,
|
|
piston.name,
|
|
chamber_endpoint.component,
|
|
pipe_endpoint.component,
|
|
contact_endpoint.component,
|
|
),
|
|
chamber_connection_port=chamber_endpoint.port,
|
|
)
|
|
)
|
|
|
|
# Port uniqueness already proves unique pistons and masses, but explicitly
|
|
# reject a custom multi-port chamber/contact/pipe shared by two roots. The
|
|
# tangent propagation assumes one geometry seed per selected state owner.
|
|
for role_index in range(5):
|
|
if len({spec.names[role_index] for spec in specs}) != len(specs):
|
|
return _failed("unsupportedPistonBranchTopology:sharedComponent")
|
|
return tuple(specs)
|
|
|
|
|
|
def compile_supported_piston_tangent_provider(
|
|
system: "GenericFluidSystem",
|
|
) -> ThreePistonTangentCompilation:
|
|
"""Compile all name-independent, topology-supported piston branches."""
|
|
|
|
discovered = _discover_supported_piston_branch_specs(system)
|
|
if isinstance(discovered, ThreePistonTangentCompilation):
|
|
return discovered
|
|
return _compile_named_piston_tangent_provider(system, discovered)
|
|
|
|
|
|
def compile_three_piston_tangent_provider(
|
|
system: "GenericFluidSystem",
|
|
) -> ThreePistonTangentCompilation:
|
|
"""Compile the committed legacy three-piston target by its stable names."""
|
|
|
|
return _compile_named_piston_tangent_provider(
|
|
system,
|
|
tuple(
|
|
_PistonBranchSpec(
|
|
names=names,
|
|
chamber_connection_port="port_3",
|
|
)
|
|
for names in _TARGET_BRANCH_NAMES
|
|
),
|
|
)
|