Expanded capabilities for heterogeneity (#902)

docs/_toc.yml

@@ -14,6 +14,7 @@ parts:
     - file: intro_concepts
     - file: advanced_concepts
     - file: wind_data_user
+    - file: heterogeneous_map
     - file: floating_wind_turbine
     - file: turbine_interaction
     - file: operation_models_user

examples/examples_heterogeneous/001_heterogeneous_inflow_single.py

@@ -17,7 +17,7 @@
 import numpy as np
 
 from floris import FlorisModel, TimeSeries
-from floris.flow_visualization import visualize_cut_plane
+from floris.flow_visualization import visualize_heterogeneous_cut_plane
 from floris.layout_visualization import plot_turbine_labels
 
 
@@ -65,15 +65,20 @@
     x_resolution=200, y_resolution=100, height=90.0
 )
 
-# Plot the horizontal plane
+# Plot the horizontal plane using the visualize_heterogeneous_cut_plane.
+# Note that this function is not very different than the standard
+# visualize_cut_plane except that it accepts the fmodel object in order to
+# visualize the boundary of the heterogeneous inflow region.
 fig, ax = plt.subplots()
-visualize_cut_plane(
+visualize_heterogeneous_cut_plane(
     horizontal_plane,
+    fmodel=fmodel,
     ax=ax,
     title="Horizontal plane at hub height",
     color_bar=True,
     label_contours=True,
 )
 plot_turbine_labels(fmodel, ax)
+ax.legend()
 
 plt.show()

examples/examples_heterogeneous/002_heterogeneous_inflow_multi.py → examples/examples_heterogeneous/002_heterogeneous_using_wind_data.py

@@ -1,16 +1,14 @@
-"""Example: Heterogeneous Inflow for multiple conditions
+"""Example: Heterogeneous Inflow using wind data
 
 When multiple cases are considered, the heterogeneous inflow conditions can be defined in two ways:
 
   1. Passing heterogeneous_inflow_config to the set method, with P points,
-        and speedups of size n_findex X P
-  2. Assigning heterogeneous_inflow_config_by_wd to the wind_data object
-        used to drive FLORIS.  This object includes
-        n_wd wind_directions, and speedups is of size n_wd X P.  When applied
-          to set, the heterogeneous_inflow_config
-        is automatically generated by using the nearest wind direction
-        defined in heterogeneous_inflow_config_by_wd
-        for each findex.
+        and speed_multipliers of size n_findex X P
+  2. More conveniently, building a HeterogeneousMap object that defines the speed_multipliers as a
+        function of wind direction and/or wind speed and passing that to a WindData object.  When
+        the WindData object is passed to the set method, the heterogeneous_inflow_config is
+        automatically generated for each findex by finding the nearest wind direction and/or wind
+        speed in the HeterogeneousMap object.
 
 This example:
 
@@ -23,7 +21,11 @@
 import matplotlib.pyplot as plt
 import numpy as np
 
-from floris import FlorisModel, TimeSeries
+from floris import (
+    FlorisModel,
+    HeterogeneousMap,
+    TimeSeries,
+)
 
 
 # Initialize FlorisModel
@@ -34,7 +36,6 @@
 
 # Define a TimeSeries object with 4 wind directions and constant wind speed
 # and turbulence intensity
-
 time_series = TimeSeries(
     wind_directions=np.array([269.0, 270.0, 271.0, 282.0]),
     wind_speeds=8.0,
@@ -51,17 +52,17 @@
 
 # Assume the speed-ups are defined such that they are the same 265-275 degrees and 275-285 degrees
 
-# If defining heterogeneous_inflow_config directly, then the speedups are of size n_findex X P
-# where the first 3 rows are identical, and the last row is different
-speed_ups = [
+# If defining heterogeneous_inflow_config directly, then the speed_multipliers are of size
+# n_findex x P, where the first 3 rows are identical and the last row is different
+speed_multipliers = [
     [1.0, 1.25, 1.0, 1.25],
     [1.0, 1.25, 1.0, 1.25],
     [1.0, 1.25, 1.0, 1.25],
     [1.0, 1.35, 1.0, 1.35],
 ]
 
 heterogeneous_inflow_config = {
-    "speed_multipliers": speed_ups,
+    "speed_multipliers": speed_multipliers,
     "x": x_locs,
     "y": y_locs,
 }
@@ -75,26 +76,45 @@
 # Get the power output of the turbines
 turbine_powers = fmodel.get_turbine_powers() / 1000.0
 
-# Now repeat using the wind_data object and heterogeneous_inflow_config_by_wd
-# First, create the speedups for the two wind directions
-speed_ups = [[1.0, 1.25, 1.0, 1.25], [1.0, 1.35, 1.0, 1.35]]
+# Now repeat using the wind_data object and HeterogeneousMap object
+# First, create the speed multipliers for the two wind directions
+speed_multipliers = [[1.0, 1.25, 1.0, 1.25], [1.0, 1.35, 1.0, 1.35]]
+
+# Now define the HeterogeneousMap object
+heterogeneous_map = HeterogeneousMap(
+    x=x_locs,
+    y=y_locs,
+    speed_multipliers=speed_multipliers,
+    wind_directions=[270.0, 280.0],
+)
 
-# Create the heterogeneous_inflow_config_by_wd dictionary
+# Now create a new TimeSeries object including the heterogeneous_inflow_config_by_wd
+time_series = TimeSeries(
+    wind_directions=np.array([269.0, 270.0, 271.0, 282.0]),
+    wind_speeds=8.0,
+    turbulence_intensities=0.06,
+    heterogeneous_map=heterogeneous_map,
+)
+
+# Note that previously, the a heterogeneous_inflow_config_by_wd, which only only
+# specification by wind direction was defined, and for backwards compatibility,
+# this is still accepted.  However, the HeterogeneousMap object is more flexible.
+# The following code produces the same results as the previous code block.
 heterogeneous_inflow_config_by_wd = {
-    "speed_multipliers": speed_ups,
+    "speed_multipliers": speed_multipliers,
     "x": x_locs,
     "y": y_locs,
     "wind_directions": [270.0, 280.0],
 }
 
-# Now create a new TimeSeries object including the heterogeneous_inflow_config_by_wd
 time_series = TimeSeries(
     wind_directions=np.array([269.0, 270.0, 271.0, 282.0]),
     wind_speeds=8.0,
     turbulence_intensities=0.06,
     heterogeneous_inflow_config_by_wd=heterogeneous_inflow_config_by_wd,
 )
 
+
 # Apply the time series to the FlorisModel
 fmodel.set(wind_data=time_series)
 

examples/examples_heterogeneous/003_heterogeneous_speedup_by_wd_and_ws.py

@@ -0,0 +1,113 @@
+"""Example: Heterogeneous Speedup by Wind Direction and Wind Speed
+
+The HeterogeneousMap object is a flexible way to define speedups as a function of wind direction
+and/or wind speed.  It also contains methods to plot the speedup map for a given wind direction
+and wind speed.
+
+This example:
+
+    1) Instantiates a HeterogeneousMap object with speedups defined for two wind directions
+        and two wind speeds
+    2) Visualizes the speedups for two particular combinations of wind direction and wind speed
+    3) Runs a FLORIS simulation using the HeterogeneousMap and visualizes the results
+
+"""
+
+
+import matplotlib.pyplot as plt
+import numpy as np
+
+from floris import (
+    FlorisModel,
+    HeterogeneousMap,
+    TimeSeries,
+)
+from floris.flow_visualization import visualize_heterogeneous_cut_plane
+
+
+# Define a HeterogeneousMap object with speedups defined for two wind directions
+# and two wind speeds.  The speedups imply no heterogeneity for the first wind direction
+# (0 degrees) with heterogeneity for the second wind direction (180 degrees) with the
+# specific speedups for this direction depending on the wind speed.
+heterogeneous_map = HeterogeneousMap(
+    x=np.array([0.0, 0.0, 250.0, 500.0, 500.0]),
+    y=np.array([0.0, 500.0, 250.0, 0.0, 500.0]),
+    speed_multipliers=np.array(
+        [
+            [1.0, 1.0, 1.0, 1.0, 1.0],
+            [1.0, 1.0, 1.0, 1.0, 1.0],
+            [1.5, 1.0, 1.25, 1.5, 1.0],
+            [1.0, 1.5, 1.25, 1.0, 1.5],
+        ]
+    ),
+    wind_directions=np.array([270.0, 270.0, 90.0, 90.0]),
+    wind_speeds=np.array([5.0, 10.0, 5.0, 10.0]),
+)
+
+# Use the HeterogeneousMap object to plot the speedup map for 3 wd/ws combinations
+fig, axarr = plt.subplots(1, 3, sharex=True, sharey=True, figsize=(15, 5))
+
+
+ax = axarr[0]
+heterogeneous_map.plot_single_speed_multiplier(
+    wind_direction=60.0, wind_speed=8.5, ax=ax, vmin=1.0, vmax=1.2
+)
+ax.set_title("Wind Direction = 60.0\nWind Speed = 8.5")
+
+ax = axarr[1]
+heterogeneous_map.plot_single_speed_multiplier(
+    wind_direction=130.0, wind_speed=4.0, ax=ax, vmin=1.0, vmax=1.2
+)
+ax.set_title("Wind Direction = 130.0\nWind Speed = 4.0")
+
+ax = axarr[2]
+heterogeneous_map.plot_single_speed_multiplier(
+    wind_direction=280.0, wind_speed=16.0, ax=ax, vmin=1.0, vmax=1.2
+)
+ax.set_title("Wind Direction = 280.0\nWind Speed = 16.0")
+fig.suptitle("Heterogeneous speedup map for several directions and wind speeds")
+
+
+# Initialize FlorisModel
+fmodel = FlorisModel("../inputs/gch.yaml")
+
+# Change the layout to a 2 turbine layout within the heterogeneous domain
+fmodel.set(layout_x=[200, 200.0], layout_y=[50, 450.0])
+
+# Define a TimeSeries object with 3 wind directions and wind speeds
+# and turbulence intensity and using the above HeterogeneousMap object
+time_series = TimeSeries(
+    wind_directions=np.array([275.0, 95.0, 75.0]),
+    wind_speeds=np.array([7.0, 6.2, 8.0]),
+    turbulence_intensities=0.06,
+    heterogeneous_map=heterogeneous_map,
+)
+
+# Apply the time series to the FlorisModel
+fmodel.set(wind_data=time_series)
+
+# Run the FLORIS simulation
+fmodel.run()
+
+# Visualize each of the findices
+fig, axarr = plt.subplots(3, 1, sharex=True, sharey=True, figsize=(10, 10))
+
+for findex in range(3):
+    ax = axarr[findex]
+
+    horizontal_plane = fmodel.calculate_horizontal_plane(
+        x_resolution=200, y_resolution=100, height=90.0, findex_for_viz=findex
+    )
+
+    visualize_heterogeneous_cut_plane(
+        cut_plane=horizontal_plane,
+        fmodel=fmodel,
+        ax=ax,
+        title=(
+            f"Wind Direction = {time_series.wind_directions[findex]}\n"
+            f"Wind Speed = {time_series.wind_speeds[findex]}"
+        ),
+    )
+
+
+plt.show()

floris/__init__.py

@@ -12,6 +12,7 @@
     visualize_cut_plane,
     visualize_quiver,
 )
+from .heterogeneous_map import HeterogeneousMap
 from .parallel_floris_model import ParallelFlorisModel
 from .uncertain_floris_model import ApproxFlorisModel, UncertainFlorisModel
 from .wind_data import (

floris/core/flow_field.py

@@ -291,16 +291,57 @@ def generate_heterogeneous_wind_map(self):
             # Linear interpolation is used for points within the user-defined area of values,
             # while the freestream wind speed is used for points outside that region
             in_region = [
-                LinearNDInterpolator(list(zip(x, y, z)), multiplier, fill_value=1.0)
+                self.interpolate_multiplier_xyz(x, y, z, multiplier, fill_value=1.0)
                 for multiplier in speed_multipliers
             ]
         else:
             # Compute the 2-dimensional interpolants for each wind direction
             # Linear interpolation is used for points within the user-defined area of values,
             # while the freestream wind speed is used for points outside that region
             in_region = [
-                LinearNDInterpolator(list(zip(x, y)), multiplier, fill_value=1.0)
+                self.interpolate_multiplier_xy(x, y, multiplier, fill_value=1.0)
                 for multiplier in speed_multipliers
             ]
 
         self.het_map = in_region
+
+    @staticmethod
+    def interpolate_multiplier_xy(x: NDArrayFloat,
+                                  y: NDArrayFloat,
+                                  multiplier: NDArrayFloat,
+                                  fill_value: float = 1.0):
+        """Return an interpolant for a 2D multiplier field.
+
+        Args:
+            x (NDArrayFloat): x locations
+            y (NDArrayFloat): y locations
+            multiplier (NDArrayFloat): multipliers
+            fill_value (float): fill value for points outside the region
+
+        Returns:
+            LinearNDInterpolator: interpolant
+        """
+
+        return LinearNDInterpolator(list(zip(x, y)), multiplier, fill_value=fill_value)
+
+
+    @staticmethod
+    def interpolate_multiplier_xyz(x: NDArrayFloat,
+                                   y: NDArrayFloat,
+                                   z: NDArrayFloat,
+                                   multiplier: NDArrayFloat,
+                                   fill_value: float = 1.0):
+        """Return an interpolant for a 3D multiplier field.
+
+        Args:
+            x (NDArrayFloat): x locations
+            y (NDArrayFloat): y locations
+            z (NDArrayFloat): z locations
+            multiplier (NDArrayFloat): multipliers
+            fill_value (float): fill value for points outside the region
+
+        Returns:
+            LinearNDInterpolator: interpolant
+        """
+
+        return LinearNDInterpolator(list(zip(x, y, z)), multiplier, fill_value=fill_value)

floris/floris_model.py

@@ -954,6 +954,22 @@ def set_for_viz(self, findex: int, solver_settings: dict) -> None:
             findex (int): The findex to set the floris object to.
             solver_settings (dict): The solver settings to use for visualization.
         """
+
+        # If not None, set the heterogeneous inflow configuration
+        if self.core.flow_field.heterogeneous_inflow_config is not None:
+            heterogeneous_inflow_config = {
+                'x': self.core.flow_field.heterogeneous_inflow_config['x'],
+                'y': self.core.flow_field.heterogeneous_inflow_config['y'],
+                'speed_multipliers':
+                    self.core.flow_field.heterogeneous_inflow_config['speed_multipliers'][findex:findex+1],
+            }
+            if 'z' in self.core.flow_field.heterogeneous_inflow_config:
+                heterogeneous_inflow_config['z'] = (
+                    self.core.flow_field.heterogeneous_inflow_config['z']
+                )
+        else:
+            heterogeneous_inflow_config = None
+
         self.set(
             wind_speeds=self.wind_speeds[findex:findex+1],
             wind_directions=self.wind_directions[findex:findex+1],
@@ -962,6 +978,7 @@ def set_for_viz(self, findex: int, solver_settings: dict) -> None:
             power_setpoints=self.core.farm.power_setpoints[findex:findex+1,:],
             awc_modes=self.core.farm.awc_modes[findex:findex+1,:],
             awc_amplitudes=self.core.farm.awc_amplitudes[findex:findex+1,:],
+            heterogeneous_inflow_config = heterogeneous_inflow_config,
             solver_settings=solver_settings,
         )