Skip to content

Commit

Permalink
AD15: add ElemInflowType for blades and tower
Browse files Browse the repository at this point in the history 
Simplification based on comment from @bjonkman
  • Loading branch information
@andrew-platt
andrew-platt committed Apr 30, 2024
1 parent 8024b98 commit 479bd62
Show file tree
Hide file tree
Showing 5 changed files with 135 additions and 168 deletions.
42 changes: 21 additions & 21 deletions modules/aerodyn/src/AeroDyn.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -657,7 +657,7 @@ subroutine Init_MiscVars(m, p, u, y, errStat, errMsg)
errStat = ErrID_None
errMsg = ""

call AllocAry( m%DisturbedInflow, 3_IntKi, p%NumBlNds, p%numBlades, 'm%DisturbedInflow', ErrStat2, ErrMsg2 ) ! must be same size as u%InflowOnBlade
call AllocAry( m%DisturbedInflow, 3_IntKi, p%NumBlNds, p%numBlades, 'm%DisturbedInflow', ErrStat2, ErrMsg2 ) ! must be same size as RotInflow%Blade(k)%InflowVel
call SetErrStat( errStat2, errMsg2, errStat, errMsg, RoutineName )
call AllocAry( m%orientationAnnulus, 3_IntKi, 3_IntKi, p%NumBlNds, p%numBlades, 'm%orientationAnnulus', ErrStat2, ErrMsg2 )
call SetErrStat( errStat2, errMsg2, errStat, errMsg, RoutineName )
Expand Down Expand Up @@ -1308,29 +1308,29 @@ subroutine Init_RotInflow( p, RotInflow, errStat, ErrMsg )
ErrMsg = ""

! Arrays for InflowWind inputs:
allocate(RotInflow%Bld(p%numBlades), stat=ErrStat2)
allocate(RotInflow%Blade(p%numBlades), stat=ErrStat2)
if (ErrStat2 /= 0) then
call SetErrStat( ErrID_Fatal, 'Error allocating RotInflow%Bld', errStat, errMsg, RoutineName )
call SetErrStat( ErrID_Fatal, 'Error allocating RotInflow%Blade', errStat, errMsg, RoutineName )
if (Failed()) return
end if

do k = 1, p%NumBlades
call AllocAry( RotInflow%Bld(k)%InflowOnBlade, 3_IntKi, p%NumBlNds, 'RotInflow%Bld(k)%InflowOnBlade', ErrStat2, ErrMsg2 )
call AllocAry( RotInflow%Blade(k)%InflowVel, 3_IntKi, p%NumBlNds, 'RotInflow%Blade(k)%InflowVel', ErrStat2, ErrMsg2 )
if (Failed()) return
RotInflow%Bld(k)%InflowOnBlade = 0.0_ReKi
RotInflow%Blade(k)%InflowVel = 0.0_ReKi

if (p%MHK > 0) then
call AllocAry( RotInflow%Bld(k)%AccelOnBlade, 3_IntKi, p%NumBlNds, 'RotInflow%Bld(k)%AccelOnBlade', ErrStat2, ErrMsg2 )
call AllocAry( RotInflow%Blade(k)%InflowAcc, 3_IntKi, p%NumBlNds, 'RotInflow%Blade(k)%InflowAcc', ErrStat2, ErrMsg2 )
if (Failed()) return
RotInflow%Bld(k)%AccelOnBlade = 0.0_ReKi
RotInflow%Blade(k)%InflowAcc = 0.0_ReKi
end if
end do

call AllocAry( RotInflow%InflowOnTower, 3_IntKi, p%NumTwrNds, 'RotInflow%InflowOnTower', ErrStat2, ErrMsg2 ) ! could be size zero
call AllocAry( RotInflow%Tower%InflowVel, 3_IntKi, p%NumTwrNds, 'RotInflow%Tower%InflowVel', ErrStat2, ErrMsg2 ) ! could be size zero
if (Failed()) return

if (p%MHK > 0) then
call AllocAry( RotInflow%AccelOnTower, 3_IntKi, p%NumTwrNds, 'RotInflow%AccelOnTower', ErrStat2, ErrMsg2 ) ! could be size zero
call AllocAry( RotInflow%Tower%InflowAcc, 3_IntKi, p%NumTwrNds, 'RotInflow%Tower%InflowAcc', ErrStat2, ErrMsg2 ) ! could be size zero
if (Failed()) return
end if

