Configure paired RM3 floating joint in Python WEC API

.github/workflows/reference-model-baselines.yml

@@ -19,10 +19,12 @@ on:
       - 'tests/test_rm3_radiation_options_parity.py'
       - 'tests/test_rm3_pto_extension_parity.py'
       - 'tests/test_rm3_mcr_parity.py'
+      - 'tests/test_rm3_public_floating_joint_parity.py'
       - 'tests/test_rm3_mcr_orchestration_parity.py'
       - 'tests/test_rm3_mcr_seastate_parity.py'
       - 'tests/test_rm3_mooring_source.py'
       - 'tests/test_rm3_mooring_parity.py'
+      - 'tests/test_rm3_public_mooring_parity.py'
       - 'tests/test_sphere_mean_drift_source.py'
       - 'tests/test_sphere_mean_drift_parity.py'
       - 'tests/test_sphere_pto_config_parity.py'
@@ -82,7 +84,7 @@ on:
       model:
         description: Run all models or a selected reference application
         type: choice
-        options: [all, SPHERE_ELEVATION_IMPORT, MORISON_FIXED, SPHERE_MOVING_MORISON, SPHERE_MOVING_MORISON_WAVE, RM3_DD_PTO, SPHERE_MPC, GBM_BARGE, SPHERE_REACTIVE_DDPTO, SPHERE_VARIABLE_MASS, ELLIPSOID_NLH_REG, ELLIPSOID_NLH_CIC, ELLIPSOID_NLH_ODE45, OSWEC_PASSIVE_YAW_IRR_CONT, OSWEC_PASSIVE_YAW_IRR, OSWEC_PASSIVE_YAW, OSWEC_FULL_DIR, OSWEC_MULTI_WAVE, RM3_B2B, RM3_MOORING_MATRIX, SPHERE_MEAN_DRIFT, RM3_MCR_SEASTATE, RM3_END_STOPS, RM3_END_STOPS_STEP, RM3_END_STOPS_STEP_FINE, RM3_END_STOPS_STEP_FINER, RM3_END_STOPS_FULL_FINE, RM3_END_STOPS_FULL_FINER, RM3_END_STOPS_FULL_PAIR]
+        options: [all, RM3_MCR_ARRAY, SPHERE_ELEVATION_IMPORT, MORISON_FIXED, SPHERE_MOVING_MORISON, SPHERE_MOVING_MORISON_WAVE, RM3_DD_PTO, SPHERE_MPC, GBM_BARGE, SPHERE_REACTIVE_DDPTO, SPHERE_VARIABLE_MASS, ELLIPSOID_NLH_REG, ELLIPSOID_NLH_CIC, ELLIPSOID_NLH_ODE45, OSWEC_PASSIVE_YAW_IRR_CONT, OSWEC_PASSIVE_YAW_IRR, OSWEC_PASSIVE_YAW, OSWEC_FULL_DIR, OSWEC_MULTI_WAVE, RM3_B2B, RM3_MOORING_MATRIX, SPHERE_MEAN_DRIFT, RM3_MCR_SEASTATE, RM3_END_STOPS, RM3_END_STOPS_STEP, RM3_END_STOPS_STEP_FINE, RM3_END_STOPS_STEP_FINER, RM3_END_STOPS_FULL_FINE, RM3_END_STOPS_FULL_FINER, RM3_END_STOPS_FULL_PAIR]
         default: all
 
 jobs:
@@ -92,7 +94,7 @@ jobs:
     strategy:
       fail-fast: false
       matrix:
