Skip to content

Commit

Permalink
AD: adding skewModRedistr to BEMT (partial)
Browse files Browse the repository at this point in the history
  • Loading branch information
ebranlard committed Dec 1, 2023
1 parent 0ada500 commit f20aaeb
Show file tree
Hide file tree
Showing 8 changed files with 43 additions and 20 deletions.
3 changes: 2 additions & 1 deletion modules/aerodyn/src/AeroDyn.f90
Original file line number Diff line number Diff line change
Expand Up @@ -3847,7 +3847,7 @@ SUBROUTINE ValidateInputData( InitInp, InputFileData, NumBl, ErrStat, ErrMsg )
if ( InputFileData%IndToler < 0.0 .or. EqualRealNos(InputFileData%IndToler, 0.0_ReKi) ) &
call SetErrStat( ErrID_Fatal, 'IndToler must be greater than 0.', ErrStat, ErrMsg, RoutineName )

if ( InputFileData%Skew_Mod /= Skew_Mod_Orthogonal .and. InputFileData%Skew_Mod /= Skew_Mod_None .and. InputFileData%Skew_Mod /= Skew_Mod_Glauert) &
if ( InputFileData%Skew_Mod /= Skew_Mod_Orthogonal .and. InputFileData%Skew_Mod /= Skew_Mod_None .and. InputFileData%Skew_Mod /= Skew_Mod_Active) &
call SetErrStat( ErrID_Fatal, 'Skew_Mod must be -1, 0, or 1.', ErrStat, ErrMsg, RoutineName )