Expand All @@ -1339,7 +1339,7 @@ subroutine Init_RotInflow( p, RotInflow, errStat, ErrMsg )
RotInflow%InflowOnNacelle = 0.0_ReKi
RotInflow%InflowOnTailFin = 0.0_ReKi
RotInflow%AvgDiskVel = 0.0_ReKi
RotInflow%InflowOnTower = 0.0_ReKi
RotInflow%Tower%InflowVel = 0.0_ReKi

contains
logical function Failed()
Expand Down Expand Up @@ -1869,7 +1869,7 @@ subroutine AD_CalcWind_Rotor(t, u, FlowField, p, RotInflow, StartNode, ErrStat,
do k = 1, p%NumBlades
call IfW_FlowField_GetVelAcc(FlowField, StartNode, t, &
real(u%BladeMotion(k)%TranslationDisp + u%BladeMotion(k)%Position, ReKi), &
RotInflow%Bld(k)%InflowOnBlade, RotInflow%Bld(k)%AccelOnBlade, ErrStat2, ErrMsg2, PosOffset=PosOffset)
RotInflow%Blade(k)%InflowVel, RotInflow%Blade(k)%InflowAcc, ErrStat2, ErrMsg2, PosOffset=PosOffset)
if(Failed()) return
StartNode = StartNode + p%NumBlNds
end do
Expand All @@ -1878,7 +1878,7 @@ subroutine AD_CalcWind_Rotor(t, u, FlowField, p, RotInflow, StartNode, ErrStat,
if (u%TowerMotion%Nnodes > 0) then
call IfW_FlowField_GetVelAcc(FlowField, StartNode, t, &
real(u%TowerMotion%TranslationDisp + u%TowerMotion%Position, ReKi), &
RotInflow%InflowOnTower, RotInflow%AccelOnTower, ErrStat2, ErrMsg2, PosOffset=PosOffset)
RotInflow%Tower%InflowVel, RotInflow%Tower%InflowAcc, ErrStat2, ErrMsg2, PosOffset=PosOffset)
if(Failed()) return
StartNode = StartNode + p%NumTwrNds
end if
Expand Down Expand Up @@ -2780,7 +2780,7 @@ subroutine SetDisturbedInflow(p, p_AD, u, RotInflow, m, errStat, errMsg)
call SetErrStat(errStat2, errMsg2, errStat, errMsg, RoutineName)
else
do k = 1, p%NumBlades
m%DisturbedInflow(:,:,k) = RotInflow%Bld(k)%InflowOnBlade
m%DisturbedInflow(:,:,k) = RotInflow%Blade(k)%InflowVel
end do
end if

