Merge branch 'dev_NewTurbines'

johnjasa · Oct 8, 2019 · 0ffa2ab · 0ffa2ab
2 parents c9ee384 + d47c8be
commit 0ffa2ab
Show file tree

Hide file tree

Showing 32 changed files with 1,163 additions and 684 deletions.
diff --git a/.gitignore b/.gitignore
@@ -14,6 +14,7 @@ _build
 
 # MATLAB files
 *.m~
+*slprj*
 
 # ide settings and files
 .idea

diff --git a/DRC_Fortran/Source/ControllerBlocks.f90 b/DRC_Fortran/Source/ControllerBlocks.f90
@@ -148,7 +148,7 @@ SUBROUTINE WindSpeedEstimator(LocalVar, CntrPar, objInst, PerfData)
         ! Extended Kalman Filter (EKF) implementation
         ELSEIF (CntrPar%WE_Mode == 2) THEN
             ! Define contant values
-            L = 2.0 * CntrPar%WE_BladeRadius
+            L = 4.0 * CntrPar%WE_BladeRadius
             Ti = 0.2
             R_m = 0.02
             H = RESHAPE((/1.0 , 0.0 , 0.0/),(/1,3/))
@@ -185,10 +185,10 @@ SUBROUTINE WindSpeedEstimator(LocalVar, CntrPar, objInst, PerfData)
                 F(2,3) = PI * v_t/(2.0*L)
 
                 ! Update process noise covariance
-                Q(1,1) = 0.0001
+                Q(1,1) = 0.00001
                 Q(2,2) =(PI * (v_m**3.0) * (Ti**2.0)) / L
                 Q(3,3) = (2.0**2.0)/600.0
-\
+
                 ! Prediction update
                 ! Tau_r = 0.5 * CntrPar%WE_RhoAir * PI *CntrPar%WE_BladeRadius**3 * Cp_op * v_h**2 * 1.0/(lambda)
                 Tau_r = AeroDynTorque(LocalVar,CntrPar,PerfData)
@@ -215,6 +215,7 @@ SUBROUTINE WindSpeedEstimator(LocalVar, CntrPar, objInst, PerfData)
                 v_h = v_t + v_m
                 LocalVar%TestType = v_m + v_t
                 LocalVar%WE_Vw = v_m + v_t
+                ! LocalVar%WE_Vw = LPFilter(v_m + v_t,LocalVar%DT,0.5,LocalVar%iStatus,.FALSE.,objInst%instLPF)
             ENDIF
 
         ELSE        
@@ -274,7 +275,7 @@ REAL FUNCTION PeakShaving(LocalVar, CntrPar, objInst)
         ! V_towertop = PIController(LocalVar%FA_Acc, 0.0, 1.0, -100.00, 100.00, LocalVar%DT, 0.0, .FALSE., objInst%instPI)
 
         Vhat = LocalVar%WE_Vw
-        Vhatf = LPFilter(Vhat,LocalVar%DT,0.04,LocalVar%iStatus,.FALSE.,objInst%instLPF)
+        Vhatf = LPFilter(Vhat,LocalVar%DT,0.2,LocalVar%iStatus,.FALSE.,objInst%instLPF)
         LocalVar%TestType = Vhatf
         ! Define minimum blade pitch angle as a function of estimated wind speed
         PeakShaving = interp1d(CntrPar%PS_WindSpeeds, CntrPar%PS_BldPitchMin, Vhatf)

diff --git a/DRC_Fortran/Source/Controllers.f90 b/DRC_Fortran/Source/Controllers.f90
@@ -118,12 +118,12 @@ SUBROUTINE VariableSpeedControl(avrSWAP, CntrPar, LocalVar, objInst)
             IF (LocalVar%VS_State >= 4) THEN
                 VS_MaxTq = CntrPar%VS_RtTq
             ELSE
