Skip to content

BeamDyn: return DTBeam from BD_Init and sub-step the driver; error checks on registry checkpoint reads/writes #3517

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
andrew-platt wants to merge 10 commits into
base: rc-5.0.1
Choose a base branch
from
Open

BeamDyn: return DTBeam from BD_Init and sub-step the driver; error checks on registry checkpoint reads/writes #3517

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
b879d65
BeamDyn: return DTBeam from BD_Init; sub-step BeamDyn in the driver
andrew-platt Nov 10, 2026
6f77541
ModReg: check iostat on registry pack/unpack reads and writes
andrew-platt Nov 10, 2026
6ecfad7
BeamDyn: reset summary file unit in BD_End after closing it
andrew-platt Nov 10, 2026
5ec44bd
FAST_Solver: index tight-coupling module check through p%iModTC
andrew-platt Nov 10, 2026
fd1b6cc
OpenFAST: allow BeamDyn sub-stepping; remove DTBeam == DT check
andrew-platt Nov 10, 2026
145e182
BeamDyn: reject nonpositive DTBeam
andrew-platt Nov 10, 2026
343117d
ModReg: check for errors in RegUnpackPointer before using the index
andrew-platt Nov 10, 2026
8fa7a2a
BeamDyn driver: do not sub-step static analysis
andrew-platt Nov 10, 2026
6db3610
BeamDyn: use inputs interpolated to t+dt for the reference frame upda…
andrew-platt Nov 10, 2026
24fe391
BeamDyn: add bd_5MW_dynamic_substep regression test
andrew-platt Nov 10, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 6 additions & 0 deletions docs/source/user/beamdyn/input_files.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -266,6 +266,12 @@ keyword “DEFAULT” sets ``N_Fact = 5``.
time-integration in seconds. The keyword “DEFAULT” may be used to
indicate that the module should employ the time increment prescribed by
the driver code (FAST/stand-alone driver program).
``DTBeam`` must be equal to or an integer divisor of the glue-code or
driver time step ``DT``. When it is smaller, BeamDyn takes
``DT``/``DTBeam`` steps per glue-code or driver step, and output is
written every ``DT``. In OpenFAST with tight coupling (``ModCoupling`` =
2 or 3), the glue code integrates the BeamDyn states at ``DT`` and
``DTBeam`` is not used.