Expand Down Expand Up @@ -3141,7 +3141,7 @@ subroutine DiskAvgValues(p, u, RotInflow, m, x_hat_disk, y_hat_disk, z_hat_disk,
do k=1,p%NumBlades
do j=1,p%NumBlNds
m%AvgDiskVelDist = m%AvgDiskVelDist + m%DisturbedInflow(:,j,k)
m%AvgDiskVel = m%AvgDiskVel + RotInflow%Bld(k)%InflowOnBlade(:,j)
m%AvgDiskVel = m%AvgDiskVel + RotInflow%Blade(k)%InflowVel(:,j)
end do
end do
m%AvgDiskVelDist = m%AvgDiskVelDist / real( p%NumBlades * p%NumBlNds, ReKi )
Expand Down Expand Up @@ -3499,7 +3499,7 @@ subroutine SetInputsForAA(p, u, RotInflow, m, errStat, errMsg)
m%AA_u%AoANoise(i,j) = m%BEMT_y%AOA(i,j)

! Set the blade element undisturbed flow
m%AA_u%Inflow(:,i,j) = RotInflow%Bld(j)%InflowonBlade(:,i)
m%AA_u%Inflow(:,i,j) = RotInflow%Blade(j)%InflowVel(:,i)
end do
end do
end subroutine SetInputsForAA
Expand Down Expand Up @@ -4738,7 +4738,7 @@ SUBROUTINE ADTwr_CalcOutput(p, u, RotInflow, m, y, ErrStat, ErrMsg )

do j=1,p%NumTwrNds

V_rel = RotInflow%InflowOnTower(:,j) - u%TowerMotion%TranslationVel(:,j) ! relative wind speed at tower node
V_rel = RotInflow%Tower%InflowVel(:,j) - u%TowerMotion%TranslationVel(:,j) ! relative wind speed at tower node

tmp = u%TowerMotion%Orientation(1,:,j)
VL(1) = dot_product( V_Rel, tmp ) ! relative local x-component of wind speed of the jth node in the tower
Expand Down Expand Up @@ -4858,9 +4858,9 @@ SUBROUTINE TwrInfl( p, u, RotInflow, m, ErrStat, ErrMsg )

if ( DisturbInflow ) then
v = CalculateTowerInfluence(p, xbar, ybar, zbar, W_tower, TwrCd, TwrTI)
m%DisturbedInflow(:,j,k) = RotInflow%Bld(k)%InflowOnBlade(:,j) + matmul( theta_tower_trans, v )
m%DisturbedInflow(:,j,k) = RotInflow%Blade(k)%InflowVel(:,j) + matmul( theta_tower_trans, v )
else
m%DisturbedInflow(:,j,k) = RotInflow%Bld(k)%InflowOnBlade(:,j)
m%DisturbedInflow(:,j,k) = RotInflow%Blade(k)%InflowVel(:,j)
end if

end do !j=NumBlNds
Expand Down Expand Up @@ -5169,8 +5169,8 @@ SUBROUTINE TwrInfl_NearestLine2Element(p, u, RotInflow, BladeNodePosition, r_Tow
found = .true.
min_dist = dist

V_rel_tower = ( RotInflow%InflowOnTower(:,n1) - u%TowerMotion%TranslationVel(:,n1) ) * elem_position2 &
+ ( RotInflow%InflowOnTower(:,n2) - u%TowerMotion%TranslationVel(:,n2) ) * elem_position
V_rel_tower = ( RotInflow%Tower%InflowVel(:,n1) - u%TowerMotion%TranslationVel(:,n1) ) * elem_position2 &
+ ( RotInflow%Tower%InflowVel(:,n2) - u%TowerMotion%TranslationVel(:,n2) ) * elem_position

TwrDiam = elem_position2*p%TwrDiam(n1) + elem_position*p%TwrDiam(n2)
TwrCd = elem_position2*p%TwrCd( n1) + elem_position*p%TwrCd( n2)
Expand Down Expand Up @@ -5285,7 +5285,7 @@ SUBROUTINE TwrInfl_NearestPoint(p, u, RotInflow, BladeNodePosition, r_TowerBlade
n1 = node_with_min_distance

r_TowerBlade = BladeNodePosition - u%TowerMotion%Position(:,n1) - u%TowerMotion%TranslationDisp(:,n1)
V_rel_tower = RotInflow%InflowOnTower(:,n1) - u%TowerMotion%TranslationVel(:,n1)
V_rel_tower = RotInflow%Tower%InflowVel(:,n1) - u%TowerMotion%TranslationVel(:,n1)
TwrDiam = p%TwrDiam(n1)
TwrCd = p%TwrCd( n1)
TwrTI = p%TwrTI( n1)
Expand Down
30 changes: 15 additions & 15 deletions modules/aerodyn/src/AeroDyn_AllBldNdOuts_IO.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -336,26 +336,26 @@ SUBROUTINE Calc_WriteAllBldNdOutput( p, p_AD, u, m, m_AD, x, y, OtherState, RotI
CASE (0 ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = 0.0_ReKi; iOut = iOut + 1; enddo;enddo

! ***** Undisturbed wind velocity in inertial, polar, local and airfoil systems*****
CASE( BldNd_VUndxi ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = RotInflow%Bld(iB)%InflowOnBlade(1,iNd); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndyi ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = RotInflow%Bld(iB)%InflowOnBlade(2,iNd); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndzi ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = RotInflow%Bld(iB)%InflowOnBlade(3,iNd); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndxi ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = RotInflow%Blade(iB)%InflowVel(1,iNd); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndyi ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = RotInflow%Blade(iB)%InflowVel(2,iNd); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndzi ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = RotInflow%Blade(iB)%InflowVel(3,iNd); iOut = iOut + 1; enddo;enddo

CASE( BldNd_VUndxp ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), R_pi(1,:,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndyp ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), R_pi(2,:,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndzp ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), R_pi(3,:,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndxp ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), R_pi(1,:,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndyp ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), R_pi(2,:,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndzp ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), R_pi(3,:,iB) ); iOut = iOut + 1; enddo;enddo

CASE( BldNd_VUndxl ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), R_li(1,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndyl ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), R_li(2,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndzl ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), R_li(3,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndxl ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), R_li(1,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndyl ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), R_li(2,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndzl ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), R_li(3,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo

CASE( BldNd_VUndxa ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), u%BladeMotion(iB)%Orientation(1,:,iNd) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndya ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), u%BladeMotion(iB)%Orientation(2,:,iNd) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndza ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), u%BladeMotion(iB)%Orientation(3,:,iNd) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndxa ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), u%BladeMotion(iB)%Orientation(1,:,iNd) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndya ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), u%BladeMotion(iB)%Orientation(2,:,iNd) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndza ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), u%BladeMotion(iB)%Orientation(3,:,iNd) ); iOut = iOut + 1; enddo;enddo

