feat: WindIO material overrides (#133) and per-mode modal mass/stiffness (#134)

CHANGELOG.md

@@ -8,7 +8,22 @@ The format follows [Keep a Changelog](https://keepachangelog.com/en/1.1.0/).
 
 ## [Unreleased]
 
-(nothing yet)
+### Added
+
+- **Material / outfitting overrides on the WindIO tower paths (#133).**
+  `Tower.from_windio` and `Tower.from_windio_with_monopile` accept optional
+  `E`, `rho`, `nu` and `outfitting_factor` keywords. Each defaults to
+  `None` (use the ontology value); pass a number to override it, so a
+  material or outfitting sensitivity sweep runs straight off a WindIO file
+  without editing the yaml. On the combined monopile+tower path the
+  override applies to both segments. Existing callers are unaffected.
+- **Per-mode generalised mass and stiffness on `ModalResult` (#134).**
+  A solve now attaches `generalized_mass` (kg) and `generalized_stiffness`
+  (N/m) arrays, one entry per mode, normalised to unit tip lateral
+  displacement so a mode's modal mass and stiffness can be read off and
+  compared against another tool's modal report. `stiffness = (2π·f)² ·
+  mass` by construction; a rigid-body platform mode (tip barely moves)
+  yields `NaN`. Not serialised (derivable from the mode shapes).
 
 ## [1.16.1] — 2026-07-06
 

src/pybmodes/io/windio.py

@@ -474,6 +474,10 @@ def read_windio_monopile_tower(
     thickness_interp: str = "linear",
     n_nodes: int | None = None,
     water_depth: float | None = None,
+    E: float | None = None,
+    rho: float | None = None,
+    nu: float | None = None,
+    outfitting_factor: float | None = None,
 ) -> WindIOMonopileTower:
     """Reduce the ``monopile`` and ``tower`` components and splice them
     into one fixed-bottom cantilever (issue #92).
@@ -509,6 +513,10 @@ def read_windio_monopile_tower(
         present; ``None`` with no ontology value keeps the monopile base
         as the clamp (correct only when the axis already starts at the
         mudline).
+    E, rho, nu, outfitting_factor : optional material / outfitting overrides
+        (issue #133), each ``None`` by default (use the ontology value).
+        When given, the override is applied to **both** the monopile and the
+        tower segment before the tube reduction.
 
     Raises
     ------
@@ -532,6 +540,19 @@ def read_windio_monopile_tower(
         yaml_path, component=component_tower, thickness_interp=thickness_interp,
     )
 
+    # Optional material / outfitting overrides (issue #133), applied to both
+    # freshly-read segments so a single value drives a whole-structure
+    # sensitivity sweep. mp / tw are local dataclasses, so mutate in place.
+    for _seg in (mp, tw):
+        if E is not None:
+            _seg.E = E
+        if rho is not None:
+            _seg.rho = rho
+        if nu is not None:
+            _seg.nu = nu
+        if outfitting_factor is not None:
+            _seg.outfitting_factor = outfitting_factor
+
     # Rigid fixed-base monopiles are clamped at the mudline, not at the
     # embedded pile tip. When the monopile reference_axis extends below the
     # seabed (e.g. IEA-15: axis -75 -> +15 with the mudline at -30), the

src/pybmodes/models/_pipeline.py

@@ -31,6 +31,55 @@
 from pybmodes.models.result import ModalResult
 
 
+def _generalised_mass_stiffness(
+    eigvecs: np.ndarray,
+    gm: np.ndarray,
+    active: np.ndarray,
+    freqs_hz: np.ndarray,
+    ref_mr: float,
+) -> tuple[np.ndarray, np.ndarray]:
+    """Physical modal (generalised) mass in kg and stiffness in N/m per mode,
+    normalised to unit tip lateral displacement (issue #134).
+
+    The non-dimensional generalised mass ``phi^T gm phi`` re-dimensionalises
+    by ``ref_mr = RM * radius``; dividing by the squared tip lateral
+    displacement DOF applies the unit-tip normalisation (the arbitrary
+    eigenvector scale cancels). Generalised stiffness follows from
+    ``omega^2 * mass``. The tip node is global DOFs 0-5 with lateral
+    displacements at ``v_disp`` (1) and ``w_disp`` (3); a mode whose tip
+    barely moves laterally (a rigid-body platform mode) yields ``NaN``.
+    """
+    active = np.asarray(active)
+
+    def _row(global_dof: int) -> int | None:
+        hits = np.nonzero(active == global_dof)[0]
+        return int(hits[0]) if hits.size else None
+
+    iv, iw = _row(1), _row(3)
+    n = eigvecs.shape[1]
+    gen_mass = np.full(n, np.nan)
+    gen_stiff = np.full(n, np.nan)
+    for i in range(n):
+        phi = eigvecs[:, i]
+        gm_form = float(phi @ gm @ phi)
+        v = float(eigvecs[iv, i]) if iv is not None else 0.0
+        w = float(eigvecs[iw, i]) if iw is not None else 0.0
+        tip = float(np.hypot(v, w))
+        # Scale-aware guard: a torsion / heave / yaw-dominated mode has only
+        # numerical lateral tip motion (tip ~ eps of the mode's peak), for
+        # which the unit-tip normalisation is undefined. An exact-zero check
+        # would let that tiny value through and divide by it, so gate on the
+        # tip motion being a real fraction of the mode's largest DOF.
+        col_scale = float(np.max(np.abs(phi)))
+        if (col_scale <= 0.0 or tip <= 1.0e-8 * col_scale
+                or not np.isfinite(gm_form) or gm_form <= 0.0):
+            continue
+        m = ref_mr * gm_form / (tip * tip)
+        gen_mass[i] = m
+        gen_stiff[i] = (2.0 * np.pi * float(freqs_hz[i])) ** 2 * m
+    return gen_mass, gen_stiff
+
+
 @dataclass
 class _SectionPropsView:
     """Minimal struct of pre-nondimensionalised section-property
@@ -268,7 +317,12 @@ def run_fem(
             and np.asarray(bmi.support.distr_m).size > 0:
         ignored.append("distributed added mass (distr_m)")
 
+    gen_mass, gen_stiff = _generalised_mass_stiffness(
+        eigvecs, gm, active, freqs_hz, nd.ref_mr
+    )
+
     return ModalResult(
         frequencies=freqs_hz, shapes=shapes, mode_labels=mode_labels,
         ignored_physics=tuple(ignored), diagnostics=diagnostics,
+        generalized_mass=gen_mass, generalized_stiffness=gen_stiff,
     )

src/pybmodes/models/result.py

@@ -79,6 +79,15 @@ class ModalResult:
         mass-matrix conditioning). Attached by :func:`run_fem`; **not**
         serialised (it describes the live solve run, not the persistent
         result data) and excluded from equality. Added 1.14.0.
+    generalized_mass, generalized_stiffness : optional ``(n_modes,)`` arrays
+        of the physical modal (generalised) mass in kg and stiffness in N/m
+        for each mode, normalised to **unit tip lateral displacement** so
+        they are directly comparable to another tool's modal output (e.g.
+        an OrcaFlex modal report). ``stiffness = (2π·f)² · mass`` by
+        construction. Most meaningful for bending modes; for a rigid-body
+        platform mode (tip displacement ≈ 0) the entry is ``NaN``. Attached
+        by :func:`run_fem` on a live solve; not serialised (derivable from
+        the shapes). Added 1.17.0 (issue #134).
     """
 
     # NOTE: field order is the positional constructor ABI for this
@@ -99,6 +108,12 @@ class ModalResult:
     # Transient solve telemetry — excluded from equality and from both
     # serialisers (it is about the run, not the persisted result).
     diagnostics: SolverDiagnostics | None = field(default=None, compare=False)
+    # Per-mode physical generalised (modal) mass and stiffness, normalised
+    # to unit tip lateral displacement (issue #134). Populated on a live
+    # solve; excluded from equality and not serialised (derivable from the
+    # shapes). ``generalized_mass`` in kg, ``generalized_stiffness`` in N/m.
+    generalized_mass: np.ndarray | None = field(default=None, compare=False)
+    generalized_stiffness: np.ndarray | None = field(default=None, compare=False)
 
     # ------------------------------------------------------------------
     # Shared integrity check

src/pybmodes/models/tower.py

@@ -326,6 +326,10 @@ def from_windio(
         n_nodes: int | None = None,
         lumped_rna_cal: bool = False,
         rna_angle_units: str = "auto",
+        E: float | None = None,
+        rho: float | None = None,
+        nu: float | None = None,
+        outfitting_factor: float | None = None,
     ) -> Tower:
         """Build a tower (or monopile) model directly from a **WindIO**
         ontology ``.yaml`` (issue #35).
@@ -378,6 +382,13 @@ def from_windio(
             the WindIO rad/deg ambiguity by magnitude; pass ``"rad"`` or
             ``"deg"`` to take the file at its word. Ignored unless
             ``lumped_rna_cal`` is set.
+        E, rho, nu, outfitting_factor : optional material / outfitting
+            overrides (issue #133). Each defaults to ``None``, meaning
+            "use the value from the ontology"; pass a number to override
+            it, e.g. for a sensitivity sweep straight off a WindIO file
+            without editing the yaml. ``E`` (Pa), ``rho`` (kg/m^3), ``nu``
+            (-) set the isotropic tower-wall material; ``outfitting_factor``
+            scales the non-structural mass (see :meth:`from_geometry`).
 
         Notes
         -----
@@ -427,8 +438,12 @@ def from_windio(
             g.outer_diameter,
             g.wall_thickness,
             flexible_length=g.flexible_length,
-            E=g.E, rho=g.rho, nu=g.nu,
-            outfitting_factor=g.outfitting_factor,
+            E=g.E if E is None else E,
+            rho=g.rho if rho is None else rho,
+            nu=g.nu if nu is None else nu,
+            outfitting_factor=(
+                g.outfitting_factor if outfitting_factor is None else outfitting_factor
+            ),
             hub_conn=hub_conn,
             tip_mass=tip_mass,
             n_nodes=n_nodes,
@@ -447,6 +462,10 @@ def from_windio_with_monopile(
         water_depth: float | None = None,
         lumped_rna_cal: bool = False,
         rna_angle_units: str = "auto",
+        E: float | None = None,
+        rho: float | None = None,
+        nu: float | None = None,
+        outfitting_factor: float | None = None,
     ) -> Tower:
         """Build a combined **monopile + tower** fixed-bottom cantilever
         from a WindIO ontology ``.yaml`` (issue #92).
@@ -494,6 +513,12 @@ def from_windio_with_monopile(
         rna_angle_units : how the auto-RNA reads ``cone_angle`` / ``uptilt``
             when ``lumped_rna_cal=True`` (``"auto"`` / ``"rad"`` / ``"deg"``);
             see :meth:`from_windio`. Ignored unless ``lumped_rna_cal`` is set.
+        E, rho, nu, outfitting_factor : optional material / outfitting
+            overrides (issue #133), each defaulting to ``None`` (use the
+            ontology value). When given, the override is applied to **both**
+            the monopile and the tower segment, so a single value drives a
+            whole-structure sensitivity sweep off a WindIO file. For
+            per-segment materials, edit the ontology instead.
 
         Notes
         -----
@@ -528,6 +553,7 @@ def from_windio_with_monopile(
             thickness_interp=thickness_interp,
             n_nodes=n_nodes,
             water_depth=water_depth,
+            E=E, rho=rho, nu=nu, outfitting_factor=outfitting_factor,
         )
         tip = _coerce_tip_mass(tip_mass)
         bmi = _build_bmi_skeleton(