Skip to content

Commit

Permalink
AD: preliminary move of BN outputs from FVW to AD
Browse files Browse the repository at this point in the history
  • Loading branch information
@ebranlard
ebranlard committed Jul 18, 2023
1 parent dadfa6a commit 7b8bf96
Show file tree
Hide file tree
Showing 5 changed files with 1,060 additions and 15 deletions.
43 changes: 37 additions & 6 deletions modules/aerodyn/src/AeroDyn.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -455,7 +455,7 @@ subroutine AD_Init( InitInp, u, p, x, xd, z, OtherState, y, m, Interval, InitOut
! many states are in the BEMT module, which were initialized in BEMT_Init()

do iR = 1, nRotors
call Init_MiscVars(m%rotors(iR), p%rotors(iR), u%rotors(iR), y%rotors(iR), errStat2, errMsg2)
call Init_MiscVars(m%rotors(iR), p%rotors(iR), u%rotors(iR), y%rotors(iR), p, errStat2, errMsg2)
if (Failed()) return;
enddo

Expand Down Expand Up @@ -583,11 +583,12 @@ subroutine AD_ReInit(p, x, xd, z, OtherState, m, Interval, ErrStat, ErrMsg )
end subroutine AD_ReInit
!----------------------------------------------------------------------------------------------------------------------------------
!> This routine initializes (allocates) the misc variables for use during the simulation.
subroutine Init_MiscVars(m, p, u, y, errStat, errMsg)
subroutine Init_MiscVars(m, p, u, y, p_AD, errStat, errMsg)
type(RotMiscVarType), intent(inout) :: m !< misc/optimization data (not defined in submodules)
type(RotParameterType), intent(in ) :: p !< Parameters
type(RotInputType), intent(inout) :: u !< input for HubMotion mesh (create sibling mesh here)
type(RotOutputType), intent(inout) :: y !< output (create mapping between output and otherstate mesh here)
type(AD_ParameterType), intent(in ) :: p_AD !< Parameters
integer(IntKi), intent( out) :: errStat !< Error status of the operation
character(*), intent( out) :: errMsg !< Error message if ErrStat /= ErrID_None

Expand Down Expand Up @@ -656,6 +657,28 @@ subroutine Init_MiscVars(m, p, u, y, errStat, errMsg)
call SetErrStat( errStat2, errMsg2, errStat, errMsg, RoutineName )
call AllocAry( m%Vind_i, 3, p%NumBlNds, p%NumBlades, 'm%Vind_i', ErrStat2, ErrMsg2 )
call SetErrStat( errStat2, errMsg2, errStat, errMsg, RoutineName )
if(Failed()) return

! Lifting line calculations for output (OLAF only for now)
if (p_AD%WakeMod == WakeMod_FVW) then
call AllocAry( m%BN_AxInd , p%NumBlNds , p%NumBlades, 'm%BN_AxInd ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_AxInd = -999999_ReKi;
call AllocAry( m%BN_TnInd , p%NumBlNds , p%NumBlades, 'm%BN_TnInd' , ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_TnInd = -999999_ReKi;
call AllocAry( m%BN_Vrel , p%NumBlNds , p%NumBlades, 'm%BN_Vrel ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_Vrel = -999999_ReKi;
call AllocAry( m%BN_alpha , p%NumBlNds , p%NumBlades, 'm%BN_alpha ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_alpha = -999999_ReKi;
call AllocAry( m%BN_phi , p%NumBlNds , p%NumBlades, 'm%BN_phi ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_phi = -999999_ReKi;
call AllocAry( m%BN_Re , p%NumBlNds , p%NumBlades, 'm%BN_Re ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_Re = -999999_ReKi;
call AllocAry( m%BN_Cl_qs , p%NumBlNds , p%NumBlades, 'm%BN_Cl_qs ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_Cl_qs = -999999_ReKi;
call AllocAry( m%BN_Cd_qs , p%NumBlNds , p%NumBlades, 'm%BN_Cd_qs ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_Cd_qs = -999999_ReKi;
call AllocAry( m%BN_Cm_qs , p%NumBlNds , p%NumBlades, 'm%BN_Cm_qs ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_Cm_qs = -999999_ReKi;
call AllocAry( m%BN_Cpmin , p%NumBlNds , p%NumBlades, 'm%BN_Cpmin ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_Cpmin = -999999_ReKi;
call AllocAry( m%BN_Cl , p%NumBlNds , p%NumBlades, 'm%BN_Cl ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_Cl = -999999_ReKi;
call AllocAry( m%BN_Cd , p%NumBlNds , p%NumBlades, 'm%BN_Cd ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_Cd = -999999_ReKi;
call AllocAry( m%BN_Cm , p%NumBlNds , p%NumBlades, 'm%BN_Cm ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_Cm = -999999_ReKi;
call AllocAry( m%BN_Cx , p%NumBlNds , p%NumBlades, 'm%BN_Cx ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_Cx = -999999_ReKi;
call AllocAry( m%BN_Cy , p%NumBlNds , p%NumBlades, 'm%BN_Cy ', ErrStat2, ErrMsg2 );call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName); m%BN_Cy = -999999_ReKi;
if(Failed()) return
endif

! mesh mapping data for integrating load over entire rotor:
allocate( m%B_L_2_H_P(p%NumBlades), Stat = ErrStat2)
if (ErrStat2 /= 0) then
Expand Down Expand Up @@ -871,6 +894,11 @@ subroutine Init_MiscVars(m, p, u, y, errStat, errMsg)

m%FirstWarn_TowerStrike = .true.

contains
logical function Failed()
call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName)
Failed = ErrStat >= AbortErrLev
end function Failed
end subroutine Init_MiscVars
!----------------------------------------------------------------------------------------------------------------------------------
!> This routine initializes (allocates) the misc variables for use during the simulation.
Expand Down Expand Up @@ -1995,8 +2023,8 @@ subroutine CalcBuoyantLoads( u, p, m, y, ErrStat, ErrMsg )
REAL(ReKi), DIMENSION(3) :: BlmomentTip !< Buoyant moment on element tip in global coordinates
REAL(ReKi), DIMENSION(3) :: TwrforceBtop !< Buoyant force on tower top in global coordinates
REAL(ReKi), DIMENSION(3) :: TwrmomentBtop !< Buoyant moment on tower top in global coordinates
REAL(ReKi), DIMENSION(p%NumBlNds,p%NumBlades,3) :: BlFBtmp !< Buoyant force at blade nodes in global coordinates
REAL(ReKi), DIMENSION(p%NumBlNds,p%NumBlades,3) :: BlMBtmp !< Buoyant moment at blade nodes in global coordinates
REAL(ReKi), DIMENSION(p%NumBlNds,p%NumBlades,3) :: BlFBtmp !< Buoyant force at blade nodes in global coordinates ! TODO REVERSE MEMORY ORDER
REAL(ReKi), DIMENSION(p%NumBlNds,p%NumBlades,3) :: BlMBtmp !< Buoyant moment at blade nodes in global coordinates ! TODO REVERSE MEMORY ORDER
REAL(ReKi), DIMENSION(p%NumTwrNds,3) :: TwrFBtmp !< Buoyant force at tower nodes in global coordinates
REAL(ReKi), DIMENSION(p%NumTwrNds,3) :: TwrMBtmp !< Buoyant moment at tower nodes in global coordinates
REAL(ReKi), DIMENSION(3) :: HubFBtmp !< Buoyant force at hub node in global coordinates, passed to m%HubFB
Expand Down Expand Up @@ -3408,8 +3436,11 @@ subroutine SetOutputsFromFVW(t, u, p, OtherState, x, xd, m, y, ErrStat, ErrMsg)
y%rotors(iR)%BladeLoad(k)%Moment(:,j) = matmul( moment, m%rotors(iR)%orientationAnnulus(:,:,j,k) ) ! moment per unit length of the jth node in the kth blade

! Save results for outputs so we don't have to recalculate them all when we write outputs
m%FVW%W(iW)%BN_AxInd(j) = AxInd
m%FVW%W(iW)%BN_TanInd(j) = TanInd
m%rotors(iR)%BN_AxInd(j,k) = AxInd
m%rotors(iR)%BN_TnInd(j,k) = TanInd
m%rotors(iR)%BN_Vrel(j,k) = Vrel
! m%FVW%W(iW)%BN_AxInd(j) = AxInd
! m%FVW%W(iW)%BN_TanInd(j) = TanInd
m%FVW%W(iW)%BN_Vrel(j) = Vrel
m%FVW%W(iW)%BN_alpha(j) = alpha
m%FVW%W(iW)%BN_phi(j) = phi
Expand Down
12 changes: 8 additions & 4 deletions modules/aerodyn/src/AeroDyn_AllBldNdOuts_IO.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -417,7 +417,8 @@ SUBROUTINE Calc_WriteAllBldNdOutput( p, p_AD, u, m, m_AD, x, y, OtherState, Indx
do iB=1,nB
iW = W2B(iB)
do iNd=1,nNd;
y%WriteOutput(iOut) = -m_AD%FVW%W(iW)%BN_UrelWind_s(1,iNd) * m_AD%FVW%W(iW)%BN_AxInd(iNd)
y%WriteOutput(iOut) = -m_AD%FVW%W(iW)%BN_UrelWind_s(1,iNd) * m%BN_AxInd(iNd,iB)
!y%WriteOutput(iOut) = -m_AD%FVW%W(iW)%BN_UrelWind_s(1,iNd) * m_AD%FVW%W(iW)%BN_AxInd(iNd)
iOut = iOut + 1
enddo
enddo
Expand All @@ -435,7 +436,8 @@ SUBROUTINE Calc_WriteAllBldNdOutput( p, p_AD, u, m, m_AD, x, y, OtherState, Indx
DO iB=1,nB
iW = W2B(iB)
DO iNd=1,nNd
y%WriteOutput(iOut) = m_AD%FVW%W(iW)%BN_UrelWind_s(2,iNd) * m_AD%FVW%W(iW)%BN_TanInd(iNd)
y%WriteOutput(iOut) = m_AD%FVW%W(iW)%BN_UrelWind_s(2,iNd) * m%BN_TnInd(iNd,iB)
!y%WriteOutput(iOut) = m_AD%FVW%W(iW)%BN_UrelWind_s(2,iNd) * m_AD%FVW%W(iW)%BN_TanInd(iNd)
iOut = iOut + 1
END DO
END DO
Expand Down Expand Up @@ -549,7 +551,8 @@ SUBROUTINE Calc_WriteAllBldNdOutput( p, p_AD, u, m, m_AD, x, y, OtherState, Indx
DO iB=1,nB
iW = W2B(iB)
DO iNd=1,nNd
y%WriteOutput(iOut) = m_AD%FVW%W(iW)%BN_AxInd(iNd)
y%WriteOutput(iOut) = m%BN_AxInd(iNd,iB)
!y%WriteOutput(iOut) = m_AD%FVW%W(iW)%BN_AxInd(iNd)
iOut = iOut + 1
END DO
END DO
Expand All @@ -567,7 +570,8 @@ SUBROUTINE Calc_WriteAllBldNdOutput( p, p_AD, u, m, m_AD, x, y, OtherState, Indx
DO iB=1,nB
iW = W2B(iB)
DO iNd=1,nNd
y%WriteOutput(iOut) = m_AD%FVW%W(iW)%BN_TanInd(iNd)
y%WriteOutput(iOut) = m%BN_TnInd(iNd,iB)
!y%WriteOutput(iOut) = m_AD%FVW%W(iW)%BN_TanInd(iNd)
iOut = iOut + 1
END DO
END DO
Expand Down
14 changes: 9 additions & 5 deletions modules/aerodyn/src/AeroDyn_IO.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -496,11 +496,15 @@ subroutine Calc_WriteOutput_FVW
m%AllOuts( BNRe( beta,k) ) = m_AD%FVW%W(iW)%BN_Re(j) / 1.0E6
m%AllOuts( BNM( beta,k) ) = m_AD%FVW%W(iW)%BN_Vrel(j) / p%SpdSound

m%AllOuts( BNVIndx(beta,k) ) = -m_AD%FVW%W(iW)%BN_UrelWind_s(1,j) * m_AD%FVW%W(iW)%BN_AxInd(j)
m%AllOuts( BNVIndy(beta,k) ) = m_AD%FVW%W(iW)%BN_UrelWind_s(2,j) * m_AD%FVW%W(iW)%BN_TanInd(j)

m%AllOuts( BNAxInd(beta,k) ) = m_AD%FVW%W(iW)%BN_AxInd(j)
m%AllOuts( BNTnInd(beta,k) ) = m_AD%FVW%W(iW)%BN_TanInd(j)
!m%AllOuts( BNVIndx(beta,k) ) = -m_AD%FVW%W(iW)%BN_UrelWind_s(1,j) * m_AD%FVW%W(iW)%BN_AxInd(j)
!m%AllOuts( BNVIndy(beta,k) ) = m_AD%FVW%W(iW)%BN_UrelWind_s(2,j) * m_AD%FVW%W(iW)%BN_TanInd(j)
m%AllOuts( BNVIndx(beta,k) ) = -m_AD%FVW%W(iW)%BN_UrelWind_s(1,j) * m%BN_AxInd(j,k)
m%AllOuts( BNVIndy(beta,k) ) = m_AD%FVW%W(iW)%BN_UrelWind_s(2,j) * m%BN_TnInd(j,k)

!m%AllOuts( BNAxInd(beta,k) ) = m_AD%FVW%W(iW)%BN_AxInd(j)
!m%AllOuts( BNTnInd(beta,k) ) = m_AD%FVW%W(iW)%BN_TanInd(j)
m%AllOuts( BNAxInd(beta,k) ) = m%BN_AxInd(j,k)
m%AllOuts( BNTnInd(beta,k) ) = m%BN_TnInd(j,k)

m%AllOuts( BNAlpha(beta,k) ) = m_AD%FVW%W(iW)%BN_alpha(j)*R2D
m%AllOuts( BNTheta(beta,k) ) = m_AD%FVW%W(iW)%PitchAndTwist(j)*R2D
Expand Down
16 changes: 16 additions & 0 deletions modules/aerodyn/src/AeroDyn_Registry.txt
View file
Original file line number Diff line number Diff line change
Expand Up @@ -307,6 +307,22 @@ typedef ^ RotMiscVarType MeshMapType T_P_2_T_L - - - "mapping data structure to
typedef ^ RotMiscVarType Logical FirstWarn_TowerStrike - - - "flag to avoid printing tower strike multiple times" -
typedef ^ RotMiscVarType ReKi AvgDiskVel {3} - - "disk-averaged U,V,W (undisturbed)" m/s
typedef ^ RotMiscVarType ReKi AvgDiskVelDist {3} - - "disk-averaged U,V,W (disturbed)" m/s
# Lifting line outputs (for FVW only for now)
typedef ^ RotMiscVarType ReKi BN_AxInd :: - - "Axial induction [size (NumBlNds,numBlades)]" -
typedef ^ RotMiscVarType ReKi BN_TnInd :: - - "Tangential induction [size (NumBlNds,numBlades)]" -
typedef ^ RotMiscVarType ReKi BN_Vrel :: - - "Relative velocity [size (NumBlNds,numBlades)]" m/s
typedef ^ RotMiscVarType ReKi BN_alpha :: - - "Angle of attack [size (NumBlNds,numBlades)]" rad
typedef ^ RotMiscVarType ReKi BN_phi :: - - "angle between the plane of rotation and the direction of the local wind [size (NumBlNds,numBlades)]" rad
typedef ^ RotMiscVarType ReKi BN_Re :: - - "Reynolds number [size (NumBlNds,numBlades)]" -
typedef ^ RotMiscVarType ReKi BN_Cl_qs :: - - "Coefficient lift, excluding unsteady aero effects" -
typedef ^ RotMiscVarType ReKi BN_Cd_qs :: - - "Coefficient drag. excluding unsteady aero effects" -
typedef ^ RotMiscVarType ReKi BN_Cm_qs :: - - "Coefficient moment, excluding unsteady aero effects" -
typedef ^ RotMiscVarType ReKi BN_Cpmin :: - - "Coefficient minimum pressure, excluding unsteady aero effects" -
typedef ^ RotMiscVarType ReKi BN_Cl :: - - "Coefficient lift, including unsteady aero effects" -
typedef ^ RotMiscVarType ReKi BN_Cd :: - - "Coefficient drag, including unsteady aero effects" -
typedef ^ RotMiscVarType ReKi BN_Cm :: - - "Coefficient moment, including unsteady aero effects" -
typedef ^ RotMiscVarType ReKi BN_Cx :: - - "normal force coefficient (normal to the plane, not chord) of the jth node in the kth blade" -
typedef ^ RotMiscVarType ReKi BN_Cy :: - - "tangential force coefficient (tangential to the plane, not chord) of the jth node in the kth blade" -
# TailFin
typedef ^ RotMiscVarType ReKi TFinAlpha - - - "Angle of attack for tailfin"
typedef ^ RotMiscVarType ReKi TFinRe - - - "Reynolds number for tailfin"
Expand Down
Loading

0 comments on commit 7b8bf96

Please sign in to comment.