``load_retries`` specifies the maximum number of load retries allowed. This option
currently works only for static analysis. For every load retry, the applied load is
Expand Down
32 changes: 28 additions & 4 deletions modules/beamdyn/src/BeamDyn.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -127,6 +127,9 @@ SUBROUTINE BD_Init( InitInp, u, p, x, xd, z, OtherState, y, MiscVar, Interval, I
! this routine sets *some* of the parameters (basically the "easy" ones)
call SetParameters(InitInp, InputFileData, p, OtherState, ErrStat2, ErrMsg2); if (Failed()) return

! Return the time step BeamDyn uses (DTBeam) so the caller can step (or sub-step) consistently
Interval = p%dt
Comment thread
andrew-platt marked this conversation as resolved.
Resolved


! Temporary GLL point intrinsic coordinates array
CALL BD_GenerateGLL(p%nodes_per_elem,GLL_nodes,ErrStat2,ErrMsg2); if (Failed()) return
Expand Down Expand Up @@ -2010,7 +2013,10 @@ SUBROUTINE BD_End( u, p, x, xd, z, OtherState, y, m, ErrStat, ErrMsg )
! Place any last minute operations or calculations here:

! Close files here:
if (m%Un_Sum > 0) CLOSE(m%Un_Sum)
if (m%Un_Sum > 0) then
CLOSE(m%Un_Sum)
m%Un_Sum = -1
endif


! Destroy the input data:
Expand Down Expand Up @@ -2060,6 +2066,7 @@ SUBROUTINE BD_UpdateStates( t, n, u, utimes, p, x, xd, z, OtherState, m, ErrStat
INTEGER(IntKi), INTENT( OUT) :: ErrStat !< Error status of the operation
CHARACTER(*), INTENT( OUT) :: ErrMsg !< Error message if ErrStat /= ErrID_None

TYPE(BD_InputType) :: u_interp ! Inputs interpolated to t+dt (used only when sub-stepping)
INTEGER(IntKi) :: ErrStat2 ! Temporary Error status
CHARACTER(ErrMsgLen) :: ErrMsg2 ! Temporary Error message

Expand All @@ -2071,9 +2078,26 @@ SUBROUTINE BD_UpdateStates( t, n, u, utimes, p, x, xd, z, OtherState, m, ErrStat
CALL BD_GA2( t, n, u, utimes, p, x, xd, z, OtherState, m, ErrStat2, ErrMsg2 )
call SetErrStat(ErrStat2,ErrMsg2,ErrStat,ErrMsg,''); if (ErrStat >= AbortErrLev) return

! change reference frame to root motion at t=T+DT (u(1)%RootMotionMesh)
call BD_UpdateGlobalRef(u(1),p,x,OtherState,ErrStat2,ErrMsg2)
call SetErrStat(ErrStat2,ErrMsg2,ErrStat,ErrMsg,'')
! change reference frame to root motion at t+dt
if (EqualRealNos(utimes(1), t + p%dt)) then
! u(1) is at t+dt (BeamDyn time step is the same as the caller's)
call BD_UpdateGlobalRef(u(1),p,x,OtherState,ErrStat2,ErrMsg2)
call SetErrStat(ErrStat2,ErrMsg2,ErrStat,ErrMsg,'')
else
! BeamDyn is being sub-stepped, so u(1) is at the end of the caller's step, not at t+dt.
! Interpolate the (global frame) inputs to t+dt so the root motion matches the states.
call BD_CopyInput(u(1), u_interp, MESH_NEWCOPY, ErrStat2, ErrMsg2)
call SetErrStat(ErrStat2,ErrMsg2,ErrStat,ErrMsg,'')
if (ErrStat < AbortErrLev) then
call BD_Input_ExtrapInterp(u, utimes, u_interp, t + p%dt, ErrStat2, ErrMsg2)
call SetErrStat(ErrStat2,ErrMsg2,ErrStat,ErrMsg,'')
end if
if (ErrStat < AbortErrLev) then
call BD_UpdateGlobalRef(u_interp,p,x,OtherState,ErrStat2,ErrMsg2)
call SetErrStat(ErrStat2,ErrMsg2,ErrStat,ErrMsg,'')
end if
call BD_DestroyInput(u_interp, ErrStat2, ErrMsg2)
end if
ELSE !IF(p%analysis_type == BD_STATIC_ANALYSIS) THEN
CALL BD_Static( t, u, utimes, p, x, OtherState, m, ErrStat, ErrMsg )
ENDIF
Expand Down
4 changes: 4 additions & 0 deletions modules/beamdyn/src/BeamDyn_IO.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1534,6 +1534,10 @@ SUBROUTINE BD_ValidateInputData( InitInp, InputFileData, ErrStat, ErrMsg )
IF(InputFileData%QuasiStaticInit .AND. (.NOT. InitInp%DynamicSolve)) &
CALL SetErrStat ( ErrID_Warn, 'QuasiStaticInit option only valid with dynamic solve. Ignoring.', ErrStat, ErrMsg, RoutineName )

! Time step must be positive for dynamic solves (not used in static solves)
IF(InitInp%DynamicSolve .AND. InputFileData%DTBeam <= 0.0_DbKi) &
CALL SetErrStat ( ErrID_Fatal, 'DTBeam must be greater than 0.', ErrStat, ErrMsg, RoutineName )

IF(InputFileData%rhoinf .LT. 0.0_BDKi .OR. InputFileData%rhoinf .GT. 1.0_BDKi) &
CALL SetErrStat ( ErrID_Fatal, 'Numerical damping parameter \rho_{inf} must be in the range of [0.0,1.0]', ErrStat, ErrMsg, RoutineName )
IF(InputFileData%quadrature .NE. GAUSS_QUADRATURE .AND. InputFileData%quadrature .NE. TRAP_QUADRATURE) &
Expand Down
87 changes: 66 additions & 21 deletions modules/beamdyn/src/Driver_Beam.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -32,6 +32,11 @@ PROGRAM BeamDyn_Driver_Program
CHARACTER(1024) :: ErrMsg ! Error message if ErrStat /= ErrID_None
REAL(DbKi) :: dt_global ! fixed/constant global time step
REAL(DbKi) :: t_global ! global-loop time marker
REAL(DbKi) :: dt_BD ! BeamDyn time step (DTBeam, returned from BD_Init)
REAL(DbKi) :: t_BD ! BeamDyn sub-step time marker
INTEGER(IntKi) :: n_substeps ! number of BeamDyn time steps per driver time step
INTEGER(IntKi) :: j_ss ! BeamDyn sub-step counter
INTEGER(IntKi) :: n_t_BD ! BeamDyn time step counter
INTEGER(IntKi) :: n_t_final ! total number of time steps
INTEGER(IntKi) :: n_t_global ! global-loop time counter
INTEGER(IntKi) :: n_t_vtk ! global vtk step counter
Expand Down Expand Up @@ -112,6 +117,9 @@ PROGRAM BeamDyn_Driver_Program
ALLOCATE(BD_Input(BD_interp_order + 1))
ALLOCATE(BD_InputTimes(BD_interp_order + 1))

! Suggested BeamDyn time step. BD_Init returns the time step BeamDyn actually uses (DTBeam)
dt_BD = dt_global

CALL BD_Init(BD_InitInput &
, BD_Input(1) &
, BD_Parameter &
Expand All @@ -121,11 +129,15 @@ PROGRAM BeamDyn_Driver_Program
, BD_OtherState &
, BD_Output &
, BD_MiscVar &
, dt_global &
, dt_BD &
, BD_InitOutput &
, ErrStat &
, ErrMsg )
CALL CheckError()

! BeamDyn may use a smaller time step (DTBeam) than the driver. Sub-step BeamDyn in that case.
call SetBDSubSteps()
CALL CheckError()

! If the Quasi-Static solve is in use, rerun the initialization with loads at t=0
! (HACK: set in the driver only because computing Jacobians with this option [as in FAST glue code] is problematic)
Expand Down Expand Up @@ -155,7 +167,7 @@ PROGRAM BeamDyn_Driver_Program
CALL CheckError()

! solve for inputs at previous time steps
BD_InputTimes(j) = DvrData%t_initial - (j - 1) * dt_global
BD_InputTimes(j) = DvrData%t_initial - (j - 1) * dt_BD
CALL BD_InputSolve( BD_InputTimes(j), BD_Input(j), DvrData, ErrStat, ErrMsg)
CALL CheckError()
END DO
Expand Down Expand Up @@ -194,27 +206,30 @@ PROGRAM BeamDyn_Driver_Program

DO n_t_global = 0, n_t_final

! Shift "window" of BD_Input
DO j = BD_interp_order, 1, -1
CALL BD_CopyInput (BD_Input(j), BD_Input(j+1), MESH_UPDATECOPY, Errstat, ErrMsg)
IF(BD_Parameter%analysis_type .EQ. BD_STATIC_ANALYSIS .AND. n_t_global > max_ld_step) EXIT

! update states from n_t_global to n_t_global + 1, sub-stepping BeamDyn at its own time step
DO j_ss = 1, n_substeps
n_t_BD = n_t_global*n_substeps + j_ss - 1
t_BD = n_t_BD * dt_BD + DvrData%t_initial
Comment thread
andrew-platt marked this conversation as resolved.
Resolved

! Shift "window" of BD_Input
DO j = BD_interp_order, 1, -1
CALL BD_CopyInput (BD_Input(j), BD_Input(j+1), MESH_UPDATECOPY, Errstat, ErrMsg)
CALL CheckError()
BD_InputTimes(j+1) = BD_InputTimes(j)
END DO

BD_InputTimes(1) = t_BD + dt_BD
CALL BD_InputSolve( BD_InputTimes(1), BD_Input(1), DvrData, ErrStat, ErrMsg)
CALL CheckError()
BD_InputTimes(j+1) = BD_InputTimes(j)
END DO

BD_InputTimes(1) = t_global + dt_global
CALL BD_InputSolve( BD_InputTimes(1), BD_Input(1), DvrData, ErrStat, ErrMsg)
CALL CheckError()


IF(BD_Parameter%analysis_type .EQ. BD_STATIC_ANALYSIS .AND. n_t_global > max_ld_step) EXIT

! update states from n_t_global to n_t_global + 1
CALL BD_UpdateStates( t_global, n_t_global, BD_Input, BD_InputTimes, BD_Parameter, &
BD_ContinuousState, &
BD_DiscreteState, BD_ConstraintState, &
BD_OtherState, BD_MiscVar, ErrStat, ErrMsg )
CALL CheckError()

CALL BD_UpdateStates( t_BD, n_t_BD, BD_Input, BD_InputTimes, BD_Parameter, &
BD_ContinuousState, &
BD_DiscreteState, BD_ConstraintState, &
BD_OtherState, BD_MiscVar, ErrStat, ErrMsg )
CALL CheckError()
END DO

! advance time
t_global = (n_t_global+1) * dt_global + DvrData%t_initial
Expand Down Expand Up @@ -247,6 +262,36 @@ PROGRAM BeamDyn_Driver_Program
CALL Dvr_End()

CONTAINS
!----------------------------------------------------------------------------------------------------------------------------------
!> Set the number of BeamDyn sub-steps per driver time step. The BeamDyn time step (DTBeam) must be
!! equal to or an integer divisor of the driver time step. Static analysis does not use DTBeam (each
!! driver step is a load step), so it is not sub-stepped.
SUBROUTINE SetBDSubSteps()
if (BD_Parameter%analysis_type == BD_STATIC_ANALYSIS) then
n_substeps = 1
dt_BD = dt_global
elseif (EqualRealNos(dt_BD, dt_global)) then
n_substeps = 1
dt_BD = dt_global
elseif (dt_BD <= 0.0_DbKi) then
ErrStat = ErrID_Fatal
ErrMsg = 'BeamDyn time step DTBeam ('//trim(Num2LStr(dt_BD))//' s) must be greater than 0.'
elseif (dt_BD > dt_global) then
ErrStat = ErrID_Fatal
ErrMsg = 'BeamDyn time step DTBeam ('//trim(Num2LStr(dt_BD))//' s) cannot be larger than the driver time step DT ('// &
trim(Num2LStr(dt_global))//' s).'
else
n_substeps = nint(dt_global / dt_BD)
if (.not. EqualRealNos(dt_global, dt_BD*n_substeps)) then
ErrStat = ErrID_Fatal
ErrMsg = 'BeamDyn time step DTBeam ('//trim(Num2LStr(dt_BD))//' s) must be an integer divisor of the driver time step DT ('// &
trim(Num2LStr(dt_global))//' s).'
else
call WrScr(' BeamDyn time step DTBeam ('//trim(Num2LStr(dt_BD))//' s) is smaller than the driver time step DT ('// &
trim(Num2LStr(dt_global))//' s). Using '//trim(Num2LStr(n_substeps))//' BeamDyn sub-steps per driver step.')
endif
endif
END SUBROUTINE SetBDSubSteps
!----------------------------------------------------------------------------------------------------------------------------------
SUBROUTINE Dvr_End()

Expand Down
49 changes: 42 additions & 7 deletions modules/nwtc-library/ModRegGen.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -261,11 +261,24 @@
type(c_ptr), intent(out) :: Ptr
integer(B8Ki), intent(out) :: Idx

! Initialize outputs
Ptr = c_null_ptr
Idx = 0

! If error, return
if (RF%ErrStat /= ErrID_None) return

! Unpack pointer index
! Unpack pointer index, return if error
call RegUnpack(RF, Idx)
if (RegCheckErr(RF, "RegUnpackPointer")) return

! Check that index is within pointer index bounds
if (Idx < 1 .or. Idx > size(RF%Pointers, kind=B8Ki)) then
RF%ErrStat = ErrID_Fatal
write(RF%ErrMsg, *) "RegUnpackPointer: pointer index", Idx, "out of range 1 to", size(RF%Pointers, kind=B8Ki)
Idx = 0
return
end if

! Get pointer from index
Ptr = RF%Pointers(Idx)
Expand Down Expand Up @@ -303,12 +316,19 @@
type(RegFile), intent(inout) :: RF
logical :: match
integer(B8Ki) :: DataNum
integer(IntKi) :: stat

! Increment the data number to be read
RF%NumData = RF%NumData + 1

! Read the data number from the file
read(RF%Unit) DataNum
! Read the data number from the file, return if error
read(RF%Unit, iostat=stat) DataNum
if (stat /= 0) then
match = .false.
RF%ErrStat = ErrID_Fatal
write(RF%ErrMsg, *) "Error reading data number", RF%NumData, "from file, iostat =", stat
return
end if

! If data number from file does not match expected number, set match false
! and create error message; otherwise, set match to true
Expand All @@ -329,16 +349,26 @@ def gen_pack(w, dt, decl, rank):
w.write(f'\n\n subroutine {name}(RF, Data)')
w.write(f'\n type(RegFile), intent(inout) :: RF')
w.write(f'\n {decl+", intent(in)":<35s} :: Data{dims}')
w.write(f'\n integer(IntKi) :: stat')
w.write(f'\n')
w.write(f'\n ! If error, return')
w.write(f'\n if (RF%ErrStat /= ErrID_None) return')
w.write(f'\n')
w.write(f'\n ! Increment data number and write to file')
w.write(f'\n RF%NumData = RF%NumData + 1')
w.write(f'\n write(RF%Unit) RF%NumData')
w.write(f'\n write(RF%Unit, iostat=stat) RF%NumData')
w.write(f'\n if (stat /= 0) then')
w.write(f'\n RF%ErrStat = ErrID_Fatal')
w.write(f'\n write(RF%ErrMsg, *) "{name}: error writing data number", RF%NumData, "to file, iostat =", stat')
w.write(f'\n return')
w.write(f'\n end if')
w.write(f'\n')
w.write(f'\n ! Write data to file')
w.write(f'\n write(RF%Unit) Data')
w.write(f'\n write(RF%Unit, iostat=stat) Data')
w.write(f'\n if (stat /= 0) then')
w.write(f'\n RF%ErrStat = ErrID_Fatal')
w.write(f'\n write(RF%ErrMsg, *) "{name}: error writing data", RF%NumData, "to file, iostat =", stat')
w.write(f'\n end if')
w.write(f'\n end subroutine')


Expand All @@ -349,15 +379,20 @@ def gen_unpack(w, dt, decl, rank):
w.write(f'\n subroutine {name}(RF, Data)')
w.write(f'\n type(RegFile), intent(inout) :: RF')
w.write(f'\n {decl+", intent(out)":<35s} :: Data{dims}')
w.write(f'\n integer(IntKi) :: stat')
w.write(f'\n')
w.write(f'\n ! If error, return')
w.write(f'\n if (RF%ErrStat /= ErrID_None) return')
w.write(f'\n')
w.write(f'\n ! Read data number, return if invalid')
w.write(f'\n if (.not. DataNumValid(RF)) return')
w.write(f'\n')
w.write(f'\n ! Read data from file')
w.write(f'\n read(RF%Unit) Data')
w.write(f'\n ! Read data from file, return if error')
w.write(f'\n read(RF%Unit, iostat=stat) Data')
w.write(f'\n if (stat /= 0) then')
w.write(f'\n RF%ErrStat = ErrID_Fatal')
w.write(f'\n write(RF%ErrMsg, *) "{name}: error reading data", RF%NumData, "from file, iostat =", stat')
w.write(f'\n end if')
Comment thread
andrew-platt marked this conversation as resolved.
Resolved
w.write(f'\n end subroutine')

def gen_pack_alloc(w, dt, decl, rank):
Expand Down
Loading
Loading