diff --git a/examples/41_test_disable_turbines.py b/examples/41_test_disable_turbines.py new file mode 100644 index 000000000..517845bad --- /dev/null +++ b/examples/41_test_disable_turbines.py @@ -0,0 +1,97 @@ +# Copyright 2023 NREL + +# Licensed under the Apache License, Version 2.0 (the "License"); you may not +# use this file except in compliance with the License. You may obtain a copy of +# the License at http://www.apache.org/licenses/LICENSE-2.0 + +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, WITHOUT +# WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the +# License for the specific language governing permissions and limitations under +# the License. + +# See https://floris.readthedocs.io for documentation + +# Example adapted from https://github.com/NREL/floris/pull/693 contributed by Elie Kadoche + + +import matplotlib.pyplot as plt +import numpy as np +import yaml + +from floris.tools import FlorisInterface + + +""" +This example demonstrates the ability of FLORIS to shut down some turbines +during a simulation. +""" + +# Initialize the FLORIS interface +fi = FlorisInterface("inputs/gch.yaml") + +# Change to the mixed model turbine +with open( + str( + fi.floris.as_dict()["farm"]["turbine_library_path"] + / (fi.floris.as_dict()["farm"]["turbine_type"][0] + ".yaml") + ) +) as t: + turbine_type = yaml.safe_load(t) +turbine_type["power_thrust_model"] = "mixed" +fi.reinitialize(turbine_type=[turbine_type]) + +# Consider a wind farm of 3 aligned wind turbines +layout = np.array([[0.0, 0.0], [500.0, 0.0], [1000.0, 0.0]]) + +# Run the computations for 2 identical wind data +# (n_findex = 2) +wind_directions = np.array([270.0, 270.0]) +wind_speeds = np.array([8.0, 8.0]) + +# Shut down the first 2 turbines for the second findex +# 2 findex x 3 turbines +disable_turbines = np.array([[False, False, False], [True, True, False]]) + +# Simulation +# ------------------------------------------ + +# Reinitialize flow field +fi.reinitialize( + layout_x=layout[:, 0], + layout_y=layout[:, 1], + wind_directions=wind_directions, + wind_speeds=wind_speeds, +) + +# # Compute wakes +fi.calculate_wake(disable_turbines=disable_turbines) + +# Results +# ------------------------------------------ + +# Get powers and effective wind speeds +turbine_powers = fi.get_turbine_powers() +turbine_powers = np.round(turbine_powers * 1e-3, decimals=2) +effective_wind_speeds = fi.turbine_average_velocities + + +# Plot the results +fig, axarr = plt.subplots(2, 1, sharex=True) + +# Plot the power +ax = axarr[0] +ax.plot(["T0", "T1", "T2"], turbine_powers[0, :], "ks-", label="All on") +ax.plot(["T0", "T1", "T2"], turbine_powers[1, :], "ro-", label="T0 & T1 disabled") +ax.set_ylabel("Power (kW)") +ax.grid(True) +ax.legend() + +ax = axarr[1] +ax.plot(["T0", "T1", "T2"], effective_wind_speeds[0, :], "ks-", label="All on") +ax.plot(["T0", "T1", "T2"], effective_wind_speeds[1, :], "ro-", label="T0 & T1 disabled") +ax.set_ylabel("Effective wind speeds (m/s)") +ax.grid(True) +ax.legend() + +plt.show() diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index 1134c7842..2a2a24812 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -31,7 +31,11 @@ ) from floris.tools.cut_plane import CutPlane from floris.tools.wind_data import WindDataBase -from floris.type_dec import floris_array_converter, NDArrayFloat +from floris.type_dec import ( + floris_array_converter, + NDArrayBool, + NDArrayFloat, +) class FlorisInterface(LoggingManager): @@ -122,6 +126,7 @@ def calculate_wake( yaw_angles: NDArrayFloat | list[float] | None = None, # tilt_angles: NDArrayFloat | list[float] | None = None, power_setpoints: NDArrayFloat | list[float] | list[float, None] | None = None, + disable_turbines: NDArrayBool | list[bool] | None = None, ) -> None: """ Wrapper to the :py:meth:`~.Farm.set_yaw_angles` and @@ -133,6 +138,9 @@ def calculate_wake( power_setpoints (NDArrayFloat | list[float] | None, optional): Turbine power setpoints. May be specified with some float values and some None values; power maximization will be assumed for any None value. Defaults to None. + disable_turbines (NDArrayBool | list[bool] | None, optional): NDArray with dimensions + n_findex x n_turbines. True values indicate the turbine is disabled at that findex + and the power setpoint at that position is set to 0. Defaults to None """ if yaw_angles is None: @@ -160,6 +168,33 @@ def calculate_wake( ] = POWER_SETPOINT_DEFAULT power_setpoints = floris_array_converter(power_setpoints) + # Check for turbines to disable + if disable_turbines is not None: + + # Force to numpy array + disable_turbines = np.array(disable_turbines) + + # Must have first dimension = n_findex + if disable_turbines.shape[0] != self.floris.flow_field.n_findex: + raise ValueError( + f"disable_turbines has a size of {disable_turbines.shape[0]} " + f"in the 0th dimension, must be equal to " + f"n_findex={self.floris.flow_field.n_findex}" + ) + + # Must have first dimension = n_turbines + if disable_turbines.shape[1] != self.floris.farm.n_turbines: + raise ValueError( + f"disable_turbines has a size of {disable_turbines.shape[1]} " + f"in the 1th dimension, must be equal to " + f"n_turbines={self.floris.farm.n_turbines}" + ) + + # Set power_setpoints and yaw_angles to 0 in all locations where + # disable_turbines is True + yaw_angles[disable_turbines] = 0.0 + power_setpoints[disable_turbines] = 0.001 # Not zero to avoid numerical problems + self.floris.farm.power_setpoints = power_setpoints # # TODO is this required? @@ -179,6 +214,8 @@ def calculate_wake( def calculate_no_wake( self, yaw_angles: NDArrayFloat | list[float] | None = None, + power_setpoints: NDArrayFloat | list[float] | list[float, None] | None = None, + disable_turbines: NDArrayBool | list[bool] | None = None, ) -> None: """ This function is similar to `calculate_wake()` except @@ -201,6 +238,51 @@ def calculate_no_wake( ) self.floris.farm.yaw_angles = yaw_angles + if power_setpoints is None: + power_setpoints = POWER_SETPOINT_DEFAULT * np.ones( + ( + self.floris.flow_field.n_findex, + self.floris.farm.n_turbines, + ) + ) + else: + power_setpoints = np.array(power_setpoints) + + # Convert any None values to the default power setpoint + power_setpoints[ + power_setpoints == np.full(power_setpoints.shape, None) + ] = POWER_SETPOINT_DEFAULT + power_setpoints = floris_array_converter(power_setpoints) + + # Check for turbines to disable + if disable_turbines is not None: + + # Force to numpy array + # disable_turbines = np.array(disable_turbines) + + # Must have first dimension = n_findex + if disable_turbines.shape[0] != self.floris.flow_field.n_findex: + raise ValueError( + f"disable_turbines has a size of {disable_turbines.shape[0]} " + f"in the 0th dimension, must be equal to " + f"n_findex={self.floris.flow_field.n_findex}" + ) + + # Must have first dimension = n_turbines + if disable_turbines.shape[1] != self.floris.farm.n_turbines: + raise ValueError( + f"disable_turbines has a size of {disable_turbines.shape[1]} " + f"in the 1th dimension, must be equal to " + f"n_turbines={self.floris.farm.n_turbines}" + ) + + # Set power_setpoints and yaw_angles to 0 in all locations where + # disable_turbines is True + yaw_angles[disable_turbines] = 0.0 + power_setpoints[disable_turbines] = 0.001 # Not zero to avoid numerical problems + + self.floris.farm.power_setpoints = power_setpoints + # Initialize solution space self.floris.initialize_domain() diff --git a/tests/floris_interface_test.py b/tests/floris_interface_test.py index 7e41fc90d..694322c7f 100644 --- a/tests/floris_interface_test.py +++ b/tests/floris_interface_test.py @@ -2,6 +2,7 @@ import numpy as np import pytest +import yaml from floris.simulation.turbine.operation_models import POWER_SETPOINT_DEFAULT from floris.tools.floris_interface import FlorisInterface @@ -15,8 +16,6 @@ def test_read_yaml(): fi = FlorisInterface(configuration=YAML_INPUT) assert isinstance(fi, FlorisInterface) - - def test_calculate_wake(): """ In FLORIS v3.2, running calculate_wake twice incorrectly set the yaw angles when the first time @@ -143,6 +142,79 @@ def test_get_farm_power(): farm_power_from_turbine = turbine_powers.sum(axis=1) np.testing.assert_almost_equal(farm_power_from_turbine, farm_powers) +def test_disable_turbines(): + + fi = FlorisInterface(configuration=YAML_INPUT) + + # Set to mixed turbine model + with open( + str( + fi.floris.as_dict()["farm"]["turbine_library_path"] + / (fi.floris.as_dict()["farm"]["turbine_type"][0] + ".yaml") + ) + ) as t: + turbine_type = yaml.safe_load(t) + turbine_type["power_thrust_model"] = "mixed" + fi.reinitialize(turbine_type=[turbine_type]) + + # Init to n-findex = 2, n_turbines = 3 + fi.reinitialize( + wind_speeds=np.array([8.,8.,]), + wind_directions=np.array([270.,270.]), + layout_x = [0,1000,2000], + layout_y=[0,0,0] + ) + + # Confirm that passing in a disable value with wrong n_findex raises error + with pytest.raises(ValueError): + fi.calculate_wake(disable_turbines=np.zeros((10, 3), dtype=bool)) + + # Confirm that passing in a disable value with wrong n_turbines raises error + with pytest.raises(ValueError): + fi.calculate_wake(disable_turbines=np.zeros((2, 10), dtype=bool)) + + # Confirm that if all turbines are disabled, power is near 0 for all turbines + fi.calculate_wake(disable_turbines=np.ones((2, 3), dtype=bool)) + turbines_powers = fi.get_turbine_powers() + np.testing.assert_allclose(turbines_powers,0,atol=0.1) + + # Confirm the same for calculate_no_wake + fi.calculate_no_wake(disable_turbines=np.ones((2, 3), dtype=bool)) + turbines_powers = fi.get_turbine_powers() + np.testing.assert_allclose(turbines_powers,0,atol=0.1) + + # Confirm that if all disabled values set to false, equivalent to running normally + fi.calculate_wake() + turbines_powers_normal = fi.get_turbine_powers() + fi.calculate_wake(disable_turbines=np.zeros((2, 3), dtype=bool)) + turbines_powers_false_disable = fi.get_turbine_powers() + np.testing.assert_allclose(turbines_powers_normal,turbines_powers_false_disable,atol=0.1) + + # Confirm the same for calculate_no_wake + fi.calculate_no_wake() + turbines_powers_normal = fi.get_turbine_powers() + fi.calculate_no_wake(disable_turbines=np.zeros((2, 3), dtype=bool)) + turbines_powers_false_disable = fi.get_turbine_powers() + np.testing.assert_allclose(turbines_powers_normal,turbines_powers_false_disable,atol=0.1) + + # Confirm the shutting off the middle turbine is like removing from the layout + # In terms of impact on third turbine + disable_turbines = np.zeros((2, 3), dtype=bool) + disable_turbines[:,1] = [True, True] + fi.calculate_wake(disable_turbines=disable_turbines) + power_with_middle_disabled = fi.get_turbine_powers() + + fi.reinitialize(layout_x = [0,2000],layout_y = [0, 0]) + fi.calculate_wake() + power_with_middle_removed = fi.get_turbine_powers() + + np.testing.assert_almost_equal(power_with_middle_disabled[0,2], power_with_middle_removed[0,1]) + np.testing.assert_almost_equal(power_with_middle_disabled[1,2], power_with_middle_removed[1,1]) + + # Check that yaw angles are correctly set when turbines are disabled + fi.reinitialize(layout_x = [0,1000,2000],layout_y = [0,0,0]) + fi.calculate_wake(disable_turbines=disable_turbines, yaw_angles=np.ones((2, 3))) + assert (fi.floris.farm.yaw_angles == np.array([[1.0, 0.0, 1.0], [1.0, 0.0, 1.0]])).all() def test_get_farm_aep(): fi = FlorisInterface(configuration=YAML_INPUT)