Merge branch 'oneROSCO2' into mult_omega

.github/workflows/CI_rosco-toolbox.yml

@@ -86,6 +86,11 @@ jobs:
         run: |
           conda install -y compilers
           conda install -y wisdem
+      - name: Add dependencies windows specific
+        if: true == contains( matrix.os, 'windows')
+        run: |
+          conda install -y m2w64-toolchain libpython 
+          conda install -y wisdem
         
 
       # Install ROSCO toolbox

Examples/example_05.py

@@ -72,10 +72,10 @@
 
 
 # Load controller library
-controller_int = ROSCO_ci.ControllerInterface(lib_name,param_filename=param_filename)
+controller_int = ROSCO_ci.ControllerInterface(lib_name,param_filename=param_filename,sim_name='sim1')
 
 # Load the simulator
-sim = ROSCO_sim.Sim(turbine,controller_int)
+sim_1 = ROSCO_sim.Sim(turbine,controller_int)
 
 # Define a wind speed history
 dt = 0.025
@@ -88,7 +88,17 @@
     ws[i] = ws[i] + t[i]//100
 
 # Run simulator and plot results
-sim.sim_ws_series(t,ws,rotor_rpm_init=4)
+sim_1.sim_ws_series(t,ws,rotor_rpm_init=4)
+
+# Load controller library again to see if we deallocated properly
+controller_int = ROSCO_ci.ControllerInterface(lib_name,param_filename=param_filename,sim_name='sim_2')
+
+# Run simulator again and plot results
+sim_2 = ROSCO_sim.Sim(turbine,controller_int)
+sim_2.sim_ws_series(t,ws,rotor_rpm_init=4)
+
+# Check if simulations are equal
+np.testing.assert_almost_equal(sim_1.gen_speed,sim_2.gen_speed)
 
 if False:
   plt.show()

Examples/example_07.py

@@ -34,7 +34,7 @@
 controller_params   = inps['controller_params']
 
 # Ensure minimum generator speed at 50 rpm (for example's sake), turn on peak shaving and cp-maximizing min pitch
-controller_params['vs_minspd'] = 50
+controller_params['vs_minspd'] = 50/97
 controller_params['PS_Mode'] = 3
 
 # Instantiate turbine, controller, and file processing classes

ROSCO/src/ControllerBlocks.f90