! TODO: deprecate this
CASE( BldNd_VUndx ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), R_wi(1,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndy ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), R_wi(2,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndz ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Bld(iB)%InflowOnBlade(:,iNd), R_wi(3,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndx ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), R_wi(1,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndy ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), R_wi(2,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo
CASE( BldNd_VUndz ); do iB=1,nB; do iNd=1,nNd; y%WriteOutput(iOut) = dot_product( RotInflow%Blade(iB)%InflowVel(:,iNd), R_wi(3,:,iNd,iB) ); iOut = iOut + 1; enddo;enddo


! ***** Disturbed wind velocity in inertial, polar, local and airfoil systems*****
Expand Down
4 changes: 2 additions & 2 deletions modules/aerodyn/src/AeroDyn_IO.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -163,7 +163,7 @@ subroutine Calc_WriteOutput_AD()
do beta=1,p%NTwOuts
j = p%TwOutNd(beta)

tmp = matmul( u%TowerMotion%Orientation(:,:,j) , RotInflow%InflowOnTower(:,j) )
tmp = matmul( u%TowerMotion%Orientation(:,:,j) , RotInflow%Tower%InflowVel(:,j) )
m%AllOuts( TwNVUnd(:,beta) ) = tmp

tmp = matmul( u%TowerMotion%Orientation(:,:,j) , u%TowerMotion%TranslationVel(:,j) )
Expand Down Expand Up @@ -221,7 +221,7 @@ subroutine Calc_WriteOutput_AD()
do beta=1,p%NBlOuts
j=p%BlOutNd(beta)

tmp = matmul( m%orientationAnnulus(:,:,j,k), RotInflow%Bld(k)%InflowOnBlade(:,j) )
tmp = matmul( m%orientationAnnulus(:,:,j,k), RotInflow%Blade(k)%InflowVel(:,j) )
m%AllOuts( BNVUndx(beta,k) ) = tmp(1)
m%AllOuts( BNVUndy(beta,k) ) = tmp(2)
m%AllOuts( BNVUndz(beta,k) ) = tmp(3)
Expand Down
9 changes: 4 additions & 5 deletions modules/aerodyn/src/AeroDyn_Registry.txt
View file
Original file line number Diff line number Diff line change
Expand Up @@ -341,11 +341,10 @@ typedef ^ MiscVarType ReKi WindVel {:}{:} - - "XYZ components of wind
typedef ^ MiscVarType ReKi WindAcc {:}{:} - - "XYZ components of wind acceleration" -

# Inflow data storage
typedef ^ BldInflowType ReKi InflowOnBlade {:}{:} - - "U,V,W at nodes on each blade (note if we change the requirement that NumNodes is the same for each blade, this will need to change)" m/s
typedef ^ BldInflowType ReKi AccelOnBlade {:}{:} - - "Wind acceleration at nodes on each blade (note if we change the requirement that NumNodes is the same for each blade, this will need to change)" m/s
typedef ^ RotInflowType BldInflowType Bld {:} - - "Blade Inputs" -
typedef ^ RotInflowType ReKi InflowOnTower {:}{:} - - "U,V,W at nodes on the tower" m/s
typedef ^ RotInflowType ReKi AccelOnTower {:}{:} - - "Wind acceleration at nodes on the tower" m/s
typedef ^ ElemInflowType ReKi InflowVel {:}{:} - - "U,V,W at nodes on element (note if we change the requirement that NumNodes is the same for each blade, this will need to change)" m/s
typedef ^ ElemInflowType ReKi InflowAcc {:}{:} - - "Wind acceleration at nodes on element (blade or tower) (note if we change the requirement that NumNodes is the same for each blade, this will need to change)" m/s
typedef ^ RotInflowType ElemInflowType Blade {:} - - "Blade wind inputs" -
typedef ^ RotInflowType ElemInflowType Tower - - - "Blade wind inputs" -
typedef ^ RotInflowType ReKi InflowOnHub {3}{1} - - "U,V,W at hub" m/s
typedef ^ RotInflowType ReKi InflowOnNacelle {3}{1} - - "U,V,W at nacelle" m/s
typedef ^ RotInflowType ReKi InflowOnTailFin {3}{1} - - "U,V,W at tailfin" m/s
Expand Down
Loading

0 comments on commit 479bd62

Please sign in to comment.