-        model: ${{ fromJSON(inputs.model == 'SPHERE_ELEVATION_IMPORT' && '["SPHERE_ELEVATION_IMPORT"]' || inputs.model == 'MORISON_FIXED' && '["MORISON_FIXED"]' || inputs.model == 'SPHERE_MOVING_MORISON' && '["SPHERE_MOVING_MORISON"]' || inputs.model == 'SPHERE_MOVING_MORISON_WAVE' && '["SPHERE_MOVING_MORISON_WAVE"]' || inputs.model == 'RM3_DD_PTO' && '["RM3_DD_PTO"]' || inputs.model == 'SPHERE_MPC' && '["SPHERE_MPC"]' || inputs.model == 'GBM_BARGE' && '["GBM_BARGE"]' || inputs.model == 'SPHERE_REACTIVE_DDPTO' && '["SPHERE_REACTIVE_DDPTO"]' || inputs.model == 'SPHERE_VARIABLE_MASS' && '["SPHERE_VARIABLE_MASS"]' || inputs.model == 'ELLIPSOID_NLH_REG' && '["ELLIPSOID_NLH_REG"]' || inputs.model == 'ELLIPSOID_NLH_CIC' && '["ELLIPSOID_NLH_CIC"]' || inputs.model == 'ELLIPSOID_NLH_ODE45' && '["ELLIPSOID_NLH_ODE45"]' || inputs.model == 'OSWEC_PASSIVE_YAW_IRR_CONT' && '["OSWEC_PASSIVE_YAW_IRR_CONT"]' || inputs.model == 'OSWEC_PASSIVE_YAW_IRR' && '["OSWEC_PASSIVE_YAW_IRR"]' || inputs.model == 'OSWEC_PASSIVE_YAW' && '["OSWEC_PASSIVE_YAW"]' || inputs.model == 'OSWEC_FULL_DIR' && '["OSWEC_FULL_DIR"]' || inputs.model == 'OSWEC_MULTI_WAVE' && '["OSWEC_MULTI_WAVE"]' || inputs.model == 'RM3_B2B' && '["RM3_B2B"]' || inputs.model == 'RM3_MOORING_MATRIX' && '["RM3_MOORING_MATRIX"]' || inputs.model == 'SPHERE_MEAN_DRIFT' && '["SPHERE_MEAN_DRIFT"]' || inputs.model == 'RM3_MCR_SEASTATE' && '["RM3_MCR_SEASTATE"]' || inputs.model == 'RM3_END_STOPS' && '["RM3_END_STOPS"]' || inputs.model == 'RM3_END_STOPS_STEP' && '["RM3_END_STOPS_STEP"]' || inputs.model == 'RM3_END_STOPS_STEP_FINE' && '["RM3_END_STOPS_STEP_FINE"]' || inputs.model == 'RM3_END_STOPS_STEP_FINER' && '["RM3_END_STOPS_STEP_FINER"]' || inputs.model == 'RM3_END_STOPS_FULL_FINE' && '["RM3_END_STOPS_FULL_FINE"]' || inputs.model == 'RM3_END_STOPS_FULL_FINER' && '["RM3_END_STOPS_FULL_FINER"]' || inputs.model == 'RM3_END_STOPS_FULL_PAIR' && '["RM3_END_STOPS_FULL_FINE", "RM3_END_STOPS_FULL_FINER"]' || '["SPHERE_ELEVATION_IMPORT", "MORISON_FIXED", "SPHERE_MOVING_MORISON", "SPHERE_MOVING_MORISON_WAVE", "RM3_DD_PTO", "RM3", "OSWEC", "OSWEC_Nonhydro", "GBM_BARGE", "SPHERE_VARIABLE_MASS", "SPHERE_REACTIVE_DDPTO", "OSWEC_FULL_DIR", "OSWEC_MULTI_WAVE", "OSWEC_PASSIVE_YAW", "Sphere", "SPHERE_MEAN_DRIFT", "RM3_B2B", "RM3_END_STOPS", "RM3_PTO_Extension", "RM3_Radiation_Options", "Sphere_Passive", "Sphere_PTO_Config", "Sphere_Reactive_PI", "Sphere_Declutching", "Sphere_Latching", "RM3_MCR", "RM3_MCR_ARRAY", "RM3_MCR_EXCEL", "RM3_MCR_MAT", "RM3_MCR_SEASTATE", "ELLIPSOID_NLH_REG", "ELLIPSOID_NLH_CIC", "ELLIPSOID_NLH_ODE45"]') }}
+        model: ${{ fromJSON(inputs.model == 'RM3_MCR_ARRAY' && '["RM3_MCR_ARRAY"]' || inputs.model == 'SPHERE_ELEVATION_IMPORT' && '["SPHERE_ELEVATION_IMPORT"]' || inputs.model == 'MORISON_FIXED' && '["MORISON_FIXED"]' || inputs.model == 'SPHERE_MOVING_MORISON' && '["SPHERE_MOVING_MORISON"]' || inputs.model == 'SPHERE_MOVING_MORISON_WAVE' && '["SPHERE_MOVING_MORISON_WAVE"]' || inputs.model == 'RM3_DD_PTO' && '["RM3_DD_PTO"]' || inputs.model == 'SPHERE_MPC' && '["SPHERE_MPC"]' || inputs.model == 'GBM_BARGE' && '["GBM_BARGE"]' || inputs.model == 'SPHERE_REACTIVE_DDPTO' && '["SPHERE_REACTIVE_DDPTO"]' || inputs.model == 'SPHERE_VARIABLE_MASS' && '["SPHERE_VARIABLE_MASS"]' || inputs.model == 'ELLIPSOID_NLH_REG' && '["ELLIPSOID_NLH_REG"]' || inputs.model == 'ELLIPSOID_NLH_CIC' && '["ELLIPSOID_NLH_CIC"]' || inputs.model == 'ELLIPSOID_NLH_ODE45' && '["ELLIPSOID_NLH_ODE45"]' || inputs.model == 'OSWEC_PASSIVE_YAW_IRR_CONT' && '["OSWEC_PASSIVE_YAW_IRR_CONT"]' || inputs.model == 'OSWEC_PASSIVE_YAW_IRR' && '["OSWEC_PASSIVE_YAW_IRR"]' || inputs.model == 'OSWEC_PASSIVE_YAW' && '["OSWEC_PASSIVE_YAW"]' || inputs.model == 'OSWEC_FULL_DIR' && '["OSWEC_FULL_DIR"]' || inputs.model == 'OSWEC_MULTI_WAVE' && '["OSWEC_MULTI_WAVE"]' || inputs.model == 'RM3_B2B' && '["RM3_B2B"]' || inputs.model == 'RM3_MOORING_MATRIX' && '["RM3_MOORING_MATRIX"]' || inputs.model == 'SPHERE_MEAN_DRIFT' && '["SPHERE_MEAN_DRIFT"]' || inputs.model == 'RM3_MCR_SEASTATE' && '["RM3_MCR_SEASTATE"]' || inputs.model == 'RM3_END_STOPS' && '["RM3_END_STOPS"]' || inputs.model == 'RM3_END_STOPS_STEP' && '["RM3_END_STOPS_STEP"]' || inputs.model == 'RM3_END_STOPS_STEP_FINE' && '["RM3_END_STOPS_STEP_FINE"]' || inputs.model == 'RM3_END_STOPS_STEP_FINER' && '["RM3_END_STOPS_STEP_FINER"]' || inputs.model == 'RM3_END_STOPS_FULL_FINE' && '["RM3_END_STOPS_FULL_FINE"]' || inputs.model == 'RM3_END_STOPS_FULL_FINER' && '["RM3_END_STOPS_FULL_FINER"]' || inputs.model == 'RM3_END_STOPS_FULL_PAIR' && '["RM3_END_STOPS_FULL_FINE", "RM3_END_STOPS_FULL_FINER"]' || '["SPHERE_ELEVATION_IMPORT", "MORISON_FIXED", "SPHERE_MOVING_MORISON", "SPHERE_MOVING_MORISON_WAVE", "RM3_DD_PTO", "RM3", "OSWEC", "OSWEC_Nonhydro", "GBM_BARGE", "SPHERE_VARIABLE_MASS", "SPHERE_REACTIVE_DDPTO", "OSWEC_FULL_DIR", "OSWEC_MULTI_WAVE", "OSWEC_PASSIVE_YAW", "Sphere", "SPHERE_MEAN_DRIFT", "RM3_B2B", "RM3_END_STOPS", "RM3_PTO_Extension", "RM3_Radiation_Options", "Sphere_Passive", "Sphere_PTO_Config", "Sphere_Reactive_PI", "Sphere_Declutching", "Sphere_Latching", "RM3_MCR", "RM3_MCR_ARRAY", "RM3_MCR_EXCEL", "RM3_MCR_MAT", "RM3_MCR_SEASTATE", "ELLIPSOID_NLH_REG", "ELLIPSOID_NLH_CIC", "ELLIPSOID_NLH_ODE45"]') }}
     name: ${{ matrix.model }}
     steps:
       - uses: actions/checkout@v7
