BeamDyn: rotating-frame-consistent linearization Jacobians when RotStates=True#3400
Open
mayankchetan wants to merge 4 commits into
Open
BeamDyn: rotating-frame-consistent linearization Jacobians when RotStates=True#3400mayankchetan wants to merge 4 commits into
mayankchetan wants to merge 4 commits into
Conversation
…RotStates=T BD's linearized states are a frozen root-aligned inertial snapshot but were flagged RotatingFrame=T, so MBC3's Omega-corrections injected spurious spin-softening (tower-mode 'dip', w^2 = w0^2 - alpha*Omega^2, alpha ~ blade mass fraction). Re-express dXdx/dXdu/dYdx in the true rotating frame (A_r = S A S^-1 + blkdiag(-wt), B_r = S B, C_r = C S^-1, with S = [[I,0],[-wt,I]] per paired (q, dqdt) state triplet and wt built from the root angular velocity in the frozen BD basis) so the exported metadata is honest. Identity at standstill (.lin numerics bit-identical). Replaces the rotation-only transform removed in 5fbd380 (which was ~identity at the OP and lacked the velocity-mixing and transport terms). The transform applies only on the glue linearization path: a new optional IsLin flag threads from ModGlue_Linearize through FAST_JacobianPInput / FAST_JacobianPContState into BD. Solver (BuildJacobianTC) and AeroMap calls are unflagged and keep the raw frozen-frame Jacobians their iterations require (the linearization-only invariant that root orientation matches the frozen reference frame does not hold mid-step, and the integrator's Newton matrix must stay consistent with its residual). Guards: fatal on state-layout or re-anchoring-invariant violation; warnings for nonzero root angular acceleration (omega-dot transport term omitted) and large WM rotation states (first-order transform). Acceptance (NREL 5-MW BD 0-12 rpm, wide-OutFmt baseline): unpatched SS branch droops 2.29%; patched flat to 0.015%, matches python frame-fix expectation within 0.02% and ElastoDyn within 0.6% at every speed.
The root-orientation-vs-GlbRot invariant check was fatal at >1e-8, which trips on healthy linearization snapshots (deviation between the per-step re-anchoring input and the lin-time input is solver-convergence noise; observed >1e-8 on IEA-15 at 4 rpm while NREL 5-MW stays below). The transform basis remains exact regardless (omega is mapped through GlbRot, the frame the states are actually anchored to); the deviation only bounds the S-basis error. Now: warn above 1e-6 (reporting the measured value), fatal only above 1e-3 (gross misalignment = states not root-aligned at all, e.g. a caller outside the linearization snapshot path). 5-MW acceptance checker re-run after the change: PASS, unchanged (SS flatness 0.015%).
…ion snapshot Gate the RootDev warning behind DoTransport (the dXdx call) so it fires once per snapshot per BD instance instead of on all three transformed blocks; the fatal gross-misalignment check remains unconditional. Also fill the IEA-15 fleet-verdict row (patched SS flatness 0.111%, max dev vs python reference 0.022%). 5-MW acceptance checker re-run: PASS unchanged (0.015%).
…till no-op, no silent return under root acceleration Two edge-case orderings found by an adversarial review of the branch: - The re-anchoring invariant guard ran before the standstill check, so a parked linearization with a grossly mis-anchored root would abort fatally even though every transform operation is an exact identity at omega = 0. The standstill early-return now precedes the guard, making the advertised standstill no-op unconditional. - The standstill early-return could skip the omega-dot warning, so a snapshot taken at exactly zero speed while the root is accelerating (e.g. LinTimes during startup) returned silently despite the omitted -skew(alpha) transport term being nonzero there. The zero-speed path now warns before returning. No change to the transform itself; the 7-speed 5-MW acceptance checker passes with identical numbers (SS flatness 0.015%). Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01CyaECKUrqvcsDCu78FAXi1
Contributor
There was a problem hiding this comment.
Pull request overview
This PR fixes a long-standing BeamDyn linearization inconsistency when RotStates=True: Jacobians were returned in a frozen inertial snapshot while being flagged as RotatingFrame=T, which causes MBC3-based postprocessors to inject spurious speed-dependent “spin-softening.” The changes introduce an explicit rotating-frame re-expression of the BeamDyn Jacobian blocks, gated so it only affects exported linearization matrices (not solver/Newton-use Jacobians).
Changes:
- Thread a new optional
IsLinflag throughModGlue_Linearize → FAST_JacobianPInput/FAST_JacobianPContState → BeamDynso the transform is applied only during linearization output generation. - Add
BD_JacRotFrameand apply it todXdu,dYdx, anddXdxin BeamDyn whenp%RotStatesandIsLinare true, including the-ω̃transport term. - Remove dead
RotateStatesdeclarations and add runtime guards/warnings around state layout and frame re-anchoring.
Reviewed changes
Copilot reviewed 3 out of 3 changed files in this pull request and generated 4 comments.
| File | Description |
|---|---|
| modules/openfast-library/src/FAST_ModGlue.f90 | Passes IsLin=.true. when building Jacobians for .lin output so BeamDyn can apply the rotating-frame transform only on the linearization path. |
| modules/openfast-library/src/FAST_Funcs.f90 | Extends Jacobian wrapper APIs with optional IsLin and forwards it to BeamDyn to gate the transform. |
| modules/beamdyn/src/BeamDyn.f90 | Implements and applies the rotating-frame Jacobian similarity/transport transform (BD_JacRotFrame) when RotStates=True and IsLin=True. |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
Comment on lines
1387
to
+1391
| call BD_JacobianPInput(ModData%Vars, ThisTime, T%BD%Input(iInput, ModData%Ins), T%BD%p(ModData%Ins), & | ||
| T%BD%x(ModData%Ins, iState), T%BD%xd(ModData%Ins, iState), & | ||
| T%BD%z(ModData%Ins, iState), T%BD%OtherSt(ModData%Ins, iState), & | ||
| T%BD%y(ModData%Ins), T%BD%m(ModData%Ins), ErrStat2, ErrMsg2, & | ||
| dYdu=dYdu, dXdu=dXdu) | ||
| dYdu=dYdu, dXdu=dXdu, IsLin=IsLin) |
Comment on lines
1497
to
+1501
| call BD_JacobianPContState(ModData%Vars, ThisTime, T%BD%Input(iInput, ModData%Ins), T%BD%p(ModData%Ins), & | ||
| T%BD%x(ModData%Ins, iState), T%BD%xd(ModData%Ins, iState), & | ||
| T%BD%z(ModData%Ins, iState), T%BD%OtherSt(ModData%Ins, iState), & | ||
| T%BD%y(ModData%Ins), T%BD%m(ModData%Ins), ErrStat2, ErrMsg2, & | ||
| dYdx=dYdx, dXdx=dXdx) | ||
| dYdx=dYdx, dXdx=dXdx, IsLin=IsLin) |
| ! silently if the root is accelerating through zero speed — the omitted omega-dot transport | ||
| ! term is nonzero there and the exported Jacobian is missing it. | ||
| if (all(abs(omega_l) < 1.0e-12_R8Ki)) then | ||
| if (DoTransport .and. norm2(u%RootMotion%RotationAcc(:,1)) > 1.0e-4_R8Ki) then |
|
|
||
| if (DoTransport) then | ||
| ! The transform assumes a steady operating point: it has no omega-dot transport term. | ||
| if (norm2(u%RootMotion%RotationAcc(:,1)) > 1.0e-4_R8Ki) then |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Feature or improvement description
Fixes #3394. Supersedes the closed draft PR #3396, which targeted
dev. This PRcarries the same commits and targets
rc-5.0.1directly.With
RotStates = True(the default), the BeamDyn linearization Jacobians are afrozen, root-aligned snapshot in the inertial frame. But the
.linrows and columns areflagged
RotatingFrame = T. MBC3-based postprocessors thus inject spurious spin-softening(
f² = f₀² − αΩ²) into speed-independent branches. Example: the 5-MW tower side-side branchdroops 2.29 % at 12 rpm, and the IEA-22 droops ~8 % at rated speed.
This PR makes the returned Jacobians frame-consistent with the flag. A new subroutine
BD_JacRotFrameinmodules/beamdyn/src/BeamDyn.f90applies, per paired(q, dq/dt)statetriplet (pairing offset Nx/2):
to
dXdx/dXdu/dYdx. Here ω is the full rootRotationVel3-vector, mapped into thefrozen
GlbRotbasis. A scalar Ω is not sufficient: precone alone puts −0.0548 rad/s in z at12 rpm on the 5-MW. At standstill the transform is exactly the identity. This is not a
revert of the PR #2063 code that
5fbd380f8removed. That code was ≈identity at the operatingpoint, and it lacked the transport (Coriolis) term. It therefore never produced rotating-frame
Jacobians.
Design choice — where the transform lives. A new optional
IsLinflag gates thetransform. The flag is threaded
ModGlue_Linearize → FAST_JacobianPInput / FAST_JacobianPContState (FAST_Funcs.f90) → BD. The transform thus applies only atlinearization snapshots. The tight-coupling solver (
BuildJacobianTC) and the AeroMap callsare unflagged. They keep the raw frozen-frame Jacobians that their Newton iterations require.
All existing callers use keyword arguments. No C-binding, Simulink, or module-driver
interfaces are affected. No registry regeneration is needed. An alternative placement is also
workable: apply the transform glue-side in
FAST_ModGlue.f90and keep the module APIflag-free. We are glad to rework to that form if the maintainers prefer it.
Runtime guards: fatal on a state-layout-invariant violation, and fatal on a gross
re-anchoring deviation (>1e-3). The code warns once per snapshot on a measurable re-anchoring
deviation (>1e-6), on ω̇ ≠ 0, and on |WM parameters| > 0.2 (the first-order
rotation-composition limit).
Known limitation (documented in #3394): floating platforms. Platform and yaw rates fold into
the root ω. There, the scalar-Ω assumption of MBC3 independently breaks down.
Related issue, if one exists
Fixes #3394. Companion defect: #3395, fixed by PR #3393 (the
.linfull-precision writer).That defect is independent, and it also corrupts BD Campbell work. The two PRs are
independent, but one investigation isolated both. The verification input files here
neutralize #3395 with
OutFmt "ES18.9E3".Impacted areas of the software
modules/beamdyn/src/BeamDyn.f90— newBD_JacRotFrame, applied inBD_JacobianPInput/BD_JacobianPContStatewhenp%RotStates .and. IsLin. The deadRotateStatesdeclarations are removed.
modules/openfast-library/src/FAST_Funcs.f90,FAST_ModGlue.f90— these thread theoptional
IsLinflag fromModGlue_Linearizeonly.Additional supporting information
This PR replaces the closed draft #3396. That draft targeted
devand could not change itsbase after it was closed. The branch itself is unchanged in content: it was rebased onto
rc-5.0.1(from dev @0160d29a0; the three changed files are identical in both bases).The acceptance gate was re-run on an
rc-5.0.1build of this branch and passes: tower-SSflatness 0.015 % over 0–12 rpm, all damping ratios positive.
Verification (double precision, dev @
0160d29a0, wide-format input files):NREL 5-MW Campbell 0–25 rpm (CompElast=2): the tower-SS flatness improves from 2.29 %
to 0.015 %. The patched results match an independent python-side transform of the
unpatched matrices to ≤0.02 %. They match the ElastoDyn reference to ≤0.6 %. All damping
ratios are positive.
IEA fleet (29 cases; IEA-15/3.4/10/22 ported to dev): patched Fortran vs the python
reference agrees to ≤0.022 % per point. The SS flatness is 0.111 % / 0.029 % / 0.026 % /
0.379 %. The IEA-22 residual equals the python reference's own 0.361 % (a second-order
WM/prebend effect, not introduced here). The raw unpatched droop there is ~8 %.
Standstill: the
.linnumbers are bit-identical, patched vs unpatched(timestamp-only diff).
ω mapping: an instrumented print reproduced the per-blade target vector, with the
−2.5° precone z-component included (the instrumentation is removed).
Report detailing investigation and implementations:
report-bd-lin-tower-freq.html
Generative AI usage
Substantial portions of this work were developed with Anthropic Claude (Claude Code). The
human author directed all design decisions. A multi-agent workflow was used: Claude Fable 5
performed the root-cause investigation, the fix design, and the implementation; Claude
Haiku agents ran the builds and the Campbell sweeps; Claude Opus performed an
independent whole-branch review before submission (no Critical or Important findings). A
second, adversarial review by OpenAI Codex returned two guard-ordering findings. Both are
fixed on the branch, and the acceptance gate was re-run after the fix. The fix is also
cross-checked against an independent python-side implementation of the same transform (see
the verification list above). Development relied on the OpenFAST Agent Kit (OAK) —
OpenFAST-specific sub-agents and MCP tooling for reading/editing/running cases and searching
the docs/code corpus.
Co-authored-by: Anthropic Claude claude@anthropic.com
Test results, if applicable
r-test impact: the baselines change for
5MW_Land_BD_Linear,5MW_Land_BD_Linear_Aero, andDamped_Beam_Rotating. This change is the fix at work: the rotating-case Jacobians change bydesign.
Damped_Beam_Rotatingruns clean on the patched binary. The standstill beam casesare identical. All
CompElast = 1cases are unaffected. The companion r-test branch(17 regenerated
.linbaselines, branchf/bd-rotstates-framefix) is rebased onto r-testrc-5.0.1and is pushed to the fork. We will open it as an r-test PR together with this one.We can regenerate baselines on request, or follow the maintainer process.
OpenFAST/r-test#184