Add turbine_average_velocities property to UncertainFlorisModel

floris/parallel_floris_model.py

@@ -9,7 +9,7 @@
 from floris.floris_model import FlorisModel
 from floris.logging_manager import LoggingManager
 from floris.optimization.yaw_optimization.yaw_optimizer_sr import YawOptimizationSR
-from floris.uncertain_floris_model import map_turbine_powers_uncertain, UncertainFlorisModel
+from floris.uncertain_floris_model import map_turbine_values_uncertain, UncertainFlorisModel
 
 
 def _get_turbine_powers_serial(fmodel_information, yaw_angles=None):
@@ -367,8 +367,8 @@ def get_turbine_powers(self, yaw_angles=None, no_wake=False):
         t2 = timerpc()
         turbine_powers = self._postprocessing(out)
         if self._is_uncertain:
-            turbine_powers = map_turbine_powers_uncertain(
-                unique_turbine_powers=turbine_powers,
+            turbine_powers = map_turbine_values_uncertain(
+                unique_turbine_values=turbine_powers,
                 map_to_expanded_inputs=self._map_to_expanded_inputs,
                 weights=self._weights,
                 n_unexpanded=self._n_unexpanded,

floris/uncertain_floris_model.py

@@ -10,7 +10,7 @@
 import numpy as np
 
 from floris import FlorisModel
-from floris.core import State
+from floris.core import average_velocity, State
 from floris.logging_manager import LoggingManager
 from floris.par_floris_model import ParFlorisModel
 from floris.type_dec import (
@@ -85,6 +85,10 @@ def __init__(
         fix_yaw_to_nominal_direction=False,
         verbose=False,
     ):
+        # Check validity of inputs
+        if wd_std <= 0:
+            raise ValueError("wd_std must be strictly greater than 0.")
+
         # Save these inputs
         self.wd_resolution = wd_resolution
         self.ws_resolution = ws_resolution
@@ -272,8 +276,8 @@ def _get_turbine_powers(self):
         """
 
         # Pass to off-class function
-        result = map_turbine_powers_uncertain(
-            unique_turbine_powers=self.fmodel_expanded._get_turbine_powers(),
+        result = map_turbine_values_uncertain(
+            unique_turbine_values=self.fmodel_expanded._get_turbine_powers(),
             map_to_expanded_inputs=self.map_to_expanded_inputs,
             weights=self.weights,
             n_unexpanded=self.n_unexpanded,
@@ -1142,6 +1146,27 @@ def core(self):
         """
         return self.fmodel_unexpanded.core
 
+    @property
+    def turbine_average_velocities(self) -> NDArrayFloat:
+        # Get the expanded velocities
+        expanded_velocities = average_velocity(
+            velocities=self.fmodel_expanded.core.flow_field.u,
+            method=self.fmodel_expanded.core.grid.average_method,
+            cubature_weights=self.fmodel_expanded.core.grid.cubature_weights,
+        )
+
+        # Pass to off-class function
+        result = map_turbine_values_uncertain(
+            unique_turbine_values=expanded_velocities,
+            map_to_expanded_inputs=self.map_to_expanded_inputs,
+            weights=self.weights,
+            n_unexpanded=self.n_unexpanded,
+            n_sample_points=self.n_sample_points,
+            n_turbines=self.fmodel_unexpanded.core.farm.n_turbines,
+        )
+
+        return result
+
 
 def map_turbine_powers_uncertain(
     unique_turbine_powers,
@@ -1151,35 +1176,57 @@ def map_turbine_powers_uncertain(
     n_sample_points,
     n_turbines,
 ):
-    """Calculates the power at each turbine in the wind farm based on uncertainty weights.
+    """
+    Alias for map_turbine_values_uncertain.
+    """
+    # Deprecation warning
+    print("map_turbine_powers_uncertain is deprecated, use map_turbine_values_uncertain instead.")
+    return map_turbine_values_uncertain(
+        unique_turbine_values=unique_turbine_powers,
+        map_to_expanded_inputs=map_to_expanded_inputs,
+        weights=weights,
+        n_unexpanded=n_unexpanded,
+        n_sample_points=n_sample_points,
+        n_turbines=n_turbines,
+    )
+
+def map_turbine_values_uncertain(
+    unique_turbine_values,
+    map_to_expanded_inputs,
+    weights,
+    n_unexpanded,
+    n_sample_points,
+    n_turbines,
+):
+    """Calculates values at each turbine in the wind farm based on uncertainty weights.
 
-    This function calculates the power at each turbine in the wind farm, considering
-    the underlying turbine powers and applying a weighted sum to handle uncertainty.
+    This function calculates the values (e.g. power, velocity) at each turbine in the wind farm,
+    considering the underlying turbine values and applying a weighted sum to handle uncertainty.
 
     Args:
-        unique_turbine_powers (NDArrayFloat): An array of unique turbine powers from the
-            underlying FlorisModel
-        map_to_expanded_inputs (NDArrayFloat): An array of indices mapping the unique powers to
-            the expanded powers
+        unique_turbine_values (NDArrayFloat): An array of unique turbine powers, velocities, etc
+            from the underlying FlorisModel
+        map_to_expanded_inputs (NDArrayFloat): An array of indices mapping the unique values to
+            the expanded values
         weights (NDArrayFloat): An array of weights for each wind direction sample point
         n_unexpanded (int): The number of unexpanded conditions
         n_sample_points (int): The number of wind direction sample points
         n_turbines (int): The number of turbines in the wind farm
 
     Returns:
-        NDArrayFloat: An array containing the powers at each turbine for each findex.
+        NDArrayFloat: An array containing the values at each turbine for each findex.
 
     """
 
     # Expand back to the expanded value
-    expanded_turbine_powers = unique_turbine_powers[map_to_expanded_inputs]
+    expanded_turbine_values = unique_turbine_values[map_to_expanded_inputs]
 
     # Reshape the weights array to make it compatible with broadcasting
     weights_reshaped = weights[:, np.newaxis]
 
-    # Reshape expanded_turbine_powers into blocks
+    # Reshape expanded_turbine_values into blocks
     blocks = np.reshape(
-        expanded_turbine_powers,
+        expanded_turbine_values,
         (n_unexpanded, n_sample_points, n_turbines),
         order="F",
     )
@@ -1221,7 +1268,7 @@ def __init__(
             yaw_resolution,
             power_setpoint_resolution,
             awc_amplitude_resolution,
-            wd_std=1.0,
+            wd_std=1.0,  # Arbitrary nonzero value, not used since only one sample point
             wd_sample_points=[0],
             fix_yaw_to_nominal_direction=False,
             verbose=verbose,

tests/uncertain_floris_model_integration_test.py

@@ -521,3 +521,104 @@ def test_copy():
     pufmodel_copy = pufmodel.copy()
     assert isinstance(pufmodel_copy, UncertainFlorisModel)
     assert isinstance(pufmodel_copy.fmodel_expanded, ParFlorisModel)
+
+def test_invalid_wd_std():
+    """
+    Test that the UncertainFlorisModel raises asn error with a wd_std of 0 or negative.
+    """
+    with pytest.raises(ValueError):
+        UncertainFlorisModel(configuration=YAML_INPUT, wd_std=0.0)
+
+    with pytest.raises(ValueError):
+        UncertainFlorisModel(configuration=YAML_INPUT, wd_std=-1.0)
+
+def test_turbine_average_velocities_shape_and_type():
+    """
+    Test that turbine_average_velocities returns the correct shape and type.
+    """
+    ufmodel = UncertainFlorisModel(configuration=YAML_INPUT)
+
+    # Set up a simple 2-turbine wind farm
+    ufmodel.set(
+        layout_x=[0, 500],
+        layout_y=[0, 0],
+        wind_speeds=[8.0, 10.0],
+        wind_directions=[270.0, 280.0],
+        turbulence_intensities=[0.06, 0.06],
+    )
+
+    ufmodel.run()
+
+    # Get turbine average velocities
+    velocities = ufmodel.turbine_average_velocities
+
+    # Check type
+    assert isinstance(velocities, np.ndarray)
+
+    # Check shape: should be (n_findex, n_turbines)
+    expected_shape = (ufmodel.n_findex, ufmodel.n_turbines)
+    assert velocities.shape == expected_shape
+
+    # Check that values are positive and reasonable
+    assert np.all(velocities > 0)
+    assert np.all(velocities < 20)  # Reasonable upper bound for wind speeds
+
+
+def test_turbine_average_velocities_free_stream():
+    """
+    Test that turbine_average_velocities returns the correct shape and type.
+    """
+    ufmodel = UncertainFlorisModel(configuration=YAML_INPUT)
+
+    # Set up a simple 2-turbine wind farm with no wake interactions
+    ufmodel.set(
+        layout_x=[0, 0],
+        layout_y=[0, 1000],
+        wind_speeds=[8.0, 10.0],
+        wind_directions=[270.0, 280.0],
+        turbulence_intensities=[0.06, 0.06],
+        wind_shear=0.0,  # Turn off shear to simplify the test
+    )
+
+    ufmodel.run()
+
+    # Get turbine average velocities
+    velocities = ufmodel.turbine_average_velocities
+
+    # Velocities should be the same as the wind speeds but n_turbines columns repeated
+    assert np.allclose(velocities, np.array([[8.0, 8.0], [10.0, 10.0]]))
+
+def test_turbine_average_velocities_uncertain_vs_certain():
+    """
+    Test that turbine_average_velocities returns the same values for uncertain and certain models.
+    """
+
+    # Set up (certain) FlorisModel
+    fmodel = FlorisModel(configuration=YAML_INPUT)
+    fmodel.set(
+        layout_x=[0, 500],
+        layout_y=[0, 0],
+        wind_speeds=[8.0, 10.0],
+        wind_directions=[270.0, 270.0],
+        turbulence_intensities=[0.06, 0.06],
+    )
+    fmodel.run()
+    velocities_certain = fmodel.turbine_average_velocities
+
+    # Create equivalent uncertain model
+    ufmodel = UncertainFlorisModel(configuration=fmodel)
+    ufmodel.run()
+    velocities_uncertain = ufmodel.turbine_average_velocities
+
+    # Check that the upstream turbine matches
+    assert np.allclose(velocities_uncertain[:, 0], velocities_certain[:, 0])
+    # Downstream turbine higher than certain when aligned
+    assert np.all(velocities_uncertain[:, 1] > velocities_certain[:, 1])
+
+    # Create a near 0-std uncertain model
+    ufmodel_zero_std = UncertainFlorisModel(configuration=fmodel, wd_std=1e-5)
+    ufmodel_zero_std.run()
+    velocities_uncertain_zero_std = ufmodel_zero_std.turbine_average_velocities
+
+    # Check that the uncertain model with 0 std matches the certain model
+    assert np.allclose(velocities_uncertain_zero_std, velocities_certain)