Fix InflowWind Grid3D cubic velocity interpolation NaN bug

[andrew-platt]

Ready to merge

Feature or improvement description

This is a bug fix targeting rc-5.0.1, not a feature. It fixes a long-standing defect in InflowWind's Grid3D flow field (TurbSim/Bladed/HAWC full-field wind files, WindType=3 and similar) that produces NaN wind velocities whenever VelInterpCubic = True is set in the InflowWind input file and the calling code queries velocity without also requesting acceleration output. FAST.Farm's AWAE module is exactly such a caller: it samples ambient wind velocity for its high-resolution grid without asking for acceleration, so any FAST.Farm case with VelInterpCubic = True returns a 100%-NaN ambient wind field and fails immediately at t=0 with:

FARM_InitialCO: ... the rotor plane has left the low-resolution domain

Root cause

IfW_FlowField_GetVelAcc (in modules/inflowwind/src/IfW_FlowField.f90) determines whether the caller wants acceleration output from whether it allocated the optional AccelUVW output array:

OutputAccel = allocated(AccelUVW)

For the Grid3D_FieldType case, this flag (combined with VelInterpCubic) selects one of four interpolation code paths. When VelInterpCubic = True and OutputAccel = False ("cubic velocity, no acceleration"), Grid3DField_GetCell is called with CalcAccel = OutputAccel = .false., which causes it to skip populating the local AccCell(8,3) array entirely — it is only filled if (CalcAccel) then ....

The velocity is then computed by Grid3DField_GetVelAccCubic, which implements a cubic Hermite spline in time. Critically, the Hermite value formula (not just the optional acceleration output) depends on AccCell, since a proper cubic-in-time interpolant needs derivative (tangent) information at both time endpoints to compute the interpolated value itself:

Velocity = C1*P(:, 1) + C2*PP(:, 1) + C3*P(:, 2) + C4*PP(:, 2)

where PP is derived from AccCell. Since AccCell is a local, stack-allocated array that is never initialized in this code path, PP is computed from garbage memory, corrupting Velocity for every single query — confirmed via debugger to contain Inf in every element for this test case, while the companion VelCell (actual wind data) was fully valid. This explains why the reported failure showed a 100%-NaN velocity array with a 100%-valid position array: the defect is purely a function of which interpolation-flag combination is selected, not of any particular spatial location.

A second, related defect was found in the same area. IfW_Grid3DField_CalcAccel (which precomputes the time-derivative array G3D%Acc needed for cubic interpolation) has:

if (G3D%NTGrids < 3) then
   G3D%Acc = 0.0_SiKi
   return
end if

G3D%NTGrids is the number of tower grid points (from an optional tower file), not the number of time steps (G3D%NSteps). For the very common case of no tower file (NTGrids = 0), this unconditionally forces G3D%Acc to zero and skips the real cubic-spline derivative calculation — regardless of how many time steps are actually available (4004 in the reported case). This doesn't itself produce NaN (zero is well-defined), but it silently defeats cubic-in-time interpolation whenever no tower file is present, which is likely most cases that use VelInterpCubic = True.

Both defects are pre-existing in rc-5.0.1 (and earlier) — reproduced and confirmed on an unmodified rc-5.0.1 checkout with identical input files, so this is not a regression from any other branch. It is essentially untested in the existing test suite because every other reg-test input file in the repository sets VelInterpCubic = False; the reported case is the only one using True.

Fix

1. IfW_FlowField_GetVelAcc: pass OutputAccel .or. FF%VelInterpCubic as the CalcAccel argument to Grid3DField_GetCell, so AccCell is always populated whenever cubic-in-time interpolation is active, independent of whether the caller also wants acceleration returned.
2. IfW_Grid3DField_CalcAccel: change the early-return guard from G3D%NTGrids < 3 to G3D%NSteps < 3, so real time derivatives are computed whenever enough time samples exist, regardless of whether a tower file is present.

Ready to merge

Related issue, if one exists

None filed yet; found during investigation of a reported FAST.Farm FARM_InitialCO fatal error ("rotor plane has left the low-resolution domain") on a 2-turbine wind-tunnel test case.

Impacted areas of the software

• modules/inflowwind/src/IfW_FlowField.f90 (IfW_FlowField_GetVelAcc, IfW_Grid3DField_CalcAccel)
• Any glue code or driver using InflowWind's Grid3D flow field type (TurbSim, Bladed, HAWC full-field wind, WindType = 3, 4, 5, 7) with VelInterpCubic = True: FAST.Farm (AWAE ambient-wind sampling and per-turbine InflowWind instances), OpenFAST, InflowWind driver, AeroDyn driver.
• No effect on the (much more common) VelInterpCubic = False linear-interpolation path.