if ( InputFileData%SectAvg) then
Expand Down Expand Up @@ -4295,6 +4295,7 @@ SUBROUTINE Init_BEMTmodule( InputFileData, RotInputFileData, u_AD, u, p, p_AD, x
InitInp%airDens = InputFileData%AirDens
InitInp%kinVisc = InputFileData%KinVisc
InitInp%skewWakeMod = InputFileData%Skew_Mod
InitInp%skewRedistrMod = InputFileData%SkewRedistrMod
InitInp%yawCorrFactor = InputFileData%SkewModFactor
InitInp%aTol = InputFileData%IndToler
InitInp%useTipLoss = InputFileData%TipLoss
Expand Down
12 changes: 6 additions & 6 deletions modules/aerodyn/src/AeroDyn_IO.f90
Original file line number Diff line number Diff line change
Expand Up @@ -771,8 +771,8 @@ SUBROUTINE ParsePrimaryFileInfo( PriPath, InitInp, InputFile, RootName, NumBlade
! Skew_Mod- Select skew model {0: No skew model at all, -1:Throw away non-normal component for linearization, 1: Glauert skew model, }
call ParseVar( FileInfo_In, CurLine, "Skew_Mod", InputFileData%Skew_Mod, ErrStat2, ErrMsg2, UnEc )
if (newInputAbsent('Skew_Mod', CurLine, errStat2, errMsg2)) then
call WrScr(' Setting Skew_Mod to 1 (skew active, Gluaert) as the input is absent (typical behavior).')
InputFileData%Skew_Mod = Skew_Mod_Glauert
call WrScr(' Setting Skew_Mod to 1 (skew active) as the input is absent (typical behavior).')
InputFileData%Skew_Mod = Skew_Mod_Active
else
if (skewModProvided) then
call LegacyAbort('Cannot have both Skew_Mod and SkewMod in the input file'); return
Expand Down Expand Up @@ -1065,7 +1065,7 @@ SUBROUTINE ParsePrimaryFileInfo( PriPath, InitInp, InputFile, RootName, NumBlade
else if (InputFileData%SkewMod==1) then
InputFileData%Skew_Mod = Skew_Mod_None
else if (InputFileData%SkewMod==2) then
InputFileData%Skew_Mod = Skew_Mod_Glauert
InputFileData%Skew_Mod = Skew_Mod_Active
else
call LegacyAbort('Legacy option SkewMod is not 0, 1,2 which is not supported.')
endif
Expand Down Expand Up @@ -1589,12 +1589,12 @@ SUBROUTINE AD_PrintSum( InputFileData, p, p_AD, u, y, ErrStat, ErrMsg )
Msg = 'orthogonal'
case (Skew_Mod_None)
Msg = 'no correction'
case (Skew_Mod_Glauert)
Msg = 'Glauert/Pitt/Peters'
case (Skew_Mod_Active)
Msg = 'active'
case default
Msg = 'unknown'
end select
WRITE (UnSu,Ec_IntFrmt) InputFileData%Skew_Mod, 'Skew_Mod', 'Type of skewed-wake correction model: '//TRIM(Msg)
WRITE (UnSu,Ec_IntFrmt) InputFileData%Skew_Mod, 'Skew_Mod', 'Skewed-wake correction model: '//TRIM(Msg)


! TipLoss
Expand Down
2 changes: 1 addition & 1 deletion modules/aerodyn/src/AeroDyn_Registry.txt
Original file line number Diff line number Diff line change
Expand Up @@ -182,7 +182,7 @@ typedef ^ AD_InputFile ReKi SpdSound - - - "Speed of sound" m/s
typedef ^ AD_InputFile IntKi SkewMod - - - "LEGACY - Skew Mod" -
typedef ^ AD_InputFile IntKi Skew_Mod - - - "Select skew model {0=No skew model at all, -1=Throw away non-normal component for linearization, 1=Glauert skew model}" -
typedef ^ AD_InputFile Logical SkewMomCorr - - - "Turn the skew momentum correction on or off [used only when SkewMod=1]" -
typedef ^ AD_InputFile IntKi SkewRedistrMod - - - "Type of skewed-wake correction model (switch) {0=no redistribution, 1=Glauert/Pitt/Peters, 2=Vortex Cylinder} [unsed only when SkewMod=1]" -
typedef ^ AD_InputFile IntKi SkewRedistrMod - - - "Type of skewed-wake redistribution model (switch) {0=no redistribution, 1=Glauert/Pitt/Peters, 2=Vortex Cylinder} [unsed only when SkewMod=1]" -
typedef ^ AD_InputFile ReKi SkewModFactor - - - "Constant used in Pitt/Peters skewed wake model (default is 15*pi/32)" -
typedef ^ AD_InputFile LOGICAL TipLoss - - - "Use the Prandtl tip-loss model? [unused when WakeMod=0]" flag
typedef ^ AD_InputFile LOGICAL HubLoss - - - "Use the Prandtl hub-loss model? [unused when WakeMod=0]" flag
Expand Down
2 changes: 1 addition & 1 deletion modules/aerodyn/src/AeroDyn_Types.f90
Original file line number Diff line number Diff line change
Expand Up @@ -211,7 +211,7 @@ MODULE AeroDyn_Types
INTEGER(IntKi) :: SkewMod = 0_IntKi !< LEGACY - Skew Mod [-]
INTEGER(IntKi) :: Skew_Mod = 0_IntKi !< Select skew model {0=No skew model at all, -1=Throw away non-normal component for linearization, 1=Glauert skew model} [-]
LOGICAL :: SkewMomCorr = .false. !< Turn the skew momentum correction on or off [used only when SkewMod=1] [-]
INTEGER(IntKi) :: SkewRedistrMod = 0_IntKi !< Type of skewed-wake correction model (switch) {0=no redistribution, 1=Glauert/Pitt/Peters, 2=Vortex Cylinder} [unsed only when SkewMod=1] [-]
INTEGER(IntKi) :: SkewRedistrMod = 0_IntKi !< Type of skewed-wake redistribution model (switch) {0=no redistribution, 1=Glauert/Pitt/Peters, 2=Vortex Cylinder} [unsed only when SkewMod=1] [-]
REAL(ReKi) :: SkewModFactor = 0.0_ReKi !< Constant used in Pitt/Peters skewed wake model (default is 15*pi/32) [-]
LOGICAL :: TipLoss = .false. !< Use the Prandtl tip-loss model? [unused when WakeMod=0] [flag]
LOGICAL :: HubLoss = .false. !< Use the Prandtl hub-loss model? [unused when WakeMod=0] [flag]
Expand Down
12 changes: 9 additions & 3 deletions modules/aerodyn/src/BEMT.f90
Original file line number Diff line number Diff line change
Expand Up @@ -187,7 +187,6 @@ subroutine BEMT_SetParameters( InitInp, p, errStat, errMsg )
p%numBlades = InitInp%numBlades
p%UA_Flag = InitInp%UA_Flag
p%DBEMT_Mod = InitInp%DBEMT_Mod
p%MomentumCorr = InitInp%MomentumCorr
p%BEM_Mod = InitInp%BEM_Mod
!call WrScr('>>>> BEM_Mod '//trim(num2lstr(p%BEM_Mod)))
if (.not.(ANY( p%BEM_Mod == (/BEMMod_2D, BEMMod_3D/)))) then
Expand Down Expand Up @@ -250,6 +249,13 @@ subroutine BEMT_SetParameters( InitInp, p, errStat, errMsg )
p%airDens = InitInp%airDens
p%kinVisc = InitInp%kinVisc
p%skewWakeMod = InitInp%skewWakeMod
if (p%skewWakeMod==Skew_Mod_Active) then
p%skewRedistrMod = InitInp%skewRedistrMod
p%MomentumCorr = InitInp%MomentumCorr
else
p%skewRedistrMod = SkewRedistrMod_None
p%MomentumCorr = .false.
endif
p%yawCorrFactor = InitInp%yawCorrFactor
p%useTipLoss = InitInp%useTipLoss
p%useHubLoss = InitInp%useHubLoss
Expand Down Expand Up @@ -1590,7 +1596,7 @@ subroutine ApplySkewedWakeCorrection_AllNodes(p, u, m, x, phi, OtherState, axInd
!............................................
! Apply skewed wake correction to the axial induction (y%axInduction)
!............................................
if ( p%skewWakeMod == Skew_Mod_Glauert ) then
if ( p%skewWakeMod == Skew_Mod_Active ) then
if (p%BEM_Mod==BEMMod_2D) then
! do nothing
else
Expand All @@ -1602,7 +1608,7 @@ subroutine ApplySkewedWakeCorrection_AllNodes(p, u, m, x, phi, OtherState, axInd
do i = 1,p%numBladeNodes ! Loop through the blade nodes / elements
if ( .not. p%FixedInductions(i,j) ) then
F = getHubTipLossCorrection(p%BEM_Mod, p%useHubLoss, p%useTipLoss, p%hubLossConst(i,j), p%tipLossConst(i,j), phi(i,j), u%cantAngle(i,j) )
call ApplySkewedWakeCorrection( p%BEM_Mod, p%skewWakeMod, p%yawCorrFactor, F, u%psi_s(j), u%psiSkewOffset, u%chi0, u%rlocal(i,j)/m%Rtip(j), axInduction(i,j), chi(i,j), m%FirstWarn_Skew )
call ApplySkewedWakeCorrection( p%BEM_Mod, p%skewRedistrMod, p%yawCorrFactor, F, u%psi_s(j), u%psiSkewOffset, u%chi0, u%rlocal(i,j)/m%Rtip(j), axInduction(i,j), chi(i,j), m%FirstWarn_Skew )
end if ! .not. p%FixedInductions (special case for tip and/or hub loss)
enddo ! I - Blade nodes / elements
enddo ! J - All blades
Expand Down
8 changes: 6 additions & 2 deletions modules/aerodyn/src/BEMTUncoupled.f90
Original file line number Diff line number Diff line change
Expand Up @@ -449,10 +449,10 @@ real(ReKi) function BEMTU_InductionWithResidual(p, u, i, j, phi, AFInfo, IsValid

end function BEMTU_InductionWithResidual
!-----------------------------------------------------------------------------------------
subroutine ApplySkewedWakeCorrection(BEM_Mod, SkewMod, yawCorrFactor, F, azimuth, azimuthOffset, chi0, tipRatio, a, chi, FirstWarn )
subroutine ApplySkewedWakeCorrection(BEM_Mod, SkewRedistrMod, yawCorrFactor, F, azimuth, azimuthOffset, chi0, tipRatio, a, chi, FirstWarn )

integer(IntKi), intent(in ) :: BEM_Mod
integer(IntKi), intent(in ) :: SkewMod
integer(IntKi), intent(in ) :: SkewRedistrMod
real(ReKi), intent(in ) :: yawCorrFactor ! set to 15*pi/32 previously; now allowed to be input (to better match data)
real(ReKi), intent(in ) :: F ! tip/hub loss factor
real(ReKi), intent(in ) :: azimuth
Expand All @@ -467,6 +467,10 @@ subroutine ApplySkewedWakeCorrection(BEM_Mod, SkewMod, yawCorrFactor, F, azimuth
real(ReKi) :: yawCorr
real(ReKi) :: yawCorr_tan ! magnitude of the tan(chi/2) correction term (with possible limits)

if (SkewRedistrMod==SkewRedistrMod_None) then
print*,'>>> SkewRedistrMod is None, Shouldnt happen'
STOP
endif

! Skewed wake correction
if(BEM_Mod==BEMMod_2D) then
Expand Down
8 changes: 5 additions & 3 deletions modules/aerodyn/src/BEMT_Registry.txt
Original file line number Diff line number Diff line change
Expand Up @@ -17,8 +17,8 @@ usefrom UnsteadyAero_Registry.txt
usefrom DBEMT_Registry.txt

param BEMT/BEMT - INTEGER Skew_Mod_Orthogonal - -1 - "Inflow orthogonal to rotor [-]" -
param BEMT/BEMT - INTEGER Skew_Mod_None - 0 - "No skew model (previously called uncoupled)" -
param BEMT/BEMT - INTEGER Skew_Mod_Glauert - 1 - "Pitt/Peters/Glauert skew model (should be 1)" -
param BEMT/BEMT - INTEGER Skew_Mod_None - 0 - "No skew model" -
param BEMT/BEMT - INTEGER Skew_Mod_Active - 1 - "Skew model active" -
param BEMT/BEMT - INTEGER Skew_Mod_PittPeters_Cont - 4 - "Pitt/Peters continuous formulation" -

param BEMT/BEMT - INTEGER SkewRedistrMod_None - 0 - "No redistribution" -
Expand All @@ -39,7 +39,8 @@ typedef BEMT/BEMT InitInputType ReKi
typedef ^ ^ INTEGER numBlades - - - "Number of blades" -
typedef ^ ^ ReKi airDens - - - "Air density" kg/m^3
typedef ^ ^ ReKi kinVisc - - - "Kinematic air viscosity" m^2/s
typedef ^ ^ INTEGER skewWakeMod - - - "Type of skewed-wake correction model [switch] {1=uncoupled, 2=Pitt/Peters, 3=coupled}" -
typedef ^ ^ INTEGER skewWakeMod - - - "Type of skewed-wake model [switch] {0=None, 1=Glauert}" -
typedef ^ ^ INTEGER skewRedistrMod - - - "Type of skewed-wake redistribution model (switch) {0=no redistribution, 1=Glauert/Pitt/Peters, 2=Vortex Cylinder} [unsed only when SkewMod=1]" -
typedef ^ ^ ReKi aTol - - - "Tolerance for the induction solution" -
typedef ^ ^ LOGICAL useTipLoss - - - "Use the Prandtl tip-loss model? [flag]" -
typedef ^ ^ LOGICAL useHubLoss - - - "Use the Prandtl hub-loss model? [flag]" -
Expand Down Expand Up @@ -144,6 +145,7 @@ typedef ^ ^ INTEGER
typedef ^ ^ ReKi airDens - - - "Air density" kg/m^3
typedef ^ ^ ReKi kinVisc - - - "Kinematic air viscosity" m^2/s
typedef ^ ^ INTEGER skewWakeMod - - - "Type of skewed-wake correction model [switch] {0=None, 1=Glauert/Pitt/Peters}" -
typedef ^ ^ INTEGER skewRedistrMod - - - "Type of skewed-wake redistribution model (switch) {0=no redistribution, 1=Glauert/Pitt/Peters, 2=Vortex Cylinder} [unsed only when SkewMod=1]" -
typedef ^ ^ ReKi aTol - - - "Tolerance for the induction solution" -
typedef ^ ^ LOGICAL useTipLoss - - - "Use the Prandtl tip-loss model? [flag]" -
typedef ^ ^ LOGICAL useHubLoss - - - "Use the Prandtl hub-loss model? [flag]" -
Expand Down
16 changes: 13 additions & 3 deletions modules/aerodyn/src/BEMT_Types.f90
Original file line number Diff line number Diff line change
Expand Up @@ -37,8 +37,8 @@ MODULE BEMT_Types
USE NWTC_Library
IMPLICIT NONE
INTEGER(IntKi), PUBLIC, PARAMETER :: Skew_Mod_Orthogonal = -1 ! Inflow orthogonal to rotor [-] [-]
INTEGER(IntKi), PUBLIC, PARAMETER :: Skew_Mod_None = 0 ! No skew model (previously called uncoupled) [-]
INTEGER(IntKi), PUBLIC, PARAMETER :: Skew_Mod_Glauert = 1 ! Pitt/Peters/Glauert skew model (should be 1) [-]
INTEGER(IntKi), PUBLIC, PARAMETER :: Skew_Mod_None = 0 ! No skew model [-]
INTEGER(IntKi), PUBLIC, PARAMETER :: Skew_Mod_Active = 1 ! Skew model active [-]
INTEGER(IntKi), PUBLIC, PARAMETER :: Skew_Mod_PittPeters_Cont = 4 ! Pitt/Peters continuous formulation [-]
INTEGER(IntKi), PUBLIC, PARAMETER :: SkewRedistrMod_None = 0 ! No redistribution [-]
INTEGER(IntKi), PUBLIC, PARAMETER :: SkewRedistrMod_PittPeters = 1 ! Pitt/Petesr/Glauert redistribution [-]
Expand All @@ -52,7 +52,8 @@ MODULE BEMT_Types
INTEGER(IntKi) :: numBlades = 0_IntKi !< Number of blades [-]
REAL(ReKi) :: airDens = 0.0_ReKi !< Air density [kg/m^3]
REAL(ReKi) :: kinVisc = 0.0_ReKi !< Kinematic air viscosity [m^2/s]
INTEGER(IntKi) :: skewWakeMod = 0_IntKi !< Type of skewed-wake correction model [switch] {1=uncoupled, 2=Pitt/Peters, 3=coupled} [-]
INTEGER(IntKi) :: skewWakeMod = 0_IntKi !< Type of skewed-wake model [switch] {0=None, 1=Glauert} [-]
INTEGER(IntKi) :: skewRedistrMod = 0_IntKi !< Type of skewed-wake redistribution model (switch) {0=no redistribution, 1=Glauert/Pitt/Peters, 2=Vortex Cylinder} [unsed only when SkewMod=1] [-]
REAL(ReKi) :: aTol = 0.0_ReKi !< Tolerance for the induction solution [-]
LOGICAL :: useTipLoss = .false. !< Use the Prandtl tip-loss model? [flag] [-]
LOGICAL :: useHubLoss = .false. !< Use the Prandtl hub-loss model? [flag] [-]
Expand Down Expand Up @@ -154,6 +155,7 @@ MODULE BEMT_Types
REAL(ReKi) :: airDens = 0.0_ReKi !< Air density [kg/m^3]
REAL(ReKi) :: kinVisc = 0.0_ReKi !< Kinematic air viscosity [m^2/s]
INTEGER(IntKi) :: skewWakeMod = 0_IntKi !< Type of skewed-wake correction model [switch] {0=None, 1=Glauert/Pitt/Peters} [-]
INTEGER(IntKi) :: skewRedistrMod = 0_IntKi !< Type of skewed-wake redistribution model (switch) {0=no redistribution, 1=Glauert/Pitt/Peters, 2=Vortex Cylinder} [unsed only when SkewMod=1] [-]
REAL(ReKi) :: aTol = 0.0_ReKi !< Tolerance for the induction solution [-]
LOGICAL :: useTipLoss = .false. !< Use the Prandtl tip-loss model? [flag] [-]
LOGICAL :: useHubLoss = .false. !< Use the Prandtl hub-loss model? [flag] [-]
Expand Down Expand Up @@ -259,6 +261,7 @@ subroutine BEMT_CopyInitInput(SrcInitInputData, DstInitInputData, CtrlCode, ErrS
DstInitInputData%airDens = SrcInitInputData%airDens
DstInitInputData%kinVisc = SrcInitInputData%kinVisc
DstInitInputData%skewWakeMod = SrcInitInputData%skewWakeMod
DstInitInputData%skewRedistrMod = SrcInitInputData%skewRedistrMod
DstInitInputData%aTol = SrcInitInputData%aTol
DstInitInputData%useTipLoss = SrcInitInputData%useTipLoss
DstInitInputData%useHubLoss = SrcInitInputData%useHubLoss
Expand Down Expand Up @@ -428,6 +431,7 @@ subroutine BEMT_PackInitInput(Buf, Indata)
call RegPack(Buf, InData%airDens)
call RegPack(Buf, InData%kinVisc)
call RegPack(Buf, InData%skewWakeMod)
call RegPack(Buf, InData%skewRedistrMod)
call RegPack(Buf, InData%aTol)
call RegPack(Buf, InData%useTipLoss)
call RegPack(Buf, InData%useHubLoss)
Expand Down Expand Up @@ -522,6 +526,8 @@ subroutine BEMT_UnPackInitInput(Buf, OutData)
if (RegCheckErr(Buf, RoutineName)) return
call RegUnpack(Buf, OutData%skewWakeMod)
if (RegCheckErr(Buf, RoutineName)) return
call RegUnpack(Buf, OutData%skewRedistrMod)
if (RegCheckErr(Buf, RoutineName)) return
call RegUnpack(Buf, OutData%aTol)
if (RegCheckErr(Buf, RoutineName)) return
call RegUnpack(Buf, OutData%useTipLoss)
Expand Down Expand Up @@ -1578,6 +1584,7 @@ subroutine BEMT_CopyParam(SrcParamData, DstParamData, CtrlCode, ErrStat, ErrMsg)
DstParamData%airDens = SrcParamData%airDens
DstParamData%kinVisc = SrcParamData%kinVisc
DstParamData%skewWakeMod = SrcParamData%skewWakeMod
DstParamData%skewRedistrMod = SrcParamData%skewRedistrMod
DstParamData%aTol = SrcParamData%aTol
DstParamData%useTipLoss = SrcParamData%useTipLoss
DstParamData%useHubLoss = SrcParamData%useHubLoss
Expand Down Expand Up @@ -1726,6 +1733,7 @@ subroutine BEMT_PackParam(Buf, Indata)
call RegPack(Buf, InData%airDens)
call RegPack(Buf, InData%kinVisc)
call RegPack(Buf, InData%skewWakeMod)
call RegPack(Buf, InData%skewRedistrMod)
call RegPack(Buf, InData%aTol)
call RegPack(Buf, InData%useTipLoss)
call RegPack(Buf, InData%useHubLoss)
Expand Down Expand Up @@ -1810,6 +1818,8 @@ subroutine BEMT_UnPackParam(Buf, OutData)
if (RegCheckErr(Buf, RoutineName)) return
call RegUnpack(Buf, OutData%skewWakeMod)
if (RegCheckErr(Buf, RoutineName)) return
call RegUnpack(Buf, OutData%skewRedistrMod)
if (RegCheckErr(Buf, RoutineName)) return
call RegUnpack(Buf, OutData%aTol)
if (RegCheckErr(Buf, RoutineName)) return
call RegUnpack(Buf, OutData%useTipLoss)
Expand Down

0 comments on commit f20aaeb

Please sign in to comment.