@@ -243,11 +243,11 @@ SUBROUTINE WindSpeedEstimator(LocalVar, CntrPar, objInst, PerfData, DebugVar, Er
             Q = RESHAPE((/0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0/),(/3,3/))
             IF (LocalVar%iStatus == 0) THEN
                 ! Initialize recurring values
-                om_r = max(LocalVar%RotSpeedF, CntrPar%VS_MinOMSpd)
+                om_r = max(LocalVar%RotSpeedF, EPSILON(1.0))
                 v_t = 0.0
                 v_m = LocalVar%HorWindV
                 v_h = LocalVar%HorWindV
-                lambda = max(LocalVar%RotSpeed, CntrPar%VS_MinOMSpd) * CntrPar%WE_BladeRadius/v_h
+                lambda = max(LocalVar%RotSpeed, EPSILON(1.0)) * CntrPar%WE_BladeRadius/v_h
                 xh = RESHAPE((/om_r, v_t, v_m/),(/3,1/))
                 P = RESHAPE((/0.01, 0.0, 0.0, 0.0, 0.01, 0.0, 0.0, 0.0, 1.0/),(/3,3/))
                 K = RESHAPE((/0.0,0.0,0.0/),(/3,1/))
@@ -258,7 +258,7 @@ SUBROUTINE WindSpeedEstimator(LocalVar, CntrPar, objInst, PerfData, DebugVar, Er
                 A_op = interp1d(CntrPar%WE_FOPoles_v,CntrPar%WE_FOPoles,v_h,ErrVar)
 
                 ! TEST INTERP2D
-                lambda = max(LocalVar%RotSpeed, CntrPar%VS_MinOMSpd) * CntrPar%WE_BladeRadius/v_h
+                lambda = max(LocalVar%RotSpeed, EPSILON(1.0)) * CntrPar%WE_BladeRadius/v_h
                 Cp_op = interp2d(PerfData%Beta_vec,PerfData%TSR_vec,PerfData%Cp_mat, LocalVar%BlPitch(1)*R2D, lambda , ErrVar)
                 Cp_op = max(0.0,Cp_op)
 
@@ -292,7 +292,7 @@ SUBROUTINE WindSpeedEstimator(LocalVar, CntrPar, objInst, PerfData, DebugVar, Er
 
 
                 ! Wind Speed Estimate
-                om_r = xh(1,1)
+                om_r = max(xh(1,1), EPSILON(1.0))
                 v_t = xh(2,1)
                 v_m = xh(3,1)
                 v_h = v_t + v_m

ROSCO/src/DISCON.F90

@@ -64,7 +64,11 @@ SUBROUTINE DISCON(avrSWAP, aviFAIL, accINFILE, avcOUTNAME, avcMSG) BIND (C, NAME
 !------------------------------------------------------------------------------------------------------------------------------
 ! Read avrSWAP array into derived types/variables
 CALL ReadAvrSWAP(avrSWAP, LocalVar)
+
+! Set Control Parameters
 CALL SetParameters(avrSWAP, accINFILE, SIZE(avcMSG), CntrPar, LocalVar, objInst, PerfData, ErrVar)
+
+! Filter signals
 CALL PreFilterMeasuredSignals(CntrPar, LocalVar, objInst)
 
 IF ((LocalVar%iStatus >= 0) .AND. (ErrVar%aviFAIL >= 0))  THEN  ! Only compute control calculations if no error has occurred and we are not on the last time step

ROSCO_toolbox/control_interface.py

@@ -12,6 +12,7 @@
 from ctypes import byref, cdll, c_int, POINTER, c_float, c_char_p, c_double, create_string_buffer, c_int32, c_void_p
 import numpy as np
 from numpy.ctypeslib import ndpointer
+import platform, ctypes
 
 # Some useful constants
 deg2rad = np.deg2rad(1)
@@ -75,6 +76,14 @@ def init_discon(self):
         self.avrSWAP[20] = 1.0 # HARD CODE initial rot speed = 1 rad/s
         self.avrSWAP[82] = 0 # HARD CODE initial nacIMU = 0
         self.avrSWAP[26] = 10 # HARD CODE initial wind speed = 10 m/s
+
+        # Blade pitch initial conditions
+        self.avrSWAP[3]     = 0 * np.deg2rad(1)
+        self.avrSWAP[32]    = 0 * np.deg2rad(1)
+        self.avrSWAP[33]    = 0 * np.deg2rad(1)
+
+        # Torque initial condition
+        self.avrSWAP[22]    = 0
 
 
         # Code this as first casll
@@ -173,14 +182,75 @@ def show_control_values(self):
 
     def kill_discon(self):
         '''
-        Unload the dylib from memory: https://stackoverflow.com/questions/359498/how-can-i-unload-a-dll-using-ctypes-in-python
-        This is unix-specific, but there seems to be a windows solution as well
+        Unload the dylib from memory: https://github.com/bwoodsend/cslug/blob/master/cslug/_stdlib.py
         '''
 
         print('Shutting down {}'.format(self.lib_name))
         handle = self.discon._handle
-        dlclose_func = self.discon.dlclose
-        dlclose_func.argtypes = [c_void_p]
+
+        # Start copy here
+        OS = platform.system()
+
+        def null_free_dll(*spam):  # pragma: no cover
+            pass
+
+        # Try to find a good runtime library which is always available and contains
+        # the standard library C functions such as malloc() or printf().
+        # XXX: Keep chosen library names in sync with the table in `cslug/stdlib.py`.
+
+        extra_libs = []
+
+        if OS == "Windows":  # pragma: Windows
+            _dlclose = ctypes.windll.kernel32.FreeLibrary
+            dlclose = lambda handle: 0 if _dlclose(handle) else 1
+            # There's some controversy as to whether this DLL is guaranteed to exist.
+            # It always has so far but isn't documented. However, MinGW assumes that it
+            # is so, should this DLL be removed, then we have much bigger problems than
+            # just this line. There is also vcruntime140.dll which isn't a standard part
+            # of the OS but is always shipped with Python so we can guarantee its
+            # presence. But vcruntime140 contains only a tiny strict-subset of msvcrt.
+            stdlib = ctypes.CDLL("msvcrt")
+
+        elif OS == "Darwin":  # pragma: Darwin
+            try:
+                try:
+                    # macOS 11 (Big Sur). Possibly also later macOS 10s.
+                    stdlib = ctypes.CDLL("libc.dylib")
+                except OSError:  # pragma: no cover
+                    stdlib = ctypes.CDLL("libSystem")
+            except OSError:  # pragma: no cover
+                # Older macOSs. Not only is the name inconsistent but it's
+                # not even in PATH.
+                _stdlib = "/usr/lib/system/libsystem_c.dylib"
+                if os.path.exists(_stdlib):
+                    stdlib = ctypes.CDLL(_stdlib)
+                else:
+                    stdlib = None
+            if stdlib is not None:  # pragma: no branch
+                dlclose = stdlib.dlclose
+            else:  # pragma: no cover
+                # I hope this never happens.
+                dlclose = null_free_dll
+
+        elif OS == "Linux":  # pragma: Linux
+            try:
+                stdlib = ctypes.CDLL("")
+            except OSError:  # pragma: no cover
+                # Either Alpine Linux or Android.
+                # Unfortunately, there doesn't seem to be any practical way
+                # to tell them apart.
+                stdlib = ctypes.CDLL("libc.so")
+
+                # Android, like FreeBSD puts its math functions
+                # in a dedicated `libm.so`.
+                # The only way to know that this is not Alpine is to check if the math
+                # functions are already available in `libc.so`.
+                if not hasattr(stdlib, "sin"):
+                    extra_libs.append(ctypes.CDLL("libm.so"))
+            dlclose = stdlib.dlclose
+
+        # End copy here
+        dlclose.argtypes = [c_void_p]
+        dlclose(handle)
 
         del self.discon
-        dlclose_func(handle)

ROSCO_toolbox/controller.py

@@ -246,9 +246,9 @@ def tune_controller(self, turbine):
         # -- Define some setpoints --
         # minimum rotor speed saturation limits
         if self.vs_minspd:
-            self.vs_minspd = np.maximum(self.vs_minspd, (turbine.TSR_operational * turbine.v_min / turbine.rotor_radius) * Ng)
+            self.vs_minspd = np.maximum(self.vs_minspd, (turbine.TSR_operational * turbine.v_min / turbine.rotor_radius))
         else: 
-            self.vs_minspd = (turbine.TSR_operational * turbine.v_min / turbine.rotor_radius) * Ng
+            self.vs_minspd = (turbine.TSR_operational * turbine.v_min / turbine.rotor_radius)
         self.pc_minspd = self.vs_minspd
 
         # max pitch angle for shutdown
@@ -483,7 +483,7 @@ def peak_shaving(self,controller, turbine):
     def min_pitch_saturation(self, controller, turbine):
 
         # Find TSR associated with minimum rotor speed
-        TSR_at_minspeed = (controller.pc_minspd/turbine.Ng) * turbine.rotor_radius / controller.v_below_rated
+        TSR_at_minspeed = (controller.pc_minspd) * turbine.rotor_radius / controller.v_below_rated
         for i in range(len(TSR_at_minspeed)):
             if TSR_at_minspeed[i] > controller.TSR_op[i]:
                 controller.TSR_op[i] = TSR_at_minspeed[i]

ROSCO_toolbox/sim.py

@@ -116,6 +116,9 @@ def sim_ws_series(self,t_array,ws_array,rotor_rpm_init=10,init_pitch=0.0, make_p
         self.t_array = t_array
         self.ws_array = ws_array
 
+        # Close shared library when finished
+        self.controller_int.kill_discon()
+
         if make_plots:
             fig, axarr = plt.subplots(4,1,sharex=True,figsize=(6,10))
 

setup.py

@@ -35,7 +35,7 @@
 from io import open
 
 # Package meta-data.
-NAME = 'rosco-toolbox'
+NAME = 'rosco'
 DESCRIPTION = 'A toolbox for development of wind turbine controllers.'
 URL = 'https://github.com/NREL/ROSCO_toolbox'
 EMAIL = 'nikhar.abbas@nrel.gov'