Additional supporting information

Diagnosed via GDB by tracing a NaN observed in AeroDyn's RotInflow%Blade%InflowVel backwards through FASTWrapper, AWAE's ambient-wind fill, and into InflowWind's Grid3D cubic interpolation, where AccCell was directly observed to contain Inf in all elements at the point of use. Confirmed the same failure reproduces on an unmodified rc-5.0.1 build with AMReX disabled, ruling out any connection to unrelated work-in-progress changes on another branch (including a separately-identified, unrelated WakeDynamics array-bounds issue, which is not part of this PR).

Generative AI usage

Root-cause investigation (GDB tracing, source review), the code fix, and the new unit tests were developed with substantial assistance from an AI coding agent.
Co-authored-by: Claude Sonnet 5 noreply@anthropic.com

Test results, if applicable

• Added modules/inflowwind/tests/test_grid3d_field.F90 with two new unit tests, registered in inflowwind_utest:
  • test_grid3d_cubic_vel_only — builds a minimal synthetic Grid3D field with a spatially-uniform, linear-in-time velocity ramp and VelInterpCubic = True; queries velocity both with and without requesting acceleration output. Fails today with "not finite (NaN/Inf)" without fix AeroDyn14 Driver #1; passes with it.
  • test_grid3d_calcaccel_no_tower — calls IfW_Grid3DField_CalcAccel directly on a field with no tower grid and checks that the computed derivative recovers the known true value. Fails today (derivative forced to zero) without fix DWM Driver #2; passes with it.
• Full inflowwind_utest suite passes with no regressions in any pre-existing test.

[Copilot]

Pull request overview

Fixes InflowWind Grid3D cubic-in-time velocity interpolation behavior so callers that request velocity-only (no acceleration output) no longer receive NaN/Inf wind velocities, and corrects the Grid3D acceleration precompute guard to use the number of time steps rather than the tower-grid count.

Changes:

• Ensure Grid3DField_GetCell populates AccCell whenever VelInterpCubic is enabled, independent of whether AccelUVW is allocated by the caller.
• Fix IfW_Grid3DField_CalcAccel early-return condition to check G3D%NSteps < 3 (time steps) instead of G3D%NTGrids < 3 (tower grid points).
• Add and register new InflowWind unit-test coverage for the cubic-velocity-only bug path.

Reviewed changes

Copilot reviewed 4 out of 4 changed files in this pull request and generated 3 comments.

File	Description
unit_tests/CMakeLists.txt	Adds the new Grid3D unit test source to the InflowWind unit test target.
modules/inflowwind/tests/test_grid3d_field.F90	Introduces a new unit test module exercising the cubic velocity interpolation path.
modules/inflowwind/tests/inflowwind_utest.F90	Registers the new Grid3D test suite in the InflowWind unit test runner.
modules/inflowwind/src/IfW_FlowField.f90	Fixes the Grid3D cubic interpolation path selection and corrects the CalcAccel guard condition.

---

💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.

[Copilot]

Pull request overview

Copilot reviewed 4 out of 4 changed files in this pull request and generated 2 comments.

[andrew-platt]

ifw_BoxExceed regression baseline comparison

The build-all-test-modules-debug CI failure on ifw_BoxExceed is an expected consequence of this fix, not a regression.

ifw_BoxExceed uses a TurbSim .bts wind file with no tower grid (NTGrids=0) and requests acceleration output (CalcAccel=true). Acceleration is always computed via the cubic-Hermite-in-time formula, which depends on IfW_Grid3DField_CalcAccel. Under the old NTGrids<3 guard, this routine zeroed the time-derivative input for any no-tower case, so the reported acceleration was computed as if dU/dt = 0 at the endpoints. With the fix (NSteps-based guard), the real cubic-spline-derived derivative is used instead, which is why the accelerations reported here differ from the current committed baseline.

Position and velocity channels are byte-for-byte identical (0.0 diff across all 300 points) — only the acceleration channels (UA, VA, WA) change:

Channel	Max |diff|	Mean |diff|	Old max magnitude	New max magnitude
T, X, Y, Z, U, V, W	0.0	0.0	—	—
UA	10.16	1.30	5.37	13.57
VA	11.97	1.39	5.57	15.50
WA	9.26	1.28	7.97	12.63

The new accelerations are ~2-3x larger in typical magnitude, consistent with the old code systematically biasing derivatives toward zero rather than a small numeric drift.