Expand functionalities of UncertaintyInterface and fix no_wake compat…

…ibility (#428)

floris/tools/uncertainty_interface.py

@@ -118,6 +118,9 @@ def __init__(
             fix_yaw_in_relative_frame=fix_yaw_in_relative_frame,
         )
 
+        # Add a _no_wake switch to keep track of calculate_wake/calculate_no_wake
+        self._no_wake = False
+
     # Private methods
 
     def _generate_pdfs_from_dict(self):
@@ -364,8 +367,24 @@ def calculate_wake(self, yaw_angles=None):
             yaw_angles: NDArrayFloat | list[float] | None = None,
         """
         self._reassign_yaw_angles(yaw_angles)
+        self._no_wake = False
+
+    def calculate_no_wake(self, yaw_angles=None):
+        """Replaces the 'calculate_no_wake' function in the FlorisInterface
+        object. Fundamentally, this function only overwrites the nominal
+        yaw angles in the FlorisInterface object. The actual wake calculations
+        are performed once 'get_turbine_powers' or 'get_farm_powers' is
+        called. However, to allow users to directly replace a FlorisInterface
+        object with this UncertaintyInterface object, this function is
+        required.
 
-    def get_turbine_powers(self, no_wake=False):
+        Args:
+            yaw_angles: NDArrayFloat | list[float] | None = None,
+        """
+        self._reassign_yaw_angles(yaw_angles)
+        self._no_wake = True
+
+    def get_turbine_powers(self):
         """Calculates the probability-weighted power production of each
         turbine in the wind farm.
 
@@ -419,7 +438,7 @@ def get_turbine_powers(self, no_wake=False):
 
         # Evaluate floris for minimal probablistic set
         self.fi.reinitialize(wind_directions=wd_array_probablistic_min)
-        if no_wake:
+        if self._no_wake:
             self.fi.calculate_no_wake(yaw_angles=yaw_angles_probablistic_min)
         else:
             self.fi.calculate_wake(yaw_angles=yaw_angles_probablistic_min)
@@ -443,7 +462,7 @@ def get_turbine_powers(self, no_wake=False):
         # Now apply probability distribution weighing to get turbine powers
         return np.sum(wd_weighing * power_probablistic, axis=0)
 
-    def get_farm_power(self, no_wake=False):
+    def get_farm_power(self):
         """Calculates the probability-weighted power production of the
         collective of all turbines in the farm, for each wind direction
         and wind speed specified.
@@ -458,9 +477,119 @@ def get_farm_power(self, no_wake=False):
             NDArrayFloat: Expectation of power production of the wind farm.
             This array has the shape (num_wind_directions, num_wind_speeds).
         """
-        turbine_powers = self.get_turbine_powers(no_wake=no_wake)
+        turbine_powers = self.get_turbine_powers()
         return np.sum(turbine_powers, axis=2)
 
+    def get_farm_AEP(
+        self,
+        freq,
+        cut_in_wind_speed=0.001,
+        cut_out_wind_speed=None,
+        yaw_angles=None,
+        no_wake=False,
+    ) -> float:
+        """
+        Estimate annual energy production (AEP) for distributions of wind speed, wind
+        direction, frequency of occurrence, and yaw offset.
+
+        Args:
+            freq (NDArrayFloat): NumPy array with shape (n_wind_directions,
+                n_wind_speeds) with the frequencies of each wind direction and
+                wind speed combination. These frequencies should typically sum
+                up to 1.0 and are used to weigh the wind farm power for every
+                condition in calculating the wind farm's AEP.
+            cut_in_wind_speed (float, optional): Wind speed in m/s below which
+                any calculations are ignored and the wind farm is known to 
+                produce 0.0 W of power. Note that to prevent problems with the
+                wake models at negative / zero wind speeds, this variable must
+                always have a positive value. Defaults to 0.001 [m/s].
+            cut_out_wind_speed (float, optional): Wind speed above which the
+                wind farm is known to produce 0.0 W of power. If None is
+                specified, will assume that the wind farm does not cut out
+                at high wind speeds. Defaults to None.
+            yaw_angles (NDArrayFloat | list[float] | None, optional):
+                The relative turbine yaw angles in degrees. If None is
+                specified, will assume that the turbine yaw angles are all
+                zero degrees for all conditions. Defaults to None.
+            no_wake: (bool, optional): When *True* updates the turbine
+                quantities without calculating the wake or adding the wake to
+                the flow field. This can be useful when quantifying the loss
+                in AEP due to wakes. Defaults to *False*.
+
+        Returns:
+            float: 
+                The Annual Energy Production (AEP) for the wind farm in
+                watt-hours.
+        """
+
+        # Verify dimensions of the variable "freq"
+        if not (
+            (np.shape(freq)[0] == self.floris.flow_field.n_wind_directions)
+            & (np.shape(freq)[1] == self.floris.flow_field.n_wind_speeds)
+            & (len(np.shape(freq)) == 2)
+        ):
+            raise UserWarning(
+                "'freq' should be a two-dimensional array with dimensions"
+                + " (n_wind_directions, n_wind_speeds)."
+            )
+
+        # Check if frequency vector sums to 1.0. If not, raise a warning
+        if np.abs(np.sum(freq) - 1.0) > 0.001:
+            self.logger.warning(
+                "WARNING: The frequency array provided to get_farm_AEP() "
+                + "does not sum to 1.0. "
+            )
+
+        # Copy the full wind speed array from the floris object and initialize
+        # the the farm_power variable as an empty array.
+        wind_speeds = np.array(self.fi.floris.flow_field.wind_speeds, copy=True)
+        farm_power = np.zeros(
+            (self.fi.floris.flow_field.n_wind_directions, len(wind_speeds))
+        )
+
+        # Determine which wind speeds we must evaluate in floris
+        conditions_to_evaluate = (wind_speeds >= cut_in_wind_speed)
+        if cut_out_wind_speed is not None:
+            conditions_to_evaluate = conditions_to_evaluate & (
+                wind_speeds < cut_out_wind_speed
+            )
+
+        # Evaluate the conditions in floris
+        if np.any(conditions_to_evaluate):
+            wind_speeds_subset = wind_speeds[conditions_to_evaluate]
+            yaw_angles_subset = None
+            if yaw_angles is not None:
+                yaw_angles_subset = yaw_angles[:, conditions_to_evaluate]
+            self.reinitialize(wind_speeds=wind_speeds_subset)
+            if no_wake:
+                self.calculate_no_wake(yaw_angles=yaw_angles_subset)
+            else:
+                self.calculate_wake(yaw_angles=yaw_angles_subset)
+            farm_power[:, conditions_to_evaluate] = self.get_farm_power()
+
+        # Finally, calculate AEP in GWh
+        aep = np.sum(np.multiply(freq, farm_power) * 365 * 24)
+
+        # Reset the FLORIS object to the full wind speed array
+        self.reinitialize(wind_speeds=wind_speeds)
+
+        return aep
+
+    def assign_hub_height_to_ref_height(self):
+        return self.fi.assign_hub_height_to_ref_height()
+
+    def get_turbine_layout(self, z=False):
+        return self.fi.get_turbine_layout(z=z)
+
+    def get_turbine_Cts(self):
+        return self.fi.get_turbine_Cts()
+
+    def get_turbine_ais(self):
+        return self.fi.get_turbine_ais()
+
+    def get_turbine_average_velocities(self):
+        return self.fi.get_turbine_average_velocities()
+
     # Define getter functions that just pass information from FlorisInterface
     @property
     def floris(self):