bugfix: IfW rotor points for disk average incorrect

modules/inflowwind/src/IfW_FlowField.f90

@@ -1591,6 +1591,7 @@ subroutine Grid4DField_GetVel(G4D, Time, Position, Velocity, ErrStat, ErrMsg)
    real(ReKi)                          :: P(3, 16)    ! Point values
    real(ReKi)                          :: tmp
    integer(IntKi)                      :: i
+   character(60)                       :: PtLoc
 
    ErrStat = ErrID_None
    ErrMsg = ""
@@ -1627,11 +1628,12 @@ subroutine Grid4DField_GetVel(G4D, Time, Position, Velocity, ErrStat, ErrMsg)
    do i = 1, 4
       if (Indx_Lo(i) <= 0) then
          Indx_Lo(i) = 1
-         call SetErrStat(ErrID_Fatal, 'Outside the grid bounds.', ErrStat, ErrMsg, RoutineName)
+         write(PtLoc,'(A1,3(f8.2,A1))') '(',Position(1),',',Position(2),',',Position(3),')'
+         call SetErrStat(ErrID_Fatal, 'Outside the grid bounds: '//trim(PtLoc), ErrStat, ErrMsg, RoutineName)
          return
       elseif (Indx_Lo(i) >= G4D%n(i)) then
-         Indx_Lo(i) = max(G4D%n(i) - 1, 1)           ! make sure it's a valid index
-         call SetErrStat(ErrID_Fatal, 'Outside the grid bounds.', ErrStat, ErrMsg, RoutineName)
+         write(PtLoc,'(A1,3(f8.2,A1))') '(',Position(1),',',Position(2),',',Position(3),')'
+         call SetErrStat(ErrID_Fatal, 'Outside the grid bounds: '//trim(PtLoc), ErrStat, ErrMsg, RoutineName)
          return
       end if
       Indx_Hi(i) = min(Indx_Lo(i) + 1, G4D%n(i))     ! make sure it's a valid index

modules/inflowwind/src/InflowWind_Subs.f90

@@ -967,6 +967,7 @@ SUBROUTINE InflowWind_SetParameters( InitInp, InputFileData, p, m, ErrStat, ErrM
       ! Temporary variables
    INTEGER(IntKi)                                     :: TmpErrStat           !< Temporary error status  for subroutine and function calls
    CHARACTER(ErrMsgLen)                               :: TmpErrMsg            !< Temporary error message for subroutine and function calls
+   integer(IntKi)                                     :: NumPtsAvg            !< Number of points to use for disk average vel (1 if no radius)
 
       ! Local variables
    INTEGER(IntKi)                                     :: I                    !< Generic counter
@@ -1029,28 +1030,34 @@ SUBROUTINE InflowWind_SetParameters( InitInp, InputFileData, p, m, ErrStat, ErrM
    CALL AllocAry( m%y_Hub%VelocityUVW, 3, 1, "Array of velocities for hub values", TmpErrStat, TmpErrMsg )
       CALL SetErrStat(TmpErrStat,TmpErrMsg,ErrStat,ErrMsg,RoutineName)
 
-   CALL AllocAry( m%u_Avg%PositionXYZ, 3, IfW_NumPtsAvg, "Array of positions for rotor-averaged values", TmpErrStat, TmpErrMsg )
+   ! Disk Velocity calculations
+   if (InitInp%RadAvg < 0.0_ReKi)  then
+      NumPtsAvg = 1_IntKi              ! Use only hub point
+   else
+      NumPtsAvg = IfW_NumPtsAvg        ! Use a field of points for disk average calculations
+   endif
+
+   CALL AllocAry( m%u_Avg%PositionXYZ, 3, NumPtsAvg, "Array of positions for rotor-averaged values", TmpErrStat, TmpErrMsg )
       CALL SetErrStat(TmpErrStat,TmpErrMsg,ErrStat,ErrMsg,RoutineName)
-   CALL AllocAry( m%y_Avg%VelocityUVW, 3, IfW_NumPtsAvg, "Array of velocities for rotor-averaged values", TmpErrStat, TmpErrMsg )
+   CALL AllocAry( m%y_Avg%VelocityUVW, 3, NumPtsAvg, "Array of velocities for rotor-averaged values", TmpErrStat, TmpErrMsg )
       CALL SetErrStat(TmpErrStat,TmpErrMsg,ErrStat,ErrMsg,RoutineName)
-   CALL AllocAry( p%PositionAvg, 3, IfW_NumPtsAvg, "Array of positions for computing average wind speed", TmpErrStat, TmpErrMsg )
+   CALL AllocAry( p%PositionAvg, 3, NumPtsAvg, "Array of positions for computing average wind speed", TmpErrStat, TmpErrMsg )
       CALL SetErrStat(TmpErrStat,TmpErrMsg,ErrStat,ErrMsg,RoutineName)
       IF ( ErrStat>= AbortErrLev ) RETURN
 
 
    if (InitInp%RadAvg < 0.0_ReKi)  then
-      R  =  max(1.0_ReKi, InputFileData%Uniform_RefLength)/2.0_ReKi    ! We'll use this as a guess for the rotor radius
+      p%PositionAvg = 0.0_ReKi         ! Use only hub
    else
-      R = InitInp%RadAvg
+      ! Calculate a ring of points at 70 rotor radius
+      R = InitInp%RadAvg * 0.7_ReKi !70% radius
+      do i=1,NumPtsAvg
+         theta = pi +(i-1)*TwoPi/NumPtsAvg
+         p%PositionAvg(1,i) = 0.0_ReKi          ! Hub X (perpindicular to rotor plane)
+         p%PositionAvg(2,i) = R*cos(theta)      ! Hub Y
+         p%PositionAvg(3,i) = R*sin(theta)      ! Hub Z (in vertical plane when azimuth=0)
+      end do
    end if
-   R = R * 0.7_ReKi !70% radius
-
-   do i=1,IfW_NumPtsAvg
-      theta = pi +(i-1)*TwoPi/IfW_NumPtsAvg
-      p%PositionAvg(1,i) = R*cos(theta)
-      p%PositionAvg(2,i) = R*sin(theta)
-      p%PositionAvg(3,i) = 0.0_ReKi
-   end do
 
    p%OutputAccel = InitInp%OutputAccel
 
@@ -1709,17 +1716,17 @@ SUBROUTINE InflowWind_GetSpatialAverage( Time, InputData, p, x, xd, z, OtherStat
    m%u_Avg%HubPosition    = InputData%HubPosition
    m%u_Avg%HubOrientation = InputData%HubOrientation
 
-   do i=1,IfW_NumPtsAvg
+   do i=1,size(m%u_Avg%PositionXYZ,DIM=2)
       m%u_Avg%PositionXYZ(:,i) = matmul(InputData%HubOrientation,p%PositionAvg(:,i)) + InputData%HubPosition
    end do
 
    CALL CalculateOutput( Time, m%u_Avg, p, x, xd, z, OtherStates, m%y_Avg, m, FillWrOut, ErrStat2, ErrMsg2 )
       CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
 
-   do i=1,IfW_NumPtsAvg
+   do i=1,size(m%u_Avg%PositionXYZ,DIM=2)
       MeanVelocity = MeanVelocity + m%y_Avg%VelocityUVW(:,i)
    end do
-   MeanVelocity = MeanVelocity / REAL(IfW_NumPtsAvg,ReKi)
+   MeanVelocity = MeanVelocity / REAL(size(m%u_Avg%PositionXYZ,DIM=2),ReKi)
 
 
    RETURN