@@ -294,6 +296,11 @@ jobs:
           WEC_SIM_MCR_LABEL: ${{ matrix.model }}
           WEC_SIM_APPLICATIONS_DIR: applications
           WEC_SIM_MATLAB_MODEL_OUTPUT_DIR: matlab-reference-model-output
+      - run: python -m pytest -q tests/test_rm3_public_floating_joint_parity.py
+        if: matrix.model == 'RM3_MCR_ARRAY'
+        env:
+          WEC_SIM_APPLICATIONS_DIR: applications
+          WEC_SIM_MATLAB_MODEL_OUTPUT_DIR: matlab-reference-model-output
       - run: python -m pytest -q tests/test_rm3_mcr_seastate_parity.py
         if: matrix.model == 'RM3_MCR_SEASTATE'
         env:
@@ -309,6 +316,11 @@ jobs:
         env:
           WEC_SIM_APPLICATIONS_DIR: applications
           WEC_SIM_MATLAB_MODEL_OUTPUT_DIR: matlab-reference-model-output
+      - run: python -m pytest -q tests/test_rm3_public_mooring_parity.py
+        if: matrix.model == 'RM3_MOORING_MATRIX'
+        env:
+          WEC_SIM_APPLICATIONS_DIR: applications
+          WEC_SIM_MATLAB_MODEL_OUTPUT_DIR: matlab-reference-model-output
       - run: python -m pytest -q tests/test_sphere_mean_drift_source.py
         if: matrix.model == 'SPHERE_MEAN_DRIFT'
         env:

PARITY.md

@@ -44,6 +44,7 @@ the production Python code. The live wave comparison passed on 6 October
 | RM3 PTO extension float and spar free decays | Pinned MATLAB Applications `RM3_PTO_Extension` inputs and BEMIO-generated RM3 HDF5 | The no-wave floating-joint solver initializes the float at +5 m or the spar at −5 m, reproducing each published case's +5 m PTO stroke. Over 30 s, maximum active-body heave differences are 47.2 mm for float and 6.85 mm for spar; velocity differences are 86.1 mm/s and 1.87 mm/s. Both initial poses and the passive body's heave agree within numerical precision; PTO force and power are zero in Python and within `1.7e-8` in MATLAB output. |
 | RM3 Multiple Condition Runs Option 1 physical sweep | Pinned MATLAB Applications RM3 input, expanded to eight scalar height × period × PTO-damping combinations | The Python `regularCIC` floating-joint solver agrees across all eight 400 s conditions within 75.1 mm surge, 0.387 mm heave, 0.000183 rad pitch, 0.413 mm PTO stroke, 0.231 mm/s PTO speed, 0.479 kN PTO force, and 0.487 kW PTO power. These scalar runs validate each physical condition; all three MCR drivers are separately checked below. |
 | RM3 Multiple Condition Runs Options 1–3 orchestration | Pinned MATLAB Applications arrays, Option 2 Excel grid, Option 3 MAT file, and three actual MATLAB R2025b `wecSimMCR` runs | The Python loaders produce the same ordered eight-case table for all three published inputs. Each driver has a separate paired gate for every body's 4,001 position and velocity samples, every PTO stroke, speed, force, and power trace, the eight mean powers, and two 2×2 power matrices with period rows and height columns. Python mean absorbed powers differ from MATLAB by at most 0.155 kW or 0.059% across eight conditions in the MAT driver. The Option 2 Excel and Option 3 MAT driver outputs agree numerically exactly; Option 1 array and Option 3 MAT body positions differ by at most `1.9e-12` m. MATLAB's PTO power sign is opposite Python's positive absorbed-power convention. Phase-seed sweeps, multiple PTOs, and other postprocessing are not yet paired. |
