Skip to content

Commit c847dbc

Browse files
andrew-plattandrew-platt
authored
Merge pull request #2975 from bjonkman/b/MD_updates
Minor numerical improvements
2 parents 3750d8f + eab1548 commit c847dbc

File tree

12 files changed

+205
-161
lines changed

12 files changed

+205
-161
lines changed

‎modules/aerodyn/src/AirfoilInfo.f90‎

Lines changed: 1 addition & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -1745,12 +1745,7 @@ subroutine AFI_ComputeAirfoilCoefs1D( AOA, p, AFI_interp, errStat, errMsg, Table
17451745

17461746

17471747
! Spline interpolation of lower table based on requested AOA
1748-
1749-
IntAFCoefs(1:s1) = CubicSplineInterpM( Alpha &
1750-
, p%Table(iTab)%Alpha &
1751-
, p%Table(iTab)%Coefs &
1752-
, p%Table(iTab)%SplineCoefs &
1753-
, ErrStat, ErrMsg )
1748+
CALL CubicSplineInterpM( Alpha, p%Table(iTab)%Alpha, p%Table(iTab)%Coefs, p%Table(iTab)%SplineCoefs, IntAFCoefs(1:s1) )
17541749
end if
17551750

17561751
AFI_interp%Cl = IntAFCoefs(p%ColCl)

‎modules/hydrodyn/src/HydroDyn.txt‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -71,7 +71,7 @@ typedef HydroDyn/HydroDyn InitInputType CHARACTER(1
7171
typedef ^ ^ LOGICAL UseInputFile - .TRUE. - "Supplied by Driver: .TRUE. if using a input file, .FALSE. if all inputs are being passed in by the caller" -
7272
typedef ^ ^ FileInfoType PassedFileData - - - "If we don't use the input file, pass everything through this" -
7373
typedef ^ ^ CHARACTER(1024) OutRootName - - - "Supplied by Driver: The name of the root file (without extension) including the full path" -
74-
typedef ^ ^ Logical Linearize - .FALSE. - "Flag that tells this module if the glue code wants to linearize." -
74+
typedef ^ ^ Logical Linearize - .FALSE. - "Flag that tells this module if the glue code wants to linearize." -
7575
typedef ^ ^ ReKi Gravity - - - "Supplied by Driver: Gravitational acceleration" "(m/s^2)"
7676
typedef ^ ^ DbKi TMax - - - "Supplied by Driver: The total simulation time" "(sec)"
7777
typedef ^ ^ logical VisMeshes - .false. - "Output visualization meshes" -

‎modules/hydrodyn/src/Morison.f90‎

Lines changed: 1 addition & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -4988,8 +4988,7 @@ SUBROUTINE GetSectionHstLds_Rec( origin, pos0, k_hat, x_hat, y_hat, Sa, Sb, dSad
49884988
s2 = 0.0_DbKi;
49894989
end if
49904990

4991-
dFdl(1:3) = dFdl(1:3) + &
4992-
-n(:,i) * ( z1*(s2-s1) + 0.5_DbKi*(z2-z1)/s*(s2*s2-s1*s1) )
4991+
dFdl(1:3) = dFdl(1:3) -n(:,i) * ( z1*(s2-s1) + 0.5_DbKi*(z2-z1)/s*(s2*s2-s1*s1) )
49934992
C(1) = (z2-z1)*(x2-x1)/3.0_DbKi/(s*s)*(s2**3-s1**3) + 0.5_DbKi*((z2-z1)*(x1-x0)/s+(x2-x1)*z1/s)*(s2*s2-s1*s1) + z1*(x1-x0)*(s2-s1)
49944993
C(2) = (z2-z1)*(y2-y1)/3.0_DbKi/(s*s)*(s2**3-s1**3) + 0.5_DbKi*((z2-z1)*(y1-y0)/s+(y2-y1)*z1/s)*(s2*s2-s1*s1) + z1*(y1-y0)*(s2-s1)
49954994
C(3) = (z2-z1)*(z2-z1)/3.0_DbKi/(s*s)*(s2**3-s1**3) + 0.5_DbKi*((z2-z1)*(z1-z0)/s+(z2-z1)*z1/s)*(s2*s2-s1*s1) + z1*(z1-z0)*(s2-s1)

‎modules/moordyn/src/MoorDyn_Driver.f90‎

Lines changed: 57 additions & 57 deletions
Original file line numberDiff line numberDiff line change
@@ -78,15 +78,14 @@ PROGRAM MoorDyn_Driver
7878
! SeaState types
7979
TYPE(SeaSt_InitInputType) :: InitInData_SeaSt ! Input data for initialization
8080
TYPE(SeaSt_InitOutputType) :: InitOutData_SeaSt ! Output data from initialization
81-
type(SeaSt_ContinuousStateType) :: x_SeaSt ! Continuous states
82-
type(SeaSt_DiscreteStateType) :: xd_SeaSt ! Discrete states
83-
type(SeaSt_ConstraintStateType) :: z_SeaSt ! Constraint states
84-
type(SeaSt_OtherStateType) :: OtherState_SeaSt ! Other states
85-
type(SeaSt_MiscVarType) :: m_SeaSt ! Misc/optimization variables
86-
type(SeaSt_ParameterType) :: p_SeaSt ! Parameters
87-
type(SeaSt_InputType) :: u_SeaSt(1) ! System inputs
88-
type(SeaSt_OutputType) :: y_SeaSt ! System outputs
89-
LOGICAL :: SeaState_Initialized = .FALSE.
81+
type(SeaSt_ContinuousStateType) :: x_SeaSt ! Continuous states
82+
type(SeaSt_DiscreteStateType) :: xd_SeaSt ! Discrete states
83+
type(SeaSt_ConstraintStateType) :: z_SeaSt ! Constraint states
84+
type(SeaSt_OtherStateType) :: OtherState_SeaSt ! Other states
85+
type(SeaSt_MiscVarType) :: m_SeaSt ! Misc/optimization variables
86+
type(SeaSt_ParameterType) :: p_SeaSt ! Parameters
87+
type(SeaSt_InputType) :: u_SeaSt(1) ! System inputs
88+
type(SeaSt_OutputType) :: y_SeaSt ! System outputs
9089

9190
! Motion file parsing
9291
type(FileInfoType) :: FileInfo_PrescribeMtn !< The derived type for holding the prescribed forces input file for parsing -- we may pass this in the future
@@ -138,6 +137,7 @@ PROGRAM MoorDyn_Driver
138137
ErrStat = ErrID_None
139138
UnEcho=-1 ! set to -1 as echo is no longer used by MD
140139
UnIn =-1
140+
141141

142142
! TODO: Sort out error handling (two sets of flags currently used)
143143

@@ -210,8 +210,6 @@ PROGRAM MoorDyn_Driver
210210

211211
! allocate Input and Output arrays; used for interpolation and extrapolation
212212
Allocate(MD_uTimes(MD_interp_order + 1))
213-
214-
! @bonnie : This is in the FAST developers glue code example, but it's probably not needed here.
215213
Allocate(MD_u(MD_interp_order + 1))
216214

217215

@@ -237,8 +235,8 @@ PROGRAM MoorDyn_Driver
237235
InitInData_SeaSt%TMax = MD_InitInp%TMax
238236
InitInData_SeaSt%Linearize = MD_InitInp%Linearize
239237

240-
CALL SeaSt_Init( InitInData_SeaSt, u_SeaSt(1), p_SeaSt, x_SeaSt, xd_SeaSt, z_SeaSt, OtherState_SeaSt, y_SeaSt, m_SeaSt, dtC, InitOutData_SeaSt, ErrStat2, ErrMsg2 ); call AbortIfFailed()
241-
SeaState_Initialized = .TRUE.
238+
CALL SeaSt_Init( InitInData_SeaSt, u_SeaSt(1), p_SeaSt, x_SeaSt, xd_SeaSt, z_SeaSt, OtherState_SeaSt, y_SeaSt, m_SeaSt, dtC, InitOutData_SeaSt, ErrStat2, ErrMsg2 )
239+
call AbortIfFailed()
242240

243241
IF ( dtC /= drvrInitInp%dtC) THEN
244242
ErrMsg2 = 'The SeaState Module attempted to change the coupling timestep, but this is not allowed. The SeaState Module must use the Driver coupling timestep.'
@@ -252,13 +250,10 @@ PROGRAM MoorDyn_Driver
252250
END IF
253251

254252
! call the initialization routine
255-
CALL MD_Init( MD_InitInp, MD_u(1), MD_p, MD_x , MD_xd, MD_xc, MD_xo, MD_y, MD_m, dtC, MD_InitOut, ErrStat2, ErrMsg2 ); call AbortIfFailed()
253+
CALL MD_Init( MD_InitInp, MD_u(1), MD_p, MD_x , MD_xd, MD_xc, MD_xo, MD_y, MD_m, dtC, MD_InitOut, ErrStat2, ErrMsg2 )
254+
call AbortIfFailed()
256255

257-
CALL MD_DestroyInitInput ( MD_InitInp , ErrStat2, ErrMsg2 ); call AbortIfFailed()
258-
CALL MD_DestroyInitOutput ( MD_InitOut , ErrStat2, ErrMsg2 ); call AbortIfFailed()
259-
260-
CALL DispNVD( MD_InitOut%Ver )
261-
256+
CALL DispNVD( MD_InitOut%Ver )
262257

263258
! determine number of input channels expected from driver input file time series (DOFs including active tensioning channels)
264259
if (allocated(MD_u(1)%DeltaL)) then
@@ -693,50 +688,47 @@ PROGRAM MoorDyn_Driver
693688
CALL RunTimes( ProgStrtTime, ProgStrtCPU, SimStrtTime, SimStrtCPU, t )
694689

695690
! Destroy all objects
696-
IF (SeaState_Initialized) THEN
697-
CALL SeaSt_End( u_SeaSt(1), p_SeaSt, x_SeaSt, xd_SeaSt, z_SeaSt, OtherState_SeaSt, y_SeaSt, m_SeaSt, ErrStat2, ErrMsg2); call AbortIfFailed()
698-
ENDIF
699-
CALL MD_End( MD_u(1), MD_p, MD_x, MD_xd, MD_xc , MD_xo, MD_y, MD_m, ErrStat2, ErrMsg2 ); call AbortIfFailed()
700-
701-
do j = 2,MD_interp_order+1
702-
call MD_DestroyInput( MD_u(j), ErrStat2, ErrMsg2)
703-
end do
704-
705-
if ( ErrStat /= ErrID_None ) THEN ! Display all errors
706-
CALL WrScr1( "Errors: " )
707-
CALL WrScr( trim(GetErrStr(ErrStat))//': '//trim(ErrMsg) )
708-
endif
709-
710-
!close (un)
711-
call CleanUp()
691+
call EndAndCleanUp()
712692
CALL NormStop()
713693

714694

715695
CONTAINS
716-
696+
!-------------------------------------------------------------------------------------------------------------------------------
717697
SUBROUTINE AbortIfFailed()
718-
719-
call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver')
720-
721-
if (ErrStat >= AbortErrLev) then
722-
if (SeaState_Initialized) then
723-
call SeaSt_End( u_SeaSt(1), p_SeaSt, x_SeaSt, xd_SeaSt, z_SeaSt, OtherState_SeaSt, y_SeaSt, m_SeaSt, ErrStat2, ErrMsg2)
724-
call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver' )
725-
end if
726-
727-
CALL SeaSt_DestroyInitOutput( InitOutData_SeaSt, ErrStat2, ErrMsg2 )
728-
call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver' )
729-
CALL SeaSt_DestroyInitInput( InitInData_SeaSt, ErrStat2, ErrMsg2 )
730-
call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver' )
731-
732-
call CleanUp()
733-
Call ProgAbort(trim(ErrMsg))
734-
elseif ( ErrStat2 /= ErrID_None ) THEN
735-
CALL WrScr1( trim(GetErrStr(ErrStat2))//': '//trim(ErrMsg2)//NewLine)
736-
end if
698+
699+
if (ErrStat >= AbortErrLev .OR. ErrStat2 >= AbortErrLev) then
700+
call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver')
701+
702+
call EndAndCleanUp()
703+
Call ProgAbort(trim(ErrMsg))
704+
elseif ( ErrStat2 /= ErrID_None ) THEN ! print messages as we get them (but don't call SetErrStat or they will be printed 2x)
705+
CALL WrScr1( trim(GetErrStr(ErrStat2))//': '//trim(ErrMsg2)//NewLine)
706+
end if
707+
737708
END SUBROUTINE AbortIfFailed
709+
!-------------------------------------------------------------------------------------------------------------------------------
710+
SUBROUTINE EndAndCleanUp()
711+
CALL SeaSt_DestroyInitOutput( InitOutData_SeaSt, ErrStat2, ErrMsg2 )
712+
call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver' )
713+
CALL SeaSt_DestroyInitInput( InitInData_SeaSt, ErrStat2, ErrMsg2 )
714+
call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver' )
715+
716+
CALL MD_DestroyInitOutput( MD_InitOut, ErrStat2, ErrMsg2 )
717+
call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver' )
718+
CALL MD_DestroyInitInput( MD_InitInp, ErrStat2, ErrMsg2 )
719+
call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver' )
720+
721+
call SeaSt_End( u_SeaSt(1), p_SeaSt, x_SeaSt, xd_SeaSt, z_SeaSt, OtherState_SeaSt, y_SeaSt, m_SeaSt, ErrStat2, ErrMsg2)
722+
call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver' )
723+
724+
CALL MD_End( MD_u(1), MD_p, MD_x, MD_xd, MD_xc , MD_xo, MD_y, MD_m, ErrStat2, ErrMsg2 )
725+
call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver' )
726+
727+
do j = 2,MD_interp_order+1
728+
call MD_DestroyInput( MD_u(j), ErrStat2, ErrMsg2)
729+
call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'MoorDyn_Driver' )
730+
end do
738731

739-
SUBROUTINE CleanUp()
740732
if(UnEcho >0) CLOSE( UnEcho )
741733
if(UnIn >0) CLOSE( UnIn )
742734

@@ -749,7 +741,15 @@ SUBROUTINE CleanUp()
749741
IF (ALLOCATED(rd_in2 )) DEALLOCATE(rd_in2 )
750742
IF (ALLOCATED(rdd_in )) DEALLOCATE(rdd_in )
751743
IF (ALLOCATED(rdd_in2 )) DEALLOCATE(rdd_in2 )
752-
END SUBROUTINE CleanUp
744+
IF (ALLOCATED(TmpRe )) DEALLOCATE(TmpRe )
745+
746+
747+
748+
if ( ErrStat /= ErrID_None ) THEN ! Display all errors
749+
CALL WrScr1( "Errors: " )
750+
CALL WrScr( trim(GetErrStr(ErrStat))//': '//trim(ErrMsg) )
751+
endif
752+
END SUBROUTINE EndAndCleanUp
753753

754754
!-------------------------------------------------------------------------------------------------------------------------------
755755
SUBROUTINE ReadDriverInputFile( inputFile, InitInp)

‎modules/moordyn/src/MoorDyn_Misc.f90‎

Lines changed: 5 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -1071,7 +1071,6 @@ SUBROUTINE setupWaterKin(WaterKinString, p, Tmax, ErrStat, ErrMsg)
10711071
REAL(SiKi) :: WaveDOmega ! frequency step
10721072
REAL(SiKi), ALLOCATABLE :: SinWaveDir(:) ! SIN( WaveDirArr(I) ) -- Each wave frequency has a unique wave direction.
10731073
REAL(SiKi), ALLOCATABLE :: CosWaveDir(:) ! COS( WaveDirArr(I) ) -- Each wave frequency has a unique wave direction.
1074-
LOGICAL :: WaveMultiDir = .FALSE. ! Indicates the waves are multidirectional -- set by WaveField pointer if enabled
10751074

10761075
REAL(SiKi), ALLOCATABLE :: TmpFFTWaveElev(:) ! Data for the FFT calculation
10771076
TYPE(FFT_DataType) :: FFT_Data ! the instance of the FFT module we're using
@@ -1083,7 +1082,6 @@ SUBROUTINE setupWaterKin(WaterKinString, p, Tmax, ErrStat, ErrMsg)
10831082
COMPLEX(SiKi), PARAMETER :: ImagNmbr = (0.0,1.0) ! The imaginary number, SQRT(-1.0)
10841083
COMPLEX(SiKi) :: ImagOmega ! = ImagNmbr*Omega (rad/s)
10851084
REAL(DbKi), ALLOCATABLE :: WaveNmbr(:) ! wave number for frequency array
1086-
REAL(SiKi), ALLOCATABLE :: WaveDirArr(:) ! Wave direction array. Each frequency has a unique direction of WaveNDir > 1 (degrees). 0's for WaveKin = 1 or if disabled in SeaState.
10871085
REAL(SiKi), ALLOCATABLE :: WaveElevC0(:,:) ! Discrete Fourier transform of the instantaneous elevation of incident waves at the ref point (meters)
10881086
COMPLEX(SiKi), ALLOCATABLE :: WaveElevC( :) ! Discrete Fourier transform of the instantaneous elevation of incident waves at the ref point (meters)
10891087
COMPLEX(SiKi), ALLOCATABLE :: WaveAccCHx(:) ! Discrete Fourier transform of the instantaneous horizontal acceleration in x-direction of incident waves before applying stretching at the zi-coordinates for points (m/s^2)
@@ -1412,7 +1410,7 @@ SUBROUTINE setupWaterKin(WaterKinString, p, Tmax, ErrStat, ErrMsg)
14121410
END IF
14131411

14141412
! Warning check to make sure SeaState and MoorDyn have the same wave dir. For now, no wave spreading. This can be updated
1415-
IF (p%WaveField%WaveDir /= WaveDir) THEN
1413+
IF (p%WaveField%WaveDir /= WaveDir) THEN !bjj: the local WaveDir doesn't appear to be used when WaveKin = 2 and p%WaveField%WaveDirArr is true, so I don't think this error message is completely valid.
14161414
IF (p%writeLog > 0) THEN
14171415
WRITE(p%UnLog, '(A)' ) " WARNING SeaState WaveDir does not match MoorDyn WaveDir. Using MoorDyn values for interpolating SeaState data to MoorDyn grid."
14181416
ENDIF
@@ -1422,7 +1420,7 @@ SUBROUTINE setupWaterKin(WaterKinString, p, Tmax, ErrStat, ErrMsg)
14221420
! Info check for if MoorDyn dtWave is non-zero. Users may set this accidentially, or could be left over in an input file
14231421
IF (p%dtWave > 0) THEN
14241422
IF (p%writeLog > 0) THEN
1425-
WRITE(p%UnLog, '(A)' ) " MoorDyn dtWave is ignored when using WaveKinMod = 2 becasue wave frequency information is supplied by SeaState"
1423+
WRITE(p%UnLog, '(A)' ) " MoorDyn dtWave is ignored when using WaveKinMod = 2 because wave frequency information is supplied by SeaState"
14261424
ENDIF
14271425
END IF
14281426

@@ -1445,15 +1443,8 @@ SUBROUTINE setupWaterKin(WaterKinString, p, Tmax, ErrStat, ErrMsg)
14451443
NStepWave = p%WaveField%NStepWave
14461444

14471445
! Pull some other things out of the WaveField pointer
1448-
WaveMultiDir = p%WaveField%WaveMultiDir
14491446
p%ntWave = NStepWave ! set ntWave to NStepWave
14501447

1451-
! Set wave spreading array if enabled in SeaState, otherwise set to zero
1452-
If (WaveMultiDir) THEN
1453-
! Note: allocations not needed here because they are already allocated in SeaState
1454-
WaveDirArr = p%WaveField%WaveDirArr
1455-
ENDIF
1456-
14571448
ELSEIF (p%WaveKin == 1) THEN ! must be a filepath therefore read wave elevations from timeseries
14581449

14591450
! NOTE: there is a decent ammount of code duplication (intentional for now) with what is in SeaState that eventually
@@ -1659,9 +1650,9 @@ SUBROUTINE setupWaterKin(WaterKinString, p, Tmax, ErrStat, ErrMsg)
16591650
ALLOCATE ( SinWaveDir(0:NStepWave2), STAT=ErrStat2); ErrMsg2 = 'Cannot allocate SinWaveDir.'; IF (Failed0()) RETURN
16601651

16611652
! Set the CosWaveDir and SinWaveDir values.
1662-
IF (WaveMultiDir) THEN ! This is only possible with WaveKinMod = 2
1663-
CosWaveDir=COS(D2R*WaveDirArr)
1664-
SinWaveDir=SIN(D2R*WaveDirArr)
1653+
IF (p%WaveKin == 2 .and. p%WaveField%WaveMultiDir) THEN ! This is only possible with WaveKinMod = 2
1654+
CosWaveDir=COS(D2R*p%WaveField%WaveDirArr)
1655+
SinWaveDir=SIN(D2R*p%WaveField%WaveDirArr)
16651656
ELSE
16661657
CosWaveDir=COS(D2R*WaveDir)
16671658
SinWaveDir=SIN(D2R*WaveDir)

‎modules/nwtc-library/src/ModMesh.f90‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -2928,7 +2928,7 @@ SUBROUTINE PackMotionMesh(M, Ary, indx_first, FieldMask, TrimOP)
29282928

29292929
if (PackForTrimSolution) then
29302930
do i=1,M%NNodes
2931-
call DCM_logMap(M%Orientation(:,:,i), logmap, ErrStat2, ErrMsg2) !NOTE: we cannot use GetSmllRotAngs because we CANNOT assume that all DCMs in the code are small.
2931+
logmap = EulerExtract(M%Orientation(:,:,i)) !NOTE: we cannot use GetSmllRotAngs because we CANNOT assume that all DCMs in the code are small.
29322932
do k=1,3
29332933
Ary(indx_first) = logmap(k)
29342934
indx_first = indx_first + 1

0 commit comments

Comments
 (0)