diff --git a/modules/aerodyn/src/AeroAcoustics/AeroAcoustics.f90 b/modules/aerodyn/src/AeroAcoustics/AeroAcoustics.f90 index 07f9f76f8..95c23c05b 100644 --- a/modules/aerodyn/src/AeroAcoustics/AeroAcoustics.f90 +++ b/modules/aerodyn/src/AeroAcoustics/AeroAcoustics.f90 @@ -443,10 +443,10 @@ subroutine Init_u( u, p, InputFileData, InitInp, errStat, errMsg ) call AllocAry(u%Vrel , p%NumBlNds, p%numBlades, 'u%Vrel' , errStat2 , errMsg2); if(Failed()) return call AllocAry(u%AeroCent_G, 3 , p%NumBlNds , p%numBlades , 'u%AeroCent_G', errStat2 , errMsg2); if(Failed()) return call AllocAry(u%Inflow , 3_IntKi , p%NumBlNds , p%numBlades , 'u%Inflow' , ErrStat2 , ErrMsg2); if(Failed()) return - call AllocAry(u%RotLtoG , 3 , 3 , p%NumBlNds , p%numBlades , 'u%RotLtoG' , errStat2 , errMsg2); if(Failed()) return + call AllocAry(u%RotGtoL , 3 , 3 , p%NumBlNds , p%numBlades , 'u%RotGtoL' , errStat2 , errMsg2); if(Failed()) return u%AoANoise = 0.0_Reki u%Vrel = 0.0_Reki - u%RotLtoG = 0.0_Reki + u%RotGtoL = 0.0_Reki u%AeroCent_G = 0.0_Reki u%Inflow = 0.0_Reki contains @@ -813,22 +813,20 @@ SUBROUTINE CalcObserve(p,m,u,xd,nt,errStat,errMsg) INTEGER(IntKi), intent( out) :: errStat !< Error status of the operation CHARACTER(*), intent( out) :: errMsg !< Error message if ErrStat /= ErrID_None ! Local variables. - REAL(ReKi) :: RLEObserve (3) ! position vector from leading edge to observer in trailing edge coordinate system - !REAL(ReKi) :: BladeObserve (3) ! position vector from leading edge to observer in trailing edge coordinate system - REAL(ReKi) :: RTEObserve (3) ! position vector from trailing edge to observer in trailing edge coordinate system - REAL(ReKi) :: RTEObservereal (3) ! location of trailing edge in global coordinate system - REAL(ReKi) :: RLEObservereal (3) ! location of leading edge in global coordinate system - REAL(ReKi) :: LocalToGlobal(3,3) ! trasnformation matrix - REAL(ReKi) :: rSLE (3) ! Distance from tower base to leading edge in trailing edge coordinate system - REAL(ReKi) :: rSTE (3) ! Distance from tower base to trailing edge in trailing edge coordinate system + REAL(ReKi) :: RLEObserve (3) ! Position vector from leading edge to observer in trailing edge coordinate system + REAL(ReKi) :: RTEObserve (3) ! Position vector from trailing edge to observer in trailing edge coordinate system + REAL(ReKi) :: RTEObserveG (3) ! Position vector from trailing edge to observer in the coordinate system located at the trailing edge and rotated as the global + REAL(ReKi) :: RLEObserveG (3) ! Position vector from leading edge to observer in the coordinate system located at the leading edge and rotated as the global + REAL(ReKi) :: RTEObservereal (3) ! Location of trailing edge in global coordinate system + REAL(ReKi) :: RLEObservereal (3) ! Location of leading edge in global coordinate system + REAL(ReKi) :: LocalToGlobal(3,3) ! Transformation matrix REAL(ReKi) :: timeLE ! Time of sound propagation from leading edge to observer REAL(ReKi) :: timeTE ! Time of sound propagation from trailing edge to observer - REAL(ReKi) :: tmpR (3) ! temporary distance vector - REAL(ReKi) :: UConvect (3) ! convection velocity of noise source in trailing edge coordinate system + REAL(ReKi) :: phi_e ! Spanwise directivity angle + REAL(ReKi) :: theta_e ! Chordwise directivity angle INTEGER(intKi) :: I ! I A generic index for DO loops. INTEGER(intKi) :: J ! J A generic index for DO loops. INTEGER(intKi) :: K ! K A generic index for DO loops. - INTEGER(intKi) :: cou ! K A generic index for DO loops. INTEGER(intKi) :: ErrStat2 CHARACTER(ErrMsgLen) :: ErrMsg2 CHARACTER(*), parameter :: RoutineName = 'CalcObserveDist' @@ -836,52 +834,65 @@ SUBROUTINE CalcObserve(p,m,u,xd,nt,errStat,errMsg) ErrStat = ErrID_None ErrMsg = "" - + ! Loop thoruhg the blades DO I = 1,p%numBlades + ! Loop through the nodes along blade span DO J = 1,p%NumBlNds - ! Transpose the matrix RotLtoG - LocalToGlobal = TRANSPOSE(u%RotLtoG(:,:,J,I)) + ! Transpose the rotational vector GlobalToLocal to obtain the rotation LocalToGlobal + LocalToGlobal = TRANSPOSE(u%RotGtoL(:,:,J,I)) + ! Rotate the coordinates of leading and trailing edge from the local reference system to the global. Then add the coordinates of the aerodynamic center in the global coordinate system + ! The global coordinate system is located on the ground, has x pointing downwind, y pointing laterally, and z pointing vertically upwards RTEObservereal = MATMUL(LocalToGlobal, p%AFTeCo(:,J,I)) + u%AeroCent_G(:,J,I) RLEObservereal = MATMUL(LocalToGlobal, p%AFLeCo(:,J,I)) + u%AeroCent_G(:,J,I) - - m%LE_Location(1,J,I) = RLEObservereal(1) ! - m%LE_Location(2,J,I) = RLEObservereal(2) ! - m%LE_Location(3,J,I) = RLEObservereal(3) ! the height of leading edge + ! Compute the coordinates of the leading edge in the global coordinate system + m%LE_Location(1,J,I) = RLEObservereal(1) + m%LE_Location(2,J,I) = RLEObservereal(2) + m%LE_Location(3,J,I) = RLEObservereal(3) + ! If the time step is set to generate AA outputs IF (nt.gt.p%Comp_AA_after) THEN IF ( (mod(nt,p%saveeach).eq.0) ) THEN - + ! Loop through the observers DO K = 1,p%NrObsLoc - ! Calculate position vector from leading and trailing edge to observer in retarded trailing edge coordinate system - RTEObserve(1)=ABS(p%Obsx(K)-RTEObservereal(1)) - RTEObserve(2)=ABS(p%Obsy(K)-RTEObservereal(2)) - RTEObserve(3)=ABS(p%Obsz(K)-RTEObservereal(3)) - - RLEObserve(1)=ABS(p%Obsx(K)-RLEObservereal(1)) - RLEObserve(2)=ABS(p%Obsy(K)-RLEObservereal(2)) - RLEObserve(3)=ABS(p%Obsz(K)-RLEObservereal(3)) - ! - ! Calculate convection velocity of noise source - ! Assumes noise source convects at some constant times the mean wind speed, approximately accounts for - ! induction velocity and change in convection velocity as noise propagates to observer (likely on the ground) - UConvect (1) = 0.8*xd%MeanVxVyVz(J,I) - UConvect (2) = 0.8*xd%MeanVxVyVz(J,I) - UConvect (3) = 0.8*xd%MeanVxVyVz(J,I) - ! Calculate time of noise propagation to observer - timeTE = SQRT (RTEObserve(1)**2+RTEObserve(2)**2+RTEObserve(3)**2)/p%SpdSound - timeLE = SQRT (RLEObserve(1)**2+RLEObserve(2)**2+RLEObserve(3)**2)/p%SpdSound - ! Calculate inputs into noise subroutines + ! Calculate the position of the observer K in a reference system located at the trailing edge and oriented as the global reference system + RTEObserveG(1)=p%Obsx(K)-RTEObservereal(1) + RTEObserveG(2)=p%Obsy(K)-RTEObservereal(2) + RTEObserveG(3)=p%Obsz(K)-RTEObservereal(3) + ! Calculate the position of the observer K in a reference system located at the leading edge and oriented as the global reference system + RLEObserveG(1)=p%Obsx(K)-RLEObservereal(1) + RLEObserveG(2)=p%Obsy(K)-RLEObservereal(2) + RLEObserveG(3)=p%Obsz(K)-RLEObservereal(3) + ! Rotate back the two reference systems from global to local. + RTEObserve = MATMUL(u%RotGtoL(:,:,J,I), RTEObserveG) + RLEObserve = MATMUL(u%RotGtoL(:,:,J,I), RLEObserveG) + + ! Calculate absolute distance between node and observer m%rTEtoObserve(K,J,I) = SQRT (RTEObserve(1)**2+RTEObserve(2)**2+RTEObserve(3)**2) m%rLEtoObserve(K,J,I) = SQRT (RLEObserve(1)**2+RLEObserve(2)**2+RLEObserve(3)**2) - m%ChordAngleTE(K,J,I) = ACOS (RTEObserve(2)/SQRT(RTEObserve(1)**2+RTEObserve(2)**2+RTEObserve(3)**2))*R2D ! theta - m%SpanAngleTE(K,J,I) = ACOS (RTEObserve(3)/SQRT(RTEObserve(1)**2+RTEObserve(3)**2))*R2D !phi - IF (m%SpanAngleTE(K,J,I)< 0) m%SpanAngleTE(K,J,I)= 180+m%SpanAngleTE(K,J,I) - IF (m%ChordAngleTE(K,J,I)< 0) m%ChordAngleTE(K,J,I)= 180+m%ChordAngleTE(K,J,I) + ! Calculate time of noise propagation to observer + timeTE = m%rTEtoObserve(K,J,I) / p%SpdSound + timeLE = m%rLEtoObserve(K,J,I) / p%SpdSound + + ! The local system has y alinged with the chord, x pointing towards the airfoil suction side, and z aligned with blade span from root towards tip + ! x ---> z_e + ! y ---> x_e + ! z ---> y_e + + ! Compute spanwise directivity angle phi for the trailing edge + phi_e = ATAN (RTEObserve(1)/RTEObserve(3)) + m%SpanAngleTE(K,J,I) = phi_e * R2D + + ! Compute chordwise directivity angle theta for the trailing edge + theta_e = ATAN ((RTEObserve(3) * COS (phi_e) + RTEObserve(1) * SIN (phi_e) ) / RTEObserve(2)) + m%ChordAngleTE(K,J,I) = theta_e * R2D + + ! Compute spanwise directivity angle phi for the leading edge (it's the same angle for the trailing edge) + phi_e = ATAN (RLEObserve(1)/RLEObserve(3)) + m%SpanAngleLE(K,J,I) = phi_e * R2D - m%ChordAngleLE(K,J,I) = ACOS (RLEObserve(2)/SQRT(RLEObserve(1)**2+RLEObserve(2)**2+RLEObserve(3)**2))*R2D - m%SpanAngleLE(K,J,I) = ACOS (RLEObserve(3)/SQRT(RLEObserve(1)**2+RLEObserve(3)**2))*R2D - IF (m%SpanAngleLE(K,J,I)< 0) m%SpanAngleLE(K,J,I)= 180+m%SpanAngleLE(K,J,I) - IF (m%ChordAngleLE(K,J,I)< 0) m%ChordAngleLE(K,J,I)= 180+m%ChordAngleLE(K,J,I) + ! Compute chordwise directivity angle theta for the leading edge + theta_e = ATAN ((RLEObserve(3) * COS (phi_e) + RLEObserve(1) * SIN (phi_e) ) / RLEObserve(2)) + m%ChordAngleLE(K,J,I) = theta_e * R2D ENDDO !K, observers ENDIF ! every Xth time step or so.. diff --git a/modules/aerodyn/src/AeroAcoustics_Registry.txt b/modules/aerodyn/src/AeroAcoustics_Registry.txt index 2044574fa..46d89ec1b 100644 --- a/modules/aerodyn/src/AeroAcoustics_Registry.txt +++ b/modules/aerodyn/src/AeroAcoustics_Registry.txt @@ -229,7 +229,7 @@ typedef ^ ParameterType ReKi AFThick # ..... Inputs .................................................................................................................... # Define inputs that are contained on the mesh here: -typedef ^ InputType ReKi RotLtoG {:}{:}{:}{:} - - "3x3 rotation matrix transform a vector from the local airfoil coordinate system to the global inertial coordinate system" - +typedef ^ InputType ReKi RotGtoL {:}{:}{:}{:} - - "3x3 rotation matrix transform a vector from the local airfoil coordinate system to the global inertial coordinate system" - typedef ^ InputType ReKi AeroCent_G {:}{:}{:} - - "location in global coordinates of the blade element aerodynamic center. 1st index = vector components, 2nd index = blade node, 3rd index = blade" - typedef ^ InputType ReKi Vrel {:}{:} - - "Vrel" - typedef ^ InputType ReKi AoANoise {:}{:} - - "Angle of attack" - diff --git a/modules/aerodyn/src/AeroAcoustics_Types.f90 b/modules/aerodyn/src/AeroAcoustics_Types.f90 index 8312a8734..21e675c8a 100644 --- a/modules/aerodyn/src/AeroAcoustics_Types.f90 +++ b/modules/aerodyn/src/AeroAcoustics_Types.f90 @@ -256,7 +256,7 @@ MODULE AeroAcoustics_Types ! ======================= ! ========= AA_InputType ======= TYPE, PUBLIC :: AA_InputType - REAL(ReKi) , DIMENSION(:,:,:,:), ALLOCATABLE :: RotLtoG !< 3x3 rotation matrix transform a vector from the local airfoil coordinate system to the global inertial coordinate system [-] + REAL(ReKi) , DIMENSION(:,:,:,:), ALLOCATABLE :: RotGtoL !< 3x3 rotation matrix transform a vector from the local airfoil coordinate system to the global inertial coordinate system [-] REAL(ReKi) , DIMENSION(:,:,:), ALLOCATABLE :: AeroCent_G !< location in global coordinates of the blade element aerodynamic center. 1st index = vector components, 2nd index = blade node, 3rd index = blade [-] REAL(ReKi) , DIMENSION(:,:), ALLOCATABLE :: Vrel !< Vrel [-] REAL(ReKi) , DIMENSION(:,:), ALLOCATABLE :: AoANoise !< Angle of attack [-] @@ -8612,23 +8612,23 @@ SUBROUTINE AA_CopyInput( SrcInputData, DstInputData, CtrlCode, ErrStat, ErrMsg ) ! ErrStat = ErrID_None ErrMsg = "" -IF (ALLOCATED(SrcInputData%RotLtoG)) THEN - i1_l = LBOUND(SrcInputData%RotLtoG,1) - i1_u = UBOUND(SrcInputData%RotLtoG,1) - i2_l = LBOUND(SrcInputData%RotLtoG,2) - i2_u = UBOUND(SrcInputData%RotLtoG,2) - i3_l = LBOUND(SrcInputData%RotLtoG,3) - i3_u = UBOUND(SrcInputData%RotLtoG,3) - i4_l = LBOUND(SrcInputData%RotLtoG,4) - i4_u = UBOUND(SrcInputData%RotLtoG,4) - IF (.NOT. ALLOCATED(DstInputData%RotLtoG)) THEN - ALLOCATE(DstInputData%RotLtoG(i1_l:i1_u,i2_l:i2_u,i3_l:i3_u,i4_l:i4_u),STAT=ErrStat2) - IF (ErrStat2 /= 0) THEN - CALL SetErrStat(ErrID_Fatal, 'Error allocating DstInputData%RotLtoG.', ErrStat, ErrMsg,RoutineName) +IF (ALLOCATED(SrcInputData%RotGtoL)) THEN + i1_l = LBOUND(SrcInputData%RotGtoL,1) + i1_u = UBOUND(SrcInputData%RotGtoL,1) + i2_l = LBOUND(SrcInputData%RotGtoL,2) + i2_u = UBOUND(SrcInputData%RotGtoL,2) + i3_l = LBOUND(SrcInputData%RotGtoL,3) + i3_u = UBOUND(SrcInputData%RotGtoL,3) + i4_l = LBOUND(SrcInputData%RotGtoL,4) + i4_u = UBOUND(SrcInputData%RotGtoL,4) + IF (.NOT. ALLOCATED(DstInputData%RotGtoL)) THEN + ALLOCATE(DstInputData%RotGtoL(i1_l:i1_u,i2_l:i2_u,i3_l:i3_u,i4_l:i4_u),STAT=ErrStat2) + IF (ErrStat2 /= 0) THEN + CALL SetErrStat(ErrID_Fatal, 'Error allocating DstInputData%RotGtoL.', ErrStat, ErrMsg,RoutineName) RETURN END IF END IF - DstInputData%RotLtoG = SrcInputData%RotLtoG + DstInputData%RotGtoL = SrcInputData%RotGtoL ENDIF IF (ALLOCATED(SrcInputData%AeroCent_G)) THEN i1_l = LBOUND(SrcInputData%AeroCent_G,1) @@ -8701,8 +8701,8 @@ SUBROUTINE AA_DestroyInput( InputData, ErrStat, ErrMsg ) ! ErrStat = ErrID_None ErrMsg = "" -IF (ALLOCATED(InputData%RotLtoG)) THEN - DEALLOCATE(InputData%RotLtoG) +IF (ALLOCATED(InputData%RotGtoL)) THEN + DEALLOCATE(InputData%RotGtoL) ENDIF IF (ALLOCATED(InputData%AeroCent_G)) THEN DEALLOCATE(InputData%AeroCent_G) @@ -8753,10 +8753,10 @@ SUBROUTINE AA_PackInput( ReKiBuf, DbKiBuf, IntKiBuf, Indata, ErrStat, ErrMsg, Si Re_BufSz = 0 Db_BufSz = 0 Int_BufSz = 0 - Int_BufSz = Int_BufSz + 1 ! RotLtoG allocated yes/no - IF ( ALLOCATED(InData%RotLtoG) ) THEN - Int_BufSz = Int_BufSz + 2*4 ! RotLtoG upper/lower bounds for each dimension - Re_BufSz = Re_BufSz + SIZE(InData%RotLtoG) ! RotLtoG + Int_BufSz = Int_BufSz + 1 ! RotGtoL allocated yes/no + IF ( ALLOCATED(InData%RotGtoL) ) THEN + Int_BufSz = Int_BufSz + 2*4 ! RotGtoL upper/lower bounds for each dimension + Re_BufSz = Re_BufSz + SIZE(InData%RotGtoL) ! RotGtoL END IF Int_BufSz = Int_BufSz + 1 ! AeroCent_G allocated yes/no IF ( ALLOCATED(InData%AeroCent_G) ) THEN @@ -8805,27 +8805,27 @@ SUBROUTINE AA_PackInput( ReKiBuf, DbKiBuf, IntKiBuf, Indata, ErrStat, ErrMsg, Si Db_Xferred = 1 Int_Xferred = 1 - IF ( .NOT. ALLOCATED(InData%RotLtoG) ) THEN + IF ( .NOT. ALLOCATED(InData%RotGtoL) ) THEN IntKiBuf( Int_Xferred ) = 0 Int_Xferred = Int_Xferred + 1 ELSE IntKiBuf( Int_Xferred ) = 1 Int_Xferred = Int_Xferred + 1 - IntKiBuf( Int_Xferred ) = LBOUND(InData%RotLtoG,1) - IntKiBuf( Int_Xferred + 1) = UBOUND(InData%RotLtoG,1) + IntKiBuf( Int_Xferred ) = LBOUND(InData%RotGtoL,1) + IntKiBuf( Int_Xferred + 1) = UBOUND(InData%RotGtoL,1) Int_Xferred = Int_Xferred + 2 - IntKiBuf( Int_Xferred ) = LBOUND(InData%RotLtoG,2) - IntKiBuf( Int_Xferred + 1) = UBOUND(InData%RotLtoG,2) + IntKiBuf( Int_Xferred ) = LBOUND(InData%RotGtoL,2) + IntKiBuf( Int_Xferred + 1) = UBOUND(InData%RotGtoL,2) Int_Xferred = Int_Xferred + 2 - IntKiBuf( Int_Xferred ) = LBOUND(InData%RotLtoG,3) - IntKiBuf( Int_Xferred + 1) = UBOUND(InData%RotLtoG,3) + IntKiBuf( Int_Xferred ) = LBOUND(InData%RotGtoL,3) + IntKiBuf( Int_Xferred + 1) = UBOUND(InData%RotGtoL,3) Int_Xferred = Int_Xferred + 2 - IntKiBuf( Int_Xferred ) = LBOUND(InData%RotLtoG,4) - IntKiBuf( Int_Xferred + 1) = UBOUND(InData%RotLtoG,4) + IntKiBuf( Int_Xferred ) = LBOUND(InData%RotGtoL,4) + IntKiBuf( Int_Xferred + 1) = UBOUND(InData%RotGtoL,4) Int_Xferred = Int_Xferred + 2 - IF (SIZE(InData%RotLtoG)>0) ReKiBuf ( Re_Xferred:Re_Xferred+(SIZE(InData%RotLtoG))-1 ) = PACK(InData%RotLtoG,.TRUE.) - Re_Xferred = Re_Xferred + SIZE(InData%RotLtoG) + IF (SIZE(InData%RotGtoL)>0) ReKiBuf ( Re_Xferred:Re_Xferred+(SIZE(InData%RotGtoL))-1 ) = PACK(InData%RotGtoL,.TRUE.) + Re_Xferred = Re_Xferred + SIZE(InData%RotGtoL) END IF IF ( .NOT. ALLOCATED(InData%AeroCent_G) ) THEN IntKiBuf( Int_Xferred ) = 0 @@ -8935,7 +8935,7 @@ SUBROUTINE AA_UnPackInput( ReKiBuf, DbKiBuf, IntKiBuf, Outdata, ErrStat, ErrMsg Re_Xferred = 1 Db_Xferred = 1 Int_Xferred = 1 - IF ( IntKiBuf( Int_Xferred ) == 0 ) THEN ! RotLtoG not allocated + IF ( IntKiBuf( Int_Xferred ) == 0 ) THEN ! RotGtoL not allocated Int_Xferred = Int_Xferred + 1 ELSE Int_Xferred = Int_Xferred + 1 @@ -8951,10 +8951,10 @@ SUBROUTINE AA_UnPackInput( ReKiBuf, DbKiBuf, IntKiBuf, Outdata, ErrStat, ErrMsg i4_l = IntKiBuf( Int_Xferred ) i4_u = IntKiBuf( Int_Xferred + 1) Int_Xferred = Int_Xferred + 2 - IF (ALLOCATED(OutData%RotLtoG)) DEALLOCATE(OutData%RotLtoG) - ALLOCATE(OutData%RotLtoG(i1_l:i1_u,i2_l:i2_u,i3_l:i3_u,i4_l:i4_u),STAT=ErrStat2) + IF (ALLOCATED(OutData%RotGtoL)) DEALLOCATE(OutData%RotGtoL) + ALLOCATE(OutData%RotGtoL(i1_l:i1_u,i2_l:i2_u,i3_l:i3_u,i4_l:i4_u),STAT=ErrStat2) IF (ErrStat2 /= 0) THEN - CALL SetErrStat(ErrID_Fatal, 'Error allocating OutData%RotLtoG.', ErrStat, ErrMsg,RoutineName) + CALL SetErrStat(ErrID_Fatal, 'Error allocating OutData%RotGtoL.', ErrStat, ErrMsg,RoutineName) RETURN END IF ALLOCATE(mask4(i1_l:i1_u,i2_l:i2_u,i3_l:i3_u,i4_l:i4_u),STAT=ErrStat2) @@ -8963,8 +8963,8 @@ SUBROUTINE AA_UnPackInput( ReKiBuf, DbKiBuf, IntKiBuf, Outdata, ErrStat, ErrMsg RETURN END IF mask4 = .TRUE. - IF (SIZE(OutData%RotLtoG)>0) OutData%RotLtoG = UNPACK(ReKiBuf( Re_Xferred:Re_Xferred+(SIZE(OutData%RotLtoG))-1 ), mask4, 0.0_ReKi ) - Re_Xferred = Re_Xferred + SIZE(OutData%RotLtoG) + IF (SIZE(OutData%RotGtoL)>0) OutData%RotGtoL = UNPACK(ReKiBuf( Re_Xferred:Re_Xferred+(SIZE(OutData%RotGtoL))-1 ), mask4, 0.0_ReKi ) + Re_Xferred = Re_Xferred + SIZE(OutData%RotGtoL) DEALLOCATE(mask4) END IF IF ( IntKiBuf( Int_Xferred ) == 0 ) THEN ! AeroCent_G not allocated @@ -10031,13 +10031,13 @@ SUBROUTINE AA_Input_ExtrapInterp1(u1, u2, tin, u_out, tin_out, ErrStat, ErrMsg ) CALL SetErrStat(ErrID_Fatal, 't(1) must not equal t(2) to avoid a division-by-zero error.', ErrStat, ErrMsg,RoutineName) RETURN END IF -IF (ALLOCATED(u_out%RotLtoG) .AND. ALLOCATED(u1%RotLtoG)) THEN - ALLOCATE(b4(SIZE(u_out%RotLtoG,1),SIZE(u_out%RotLtoG,2), & - SIZE(u_out%RotLtoG,3),SIZE(u_out%RotLtoG,4) )) - ALLOCATE(c4(SIZE(u_out%RotLtoG,1),SIZE(u_out%RotLtoG,2), & - SIZE(u_out%RotLtoG,3),SIZE(u_out%RotLtoG,4) )) - b4 = -(u1%RotLtoG - u2%RotLtoG)/t(2) - u_out%RotLtoG = u1%RotLtoG + b4 * t_out +IF (ALLOCATED(u_out%RotGtoL) .AND. ALLOCATED(u1%RotGtoL)) THEN + ALLOCATE(b4(SIZE(u_out%RotGtoL,1),SIZE(u_out%RotGtoL,2), & + SIZE(u_out%RotGtoL,3),SIZE(u_out%RotGtoL,4) )) + ALLOCATE(c4(SIZE(u_out%RotGtoL,1),SIZE(u_out%RotGtoL,2), & + SIZE(u_out%RotGtoL,3),SIZE(u_out%RotGtoL,4) )) + b4 = -(u1%RotGtoL - u2%RotGtoL)/t(2) + u_out%RotGtoL = u1%RotGtoL + b4 * t_out DEALLOCATE(b4) DEALLOCATE(c4) END IF ! check if allocated @@ -10137,14 +10137,14 @@ SUBROUTINE AA_Input_ExtrapInterp2(u1, u2, u3, tin, u_out, tin_out, ErrStat, ErrM CALL SetErrStat(ErrID_Fatal, 't(1) must not equal t(3) to avoid a division-by-zero error.', ErrStat, ErrMsg,RoutineName) RETURN END IF -IF (ALLOCATED(u_out%RotLtoG) .AND. ALLOCATED(u1%RotLtoG)) THEN - ALLOCATE(b4(SIZE(u_out%RotLtoG,1),SIZE(u_out%RotLtoG,2), & - SIZE(u_out%RotLtoG,3),SIZE(u_out%RotLtoG,4) )) - ALLOCATE(c4(SIZE(u_out%RotLtoG,1),SIZE(u_out%RotLtoG,2), & - SIZE(u_out%RotLtoG,3),SIZE(u_out%RotLtoG,4) )) - b4 = (t(3)**2*(u1%RotLtoG - u2%RotLtoG) + t(2)**2*(-u1%RotLtoG + u3%RotLtoG))/(t(2)*t(3)*(t(2) - t(3))) - c4 = ( (t(2)-t(3))*u1%RotLtoG + t(3)*u2%RotLtoG - t(2)*u3%RotLtoG ) / (t(2)*t(3)*(t(2) - t(3))) - u_out%RotLtoG = u1%RotLtoG + b4 * t_out + c4 * t_out**2 +IF (ALLOCATED(u_out%RotGtoL) .AND. ALLOCATED(u1%RotGtoL)) THEN + ALLOCATE(b4(SIZE(u_out%RotGtoL,1),SIZE(u_out%RotGtoL,2), & + SIZE(u_out%RotGtoL,3),SIZE(u_out%RotGtoL,4) )) + ALLOCATE(c4(SIZE(u_out%RotGtoL,1),SIZE(u_out%RotGtoL,2), & + SIZE(u_out%RotGtoL,3),SIZE(u_out%RotGtoL,4) )) + b4 = (t(3)**2*(u1%RotGtoL - u2%RotGtoL) + t(2)**2*(-u1%RotGtoL + u3%RotGtoL))/(t(2)*t(3)*(t(2) - t(3))) + c4 = ( (t(2)-t(3))*u1%RotGtoL + t(3)*u2%RotGtoL - t(2)*u3%RotGtoL ) / (t(2)*t(3)*(t(2) - t(3))) + u_out%RotGtoL = u1%RotGtoL + b4 * t_out + c4 * t_out**2 DEALLOCATE(b4) DEALLOCATE(c4) END IF ! check if allocated diff --git a/modules/aerodyn/src/AeroDyn.f90 b/modules/aerodyn/src/AeroDyn.f90 index 5a0da3578..d4c245ba3 100644 --- a/modules/aerodyn/src/AeroDyn.f90 +++ b/modules/aerodyn/src/AeroDyn.f90 @@ -1495,8 +1495,7 @@ subroutine SetInputsForAA(p, u, m, errStat, errMsg) do j=1,p%NumBlades do i = 1,p%NumBlNds ! Get local orientation matrix to transform from blade element coordinates to global coordinates - !m%AA_u%RotLtoG(:,:,i,j) = m%WithoutSweepPitchTwist(:,:,i,j) - m%AA_u%RotLtoG(:,:,i,j) = u%BladeMotion(j)%Orientation(:,:,i) + m%AA_u%RotGtoL(:,:,i,j) = u%BladeMotion(j)%Orientation(:,:,i) ! Get blade element aerodynamic center in global coordinates m%AA_u%AeroCent_G(:,i,j) = u%BladeMotion(j)%Position(:,i) + u%BladeMotion(j)%TranslationDisp(:,i)