+| Public RM3 floating-joint configuration | Pinned MATLAB MCR Option 1 array condition 1 and published MooringMatrix application; fresh paired [MCR](https://github.com/cmudrc/wec-sim-python/actions/runs/37787926107) and [mooring](https://github.com/cmudrc/wec-sim-python/actions/runs/37787936962) runs | `WEC.floating_joint` routes a Python-configured two-body device through the already paired pitched-slider solver and returns named body, four-coordinate, and PTO histories. In the 400 s MCR case, maximum differences are 35.1 mm surge, 0.228 mm heave, `3.46e-6` rad pitch, 0.243 mm PTO stroke, 163 N PTO force, and 118 W absorbed power; explicit gates cover positions, speeds, stroke, force, and power at all 4,001 samples. A second public configuration reproduces the 400 s imported-elevation MooringMatrix run, including its initial spar offset and joint surge spring, under the existing 40,001-sample body/PTO/mooring gates. Both fresh MATLAB jobs pass. The builder's implicit added-mass default remains physical; source-compatible delayed feedback is explicit for the MCR comparison. This layout does not represent arbitrary Simscape joints or off-axis PTO endpoints. |
 | RM3 imported-spectrum sea-state MCR | Pinned MATLAB Applications `RM3_MCROPT3_SeaState` and three actual MATLAB MCR runs | All three imported component tables agree to `4.9e-15`, wave elevations to `1.0e-13` m, and active excitation forces and moments to `3.8e-8` N or N m. The public `WEC.run(ImportedSpectrumWave(...))` path now replays all 4,001 incident-wave samples in each 400 s case within `1.0e-13` m; its general WEC dynamics are not claimed to reproduce the published four-coordinate RM3 motion. On MATLAB body velocities, the 60 s radiation convolution matches every logged force component within `5.4e-9` N or N m. The applied added-mass force reconstructs from the exported matrix and delayed acceleration within `7e-7` N or N m. The exact pitched-slider geometry and source-compatible mass feedback bring all three 400 s paired trajectories within 14.4 mm surge, 1.03 mm heave, and 0.000638 rad pitch; velocity differences are below 3.26 mm/s surge, 0.883 mm/s heave, and 0.000256 rad/s pitch. PTO stroke, speed, force, and absorbed power differ by at most 2.44 mm, 1.41 mm/s, 1.69 kN, and 1.19 kW. Mean absorbed power differs by at most 29 W. Explicit paired gates passed for the specialized MCR solver; PR #18 was merged. |
 | RM3 PTO-Sim direct linear generator | Pinned MATLAB Applications `PTO-Sim/RM3/RM3_DD_PTO`, its BEMIO-generated RM3 HDF5, and a fresh 400 s `ode4` run at 0.0005 s | Both the focused Python two-heave runner and a public `WEC.pto(linear_generator=...)` configuration integrate the source block's d/q flux states and electrical angle with float and spar motion. The public configuration uses the existing body-local endpoint and coordinate maps. Source-only checks verify PTO relative motion, friction, load voltage, and mechanical/electrical power identities. Across all 40,001 saved samples, both paired gates require body heaves/velocities and PTO stroke/speed within 10 µm or µm/s, generator force and powers within 0.005 N or W, and all three phase currents/voltages within 0.0001 A/0.01 V. The [fresh source run](https://github.com/cmudrc/wec-sim-python/actions/runs/37749845161) and local paired tests pass. This validates the published vertical two-body regular-wave layout, not general PTO-Sim electrical networks, radiation-memory coupling, or arbitrary attachment geometry. |
 | Sphere passive controller | Pinned MATLAB Applications `Controls/Passive (P)` and MATLAB-generated Sphere HDF5 | The Python heave model represents the published proportional controller as an 860,870 N s/m damper. Maximum differences over 400 s are 0.063 mm position, 0.065 mm/s velocity, 56 N controller force, and 32 W controller power after accounting for MATLAB's opposite force and power signs. |
@@ -245,8 +246,10 @@ that body only heaves, the shifted point does not affect its trajectory;
 attachment-location dynamics remain unpaired.
 `wecsim.WEC` provides a Python builder for this same validated path and
 returns named NumPy body, coordinate, and PTO histories. The JSON case runner
-remains available for saved cases; the Python builder currently covers the
-`linear_subspace` layout only.
+remains available for saved cases; the Python builder covers the mapped
+`linear_subspace`, one-body `floating_gbm`, and two-body RM3-style
+`floating_joint` layouts. The latter uses a relative-heave PTO and does not
+accept arbitrary attachment points.
 `python -m wecsim CASE.json --output motion.csv` runs a
 supported dynamics configuration. The case declares wave, body, constraint,
 PTO, and time settings; the result includes all six body position and

README.md

@@ -71,10 +71,36 @@ two-body device with a shared pitch pivot and body-local PTO points. Run it
 with `python -m examples.configurable_rm3_pto`. Relative HDF5 paths in
 the Python API resolve from the current directory unless `base_dir` is passed
 to `wec.run`. Body order must match the HDF5 hydrodynamic body order. The
-Python builder covers `linear_subspace` with supported waves and the
-single-body `floating_gbm` regular-wave layout. Its named coordinates use
-small-motion kinematics and fixed-axis PTOs. Other validated reference
-layouts remain available through the case runner and focused solver functions.
+Python builder covers `linear_subspace`, the single-body `floating_gbm`
+regular-wave layout, and the paired two-body `floating_joint` layout. Mapped
+coordinates use small-motion kinematics and fixed-axis PTOs. The floating
+joint uses its own pitched-slider geometry and a relative-heave PTO; arbitrary
+PTO attachment points are not part of that reduced layout.
+
+For the published RM3 floating joint, configure the two bodies in HDF5 order:
+
+```python
+from wecsim import RegularCICWave, WEC
+
+wec = WEC("RM3 floating joint")
+float_body = wec.body("float", "path/to/rm3.h5",
+                      inertia=(0, 21_306_090.66, 0))
+spar = wec.body("spar", "path/to/rm3.h5",
+                inertia=(0, 94_407_091.24, 0))
+wec.floating_joint(float_body, spar, damping=1_200_000)
+result = wec.run(RegularCICWave(2.5, 8), dt=0.1, end_time=400,
+                 ramp_time=100, radiation_memory=60)
+stroke = result.ptos["relative_heave"].stroke
+```
+
+The default uses implicit added mass. For comparisons with the pinned MATLAB
+Simulink block, `added_mass_scheme="simulink_delay"` selects its numerical
+feedback setting explicitly. `floating_joint` also accepts a joint surge
+spring, PTO stiffness and equilibrium, hard stops, and convolution/FIR
+radiation where the underlying paired solver supports them. The returned
+`absorbed_power` counts the linear damper; spring energy exchange can be
+computed from `-force * velocity`. The returned stroke is float heave minus
+spar heave, so an initial body offset appears in the first stroke sample.
 
 For the published generalized-body-mode barge, generate its HDF5 with the
 Applications `Generalized_Body_Modes/hydroData/bemio.m`, then run:

tests/test_python_api.py

@@ -164,6 +164,32 @@ def test_python_builder_runs_imported_elevation_with_relative_mat_path():
     assert dict(response.raw.extra_outputs)["body1_excitation_force"].shape == (3, 6)
 
 
+def test_python_builder_configures_paired_floating_joint_and_pto():
+    wec = WEC("RM3 floating joint")
+    float_body = wec.body("float", HYDRO, inertia=(0, 21_306_090.66, 0))
+    spar = wec.body("spar", HYDRO, inertia=(0, 94_407_091.24, 0))
+    wec.floating_joint(
+        float_body, spar, damping=1_200_000, stiffness=100_000,
+        equilibrium_position=0.1, pto_name="main",
+    )
+    result = wec.run(RegularWave(2.5, 8), dt=0.1, end_time=0.2,
+                     ramp_time=1)
+    assert set(result.coordinates) == {"surge", "float_heave",
+                                       "spar_heave", "pitch"}
+    pto = result.ptos["main"]
+    np.testing.assert_allclose(
+        pto.force,
+        -1_200_000 * pto.velocity - 100_000 * (pto.stroke - 0.1),
+    )
+    np.testing.assert_allclose(pto.absorbed_power,
+                               1_200_000 * pto.velocity**2)
+    np.testing.assert_allclose(
+        pto.stroke,
+        result.coordinates["float_heave"].position
+        - result.coordinates["spar_heave"].position,
+    )
+
+
 def test_python_builder_rejects_foreign_body_attachment():
     first = WEC("first")
     second = WEC("second")