-                ! VS_MaxTq = CntrPar%VS_MaxTq
+                ! VS_MaxTq = CntrPar%VS_MaxTq           ! NJA: May want to boost max torque
                 VS_MaxTq = CntrPar%VS_RtTq
             ENDIF
             LocalVar%GenTq = PIController(LocalVar%VS_SpdErr, CntrPar%VS_KP(1), CntrPar%VS_KI(1), CntrPar%VS_MinTq, VS_MaxTq, LocalVar%DT, CntrPar%VS_MaxOMTq, .FALSE., objInst%instPI)
 
-            ! K*Omega^2 control law with PI torque control in transition regions
+        ! K*Omega^2 control law with PI torque control in transition regions
         ELSE
             ! Update PI loops for region 1.5 and 2.5 PI control
             ! LocalVar%GenArTq = PIController(LocalVar%VS_SpdErrAr, CntrPar%VS_KP(1), CntrPar%VS_KI(1), CntrPar%VS_MaxOMTq, CntrPar%VS_ArSatTq, LocalVar%DT, CntrPar%VS_RtTq, .TRUE., objInst%instPI)

diff --git a/DRC_Fortran/Source/DISCON.f90 b/DRC_Fortran/Source/DISCON.f90
@@ -39,6 +39,7 @@ SUBROUTINE DISCON(avrSWAP, aviFAIL, accINFILE, avcOUTNAME, avcMSG) BIND (C, NAME
 TYPE(LocalVariables), SAVE            :: LocalVar
 TYPE(ObjectInstances), SAVE           :: objInst
 TYPE(PerformanceData), SAVE           :: PerfData
+
 !------------------------------------------------------------------------------------------------------------------------------
 ! Main control calculations
 !------------------------------------------------------------------------------------------------------------------------------
@@ -48,7 +49,7 @@ SUBROUTINE DISCON(avrSWAP, aviFAIL, accINFILE, avcOUTNAME, avcMSG) BIND (C, NAME
 CALL PreFilterMeasuredSignals(CntrPar, LocalVar, objInst)
 
 IF ((LocalVar%iStatus >= 0) .AND. (aviFAIL >= 0))  THEN  ! Only compute control calculations if no error has occurred and we are not on the last time step
-    CALL ComputeVariablesSetpoints(CntrPar, LocalVar)
+    CALL ComputeVariablesSetpoints(CntrPar, LocalVar, objInst)
 
     CALL StateMachine(CntrPar, LocalVar)
     CALL WindSpeedEstimator(LocalVar, CntrPar, objInst, PerfData)

diff --git a/DRC_Fortran/Source/DRC_Types.f90 b/DRC_Fortran/Source/DRC_Types.f90
@@ -52,6 +52,7 @@ MODULE DRC_Types
     INTEGER(4)                          :: VS_n                         ! Number of controller gains
     REAL(4), DIMENSION(:), ALLOCATABLE  :: VS_KP                        ! Proportional gain for generator PI torque controller, used in the transitional 2.5 region
     REAL(4), DIMENSION(:), ALLOCATABLE  :: VS_KI                        ! Integral gain for generator PI torque controller, used in the transitional 2.5 region
+    REAL(4)                             :: VS_TSRopt                    ! Power-maximizing region 2 tip-speed ratio [rad]
 
     INTEGER(4)                          :: SS_Mode                      ! Setpoint Smoother mode {0: no setpoint smoothing, 1: introduce setpoint smoothing}
     REAL(4)                             :: SS_VSGainBias                ! Variable speed torque controller gain bias, [(rad/s)/rad].

diff --git a/DRC_Fortran/Source/ReadSetParameters.f90 b/DRC_Fortran/Source/ReadSetParameters.f90
@@ -104,6 +104,7 @@ SUBROUTINE ReadControlParameterFileSub(CntrPar)
         READ(UnControllerParameters,*) CntrPar%VS_KP
         ALLOCATE(CntrPar%VS_KI(CntrPar%VS_n))
         READ(UnControllerParameters,*) CntrPar%VS_KI
+        READ(UnControllerParameters,*) CntrPar%VS_TSRopt
         READ(UnControllerParameters, *)
 
         !------- Setpoint Smoother --------------------------------
@@ -179,17 +180,20 @@ SUBROUTINE ReadControlParameterFileSub(CntrPar)
     END SUBROUTINE ReadControlParameterFileSub
     ! -----------------------------------------------------------------------------------
     ! Calculate setpoints for primary control actions    
-    SUBROUTINE ComputeVariablesSetpoints(CntrPar, LocalVar)
-        USE DRC_Types, ONLY : ControlParameters, LocalVariables
+    SUBROUTINE ComputeVariablesSetpoints(CntrPar, LocalVar, objInst)
+        USE DRC_Types, ONLY : ControlParameters, LocalVariables, ObjectInstances
 
         ! Allocate variables
         TYPE(ControlParameters), INTENT(INOUT)  :: CntrPar
         TYPE(LocalVariables), INTENT(INOUT)     :: LocalVar
+        TYPE(ObjectInstances), INTENT(INOUT)    :: objInst
+
         REAL(4)                                 :: VS_RefSpd        ! Referece speed for variable speed torque controller, [rad/s] 
         REAL(4)                                 :: PC_RefSpd        ! Referece speed for pitch controller, [rad/s] 
         REAL(4)                                 :: Omega_op         ! Optimal TSR-tracking generator speed, [rad/s]
+        REAL(4)                                 :: WE_Vw_f          ! Filtered Wind Speed Estimate
         ! temp
-        REAL(4)                                 :: VS_TSRop = 7.5
+        ! REAL(4)                                 :: VS_TSRop = 7.5
 
         ! ----- Calculate yaw misalignment error -----
         LocalVar%Y_MErr = LocalVar%Y_M + CntrPar%Y_MErrSet ! Yaw-alignment error
@@ -207,7 +211,8 @@ SUBROUTINE ComputeVariablesSetpoints(CntrPar, LocalVar)
         ! ----- Torque controller reference errors -----
         ! Define VS reference generator speed [rad/s]
         IF (CntrPar%VS_ControlMode == 2) THEN
-            VS_RefSpd = (VS_TSRop * LocalVar%WE_Vw / CntrPar%WE_BladeRadius) * CntrPar%WE_GearboxRatio
+            WE_Vw_f = LPFilter(LocalVar%We_Vw, LocalVar%DT, 0.1, LocalVar%iStatus, .FALSE., objInst%instLPF)
+            VS_RefSpd = (CntrPar%VS_TSRopt * WE_Vw_f / CntrPar%WE_BladeRadius) * CntrPar%WE_GearboxRatio
             VS_RefSpd = saturate(VS_RefSpd,CntrPar%VS_MinOMSpd, CntrPar%PC_RefSpd)
         ELSE
             VS_RefSpd = CntrPar%VS_RefSpd
@@ -366,36 +371,38 @@ SUBROUTINE Assert(LocalVar, CntrPar, avrSWAP, aviFAIL, ErrMsg, size_avcMSG)
             ErrMsg  = 'IPC_KI(2) must be zero or greater than zero.'
         ENDIF
 
-        IF (CntrPar%Y_IPC_omegaLP <= 0.0)  THEN
-            aviFAIL = -1
-            ErrMsg  = 'Y_IPC_omegaLP must be greater than zero.'
-        ENDIF
-
-        IF (CntrPar%Y_IPC_zetaLP <= 0.0)  THEN
-            aviFAIL = -1
-            ErrMsg  = 'Y_IPC_zetaLP must be greater than zero.'
-        ENDIF
-
-        IF (CntrPar%Y_ErrThresh <= 0.0)  THEN
-            aviFAIL = -1
-            ErrMsg  = 'Y_ErrThresh must be greater than zero.'
-        ENDIF
-
-        IF (CntrPar%Y_Rate <= 0.0)  THEN
-            aviFAIL = -1
-            ErrMsg  = 'CntrPar%Y_Rate must be greater than zero.'
-        ENDIF
-
-        IF (CntrPar%Y_omegaLPFast <= 0.0)  THEN
-            aviFAIL = -1
-            ErrMsg  = 'Y_omegaLPFast must be greater than zero.'
-        ENDIF
-
-        IF (CntrPar%Y_omegaLPSlow <= 0.0)  THEN
-            aviFAIL = -1
-            ErrMsg  = 'Y_omegaLPSlow must be greater than zero.'
+        ! ---- Yaw Control ----
+        IF (CntrPar%Y_ControlMode > 0) THEN
+            IF (CntrPar%Y_IPC_omegaLP <= 0.0)  THEN
+                aviFAIL = -1
+                ErrMsg  = 'Y_IPC_omegaLP must be greater than zero.'
+            ENDIF
+
+            IF (CntrPar%Y_IPC_zetaLP <= 0.0)  THEN
+                aviFAIL = -1
+                ErrMsg  = 'Y_IPC_zetaLP must be greater than zero.'
+            ENDIF
+
+            IF (CntrPar%Y_ErrThresh <= 0.0)  THEN
+                aviFAIL = -1
+                ErrMsg  = 'Y_ErrThresh must be greater than zero.'
+            ENDIF
+
+            IF (CntrPar%Y_Rate <= 0.0)  THEN
+                aviFAIL = -1
+                ErrMsg  = 'CntrPar%Y_Rate must be greater than zero.'
+            ENDIF
+
+            IF (CntrPar%Y_omegaLPFast <= 0.0)  THEN
+                aviFAIL = -1
+                ErrMsg  = 'Y_omegaLPFast must be greater than zero.'
+            ENDIF
+
+            IF (CntrPar%Y_omegaLPSlow <= 0.0)  THEN
+                aviFAIL = -1
+                ErrMsg  = 'Y_omegaLPSlow must be greater than zero.'
+            ENDIF
         ENDIF
-
         ! Abort if the user has not requested a pitch angle actuator (See Appendix A
         ! of Bladed User's Guide):
         IF (NINT(avrSWAP(10)) /= 0)  THEN ! .TRUE. if a pitch angle actuator hasn't been requested

diff --git a/Matlab_Toolbox/Post_TestCases.m b/Matlab_Toolbox/Post_TestCases.m
@@ -44,8 +44,13 @@
  % Will want to (un)comment desired cases to plot
 
  % Usable Plot types: 
+<<<<<<< HEAD
  %  - Step, Below Rated, Near Rated, Above Rated, Floating Steady, Floating Near Rated, Floating Mexican Hat
 plottype = 'Floating Above Rated';
+=======
+ %  - Step, Below Rated, Near Rated, Above Rated, Floating Steady, Floating Near Rated, Floating Above Rated, Floating Mexican Hat
+plottype = 'Near Rated';
+>>>>>>> dev_NewTurbines
 
 close all
 % Below Rated
@@ -62,7 +67,11 @@
         Pl_FastPlots(fo.OC4_ARsteady_Legacy, fo.OC4_ARsteady_Baseline)
     case 'Floating Near Rated'
         Pl_FastPlots(fo.OC4_NR_Legacy, fo.OC4_NR_Baseline)
+<<<<<<< HEAD
    case 'Floating Above Rated'
+=======
+    case 'Floating Above Rated'
+>>>>>>> dev_NewTurbines
         Pl_FastPlots(fo.OC4_AR_Legacy, fo.OC4_AR_Baseline)
     case 'Floating Mexican Hat'
         Pl_FastPlots(fo.OC4_MH_Legacy, fo.OC4_MH_Baseline)