From 59f3a603911cc287bc67e71cb3421a8cda386053 Mon Sep 17 00:00:00 2001 From: Bart Doekemeijer Date: Wed, 15 Feb 2023 17:50:20 +0100 Subject: [PATCH 01/36] Add coordinates from inertial frame to grid and add visualization capabilities --- floris/simulation/grid.py | 75 +++++++++++++++++++++++--- floris/tools/floris_interface.py | 15 ++++-- floris/tools/visualization.py | 90 ++++++++++++++++++-------------- floris/utilities.py | 45 ++++++++++++++-- 4 files changed, 174 insertions(+), 51 deletions(-) diff --git a/floris/simulation/grid.py b/floris/simulation/grid.py index 9f8fdc3ed..9c72ae59a 100644 --- a/floris/simulation/grid.py +++ b/floris/simulation/grid.py @@ -28,8 +28,11 @@ NDArrayFloat, NDArrayInt, ) -from floris.utilities import rotate_coordinates_rel_west, Vec3 - +from floris.utilities import ( + rotate_coordinates_rel_west, + reverse_rotate_coordinates_rel_west, + Vec3 +) @define class Grid(ABC): @@ -69,6 +72,8 @@ class Grid(ABC): wind_speeds: NDArrayFloat = field(converter=floris_array_converter) time_series: bool = field() + _xc_rot: float = field(init=False) + _yc_rot: float = field(init=False) n_turbines: int = field(init=False) n_wind_speeds: int = field(init=False) n_wind_directions: int = field(init=False) @@ -79,6 +84,12 @@ class Grid(ABC): x_sorted: NDArrayFloat = field(init=False) y_sorted: NDArrayFloat = field(init=False) z_sorted: NDArrayFloat = field(init=False) + x_inertial_frame: NDArrayFloat = field(init=False) + y_inertial_frame: NDArrayFloat = field(init=False) + z_inertial_frame: NDArrayFloat = field(init=False) + x_sorted_inertial_frame: NDArrayFloat = field(init=False) + y_sorted_inertial_frame: NDArrayFloat = field(init=False) + z_sorted_inertial_frame: NDArrayFloat = field(init=False) def __attrs_post_init__(self) -> None: self.turbine_coordinates_array = np.array([c.elements for c in self.turbine_coordinates]) @@ -202,7 +213,7 @@ def set_grid(self) -> None: # the foot of the turbine where the tower meets the ground. # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z = rotate_coordinates_rel_west(self.wind_directions, self.turbine_coordinates_array) + x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west(self.wind_directions, self.turbine_coordinates_array) # - **rloc** (*float, optional): A value, from 0 to 1, that determines # the width/height of the grid of points on the rotor as a ratio of @@ -265,6 +276,33 @@ def set_grid(self) -> None: self.y = np.take_along_axis(self.y_sorted, self.unsorted_indices, axis=2) self.z = np.take_along_axis(self.z_sorted, self.unsorted_indices, axis=2) + # Now calculate grid coordinates in original frame (from 270 deg perspective) + self.x_inertial_frame, self.y_inertial_frame, self.z_inertial_frame = \ + reverse_rotate_coordinates_rel_west( + wind_directions=self.wind_directions, + grid_x=self.x, + grid_y=self.y, + grid_z=self.z, + x_center_of_rotation=xc_rot, + y_center_of_rotation=yc_rot, + ) + + # Now calculate grid coordinates in original frame (from 270 deg perspective) + self.x_sorted_inertial_frame, self.y_sorted_inertial_frame, self.z_sorted_inertial_frame = \ + reverse_rotate_coordinates_rel_west( + wind_directions=self.wind_directions, + grid_x=self.x_sorted, + grid_y=self.y_sorted, + grid_z=self.z_sorted, + x_center_of_rotation=xc_rot, + y_center_of_rotation=yc_rot, + ) + + # Save center of rotation to self + self._xc_rot = xc_rot + self._yc_rot = yc_rot + + @define class FlowFieldGrid(Grid): """ @@ -296,7 +334,7 @@ def set_grid(self) -> None: """ # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z = rotate_coordinates_rel_west(self.wind_directions, self.turbine_coordinates_array) + x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west(self.wind_directions, self.turbine_coordinates_array) # Construct the arrays storing the grid points eps = 0.01 @@ -318,6 +356,21 @@ def set_grid(self) -> None: self.y_sorted = y_points[None, None, :, :, :] self.z_sorted = z_points[None, None, :, :, :] + # Now calculate grid coordinates in original frame (from 270 deg perspective) + self.x_sorted_inertial_frame, self.y_sorted_inertial_frame, self.z_sorted_inertial_frame = \ + reverse_rotate_coordinates_rel_west( + wind_directions=self.wind_directions, + grid_x=self.x_sorted, + grid_y=self.y_sorted, + grid_z=self.z_sorted, + x_center_of_rotation=xc_rot, + y_center_of_rotation=yc_rot, + ) + + # Save center of rotation to self + self._xc_rot = xc_rot + self._yc_rot = yc_rot + @define class FlowFieldPlanarGrid(Grid): """ @@ -355,8 +408,7 @@ def set_grid(self) -> None: Then, create the grid based on this wind-from-left orientation """ # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z = rotate_coordinates_rel_west(self.wind_directions, self.turbine_coordinates_array) - + x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west(self.wind_directions, self.turbine_coordinates_array) max_diameter = np.max(self.reference_turbine_diameter) if self.normal_vector == "z": # Rules of thumb for horizontal plane @@ -418,6 +470,17 @@ def set_grid(self) -> None: self.y_sorted = y_points[None, None, :, :, :] self.z_sorted = z_points[None, None, :, :, :] + # Now calculate grid coordinates in original frame (from 270 deg perspective) + self.x_sorted_inertial_frame, self.y_sorted_inertial_frame, self.z_sorted_inertial_frame = \ + reverse_rotate_coordinates_rel_west( + wind_directions=self.wind_directions, + grid_x=self.x_sorted, + grid_y=self.y_sorted, + grid_z=self.z_sorted, + x_center_of_rotation=xc_rot, + y_center_of_rotation=yc_rot, + ) + # self.sorted_indices = self.x.argsort(axis=2) # self.unsorted_indices = self.sorted_indices.argsort(axis=2) diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index 0c90724fc..86aed21e0 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -264,6 +264,7 @@ def get_plane_of_points( self, normal_vector="z", planar_coordinate=None, + rotate_to_inertial_frame=False, ): """ Calculates velocity values through the @@ -281,9 +282,15 @@ def get_plane_of_points( :py:class:`pandas.DataFrame`: containing values of x1, x2, u, v, w """ # Get results vectors - x_flat = self.floris.grid.x_sorted[0, 0].flatten() - y_flat = self.floris.grid.y_sorted[0, 0].flatten() - z_flat = self.floris.grid.z_sorted[0, 0].flatten() + rotate_to_inertial_frame=True, + if (normal_vector == "z") and rotate_to_inertial_frame: + x_flat = self.floris.grid.x_sorted_inertial_frame[0, 0].flatten() + y_flat = self.floris.grid.y_sorted_inertial_frame[0, 0].flatten() + z_flat = self.floris.grid.z_sorted_inertial_frame[0, 0].flatten() + else: + x_flat = self.floris.grid.x_sorted[0, 0].flatten() + y_flat = self.floris.grid.y_sorted[0, 0].flatten() + z_flat = self.floris.grid.z_sorted[0, 0].flatten() u_flat = self.floris.flow_field.u_sorted[0, 0].flatten() v_flat = self.floris.flow_field.v_sorted[0, 0].flatten() w_flat = self.floris.flow_field.w_sorted[0, 0].flatten() @@ -347,6 +354,7 @@ def calculate_horizontal_plane( wd=None, ws=None, yaw_angles=None, + rotate_to_inertial_frame=False, ): """ Shortcut method to instantiate a :py:class:`~.tools.cut_plane.CutPlane` @@ -402,6 +410,7 @@ def calculate_horizontal_plane( df = self.get_plane_of_points( normal_vector="z", planar_coordinate=height, + rotate_to_inertial_frame=rotate_to_inertial_frame, ) # Compute the cutplane diff --git a/floris/tools/visualization.py b/floris/tools/visualization.py index 178ad17e5..51e23b9fb 100644 --- a/floris/tools/visualization.py +++ b/floris/tools/visualization.py @@ -39,7 +39,8 @@ def plot_turbines( yaw_angles, rotor_diameters, color: str | None = None, - wind_direction: float = 270.0 + wind_direction: float = 270.0, + rotate_to_inertial_frame=False, ): """ Plot wind plant layout from turbine locations. @@ -55,12 +56,16 @@ def plot_turbines( """ if color is None: color = "k" - - coordinates_array = np.array([[x, y, 0.0] for x, y in list(zip(layout_x, layout_y))]) - layout_x, layout_y, _ = rotate_coordinates_rel_west( - np.array([wind_direction]), - coordinates_array - ) + + if not rotate_to_inertial_frame: + coordinates_array = np.array([[x, y, 0.0] for x, y in list(zip(layout_x, layout_y))]) + layout_x, layout_y, _, _, _ = rotate_coordinates_rel_west( + np.array([wind_direction]), + coordinates_array + ) + else: + layout_x = layout_x[None, None, :] + layout_y = layout_y[None, None, :] for x, y, yaw, d in zip(layout_x[0,0], layout_y[0,0], yaw_angles, rotor_diameters): R = d / 2.0 @@ -71,7 +76,7 @@ def plot_turbines( ax.plot([x_0, x_1], [y_0, y_1], color=color) -def plot_turbines_with_fi(fi: FlorisInterface, ax=None, color=None, yaw_angles=None): +def plot_turbines_with_fi(fi: FlorisInterface, ax=None, color=None, yaw_angles=None, rotate_to_inertial_frame=False): """ Wrapper function to plot turbines which extracts the data from a FLORIS interface object @@ -86,6 +91,10 @@ def plot_turbines_with_fi(fi: FlorisInterface, ax=None, color=None, yaw_angles=N fig, ax = plt.subplots() if yaw_angles is None: yaw_angles = fi.floris.farm.yaw_angles + + if rotate_to_inertial_frame: + yaw_angles = yaw_angles + (270.0 - fi.floris.flow_field.wind_directions[0]) + plot_turbines( ax, fi.layout_x, @@ -94,10 +103,11 @@ def plot_turbines_with_fi(fi: FlorisInterface, ax=None, color=None, yaw_angles=N fi.floris.farm.rotor_diameters[0, 0], color=color, wind_direction=fi.floris.flow_field.wind_directions[0], + rotate_to_inertial_frame=rotate_to_inertial_frame, ) -def add_turbine_id_labels(fi: FlorisInterface, ax: plt.Axes, **kwargs): +def add_turbine_id_labels(fi: FlorisInterface, ax: plt.Axes, rotate_to_inertial_frame=False, **kwargs): """ Adds index labels to a plot based on the given FlorisInterface. See the pyplot.annotate docs for more info: @@ -109,15 +119,20 @@ def add_turbine_id_labels(fi: FlorisInterface, ax: plt.Axes, **kwargs): fi (FlorisInterface): Simulation object to get the layout and index information. ax (plt.Axes): Axes object to add the labels. """ - coordinates_array = np.array([ - [x, y, 0.0] - for x, y in list(zip(fi.layout_x, fi.layout_y)) - ]) - wind_direction = fi.floris.flow_field.wind_directions[0] - layout_x, layout_y, _ = rotate_coordinates_rel_west( - np.array([wind_direction]), - coordinates_array - ) + + if rotate_to_inertial_frame: + coordinates_array = np.array([ + [x, y, 0.0] + for x, y in list(zip(fi.layout_x, fi.layout_y)) + ]) + wind_direction = fi.floris.flow_field.wind_directions[0] + layout_x, layout_y, _, _, _ = rotate_coordinates_rel_west( + np.array([wind_direction]), + coordinates_array + ) + else: + layout_x = fi.layout_x[None, None, :] + layout_y = fi.layout_y[None, None, :] for i in range(fi.floris.farm.n_turbines): ax.annotate( @@ -148,15 +163,18 @@ def line_contour_cut_plane(cut_plane, ax=None, levels=None, colors=None, **kwarg if not ax: fig, ax = plt.subplots() - # Reshape UMesh internally - x1_mesh = cut_plane.df.x1.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) - x2_mesh = cut_plane.df.x2.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) - u_mesh = cut_plane.df.u.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) - Zm = np.ma.masked_where(np.isnan(u_mesh), u_mesh) rcParams["contour.negative_linestyle"] = "solid" # Plot the cut-through - contours = ax.contour(x1_mesh, x2_mesh, Zm, levels=levels, colors=colors, **kwargs) + contours = ax.tricontour( + cut_plane.df.x1, + cut_plane.df.x2, + cut_plane.df.u, + levels=levels, + colors=colors, + extend="both", + **kwargs, + ) ax.clabel(contours, contours.levels, inline=True, fontsize=10, colors="black") @@ -198,38 +216,34 @@ def visualize_cut_plane( if not ax: fig, ax = plt.subplots() if vel_component=='u': - vel_mesh = cut_plane.df.u.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) + # vel_mesh = cut_plane.df.u.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) if min_speed is None: min_speed = cut_plane.df.u.min() if max_speed is None: max_speed = cut_plane.df.u.max() elif vel_component=='v': - vel_mesh = cut_plane.df.v.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) + # vel_mesh = cut_plane.df.v.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) if min_speed is None: min_speed = cut_plane.df.v.min() if max_speed is None: max_speed = cut_plane.df.v.max() elif vel_component=='w': - vel_mesh = cut_plane.df.w.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) + # vel_mesh = cut_plane.df.w.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) if min_speed is None: min_speed = cut_plane.df.w.min() if max_speed is None: max_speed = cut_plane.df.w.max() - # Reshape to 2d for plotting - x1_mesh = cut_plane.df.x1.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) - x2_mesh = cut_plane.df.x2.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) - Zm = np.ma.masked_where(np.isnan(vel_mesh), vel_mesh) - # Plot the cut-through - im = ax.pcolormesh( - x1_mesh, - x2_mesh, - Zm, - cmap=cmap, + fig, ax = plt.subplots() + im = ax.tricontourf( + cut_plane.df.x1, + cut_plane.df.x2, + cut_plane.df.u, vmin=min_speed, vmax=max_speed, - shading="nearest" + cmap=cmap, + extend="both", ) # Add line contour diff --git a/floris/utilities.py b/floris/utilities.py index 276bf26b6..7adcddca4 100644 --- a/floris/utilities.py +++ b/floris/utilities.py @@ -195,7 +195,7 @@ def wind_delta(wind_directions): return ((wind_directions - 270) % 360 + 360) % 360 -def rotate_coordinates_rel_west(wind_directions, coordinates): +def rotate_coordinates_rel_west(wind_directions, coordinates, x_center_of_rotation=None, y_center_of_rotation=None): # Calculate the difference in given wind direction from 270 / West wind_deviation_from_west = wind_delta(wind_directions) wind_deviation_from_west = np.reshape(wind_deviation_from_west, (len(wind_directions), 1, 1)) @@ -204,8 +204,10 @@ def rotate_coordinates_rel_west(wind_directions, coordinates): x_coordinates, y_coordinates, z_coordinates = coordinates.T # Find center of rotation - this is the center of box bounding all of the turbines - x_center_of_rotation = (np.min(x_coordinates) + np.max(x_coordinates)) / 2 - y_center_of_rotation = (np.min(y_coordinates) + np.max(y_coordinates)) / 2 + if x_center_of_rotation is None: + x_center_of_rotation = (np.min(x_coordinates) + np.max(x_coordinates)) / 2 + if y_center_of_rotation is None: + y_center_of_rotation = (np.min(y_coordinates) + np.max(y_coordinates)) / 2 # Rotate turbine coordinates about the center x_coord_offset = x_coordinates - x_center_of_rotation @@ -221,7 +223,42 @@ def rotate_coordinates_rel_west(wind_directions, coordinates): + y_center_of_rotation ) z_coord_rotated = np.ones_like(wind_deviation_from_west) * z_coordinates - return x_coord_rotated, y_coord_rotated, z_coord_rotated + return x_coord_rotated, y_coord_rotated, z_coord_rotated, x_center_of_rotation, y_center_of_rotation + + +def reverse_rotate_coordinates_rel_west(wind_directions, grid_x, grid_y, grid_z, x_center_of_rotation, y_center_of_rotation): + # Calculate the difference in given wind direction from 270 / West + wind_deviation_from_west = -1.0 * wind_delta(wind_directions) # We are rotating in the other direction + + # Construct the arrays storing the turbine locations + grid_x_reversed = np.zeros_like(grid_x) + grid_y_reversed = np.zeros_like(grid_x) + grid_z_reversed = np.zeros_like(grid_x) + for wii, angle_rotation in enumerate(wind_deviation_from_west): + x_rot = grid_x[wii, :, :, :, :] + y_rot = grid_y[wii, :, :, :, :] + z_rot = grid_z[wii, :, :, :, :] + + # Rotate turbine coordinates about the center + x_rot_offset = x_rot - x_center_of_rotation + y_rot_offset = y_rot - y_center_of_rotation + x = ( + x_rot_offset * cosd(angle_rotation) + - y_rot_offset * sind(angle_rotation) + + x_center_of_rotation + ) + y = ( + x_rot_offset * sind(angle_rotation) + + y_rot_offset * cosd(angle_rotation) + + y_center_of_rotation + ) + z = z_rot # Nothing changed in this rotation + + grid_x_reversed[wii, :, :, :, :] = x + grid_y_reversed[wii, :, :, :, :] = y + grid_z_reversed[wii, :, :, :, :] = z + + return grid_x_reversed, grid_y_reversed, grid_z_reversed class Loader(yaml.SafeLoader): From 996680e31870a4a5a951ec3dc3e09172d1147fb9 Mon Sep 17 00:00:00 2001 From: Bart Doekemeijer Date: Thu, 16 Feb 2023 07:49:27 +0100 Subject: [PATCH 02/36] Update unit test for utilities, bug fix --- tests/utilities_unit_test.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/tests/utilities_unit_test.py b/tests/utilities_unit_test.py index eeefc6f93..e9f952928 100644 --- a/tests/utilities_unit_test.py +++ b/tests/utilities_unit_test.py @@ -89,7 +89,7 @@ def test_rotate_coordinates_rel_west(): # For 270, the coordinates should not change. wind_directions = np.array([270.0]) - x_rotated, y_rotated, z_rotated = rotate_coordinates_rel_west(wind_directions, coordinates) + x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west(wind_directions, coordinates) np.testing.assert_array_equal( X_COORDS, x_rotated[0,0] ) np.testing.assert_array_equal( Y_COORDS, y_rotated[0,0] ) @@ -105,7 +105,7 @@ def test_rotate_coordinates_rel_west(): # NOTE: These adjustments are not general and will fail if the coordinates in # conftest change. wind_directions = np.array([360.0]) - x_rotated, y_rotated, z_rotated = rotate_coordinates_rel_west(wind_directions, coordinates) + x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west(wind_directions, coordinates) np.testing.assert_almost_equal( Y_COORDS, x_rotated[0,0] - np.min(x_rotated[0,0])) np.testing.assert_almost_equal( X_COORDS, y_rotated[0,0] - np.min(y_rotated[0,0])) np.testing.assert_almost_equal( @@ -114,7 +114,7 @@ def test_rotate_coordinates_rel_west(): ) wind_directions = np.array([90.0]) - x_rotated, y_rotated, z_rotated = rotate_coordinates_rel_west(wind_directions, coordinates) + x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west(wind_directions, coordinates) np.testing.assert_almost_equal( X_COORDS[-1:-4:-1], x_rotated[0,0] ) np.testing.assert_almost_equal( Y_COORDS, y_rotated[0,0] ) np.testing.assert_almost_equal( Z_COORDS, z_rotated[0,0] ) From 01ce5c209507e55fe275b65cc923c057c09ae4e6 Mon Sep 17 00:00:00 2001 From: Rafael M Mudafort Date: Wed, 15 Mar 2023 10:48:22 -0500 Subject: [PATCH 03/36] Remove unused variables --- floris/simulation/grid.py | 11 ----------- 1 file changed, 11 deletions(-) diff --git a/floris/simulation/grid.py b/floris/simulation/grid.py index 9c72ae59a..9433cc9ae 100644 --- a/floris/simulation/grid.py +++ b/floris/simulation/grid.py @@ -72,8 +72,6 @@ class Grid(ABC): wind_speeds: NDArrayFloat = field(converter=floris_array_converter) time_series: bool = field() - _xc_rot: float = field(init=False) - _yc_rot: float = field(init=False) n_turbines: int = field(init=False) n_wind_speeds: int = field(init=False) n_wind_directions: int = field(init=False) @@ -298,11 +296,6 @@ def set_grid(self) -> None: y_center_of_rotation=yc_rot, ) - # Save center of rotation to self - self._xc_rot = xc_rot - self._yc_rot = yc_rot - - @define class FlowFieldGrid(Grid): """ @@ -367,10 +360,6 @@ def set_grid(self) -> None: y_center_of_rotation=yc_rot, ) - # Save center of rotation to self - self._xc_rot = xc_rot - self._yc_rot = yc_rot - @define class FlowFieldPlanarGrid(Grid): """ From abf803038f64d53eab81f16ad032273ed0419705 Mon Sep 17 00:00:00 2001 From: Bart Doekemeijer Date: Tue, 2 May 2023 14:18:19 +0200 Subject: [PATCH 04/36] Tiny bug fixes in visualization with inertial grid rotation functionality --- floris/tools/visualization.py | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/floris/tools/visualization.py b/floris/tools/visualization.py index d7083ab6c..d85d66de4 100644 --- a/floris/tools/visualization.py +++ b/floris/tools/visualization.py @@ -99,8 +99,8 @@ def plot_turbines_with_fi(fi: FlorisInterface, ax=None, color=None, yaw_angles=N ax, fi.layout_x, fi.layout_y, - yaw_angles[0, 0], - fi.floris.farm.rotor_diameters[0, 0], + yaw_angles.flatten(), + fi.floris.farm.rotor_diameters.flatten(), color=color, wind_direction=fi.floris.flow_field.wind_directions[0], rotate_to_inertial_frame=rotate_to_inertial_frame, @@ -190,6 +190,7 @@ def visualize_cut_plane( max_speed=None, cmap="coolwarm", levels=None, + clevels=None, color_bar=False, title="", **kwargs @@ -215,6 +216,7 @@ def visualize_cut_plane( if not ax: fig, ax = plt.subplots() + if vel_component=='u': # vel_mesh = cut_plane.df.u.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) if min_speed is None: @@ -234,14 +236,18 @@ def visualize_cut_plane( if max_speed is None: max_speed = cut_plane.df.w.max() + # Allow separate number of levels for tricontourf and for line_contour + if clevels is None: + clevels = levels + # Plot the cut-through - fig, ax = plt.subplots() im = ax.tricontourf( cut_plane.df.x1, cut_plane.df.x2, cut_plane.df.u, vmin=min_speed, vmax=max_speed, + levels=clevels, cmap=cmap, extend="both", ) From 624ca2e489a6fb91fa554a0b0592e9aced77064a Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 21:33:23 -0600 Subject: [PATCH 05/36] remove unnecessary calculations --- floris/simulation/grid.py | 37 ++++++++++++++++--------------------- 1 file changed, 16 insertions(+), 21 deletions(-) diff --git a/floris/simulation/grid.py b/floris/simulation/grid.py index 0d4f06a1c..59be6dc60 100644 --- a/floris/simulation/grid.py +++ b/floris/simulation/grid.py @@ -29,11 +29,12 @@ NDArrayInt, ) from floris.utilities import ( - rotate_coordinates_rel_west, reverse_rotate_coordinates_rel_west, - Vec3 + rotate_coordinates_rel_west, + Vec3, ) + @define class Grid(ABC): """ @@ -215,7 +216,10 @@ def set_grid(self) -> None: # the foot of the turbine where the tower meets the ground. # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west(self.wind_directions, self.turbine_coordinates_array) + x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west( + self.wind_directions, + self.turbine_coordinates_array, + ) # - **rloc** (*float, optional): A value, from 0 to 1, that determines # the width/height of the grid of points on the rotor as a ratio of @@ -247,7 +251,7 @@ def set_grid(self) -> None: disc_area_radius, self.grid_resolution, dtype=floris_float_type, - axis=1 + axis=1, ) # Construct the turbine grids # Here, they are already rotated to the correct orientation for each wind direction @@ -276,21 +280,6 @@ def set_grid(self) -> None: self.y_sorted = np.take_along_axis(_y, self.sorted_indices, axis=2) self.z_sorted = np.take_along_axis(_z, self.sorted_indices, axis=2) - self.x = np.take_along_axis(self.x_sorted, self.unsorted_indices, axis=2) - self.y = np.take_along_axis(self.y_sorted, self.unsorted_indices, axis=2) - self.z = np.take_along_axis(self.z_sorted, self.unsorted_indices, axis=2) - - # Now calculate grid coordinates in original frame (from 270 deg perspective) - self.x_inertial_frame, self.y_inertial_frame, self.z_inertial_frame = \ - reverse_rotate_coordinates_rel_west( - wind_directions=self.wind_directions, - grid_x=self.x, - grid_y=self.y, - grid_z=self.z, - x_center_of_rotation=xc_rot, - y_center_of_rotation=yc_rot, - ) - # Now calculate grid coordinates in original frame (from 270 deg perspective) self.x_sorted_inertial_frame, self.y_sorted_inertial_frame, self.z_sorted_inertial_frame = \ reverse_rotate_coordinates_rel_west( @@ -333,7 +322,10 @@ def set_grid(self) -> None: """ # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west(self.wind_directions, self.turbine_coordinates_array) + x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west( + self.wind_directions, + self.turbine_coordinates_array, + ) # Construct the arrays storing the grid points eps = 0.01 @@ -403,7 +395,10 @@ def set_grid(self) -> None: Then, create the grid based on this wind-from-left orientation """ # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west(self.wind_directions, self.turbine_coordinates_array) + x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west( + self.wind_directions, + self.turbine_coordinates_array, + ) max_diameter = np.max(self.reference_turbine_diameter) if self.normal_vector == "z": # Rules of thumb for horizontal plane From 9ac6e210707c496c997e05c81507c877c7e91e8a Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 21:36:35 -0600 Subject: [PATCH 06/36] remove hardcoded flag --- floris/tools/floris_interface.py | 1 - 1 file changed, 1 deletion(-) diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index f19ced717..1bb38cdbf 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -300,7 +300,6 @@ def get_plane_of_points( :py:class:`pandas.DataFrame`: containing values of x1, x2, x3, u, v, w """ # Get results vectors - rotate_to_inertial_frame=True, if (normal_vector == "z") and rotate_to_inertial_frame: x_flat = self.floris.grid.x_sorted_inertial_frame[0, 0].flatten() y_flat = self.floris.grid.y_sorted_inertial_frame[0, 0].flatten() From 1a92904707141c8b0f6e529ddbbeba91675757f5 Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 21:37:42 -0600 Subject: [PATCH 07/36] reconfigure heterogeneous map input and usage --- floris/tools/floris_interface.py | 81 ++++++++++++++++++++------------ 1 file changed, 51 insertions(+), 30 deletions(-) diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index 1bb38cdbf..b591f2ac9 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -50,7 +50,7 @@ class FlorisInterface(LoggerBase): - **logging**: See `floris.simulation.floris.Floris` for more details. """ - def __init__(self, configuration: dict | str | Path, het_map=None): + def __init__(self, configuration: dict | str | Path, het_config: dict | None = None): self.configuration = configuration if isinstance(self.configuration, (str, Path)): @@ -62,7 +62,28 @@ def __init__(self, configuration: dict | str | Path, het_map=None): else: raise TypeError("The Floris `configuration` must be of type 'dict', 'str', or 'Path'.") - # Store the heterogeneous map for use after reinitailization + # Create the heterogeneous wind map interpolant if het_config is supplied + # Check that the correct keys are supplied for the het_config dict + if het_config is not None: + for k in ["speed_ups", "x_locs", "y_locs"]: + if k not in het_config.keys(): + raise ValueError( + "het_config must contain entries for 'speeds_ups', 'x_locs', and 'y_locs'," + f" with 'z_locs' as optional. Missing '{k}'." + ) + if "z_locs" not in het_config: + het_config["z_locs"] = None + # Create the heterogeneous wind map interpolant + het_map = self.generate_heterogeneous_wind_map( + speed_ups=het_config['speed_ups'], + x=het_config['x_locs'], + y=het_config['y_locs'], + z=het_config['z_locs'], + ) + else: + het_map = None + # Store the heterogeneous info and map for use after reinitailization + self.het_config = het_config self.het_map = het_map # Assign the heterogeneous map to the flow field # Needed for a direct call to fi.calculate_wake without fi.reinitialize @@ -960,6 +981,34 @@ def get_farm_AEP_wind_rose_class( return aep + def generate_heterogeneous_wind_map(self, speed_ups, x, y, z=None): + if z is not None: + # Compute the 3-dimensional interpolants for each wind diretion + # Linear interpolation is used for points within the user-defined area of values, + # while a nearest-neighbor interpolant is used for points outside that region + in_region = [ + LinearNDInterpolator(list(zip(x, y, z)), speed_up, fill_value=1.0) + for speed_up in speed_ups + ] + out_region = [ + NearestNDInterpolator(list(zip(x, y, z)), speed_up) + for speed_up in speed_ups + ] + else: + # Compute the 2-dimensional interpolants for each wind diretion + # Linear interpolation is used for points within the user-defined area of values, + # while a nearest-neighbor interpolant is used for points outside that region + in_region = [ + LinearNDInterpolator(list(zip(x, y)), speed_up, fill_value=1.0) + for speed_up in speed_ups + ] + out_region = [ + NearestNDInterpolator(list(zip(x, y)), speed_up) + for speed_up in speed_ups + ] + + return [in_region, out_region] + @property def layout_x(self): @@ -1000,34 +1049,6 @@ def get_turbine_layout(self, z=False): return xcoords, ycoords -def generate_heterogeneous_wind_map(speed_ups, x, y, z=None): - if z is not None: - # Compute the 3-dimensional interpolants for each wind diretion - # Linear interpolation is used for points within the user-defined area of values, - # while a nearest-neighbor interpolant is used for points outside that region - in_region = [ - LinearNDInterpolator(list(zip(x, y, z)), speed_up, fill_value=np.nan) - for speed_up in speed_ups - ] - out_region = [ - NearestNDInterpolator(list(zip(x, y, z)), speed_up) - for speed_up in speed_ups - ] - else: - # Compute the 2-dimensional interpolants for each wind diretion - # Linear interpolation is used for points within the user-defined area of values, - # while a nearest-neighbor interpolant is used for points outside that region - in_region = [ - LinearNDInterpolator(list(zip(x, y)), speed_up, fill_value=np.nan) - for speed_up in speed_ups - ] - out_region = [ - NearestNDInterpolator(list(zip(x, y)), speed_up) - for speed_up in speed_ups - ] - - return [in_region, out_region] - ## Functionality removed in v3 def set_rotor_diameter(self, rotor_diameter): From e48237bcea22b37427d227054d240f1f9afd4c53 Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 21:38:45 -0600 Subject: [PATCH 08/36] update heterogeneous calculation to use rotated intertial_frame coordinates --- floris/simulation/flow_field.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/floris/simulation/flow_field.py b/floris/simulation/flow_field.py index 2e6011a12..04c4d8022 100644 --- a/floris/simulation/flow_field.py +++ b/floris/simulation/flow_field.py @@ -98,14 +98,14 @@ def initialize_velocity_field(self, grid: Grid) -> None: if len(self.het_map[0][0].points[0]) == 2: speed_ups = self.calculate_speed_ups( self.het_map, - grid.x_sorted, - grid.y_sorted + grid.x_sorted_inertial_frame, + grid.y_sorted_inertial_frame ) elif len(self.het_map[0][0].points[0]) == 3: speed_ups = self.calculate_speed_ups( self.het_map, - grid.x_sorted, - grid.y_sorted, + grid.x_sorted_inertial_frame, + grid.y_sorted_inertial_frame, grid.z_sorted ) From 57607216d58d11ffe962de901b88846a09c22b05 Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 21:47:54 -0600 Subject: [PATCH 09/36] adding specific plotting method for heterogeneous flows --- floris/tools/visualization.py | 125 ++++++++++++++++++++++++++++++++-- 1 file changed, 121 insertions(+), 4 deletions(-) diff --git a/floris/tools/visualization.py b/floris/tools/visualization.py index d85d66de4..dc8354fb6 100644 --- a/floris/tools/visualization.py +++ b/floris/tools/visualization.py @@ -22,6 +22,7 @@ import numpy as np import pandas as pd from matplotlib import rcParams +from scipy.spatial import ConvexHull from floris.simulation import Floris from floris.tools.cut_plane import CutPlane @@ -56,7 +57,7 @@ def plot_turbines( """ if color is None: color = "k" - + if not rotate_to_inertial_frame: coordinates_array = np.array([[x, y, 0.0] for x, y in list(zip(layout_x, layout_y))]) layout_x, layout_y, _, _, _ = rotate_coordinates_rel_west( @@ -76,7 +77,13 @@ def plot_turbines( ax.plot([x_0, x_1], [y_0, y_1], color=color) -def plot_turbines_with_fi(fi: FlorisInterface, ax=None, color=None, yaw_angles=None, rotate_to_inertial_frame=False): +def plot_turbines_with_fi( + fi: FlorisInterface, + ax=None, + color=None, + yaw_angles=None, + rotate_to_inertial_frame=False, +): """ Wrapper function to plot turbines which extracts the data from a FLORIS interface object @@ -91,7 +98,7 @@ def plot_turbines_with_fi(fi: FlorisInterface, ax=None, color=None, yaw_angles=N fig, ax = plt.subplots() if yaw_angles is None: yaw_angles = fi.floris.farm.yaw_angles - + if rotate_to_inertial_frame: yaw_angles = yaw_angles + (270.0 - fi.floris.flow_field.wind_directions[0]) @@ -107,7 +114,12 @@ def plot_turbines_with_fi(fi: FlorisInterface, ax=None, color=None, yaw_angles=N ) -def add_turbine_id_labels(fi: FlorisInterface, ax: plt.Axes, rotate_to_inertial_frame=False, **kwargs): +def add_turbine_id_labels( + fi: FlorisInterface, + ax: plt.Axes, + rotate_to_inertial_frame=False, + **kwargs, +): """ Adds index labels to a plot based on the given FlorisInterface. See the pyplot.annotate docs for more info: @@ -279,6 +291,111 @@ def visualize_cut_plane( return im +def visualize_heterogeneous_cut_plane( + cut_plane, + fi, + ax=None, + vel_component='u', + min_speed=None, + max_speed=None, + cmap="coolwarm", + levels=None, + color_bar=False, + title="", + plot_het_bounds=True, + **kwargs +): + """ + Generate pseudocolor mesh plot of the cut_plane. + + Args: + cut_plane (:py:class:`~.tools.cut_plane.CutPlane`): 2D + plane through wind plant. + ax (:py:class:`matplotlib.pyplot.axes`): Figure axes. Defaults + to None. + min_speed (float, optional): Minimum value of wind speed for + contours. Defaults to None. + max_speed (float, optional): Maximum value of wind speed for + contours. Defaults to None. + cmap (str, optional): Colormap specifier. Defaults to + 'coolwarm'. + + Returns: + im (:py:class:`matplotlib.plt.pcolormesh`): Image handle. + """ + + if not ax: + fig, ax = plt.subplots() + if vel_component=='u': + # vel_mesh = cut_plane.df.u.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) + if min_speed is None: + min_speed = cut_plane.df.u.min() + if max_speed is None: + max_speed = cut_plane.df.u.max() + elif vel_component=='v': + # vel_mesh = cut_plane.df.v.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) + if min_speed is None: + min_speed = cut_plane.df.v.min() + if max_speed is None: + max_speed = cut_plane.df.v.max() + elif vel_component=='w': + # vel_mesh = cut_plane.df.w.values.reshape(cut_plane.resolution[1], cut_plane.resolution[0]) + if min_speed is None: + min_speed = cut_plane.df.w.min() + if max_speed is None: + max_speed = cut_plane.df.w.max() + + # Plot the cut-through + im = ax.tricontourf( + cut_plane.df.x1, + cut_plane.df.x2, + cut_plane.df.u, + vmin=min_speed, + vmax=max_speed, + cmap=cmap, + extend="both", + ) + + # Add line contour + line_contour_cut_plane( + cut_plane, + ax=ax, + levels=levels, + colors="b", + linewidths=0.8, + alpha=0.3, + **kwargs + ) + + # Plot the user-defined heterogeneous flow area + if plot_het_bounds: + points = np.array(list(zip(fi.het_config['x_locs'], fi.het_config['y_locs']))) + hull = ConvexHull(points) + h = ax.plot( + points[np.append(hull.vertices, hull.vertices[0]),0], + points[np.append(hull.vertices, hull.vertices[0]), 1], + 'k--', + lw=2, + ) + ax.plot(points[hull.vertices,0], points[hull.vertices,1], 'ko') + ax.legend(h, ["defined heterogeneous bounds"], loc=1) + + if cut_plane.normal_vector == "x": + ax.invert_xaxis() + + if color_bar: + cbar = plt.colorbar(im, ax=ax) + cbar.set_label('m/s') + + # Set the title + ax.set_title(title) + + # Make equal axis + ax.set_aspect("equal") + + return im + + def visualize_quiver(cut_plane, ax=None, min_speed=None, max_speed=None, downSamp=1, **kwargs): """ Visualize the in-plane flows in a cut_plane using quiver. From 702b42be415aa3a4363c919a29a694160ea64a87 Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 22:16:54 -0600 Subject: [PATCH 10/36] adding wd input to plot_turbines_with_fi to properly plot turbine yaw setting relative to wd --- floris/tools/visualization.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/floris/tools/visualization.py b/floris/tools/visualization.py index dc8354fb6..509e70be9 100644 --- a/floris/tools/visualization.py +++ b/floris/tools/visualization.py @@ -27,7 +27,7 @@ from floris.simulation import Floris from floris.tools.cut_plane import CutPlane from floris.tools.floris_interface import FlorisInterface -from floris.utilities import rotate_coordinates_rel_west +from floris.utilities import rotate_coordinates_rel_west, wind_delta def show_plots(): @@ -81,6 +81,7 @@ def plot_turbines_with_fi( fi: FlorisInterface, ax=None, color=None, + wd=None, yaw_angles=None, rotate_to_inertial_frame=False, ): @@ -98,9 +99,11 @@ def plot_turbines_with_fi( fig, ax = plt.subplots() if yaw_angles is None: yaw_angles = fi.floris.farm.yaw_angles + if wd is None: + wd = fi.floris.flow_field.wind_directions[0] - if rotate_to_inertial_frame: - yaw_angles = yaw_angles + (270.0 - fi.floris.flow_field.wind_directions[0]) + # Rotate yaw angles to inertial frame for plotting turbines relative to wind direction + yaw_angles = yaw_angles - wind_delta(np.array(wd)) plot_turbines( ax, From fd34ee144ad3156a2dd177e079f5fe90fd956cc0 Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 22:17:42 -0600 Subject: [PATCH 11/36] making rotating the wind turbines in plots relative to the wind direction the default behavior --- floris/tools/floris_interface.py | 5 +--- floris/tools/visualization.py | 48 +++++++++++--------------------- 2 files changed, 18 insertions(+), 35 deletions(-) diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index b591f2ac9..6cbf399bc 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -304,7 +304,6 @@ def get_plane_of_points( self, normal_vector="z", planar_coordinate=None, - rotate_to_inertial_frame=False, ): """ Calculates velocity values through the @@ -321,7 +320,7 @@ def get_plane_of_points( :py:class:`pandas.DataFrame`: containing values of x1, x2, x3, u, v, w """ # Get results vectors - if (normal_vector == "z") and rotate_to_inertial_frame: + if (normal_vector == "z"): x_flat = self.floris.grid.x_sorted_inertial_frame[0, 0].flatten() y_flat = self.floris.grid.y_sorted_inertial_frame[0, 0].flatten() z_flat = self.floris.grid.z_sorted_inertial_frame[0, 0].flatten() @@ -392,7 +391,6 @@ def calculate_horizontal_plane( wd=None, ws=None, yaw_angles=None, - rotate_to_inertial_frame=False, ): """ Shortcut method to instantiate a :py:class:`~.tools.cut_plane.CutPlane` @@ -448,7 +446,6 @@ def calculate_horizontal_plane( df = self.get_plane_of_points( normal_vector="z", planar_coordinate=height, - rotate_to_inertial_frame=rotate_to_inertial_frame, ) # Compute the cutplane diff --git a/floris/tools/visualization.py b/floris/tools/visualization.py index 509e70be9..2e1fe7f15 100644 --- a/floris/tools/visualization.py +++ b/floris/tools/visualization.py @@ -41,7 +41,6 @@ def plot_turbines( rotor_diameters, color: str | None = None, wind_direction: float = 270.0, - rotate_to_inertial_frame=False, ): """ Plot wind plant layout from turbine locations. @@ -58,15 +57,12 @@ def plot_turbines( if color is None: color = "k" - if not rotate_to_inertial_frame: - coordinates_array = np.array([[x, y, 0.0] for x, y in list(zip(layout_x, layout_y))]) - layout_x, layout_y, _, _, _ = rotate_coordinates_rel_west( - np.array([wind_direction]), - coordinates_array - ) - else: - layout_x = layout_x[None, None, :] - layout_y = layout_y[None, None, :] + # Rotate layout to inertial frame for plotting turbines relative to wind direction + coordinates_array = np.array([[x, y, 0.0] for x, y in list(zip(layout_x, layout_y))]) + layout_x, layout_y, _, _, _ = rotate_coordinates_rel_west( + np.array([wind_direction]), + coordinates_array + ) for x, y, yaw, d in zip(layout_x[0,0], layout_y[0,0], yaw_angles, rotor_diameters): R = d / 2.0 @@ -83,7 +79,6 @@ def plot_turbines_with_fi( color=None, wd=None, yaw_angles=None, - rotate_to_inertial_frame=False, ): """ Wrapper function to plot turbines which extracts the data @@ -113,16 +108,10 @@ def plot_turbines_with_fi( fi.floris.farm.rotor_diameters.flatten(), color=color, wind_direction=fi.floris.flow_field.wind_directions[0], - rotate_to_inertial_frame=rotate_to_inertial_frame, ) -def add_turbine_id_labels( - fi: FlorisInterface, - ax: plt.Axes, - rotate_to_inertial_frame=False, - **kwargs, -): +def add_turbine_id_labels(fi: FlorisInterface, ax: plt.Axes, **kwargs): """ Adds index labels to a plot based on the given FlorisInterface. See the pyplot.annotate docs for more info: @@ -135,19 +124,16 @@ def add_turbine_id_labels( ax (plt.Axes): Axes object to add the labels. """ - if rotate_to_inertial_frame: - coordinates_array = np.array([ - [x, y, 0.0] - for x, y in list(zip(fi.layout_x, fi.layout_y)) - ]) - wind_direction = fi.floris.flow_field.wind_directions[0] - layout_x, layout_y, _, _, _ = rotate_coordinates_rel_west( - np.array([wind_direction]), - coordinates_array - ) - else: - layout_x = fi.layout_x[None, None, :] - layout_y = fi.layout_y[None, None, :] + # Rotate layout to intertial frame for plotting turbines relative to wind direction + coordinates_array = np.array([ + [x, y, 0.0] + for x, y in list(zip(fi.layout_x, fi.layout_y)) + ]) + wind_direction = fi.floris.flow_field.wind_directions[0] + layout_x, layout_y, _, _, _ = rotate_coordinates_rel_west( + np.array([wind_direction]), + coordinates_array + ) for i in range(fi.floris.farm.n_turbines): ax.annotate( From c746201be27b22fca15c43598021da00a45e2944 Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 22:25:52 -0600 Subject: [PATCH 12/36] removing unnecessary parameters from grid class definition --- floris/simulation/grid.py | 6 ------ 1 file changed, 6 deletions(-) diff --git a/floris/simulation/grid.py b/floris/simulation/grid.py index 59be6dc60..e30e987dc 100644 --- a/floris/simulation/grid.py +++ b/floris/simulation/grid.py @@ -77,15 +77,9 @@ class Grid(ABC): n_wind_speeds: int = field(init=False) n_wind_directions: int = field(init=False) turbine_coordinates_array: NDArrayFloat = field(init=False) - x: NDArrayFloat = field(init=False, default=[]) - y: NDArrayFloat = field(init=False, default=[]) - z: NDArrayFloat = field(init=False, default=[]) x_sorted: NDArrayFloat = field(init=False) y_sorted: NDArrayFloat = field(init=False) z_sorted: NDArrayFloat = field(init=False) - x_inertial_frame: NDArrayFloat = field(init=False) - y_inertial_frame: NDArrayFloat = field(init=False) - z_inertial_frame: NDArrayFloat = field(init=False) x_sorted_inertial_frame: NDArrayFloat = field(init=False) y_sorted_inertial_frame: NDArrayFloat = field(init=False) z_sorted_inertial_frame: NDArrayFloat = field(init=False) From 26312a816a62c1f7c5f680fa397590449c73da4f Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 23:15:08 -0600 Subject: [PATCH 13/36] refactored het map generation and added docstrings --- floris/tools/floris_interface.py | 88 ++++++++++++++++++++++++-------- 1 file changed, 66 insertions(+), 22 deletions(-) diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index 6cbf399bc..99674adfd 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -48,6 +48,17 @@ class FlorisInterface(LoggerBase): - **turbine**: See `floris.simulation.turbine.Turbine` for more details. - **wake**: See `floris.simulation.wake.WakeManager` for more details. - **logging**: See `floris.simulation.floris.Floris` for more details. + het_config (:py:obj:`dict`, optional): The heterogeneous inflow configuration dictionary. + The configuration should have the following inputs specified. + - **speed_ups**: See `floris.tools.floris_interface.generate_heterogeneous_wind_map' + for more details. + - **x_locs**: See `floris.tools.floris_interface.generate_heterogeneous_wind_map' + for more details. + - **y_locs**: See `floris.tools.floris_interface.generate_heterogeneous_wind_map' + for more details. + - **z_locs** (optional): See + `floris.tools.floris_interface.generate_heterogeneous_wind_map' + for more details. """ def __init__(self, configuration: dict | str | Path, het_config: dict | None = None): @@ -63,23 +74,8 @@ def __init__(self, configuration: dict | str | Path, het_config: dict | None = N raise TypeError("The Floris `configuration` must be of type 'dict', 'str', or 'Path'.") # Create the heterogeneous wind map interpolant if het_config is supplied - # Check that the correct keys are supplied for the het_config dict if het_config is not None: - for k in ["speed_ups", "x_locs", "y_locs"]: - if k not in het_config.keys(): - raise ValueError( - "het_config must contain entries for 'speeds_ups', 'x_locs', and 'y_locs'," - f" with 'z_locs' as optional. Missing '{k}'." - ) - if "z_locs" not in het_config: - het_config["z_locs"] = None - # Create the heterogeneous wind map interpolant - het_map = self.generate_heterogeneous_wind_map( - speed_ups=het_config['speed_ups'], - x=het_config['x_locs'], - y=het_config['y_locs'], - z=het_config['z_locs'], - ) + het_map = self.generate_heterogeneous_wind_map(het_config=het_config) else: het_map = None # Store the heterogeneous info and map for use after reinitailization @@ -132,7 +128,7 @@ def assign_hub_height_to_ref_height(self): def copy(self): """Create an independent copy of the current FlorisInterface object""" - return FlorisInterface(self.floris.as_dict(), het_map=self.het_map) + return FlorisInterface(self.floris.as_dict(), het_config=self.het_config) def calculate_wake( self, @@ -236,7 +232,7 @@ def reinitialize( solver_settings: dict | None = None, time_series: bool = False, layout: tuple[list[float], list[float]] | tuple[NDArrayFloat, NDArrayFloat] | None = None, - het_map=None, + het_config=None, ): # Export the floris object recursively as a dictionary floris_dict = self.floris.as_dict() @@ -260,8 +256,8 @@ def reinitialize( flow_field_dict["turbulence_intensity"] = turbulence_intensity if air_density is not None: flow_field_dict["air_density"] = air_density - if het_map is not None: - self.het_map = het_map + if het_config is not None: + self.het_config = het_config ## Farm if layout is not None: @@ -298,7 +294,12 @@ def reinitialize( # Create a new instance of floris and attach to self self.floris = Floris.from_dict(floris_dict) # Re-assign the hetergeneous inflow map to flow field - self.floris.flow_field.het_map = self.het_map + if self.het_config is not None: + het_map = self.generate_heterogeneous_wind_map(het_config=self.het_config) + else: + het_map = None + self.het_map = het_map + self.floris.flow_field.het_map = het_map def get_plane_of_points( self, @@ -978,7 +979,50 @@ def get_farm_AEP_wind_rose_class( return aep - def generate_heterogeneous_wind_map(self, speed_ups, x, y, z=None): + def generate_heterogeneous_wind_map(self, het_config): + """This function creates the heterogeneous interpolant used to calculate heterogeneous + inflows. + + Args: + het_config (dict): The heterogeneous inflow configuration dictionary. + The configuration should have the following inputs specified. + - **speed_ups** (list): A list of speed up factors that will multiply the specified + freestream wind speed. This 2-dimensional array should have an array of + multiplicative factors defined for each wind direction. + - **x_locs** (list): A list of x locations at which the speed up factors are + defined. + - **y_locs**: A list of y locations at which the speed up factors are + defined. + - **z_locs** (optional): A list of z locations at which the speed up factors are + defined. + + Raises: + ValueError: ("het_config must contain entries for 'speeds_ups', 'x_locs', and 'y_locs'," + f" with 'z_locs' as optional. Missing '{k}'.") raised if the configuration + dictionary is missing one of the required key-value pairs. + + Returns: + LinearNDInterpolator: Returns a linear interpolant for computing wind speed based on an + x and y location in the flow field. This is computed using Scipy's + LinearNDInterpolator and uses a fill value equal to the freestream for interpolated + values outside of the user-defined heterogeneous map bounds. + """ + # Check that the correct keys are supplied for the het_config dict + for k in ["speed_ups", "x_locs", "y_locs"]: + if k not in het_config.keys(): + raise ValueError( + "het_config must contain entries for 'speeds_ups', 'x_locs', and 'y_locs'," + f" with 'z_locs' as optional. Missing '{k}'." + ) + if "z_locs" not in het_config: + # If only a 2D case, add "None" for the z locations + het_config["z_locs"] = None + + speed_ups = het_config['speed_ups'] + x = het_config['x_locs'] + y = het_config['y_locs'] + z = het_config['z_locs'] + if z is not None: # Compute the 3-dimensional interpolants for each wind diretion # Linear interpolation is used for points within the user-defined area of values, From 307d4e4815f731194c8eab5f32969f34b1827a6c Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 23:21:43 -0600 Subject: [PATCH 14/36] removing out_region calculation for het_map as a freestream fill value is now used --- floris/simulation/flow_field.py | 32 +++++--------------------------- floris/tools/floris_interface.py | 14 +++----------- 2 files changed, 8 insertions(+), 38 deletions(-) diff --git a/floris/simulation/flow_field.py b/floris/simulation/flow_field.py index 04c4d8022..b78706982 100644 --- a/floris/simulation/flow_field.py +++ b/floris/simulation/flow_field.py @@ -95,13 +95,13 @@ def initialize_velocity_field(self, grid: Grid) -> None: # If heterogeneous flow data is given, the speed ups at the defined # grid locations are determined in either 2 or 3 dimensions. else: - if len(self.het_map[0][0].points[0]) == 2: + if len(self.het_map[0].points[0]) == 2: speed_ups = self.calculate_speed_ups( self.het_map, grid.x_sorted_inertial_frame, grid.y_sorted_inertial_frame ) - elif len(self.het_map[0][0].points[0]) == 3: + elif len(self.het_map[0].points[0]) == 3: speed_ups = self.calculate_speed_ups( self.het_map, grid.x_sorted_inertial_frame, @@ -174,7 +174,7 @@ def finalize(self, unsorted_indices): def calculate_speed_ups(self, het_map, x, y, z=None): # Check that the het maps wd dimension matches - if self.n_wind_directions!= np.array(het_map).shape[1]: + if self.n_wind_directions!= np.array(het_map).shape[0]: raise ValueError( "het_map's wind direction dimension not equal to number of wind directions" ) @@ -183,38 +183,16 @@ def calculate_speed_ups(self, het_map, x, y, z=None): # Calculate the 3-dimensional speed ups; reshape is needed as the generator # adds an extra dimension speed_ups = np.reshape( - [het_map[0][i](x[i:i+1], y[i:i+1], z[i:i+1]) for i in range( len(het_map[0]))], + [het_map[i](x[i:i+1], y[i:i+1], z[i:i+1]) for i in range( len(het_map))], np.shape(x) ) - # If there are any points requested outside the user-defined area, use the - # nearest-neighbor interplonat to determine those speed up values - if np.isnan(speed_ups).any(): - idx_nan = np.where(np.isnan(speed_ups)) - speed_ups_out_of_region = np.reshape( - [het_map[1][i](x[i:i+1], y[i:i+1], z[i:i+1]) for i in range(len(het_map[1]))], - np.shape(x) - ) - - speed_ups[idx_nan] = speed_ups_out_of_region[idx_nan] - else: # Calculate the 2-dimensional speed ups; reshape is needed as the generator # adds an extra dimension speed_ups = np.reshape( - [het_map[0][i](x[i:i+1], y[i:i+1]) for i in range(len(het_map[0]))], + [het_map[i](x[i:i+1], y[i:i+1]) for i in range(len(het_map))], np.shape(x) ) - # If there are any points requested outside the user-defined area, use the - # nearest-neighbor interplonat to determine those speed up values - if np.isnan(speed_ups).any(): - idx_nan = np.where(np.isnan(speed_ups)) - speed_ups_out_of_region = np.reshape( - [het_map[1][i](x[i:i+1], y[i:i+1]) for i in range(len(het_map[1]))], - np.shape(x) - ) - - speed_ups[idx_nan] = speed_ups_out_of_region[idx_nan] - return speed_ups diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index 99674adfd..a845e210e 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -1026,29 +1026,21 @@ def generate_heterogeneous_wind_map(self, het_config): if z is not None: # Compute the 3-dimensional interpolants for each wind diretion # Linear interpolation is used for points within the user-defined area of values, - # while a nearest-neighbor interpolant is used for points outside that region + # while the freestream wind speed is used for points outside that region in_region = [ LinearNDInterpolator(list(zip(x, y, z)), speed_up, fill_value=1.0) for speed_up in speed_ups ] - out_region = [ - NearestNDInterpolator(list(zip(x, y, z)), speed_up) - for speed_up in speed_ups - ] else: # Compute the 2-dimensional interpolants for each wind diretion # Linear interpolation is used for points within the user-defined area of values, - # while a nearest-neighbor interpolant is used for points outside that region + # while the freestream wind speed is used for points outside that region in_region = [ LinearNDInterpolator(list(zip(x, y)), speed_up, fill_value=1.0) for speed_up in speed_ups ] - out_region = [ - NearestNDInterpolator(list(zip(x, y)), speed_up) - for speed_up in speed_ups - ] - return [in_region, out_region] + return in_region @property From f257f42439fd7789d8fda7e715d8ac82103239b1 Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 23:22:02 -0600 Subject: [PATCH 15/36] updating uncertainty_interface to accept het_config --- floris/tools/uncertainty_interface.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/floris/tools/uncertainty_interface.py b/floris/tools/uncertainty_interface.py index 82db14d0e..877aa56f3 100644 --- a/floris/tools/uncertainty_interface.py +++ b/floris/tools/uncertainty_interface.py @@ -26,7 +26,7 @@ class UncertaintyInterface(LoggerBase): def __init__( self, configuration, - het_map=None, + het_conifg=None, unc_options=None, unc_pmfs=None, fix_yaw_in_relative_frame=False, @@ -110,7 +110,7 @@ def __init__( if isinstance(configuration, FlorisInterface): self.fi = configuration else: - self.fi = FlorisInterface(configuration, het_map=het_map) + self.fi = FlorisInterface(configuration, het_config=het_conifg) self.reinitialize_uncertainty( unc_options=unc_options, From 63e8f30a289fa1ca7b4ca5b63adbb7ead6eb0b4d Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 23:22:25 -0600 Subject: [PATCH 16/36] updating heterogeneous examples to use new het_config --- examples/16_heterogeneous_inflow.py | 22 +++++++++++------ examples/16b_heterogenaity_multiple_ws_wd.py | 25 +++++++++++++------- 2 files changed, 31 insertions(+), 16 deletions(-) diff --git a/examples/16_heterogeneous_inflow.py b/examples/16_heterogeneous_inflow.py index e0745902c..1f84d7efe 100644 --- a/examples/16_heterogeneous_inflow.py +++ b/examples/16_heterogeneous_inflow.py @@ -16,7 +16,6 @@ import matplotlib.pyplot as plt from floris.tools import FlorisInterface -from floris.tools.floris_interface import generate_heterogeneous_wind_map from floris.tools.visualization import visualize_cut_plane @@ -44,12 +43,16 @@ x_locs = [-300.0, -300.0, 2600.0, 2600.0] y_locs = [ -300.0, 300.0, -300.0, 300.0] -# Generate the linear interpolation to be used for the heterogeneous inflow. -het_map_2d = generate_heterogeneous_wind_map(speed_ups, x_locs, y_locs) +# Create the configuration dictionary to be used for the heterogeneous inflow. +het_config_2d = { + 'speed_ups': speed_ups, + 'x_locs': x_locs, + 'y_locs': y_locs, +} # Initialize FLORIS with the given input file via FlorisInterface. # Also, pass the heterogeneous map into the FlorisInterface. -fi_2d = FlorisInterface("inputs/gch.yaml", het_map=het_map_2d) +fi_2d = FlorisInterface("inputs/gch.yaml", het_config=het_config_2d) # Set shear to 0.0 to highlight the heterogeneous inflow fi_2d.reinitialize(wind_shear=0.0) @@ -92,12 +95,17 @@ y_locs = [-300.0, 300.0, -300.0, 300.0, -300.0, 300.0, -300.0, 300.0] z_locs = [540.0, 540.0, 0.0, 0.0, 540.0, 540.0, 0.0, 0.0] -# Generate the linear interpolation to be used for the heterogeneous inflow. -het_map_3d = generate_heterogeneous_wind_map(speed_ups, x_locs, y_locs, z_locs) +# Create the configuration dictionary to be used for the heterogeneous inflow. +het_config_3d = { + 'speed_ups': speed_ups, + 'x_locs': x_locs, + 'y_locs': y_locs, + 'z_locs': z_locs, +} # Initialize FLORIS with the given input file via FlorisInterface. # Also, pass the heterogeneous map into the FlorisInterface. -fi_3d = FlorisInterface("inputs/gch.yaml", het_map=het_map_3d) +fi_3d = FlorisInterface("inputs/gch.yaml", het_config=het_config_3d) # Set shear to 0.0 to highlight the heterogeneous inflow fi_3d.reinitialize(wind_shear=0.0) diff --git a/examples/16b_heterogenaity_multiple_ws_wd.py b/examples/16b_heterogenaity_multiple_ws_wd.py index d456d9f7f..a35a73aa7 100644 --- a/examples/16b_heterogenaity_multiple_ws_wd.py +++ b/examples/16b_heterogenaity_multiple_ws_wd.py @@ -17,7 +17,6 @@ import numpy as np from floris.tools import FlorisInterface -from floris.tools.floris_interface import generate_heterogeneous_wind_map from floris.tools.visualization import visualize_cut_plane @@ -34,19 +33,23 @@ x_locs = [-300.0, -300.0, 2600.0, 2600.0] y_locs = [ -300.0, 300.0, -300.0, 300.0] -# Generate the linear interpolation to be used for the heterogeneous inflow. -het_map_2d = generate_heterogeneous_wind_map(speed_ups, x_locs, y_locs) +# Create the configuration dictionary to be used for the heterogeneous inflow. +het_config_2d = { + 'speed_ups': speed_ups, + 'x_locs': x_locs, + 'y_locs': y_locs, +} # Initialize FLORIS with the given input file via FlorisInterface. # Also, pass the heterogeneous map into the FlorisInterface. -fi = FlorisInterface("inputs/gch.yaml", het_map=het_map_2d) +fi = FlorisInterface("inputs/gch.yaml", het_config=het_config_2d) # Set shear to 0.0 to highlight the heterogeneous inflow fi.reinitialize( wind_shear=0.0, wind_speeds=[8.0], wind_directions=[270.], - layout_x=[0,0], + layout_x=[0, 0], layout_y=[-300, 300], ) fi.calculate_wake() @@ -54,7 +57,7 @@ # Show the initial results print('------------------------------------------') -print('Given the speedsups and turbine locations, ') +print('Given the speedups and turbine locations, ') print(' the first turbine has an inflow wind speed') print(' twice that of the second') print(' Wind Speed = 8., Wind Direction = 270.') @@ -74,9 +77,13 @@ # To change the number of wind directions however it is necessary to make a matching # change to the dimensions of the het map -speed_ups = [[2.0, 1.0, 2.0, 1.0], [2.0, 1.0, 2.0, 1.0] ] # Expand to two wind directions -het_map_2d = generate_heterogeneous_wind_map(speed_ups, x_locs, y_locs) -fi.reinitialize(wind_directions=[270., 275.], wind_speeds=[8.], het_map=het_map_2d) +speed_ups = [[2.0, 1.0, 2.0, 1.0], [2.0, 1.0, 2.0, 1.0]] # Expand to two wind directions +het_config_2d = { + 'speed_ups': speed_ups, + 'x_locs': x_locs, + 'y_locs': y_locs, +} +fi.reinitialize(wind_directions=[270., 275.], wind_speeds=[8.], het_config=het_config_2d) fi.calculate_wake() turbine_powers = np.round(fi.get_turbine_powers() / 1000.) print('With wind directions now set to 270 and 275 deg') From ab9dbc0b706137e45d73b7cdf854b980c8d02f85 Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 23:26:16 -0600 Subject: [PATCH 17/36] formatting upates --- floris/utilities.py | 21 +++++++++++++++++---- tests/utilities_unit_test.py | 15 ++++++++++++--- 2 files changed, 29 insertions(+), 7 deletions(-) diff --git a/floris/utilities.py b/floris/utilities.py index 84a86c046..ce41648c5 100644 --- a/floris/utilities.py +++ b/floris/utilities.py @@ -218,7 +218,12 @@ def wind_delta(wind_directions: NDArrayFloat | float): return (wind_directions - 270) % 360 -def rotate_coordinates_rel_west(wind_directions, coordinates, x_center_of_rotation=None, y_center_of_rotation=None): +def rotate_coordinates_rel_west( + wind_directions, + coordinates, + x_center_of_rotation=None, + y_center_of_rotation=None +): """ This function rotates the given coordinates so that they are aligned with West (270) rather than North (0). The rotation happens about the centroid of the coordinates. @@ -256,12 +261,20 @@ def rotate_coordinates_rel_west(wind_directions, coordinates, x_center_of_rotati + y_center_of_rotation ) z_coord_rotated = np.ones_like(wind_deviation_from_west) * z_coordinates - return x_coord_rotated, y_coord_rotated, z_coord_rotated, x_center_of_rotation, y_center_of_rotation + return x_coord_rotated, y_coord_rotated, z_coord_rotated, x_center_of_rotation, \ + y_center_of_rotation -def reverse_rotate_coordinates_rel_west(wind_directions, grid_x, grid_y, grid_z, x_center_of_rotation, y_center_of_rotation): +def reverse_rotate_coordinates_rel_west( + wind_directions, + grid_x, grid_y, + grid_z, + x_center_of_rotation, + y_center_of_rotation +): # Calculate the difference in given wind direction from 270 / West - wind_deviation_from_west = -1.0 * wind_delta(wind_directions) # We are rotating in the other direction + # We are rotating in the other direction + wind_deviation_from_west = -1.0 * wind_delta(wind_directions) # Construct the arrays storing the turbine locations grid_x_reversed = np.zeros_like(grid_x) diff --git a/tests/utilities_unit_test.py b/tests/utilities_unit_test.py index 55e9fc5d3..40ea594a8 100644 --- a/tests/utilities_unit_test.py +++ b/tests/utilities_unit_test.py @@ -91,7 +91,10 @@ def test_rotate_coordinates_rel_west(): # For 270, the coordinates should not change. wind_directions = np.array([270.0]) - x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west(wind_directions, coordinates) + x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west( + wind_directions, + coordinates, + ) np.testing.assert_array_equal( X_COORDS, x_rotated[0,0] ) np.testing.assert_array_equal( Y_COORDS, y_rotated[0,0] ) @@ -107,7 +110,10 @@ def test_rotate_coordinates_rel_west(): # NOTE: These adjustments are not general and will fail if the coordinates in # conftest change. wind_directions = np.array([360.0]) - x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west(wind_directions, coordinates) + x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west( + wind_directions, + coordinates, + ) np.testing.assert_almost_equal( Y_COORDS, x_rotated[0,0] - np.min(x_rotated[0,0])) np.testing.assert_almost_equal( X_COORDS, y_rotated[0,0] - np.min(y_rotated[0,0])) np.testing.assert_almost_equal( @@ -116,7 +122,10 @@ def test_rotate_coordinates_rel_west(): ) wind_directions = np.array([90.0]) - x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west(wind_directions, coordinates) + x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west( + wind_directions, + coordinates, + ) np.testing.assert_almost_equal( X_COORDS[-1:-4:-1], x_rotated[0,0] ) np.testing.assert_almost_equal( Y_COORDS, y_rotated[0,0] ) np.testing.assert_almost_equal( Z_COORDS, z_rotated[0,0] ) From 09b4aa1e389d7708345e01366bf9c36922edfdfe Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 23:30:22 -0600 Subject: [PATCH 18/36] updating parallel_computing_interface to use het_config --- floris/tools/parallel_computing_interface.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/floris/tools/parallel_computing_interface.py b/floris/tools/parallel_computing_interface.py index ffbfaade9..c9f9923c4 100644 --- a/floris/tools/parallel_computing_interface.py +++ b/floris/tools/parallel_computing_interface.py @@ -12,17 +12,17 @@ def _load_local_floris_object( fi_dict, - het_map=None, + het_config=None, unc_pmfs=None, fix_yaw_in_relative_frame=False ): # Load local FLORIS object if unc_pmfs is None: - fi = FlorisInterface(fi_dict, het_map=het_map) + fi = FlorisInterface(fi_dict, het_config=het_config) else: fi = UncertaintyInterface( fi_dict, - het_map=het_map, + het_config=het_config, unc_pmfs=unc_pmfs, fix_yaw_in_relative_frame=fix_yaw_in_relative_frame, ) From 0b461ff37e8a0f092ba39aab52fced785f1273fb Mon Sep 17 00:00:00 2001 From: bayc Date: Wed, 3 May 2023 23:40:39 -0600 Subject: [PATCH 19/36] updating documentation --- floris/tools/visualization.py | 17 +++++++++++++---- 1 file changed, 13 insertions(+), 4 deletions(-) diff --git a/floris/tools/visualization.py b/floris/tools/visualization.py index 2e1fe7f15..f3fdaa93e 100644 --- a/floris/tools/visualization.py +++ b/floris/tools/visualization.py @@ -85,10 +85,11 @@ def plot_turbines_with_fi( from a FLORIS interface object Args: - fi (:py:class:`floris.tools.flow_data.FlowData`): - FlowData object. - ax (:py:class:`matplotlib.pyplot.axes`): figure axes. - color (str, optional): Color to plot turbines + fi (:py:class:`floris.tools.floris_interface.FlorisIntreface`): FlorisIntreface object. + ax (:py:class:`matplotlib.pyplot.axes`): Figure axes. Defaults to None. + color (str, optional): Color to plot turbines. Defaults to None. + wd (list, optional): The wind direction to plot the turbines relative to. Defaults to None. + yaw_angles (NDArray, optional): The yaw angles for the turbines. Defaults to None. """ if not ax: fig, ax = plt.subplots() @@ -210,6 +211,14 @@ def visualize_cut_plane( contours. Defaults to None. cmap (str, optional): Colormap specifier. Defaults to 'coolwarm'. + levels (np.array, optional): Contour levels for line contour plot. + Defaults to None. + clevels (np.array, optional): Contour levels for tricontourf plot. + Defaults to None. + color_bar (Boolean, optional): Flag to include a color bar on the plot. + Defaults to False. + title (str, optional): User-supplied title for the plot. Defaults to "". + **kwargs: Additional parameters to pass to line contour plot. Returns: im (:py:class:`matplotlib.plt.pcolormesh`): Image handle. From 929838f0467b1bb8dfef33fb0553f0d38895ed8c Mon Sep 17 00:00:00 2001 From: bayc Date: Thu, 4 May 2023 08:53:39 -0600 Subject: [PATCH 20/36] add separate levels parameter for contour plot for heteregeneous flows --- floris/tools/visualization.py | 25 +++++++++++++++++++++++-- 1 file changed, 23 insertions(+), 2 deletions(-) diff --git a/floris/tools/visualization.py b/floris/tools/visualization.py index f3fdaa93e..14c72c153 100644 --- a/floris/tools/visualization.py +++ b/floris/tools/visualization.py @@ -203,8 +203,10 @@ def visualize_cut_plane( Args: cut_plane (:py:class:`~.tools.cut_plane.CutPlane`): 2D plane through wind plant. - ax (:py:class:`matplotlib.pyplot.axes`): Figure axes. Defaults + ax (:py:class:`matplotlib.pyplot.axes`, optional): Figure axes. Defaults to None. + vel_component (str, optional): The velocity component that the cut plane is + perpendicular to. min_speed (float, optional): Minimum value of wind speed for contours. Defaults to None. max_speed (float, optional): Maximum value of wind speed for @@ -298,25 +300,39 @@ def visualize_heterogeneous_cut_plane( max_speed=None, cmap="coolwarm", levels=None, + clevels=None, color_bar=False, title="", plot_het_bounds=True, **kwargs ): """ - Generate pseudocolor mesh plot of the cut_plane. + Generate pseudocolor mesh plot of the heterogeneous cut_plane. Args: cut_plane (:py:class:`~.tools.cut_plane.CutPlane`): 2D plane through wind plant. + fi (:py:class:`~.tools.floris_interface.FlorisInterface`): FlorisInterface object. ax (:py:class:`matplotlib.pyplot.axes`): Figure axes. Defaults to None. + vel_component (str, optional): The velocity component that the cut plane is + perpendicular to. min_speed (float, optional): Minimum value of wind speed for contours. Defaults to None. max_speed (float, optional): Maximum value of wind speed for contours. Defaults to None. cmap (str, optional): Colormap specifier. Defaults to 'coolwarm'. + levels (np.array, optional): Contour levels for line contour plot. + Defaults to None. + clevels (np.array, optional): Contour levels for tricontourf plot. + Defaults to None. + color_bar (Boolean, optional): Flag to include a color bar on the plot. + Defaults to False. + title (str, optional): User-supplied title for the plot. Defaults to "". + plot_het_bonds (boolean, optional): Flag to include the user-defined bounds of the + heterogeneous wind speed area. Defaults to True. + **kwargs: Additional parameters to pass to line contour plot. Returns: im (:py:class:`matplotlib.plt.pcolormesh`): Image handle. @@ -343,6 +359,10 @@ def visualize_heterogeneous_cut_plane( if max_speed is None: max_speed = cut_plane.df.w.max() + # Allow separate number of levels for tricontourf and for line_contour + if clevels is None: + clevels = levels + # Plot the cut-through im = ax.tricontourf( cut_plane.df.x1, @@ -350,6 +370,7 @@ def visualize_heterogeneous_cut_plane( cut_plane.df.u, vmin=min_speed, vmax=max_speed, + levels=clevels, cmap=cmap, extend="both", ) From ec83fc590311174046d86a542b7353e697c89ea5 Mon Sep 17 00:00:00 2001 From: bayc Date: Thu, 4 May 2023 09:49:11 -0600 Subject: [PATCH 21/36] fixing typo in het_config parameter --- floris/tools/uncertainty_interface.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/floris/tools/uncertainty_interface.py b/floris/tools/uncertainty_interface.py index 877aa56f3..ad892f51d 100644 --- a/floris/tools/uncertainty_interface.py +++ b/floris/tools/uncertainty_interface.py @@ -26,7 +26,7 @@ class UncertaintyInterface(LoggerBase): def __init__( self, configuration, - het_conifg=None, + het_config=None, unc_options=None, unc_pmfs=None, fix_yaw_in_relative_frame=False, @@ -110,7 +110,7 @@ def __init__( if isinstance(configuration, FlorisInterface): self.fi = configuration else: - self.fi = FlorisInterface(configuration, het_config=het_conifg) + self.fi = FlorisInterface(configuration, het_config=het_config) self.reinitialize_uncertainty( unc_options=unc_options, From cb36a34e1d4e1d9acba2f2c3e9b089c6a46bb4b1 Mon Sep 17 00:00:00 2001 From: bayc Date: Thu, 4 May 2023 10:30:45 -0600 Subject: [PATCH 22/36] store the center of rotation in the grid classes --- floris/simulation/grid.py | 20 +++++++++++--------- 1 file changed, 11 insertions(+), 9 deletions(-) diff --git a/floris/simulation/grid.py b/floris/simulation/grid.py index e30e987dc..472cd0f6b 100644 --- a/floris/simulation/grid.py +++ b/floris/simulation/grid.py @@ -83,6 +83,8 @@ class Grid(ABC): x_sorted_inertial_frame: NDArrayFloat = field(init=False) y_sorted_inertial_frame: NDArrayFloat = field(init=False) z_sorted_inertial_frame: NDArrayFloat = field(init=False) + xc_rot: NDArrayFloat = field(init=False) + yc_rot: NDArrayFloat = field(init=False) def __attrs_post_init__(self) -> None: self.turbine_coordinates_array = np.array([c.elements for c in self.turbine_coordinates]) @@ -210,7 +212,7 @@ def set_grid(self) -> None: # the foot of the turbine where the tower meets the ground. # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west( + x, y, z, self.xc_rot, self.yc_rot = rotate_coordinates_rel_west( self.wind_directions, self.turbine_coordinates_array, ) @@ -281,8 +283,8 @@ def set_grid(self) -> None: grid_x=self.x_sorted, grid_y=self.y_sorted, grid_z=self.z_sorted, - x_center_of_rotation=xc_rot, - y_center_of_rotation=yc_rot, + x_center_of_rotation=self.xc_rot, + y_center_of_rotation=self.yc_rot, ) @define @@ -316,7 +318,7 @@ def set_grid(self) -> None: """ # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west( + x, y, z, self.xc_rot, self.yc_rot = rotate_coordinates_rel_west( self.wind_directions, self.turbine_coordinates_array, ) @@ -348,8 +350,8 @@ def set_grid(self) -> None: grid_x=self.x_sorted, grid_y=self.y_sorted, grid_z=self.z_sorted, - x_center_of_rotation=xc_rot, - y_center_of_rotation=yc_rot, + x_center_of_rotation=self.xc_rot, + y_center_of_rotation=self.yc_rot, ) @define @@ -389,7 +391,7 @@ def set_grid(self) -> None: Then, create the grid based on this wind-from-left orientation """ # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z, xc_rot, yc_rot = rotate_coordinates_rel_west( + x, y, z, self.xc_rot, self.yc_rot = rotate_coordinates_rel_west( self.wind_directions, self.turbine_coordinates_array, ) @@ -461,8 +463,8 @@ def set_grid(self) -> None: grid_x=self.x_sorted, grid_y=self.y_sorted, grid_z=self.z_sorted, - x_center_of_rotation=xc_rot, - y_center_of_rotation=yc_rot, + x_center_of_rotation=self.xc_rot, + y_center_of_rotation=self.yc_rot, ) # self.sorted_indices = self.x.argsort(axis=2) From 2d12896c151c51a3146aff35f0efed0ce21a31f0 Mon Sep 17 00:00:00 2001 From: bayc Date: Thu, 4 May 2023 15:33:23 -0600 Subject: [PATCH 23/36] add warning for users if interpolated points are calculated outside user-defined heterogeneous bounds --- floris/simulation/flow_field.py | 27 ++++++++++++++++++++++++--- 1 file changed, 24 insertions(+), 3 deletions(-) diff --git a/floris/simulation/flow_field.py b/floris/simulation/flow_field.py index b78706982..e9697cded 100644 --- a/floris/simulation/flow_field.py +++ b/floris/simulation/flow_field.py @@ -15,19 +15,22 @@ from __future__ import annotations import attrs +import matplotlib.path as mpltPath import numpy as np from attrs import define, field -from floris.simulation import Grid +from floris.simulation import ( + BaseClass, + Grid, +) from floris.type_dec import ( floris_array_converter, - FromDictMixin, NDArrayFloat, ) @define -class FlowField(FromDictMixin): +class FlowField(BaseClass): wind_speeds: NDArrayFloat = field(converter=floris_array_converter) wind_directions: NDArrayFloat = field(converter=floris_array_converter) wind_veer: float = field(converter=float) @@ -95,6 +98,24 @@ def initialize_velocity_field(self, grid: Grid) -> None: # If heterogeneous flow data is given, the speed ups at the defined # grid locations are determined in either 2 or 3 dimensions. else: + + path = mpltPath.Path(self.het_map[0].points) + points = np.column_stack( + ( + grid.x_sorted_inertial_frame.flatten(), + grid.y_sorted_inertial_frame.flatten(), + ) + ) + inside = path.contains_points(points) + if not np.all(inside): + self.logger.warning( + "The calculated flow field contains points outside of the the user-defined " + "heterogeneous inflow bounds. For these points, the interpolated value has " + "been filled with the freestream wind speed. If this is not the desired " + "behavior, the user will need to expand the heterogeneous inflow bounds to " + "fully cover the calculated flow field area." + ) + if len(self.het_map[0].points[0]) == 2: speed_ups = self.calculate_speed_ups( self.het_map, From 148f3db9241fa1ed8410818cb0f8fb990d885aa4 Mon Sep 17 00:00:00 2001 From: bayc Date: Thu, 4 May 2023 21:41:39 -0600 Subject: [PATCH 24/36] moving het_config from FlorisInterface to FlowField --- floris/simulation/flow_field.py | 72 ++++++++++++++++++++++++++------ floris/tools/floris_interface.py | 37 +--------------- floris/tools/visualization.py | 9 +++- 3 files changed, 69 insertions(+), 49 deletions(-) diff --git a/floris/simulation/flow_field.py b/floris/simulation/flow_field.py index e9697cded..e3796f5ec 100644 --- a/floris/simulation/flow_field.py +++ b/floris/simulation/flow_field.py @@ -18,6 +18,7 @@ import matplotlib.path as mpltPath import numpy as np from attrs import define, field +from scipy.interpolate import LinearNDInterpolator from floris.simulation import ( BaseClass, @@ -39,6 +40,7 @@ class FlowField(BaseClass): turbulence_intensity: float = field(converter=float) reference_wind_height: float = field(converter=float) time_series : bool = field(default=False) + het_config: dict = field(default=None) n_wind_speeds: int = field(init=False) n_wind_directions: int = field(init=False) @@ -72,6 +74,13 @@ def wind_directions_validator(self, instance: attrs.Attribute, value: NDArrayFlo """Using the validator method to keep the `n_wind_directions` attribute up to date.""" self.n_wind_directions = value.size + + + def __attrs_post_init__(self) -> None: + if self.het_config is not None: + self.generate_heterogeneous_wind_map() + + def initialize_velocity_field(self, grid: Grid) -> None: # Create an initial wind profile as a function of height. The values here will @@ -193,27 +202,64 @@ def finalize(self, unsorted_indices): ) def calculate_speed_ups(self, het_map, x, y, z=None): - - # Check that the het maps wd dimension matches - if self.n_wind_directions!= np.array(het_map).shape[0]: - raise ValueError( - "het_map's wind direction dimension not equal to number of wind directions" - ) - if z is not None: - # Calculate the 3-dimensional speed ups; reshape is needed as the generator + # Calculate the 3-dimensional speed ups; squeeze is needed as the generator # adds an extra dimension - speed_ups = np.reshape( + speed_ups = np.squeeze( [het_map[i](x[i:i+1], y[i:i+1], z[i:i+1]) for i in range( len(het_map))], - np.shape(x) + axis=1, ) else: - # Calculate the 2-dimensional speed ups; reshape is needed as the generator + # Calculate the 2-dimensional speed ups; squeeze is needed as the generator # adds an extra dimension - speed_ups = np.reshape( + speed_ups = np.squeeze( [het_map[i](x[i:i+1], y[i:i+1]) for i in range(len(het_map))], - np.shape(x) + axis=1, ) return speed_ups + + def generate_heterogeneous_wind_map(self): + """This function creates the heterogeneous interpolant used to calculate heterogeneous + inflows. The interpolant is for computing wind speed based on an x and y location in the + flow field. This is computed using Scipy's LinearNDInterpolator and uses a fill value + equal to the freestream for interpolated values outside of the user-defined heterogeneous + map bounds. + + Args: + het_config (dict): The heterogeneous inflow configuration dictionary. + The configuration should have the following inputs specified. + - **speed_ups** (list): A list of speed up factors that will multiply the specified + freestream wind speed. This 2-dimensional array should have an array of + multiplicative factors defined for each wind direction. + - **x_locs** (list): A list of x locations at which the speed up factors are + defined. + - **y_locs**: A list of y locations at which the speed up factors are + defined. + - **z_locs** (optional): A list of z locations at which the speed up factors are + defined. + """ + speed_ups = self.het_config['speed_ups'] + x = self.het_config['x_locs'] + y = self.het_config['y_locs'] + z = self.het_config['z_locs'] + + if z is not None: + # Compute the 3-dimensional interpolants for each wind diretion + # 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)), speed_up, fill_value=1.0) + for speed_up in speed_ups + ] + else: + # Compute the 2-dimensional interpolants for each wind diretion + # 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)), speed_up, fill_value=1.0) + for speed_up in speed_ups + ] + + self.het_map = in_region diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index 4627e12c7..716bb96fc 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -48,20 +48,9 @@ class FlorisInterface(LoggerBase): - **turbine**: See `floris.simulation.turbine.Turbine` for more details. - **wake**: See `floris.simulation.wake.WakeManager` for more details. - **logging**: See `floris.simulation.floris.Floris` for more details. - het_config (:py:obj:`dict`, optional): The heterogeneous inflow configuration dictionary. - The configuration should have the following inputs specified. - - **speed_ups**: See `floris.tools.floris_interface.generate_heterogeneous_wind_map' - for more details. - - **x_locs**: See `floris.tools.floris_interface.generate_heterogeneous_wind_map' - for more details. - - **y_locs**: See `floris.tools.floris_interface.generate_heterogeneous_wind_map' - for more details. - - **z_locs** (optional): See - `floris.tools.floris_interface.generate_heterogeneous_wind_map' - for more details. """ - def __init__(self, configuration: dict | str | Path, het_config: dict | None = None): + def __init__(self, configuration: dict | str | Path): self.configuration = configuration if isinstance(self.configuration, (str, Path)): @@ -73,18 +62,6 @@ def __init__(self, configuration: dict | str | Path, het_config: dict | None = N else: raise TypeError("The Floris `configuration` must be of type 'dict', 'str', or 'Path'.") - # Create the heterogeneous wind map interpolant if het_config is supplied - if het_config is not None: - het_map = self.generate_heterogeneous_wind_map(het_config=het_config) - else: - het_map = None - # Store the heterogeneous info and map for use after reinitailization - self.het_config = het_config - self.het_map = het_map - # Assign the heterogeneous map to the flow field - # Needed for a direct call to fi.calculate_wake without fi.reinitialize - self.floris.flow_field.het_map = het_map - # If ref height is -1, assign the hub height if np.abs(self.floris.flow_field.reference_wind_height + 1.0) < 1.0e-6: self.assign_hub_height_to_ref_height() @@ -256,7 +233,7 @@ def reinitialize( if air_density is not None: flow_field_dict["air_density"] = air_density if het_config is not None: - self.het_config = het_config + flow_field_dict["het_config"] = het_config ## Farm if layout_x is not None: @@ -285,13 +262,6 @@ def reinitialize( # Create a new instance of floris and attach to self self.floris = Floris.from_dict(floris_dict) - # Re-assign the hetergeneous inflow map to flow field - if self.het_config is not None: - het_map = self.generate_heterogeneous_wind_map(het_config=self.het_config) - else: - het_map = None - self.het_map = het_map - self.floris.flow_field.het_map = het_map def get_plane_of_points( self, @@ -451,7 +421,6 @@ def calculate_horizontal_plane( # Reset the fi object back to the turbine grid configuration self.floris = Floris.from_dict(floris_dict) - self.floris.flow_field.het_map = self.het_map # Run the simulation again for futher postprocessing (i.e. now we can get farm power) self.calculate_wake(yaw_angles=current_yaw_angles) @@ -530,7 +499,6 @@ def calculate_cross_plane( # Reset the fi object back to the turbine grid configuration self.floris = Floris.from_dict(floris_dict) - self.floris.flow_field.het_map = self.het_map # Run the simulation again for futher postprocessing (i.e. now we can get farm power) self.calculate_wake(yaw_angles=current_yaw_angles) @@ -609,7 +577,6 @@ def calculate_y_plane( # Reset the fi object back to the turbine grid configuration self.floris = Floris.from_dict(floris_dict) - self.floris.flow_field.het_map = self.het_map # Run the simulation again for futher postprocessing (i.e. now we can get farm power) self.calculate_wake(yaw_angles=current_yaw_angles) diff --git a/floris/tools/visualization.py b/floris/tools/visualization.py index 14c72c153..83aa8892d 100644 --- a/floris/tools/visualization.py +++ b/floris/tools/visualization.py @@ -388,7 +388,14 @@ def visualize_heterogeneous_cut_plane( # Plot the user-defined heterogeneous flow area if plot_het_bounds: - points = np.array(list(zip(fi.het_config['x_locs'], fi.het_config['y_locs']))) + points = np.array( + list( + zip( + fi.floris.flow_field.het_config['x_locs'], + fi.floris.flow_field.het_config['y_locs'], + ) + ) + ) hull = ConvexHull(points) h = ax.plot( points[np.append(hull.vertices, hull.vertices[0]),0], From b8646784bd63caccc139083947306997b88392e9 Mon Sep 17 00:00:00 2001 From: bayc Date: Thu, 4 May 2023 21:44:40 -0600 Subject: [PATCH 25/36] moving some input checking to validators --- floris/simulation/flow_field.py | 50 ++++++++++++++++++++++++++++++-- floris/tools/floris_interface.py | 11 ------- 2 files changed, 48 insertions(+), 13 deletions(-) diff --git a/floris/simulation/flow_field.py b/floris/simulation/flow_field.py index e3796f5ec..7319c0a07 100644 --- a/floris/simulation/flow_field.py +++ b/floris/simulation/flow_field.py @@ -19,6 +19,8 @@ import numpy as np from attrs import define, field from scipy.interpolate import LinearNDInterpolator +from scipy.spatial import ConvexHull +from shapely.geometry import Polygon from floris.simulation import ( BaseClass, @@ -74,6 +76,41 @@ def wind_directions_validator(self, instance: attrs.Attribute, value: NDArrayFlo """Using the validator method to keep the `n_wind_directions` attribute up to date.""" self.n_wind_directions = value.size + @het_config.validator + def het_config_validator(self, instance: attrs.Attribute, value: dict | None) -> None: + """Using the validator method to check that the het_config dictionary has the correct + key-value pairs. + """ + if value is not None: + # Check that the correct keys are supplied for the het_config dict + for k in ["speed_ups", "x_locs", "y_locs"]: + if k not in value.keys(): + raise ValueError( + "het_config must contain entries for 'speeds_ups', 'x_locs', and 'y_locs'," + f" with 'z_locs' as optional. Missing '{k}'." + ) + if "z_locs" not in value: + # If only a 2D case, add "None" for the z locations + value["z_locs"] = None + + @het_map.validator + def het_map_validator(self, instance: attrs.Attribute, value: list | None) -> None: + """Using this validator to make sure that the het_map has an interpolant defined for + each wind direction. + """ + if value is not None: + if self.n_wind_directions!= np.array(value).shape[0]: + if np.array(value).shape[0]: + self.logger.warning( + "Multiple wind directions defined with just one hetergeneous speed up " + "map. The same speed up values are being used for all the wind directions." + ) + else: + raise ValueError( + "The het_map's wind direction dimension not equal to number of wind " + "directions. Either the dimensions must be equal, or just one set up " + "of speed ups must be defined to be used across all wind directions." + ) def __attrs_post_init__(self) -> None: @@ -107,8 +144,17 @@ def initialize_velocity_field(self, grid: Grid) -> None: # If heterogeneous flow data is given, the speed ups at the defined # grid locations are determined in either 2 or 3 dimensions. else: - - path = mpltPath.Path(self.het_map[0].points) + bounds = np.array( + list( + zip( + self.het_config['x_locs'], + self.het_config['y_locs'], + ) + ) + ) + hull = ConvexHull(bounds) + polygon = Polygon(bounds[hull.vertices]) + path = mpltPath.Path(polygon.boundary.coords) points = np.column_stack( ( grid.x_sorted_inertial_frame.flatten(), diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index 716bb96fc..2f215e1bc 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -966,17 +966,6 @@ def generate_heterogeneous_wind_map(self, het_config): LinearNDInterpolator and uses a fill value equal to the freestream for interpolated values outside of the user-defined heterogeneous map bounds. """ - # Check that the correct keys are supplied for the het_config dict - for k in ["speed_ups", "x_locs", "y_locs"]: - if k not in het_config.keys(): - raise ValueError( - "het_config must contain entries for 'speeds_ups', 'x_locs', and 'y_locs'," - f" with 'z_locs' as optional. Missing '{k}'." - ) - if "z_locs" not in het_config: - # If only a 2D case, add "None" for the z locations - het_config["z_locs"] = None - speed_ups = het_config['speed_ups'] x = het_config['x_locs'] y = het_config['y_locs'] From afb5b23d02c4a0f73d452dffb3c4e95724f8249e Mon Sep 17 00:00:00 2001 From: bayc Date: Thu, 4 May 2023 22:17:44 -0600 Subject: [PATCH 26/36] removing attempt at allowing for broadcasting of a single speed up array across multiple wind directions --- floris/simulation/flow_field.py | 14 +++----------- 1 file changed, 3 insertions(+), 11 deletions(-) diff --git a/floris/simulation/flow_field.py b/floris/simulation/flow_field.py index 7319c0a07..19f60ee05 100644 --- a/floris/simulation/flow_field.py +++ b/floris/simulation/flow_field.py @@ -100,17 +100,9 @@ def het_map_validator(self, instance: attrs.Attribute, value: list | None) -> No """ if value is not None: if self.n_wind_directions!= np.array(value).shape[0]: - if np.array(value).shape[0]: - self.logger.warning( - "Multiple wind directions defined with just one hetergeneous speed up " - "map. The same speed up values are being used for all the wind directions." - ) - else: - raise ValueError( - "The het_map's wind direction dimension not equal to number of wind " - "directions. Either the dimensions must be equal, or just one set up " - "of speed ups must be defined to be used across all wind directions." - ) + raise ValueError( + "The het_map's wind direction dimension not equal to number of wind directions." + ) def __attrs_post_init__(self) -> None: From 2323469782e6274f16b7d731f2f0d2e0f4056040 Mon Sep 17 00:00:00 2001 From: bayc Date: Thu, 4 May 2023 22:18:16 -0600 Subject: [PATCH 27/36] removing old het_config definitions --- floris/tools/floris_interface.py | 56 +------------------- floris/tools/parallel_computing_interface.py | 4 +- floris/tools/uncertainty_interface.py | 3 +- 3 files changed, 3 insertions(+), 60 deletions(-) diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index 2f215e1bc..c829b8358 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -105,7 +105,7 @@ def assign_hub_height_to_ref_height(self): def copy(self): """Create an independent copy of the current FlorisInterface object""" - return FlorisInterface(self.floris.as_dict(), het_config=self.het_config) + return FlorisInterface(self.floris.as_dict()) def calculate_wake( self, @@ -937,60 +937,6 @@ def get_farm_AEP_wind_rose_class( return aep - - def generate_heterogeneous_wind_map(self, het_config): - """This function creates the heterogeneous interpolant used to calculate heterogeneous - inflows. - - Args: - het_config (dict): The heterogeneous inflow configuration dictionary. - The configuration should have the following inputs specified. - - **speed_ups** (list): A list of speed up factors that will multiply the specified - freestream wind speed. This 2-dimensional array should have an array of - multiplicative factors defined for each wind direction. - - **x_locs** (list): A list of x locations at which the speed up factors are - defined. - - **y_locs**: A list of y locations at which the speed up factors are - defined. - - **z_locs** (optional): A list of z locations at which the speed up factors are - defined. - - Raises: - ValueError: ("het_config must contain entries for 'speeds_ups', 'x_locs', and 'y_locs'," - f" with 'z_locs' as optional. Missing '{k}'.") raised if the configuration - dictionary is missing one of the required key-value pairs. - - Returns: - LinearNDInterpolator: Returns a linear interpolant for computing wind speed based on an - x and y location in the flow field. This is computed using Scipy's - LinearNDInterpolator and uses a fill value equal to the freestream for interpolated - values outside of the user-defined heterogeneous map bounds. - """ - speed_ups = het_config['speed_ups'] - x = het_config['x_locs'] - y = het_config['y_locs'] - z = het_config['z_locs'] - - if z is not None: - # Compute the 3-dimensional interpolants for each wind diretion - # 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)), speed_up, fill_value=1.0) - for speed_up in speed_ups - ] - else: - # Compute the 2-dimensional interpolants for each wind diretion - # 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)), speed_up, fill_value=1.0) - for speed_up in speed_ups - ] - - return in_region - - @property def layout_x(self): """ diff --git a/floris/tools/parallel_computing_interface.py b/floris/tools/parallel_computing_interface.py index c9f9923c4..6dd7684a3 100644 --- a/floris/tools/parallel_computing_interface.py +++ b/floris/tools/parallel_computing_interface.py @@ -12,17 +12,15 @@ def _load_local_floris_object( fi_dict, - het_config=None, unc_pmfs=None, fix_yaw_in_relative_frame=False ): # Load local FLORIS object if unc_pmfs is None: - fi = FlorisInterface(fi_dict, het_config=het_config) + fi = FlorisInterface(fi_dict) else: fi = UncertaintyInterface( fi_dict, - het_config=het_config, unc_pmfs=unc_pmfs, fix_yaw_in_relative_frame=fix_yaw_in_relative_frame, ) diff --git a/floris/tools/uncertainty_interface.py b/floris/tools/uncertainty_interface.py index 2813f5df6..747ac16ce 100644 --- a/floris/tools/uncertainty_interface.py +++ b/floris/tools/uncertainty_interface.py @@ -26,7 +26,6 @@ class UncertaintyInterface(LoggerBase): def __init__( self, configuration, - het_config=None, unc_options=None, unc_pmfs=None, fix_yaw_in_relative_frame=False, @@ -110,7 +109,7 @@ def __init__( if isinstance(configuration, FlorisInterface): self.fi = configuration else: - self.fi = FlorisInterface(configuration, het_config=het_config) + self.fi = FlorisInterface(configuration) self.reinitialize_uncertainty( unc_options=unc_options, From 8cfde459cc8b8b9f3e3f7566b072665fd1071307 Mon Sep 17 00:00:00 2001 From: bayc Date: Thu, 4 May 2023 22:18:37 -0600 Subject: [PATCH 28/36] updating the heterogeneous examples --- examples/16_heterogeneous_inflow.py | 33 +++--- examples/16b_heterogenaity_multiple_ws_wd.py | 12 +-- examples/inputs/gch_heterogeneous_inflow.yaml | 102 ++++++++++++++++++ 3 files changed, 123 insertions(+), 24 deletions(-) create mode 100644 examples/inputs/gch_heterogeneous_inflow.yaml diff --git a/examples/16_heterogeneous_inflow.py b/examples/16_heterogeneous_inflow.py index 1f84d7efe..54f7facd4 100644 --- a/examples/16_heterogeneous_inflow.py +++ b/examples/16_heterogeneous_inflow.py @@ -36,23 +36,18 @@ """ -# Define the speed ups of the heterogeneous inflow, and their locations. -# For the 2-dimensional case, this requires x and y locations. -# The speed ups are multipliers of the ambient wind speed. -speed_ups = [[2.0, 2.0, 1.0, 1.0]] -x_locs = [-300.0, -300.0, 2600.0, 2600.0] -y_locs = [ -300.0, 300.0, -300.0, 300.0] - -# Create the configuration dictionary to be used for the heterogeneous inflow. -het_config_2d = { - 'speed_ups': speed_ups, - 'x_locs': x_locs, - 'y_locs': y_locs, -} - # Initialize FLORIS with the given input file via FlorisInterface. -# Also, pass the heterogeneous map into the FlorisInterface. -fi_2d = FlorisInterface("inputs/gch.yaml", het_config=het_config_2d) +# Note that the hetergeneous flow is defined in the input file. The het_config +# dictionary is defined as below. The speed ups are multipliers of the ambient wind speed, +# and the x_locs and y_locs are the locations of the speed ups. +# +# het_config = { +# 'speed_ups': [[2.0, 1.0, 2.0, 1.0]], +# 'x_locs': [-300.0, -300.0, 2600.0, 2600.0], +# 'y_locs': [ -300.0, 300.0, -300.0, 300.0], +# } + +fi_2d = FlorisInterface("inputs/gch_heterogeneous_inflow.yaml") # Set shear to 0.0 to highlight the heterogeneous inflow fi_2d.reinitialize(wind_shear=0.0) @@ -104,8 +99,10 @@ } # Initialize FLORIS with the given input file via FlorisInterface. -# Also, pass the heterogeneous map into the FlorisInterface. -fi_3d = FlorisInterface("inputs/gch.yaml", het_config=het_config_3d) +# Note that we initialize FLORIS with a homogenous flow input file, but +# then configure the heterogeneous inflow via the reinitialize method. +fi_3d = FlorisInterface("inputs/gch.yaml") +fi_3d.reinitialize(het_config=het_config_3d) # Set shear to 0.0 to highlight the heterogeneous inflow fi_3d.reinitialize(wind_shear=0.0) diff --git a/examples/16b_heterogenaity_multiple_ws_wd.py b/examples/16b_heterogenaity_multiple_ws_wd.py index a35a73aa7..c4822df72 100644 --- a/examples/16b_heterogenaity_multiple_ws_wd.py +++ b/examples/16b_heterogenaity_multiple_ws_wd.py @@ -35,14 +35,14 @@ # Create the configuration dictionary to be used for the heterogeneous inflow. het_config_2d = { - 'speed_ups': speed_ups, - 'x_locs': x_locs, - 'y_locs': y_locs, + 'speed_ups': [[2.0, 1.0, 2.0, 1.0]], + 'x_locs': [-300.0, -300.0, 2600.0, 2600.0], + 'y_locs': [ -300.0, 300.0, -300.0, 300.0], } # Initialize FLORIS with the given input file via FlorisInterface. -# Also, pass the heterogeneous map into the FlorisInterface. -fi = FlorisInterface("inputs/gch.yaml", het_config=het_config_2d) +# Note the heterogeneous inflow is defined in the input file. +fi = FlorisInterface("inputs/gch_heterogeneous_inflow.yaml") # Set shear to 0.0 to highlight the heterogeneous inflow fi.reinitialize( @@ -50,7 +50,7 @@ wind_speeds=[8.0], wind_directions=[270.], layout_x=[0, 0], - layout_y=[-300, 300], + layout_y=[-299., 299.], ) fi.calculate_wake() turbine_powers = fi.get_turbine_powers().flatten() / 1000. diff --git a/examples/inputs/gch_heterogeneous_inflow.yaml b/examples/inputs/gch_heterogeneous_inflow.yaml new file mode 100644 index 000000000..c0458c851 --- /dev/null +++ b/examples/inputs/gch_heterogeneous_inflow.yaml @@ -0,0 +1,102 @@ +name: GCH +description: Three turbines using Gauss Curl Hybrid model +floris_version: v3.0.0 + +logging: + console: + enable: true + level: WARNING + file: + enable: false + level: WARNING + +solver: + type: turbine_grid + turbine_grid_points: 1 +farm: + layout_x: + - 0.0 + - 630.0 + - 1260.0 + layout_y: + - 0.0 + - 0.0 + - 0.0 + turbine_type: + - nrel_5MW + +flow_field: + air_density: 1.225 + het_config: + speed_ups: + - - 2.0 + - 1.0 + - 2.0 + - 1.0 + x_locs: + - -300. + - -300. + - 2600. + - 2600. + y_locs: + - -300. + - 300. + - -300. + - 300. + reference_wind_height: -1 + turbulence_intensity: 0.06 + wind_directions: + - 270.0 + wind_shear: 0.12 + wind_speeds: + - 8.0 + wind_veer: 0.0 + +wake: + model_strings: + combination_model: sosfs + deflection_model: gauss + turbulence_model: crespo_hernandez + velocity_model: gauss + enable_secondary_steering: true + enable_yaw_added_recovery: true + enable_transverse_velocities: true + + wake_deflection_parameters: + gauss: + ad: 0.0 + alpha: 0.58 + bd: 0.0 + beta: 0.077 + dm: 1.0 + ka: 0.38 + kb: 0.004 + jimenez: + ad: 0.0 + bd: 0.0 + kd: 0.05 + + wake_velocity_parameters: + cc: + a_s: 0.179367259 + b_s: 0.0118889215 + c_s1: 0.0563691592 + c_s2: 0.13290157 + a_f: 3.11 + b_f: -0.68 + c_f: 2.41 + alpha_mod: 1.0 + gauss: + alpha: 0.58 + beta: 0.077 + ka: 0.38 + kb: 0.004 + jensen: + we: 0.05 + + wake_turbulence_parameters: + crespo_hernandez: + initial: 0.1 + constant: 0.5 + ai: 0.8 + downstream: -0.32 From cbd1de5c570499e863daf4f67417a0416da73a2f Mon Sep 17 00:00:00 2001 From: bayc Date: Thu, 4 May 2023 22:35:25 -0600 Subject: [PATCH 29/36] removing one more outdated het_map reference --- floris/tools/parallel_computing_interface.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/floris/tools/parallel_computing_interface.py b/floris/tools/parallel_computing_interface.py index 6dd7684a3..600228a87 100644 --- a/floris/tools/parallel_computing_interface.py +++ b/floris/tools/parallel_computing_interface.py @@ -215,7 +215,7 @@ def _preprocessing(self, yaw_angles=None): # Prepare lightweight data to pass along if isinstance(self.fi, FlorisInterface): - fi_information = (fi_dict_split, self.fi.het_map, None, None) + fi_information = (fi_dict_split, None, None) else: fi_information = ( fi_dict_split, From c5e736e275ef3294e1adf377d4df29e5440de8f7 Mon Sep 17 00:00:00 2001 From: bayc Date: Fri, 5 May 2023 13:28:20 -0600 Subject: [PATCH 30/36] updating doctrings and comments --- floris/tools/visualization.py | 2 +- floris/utilities.py | 27 ++++++++++++++++++++++++--- 2 files changed, 25 insertions(+), 4 deletions(-) diff --git a/floris/tools/visualization.py b/floris/tools/visualization.py index 83aa8892d..47fbcf72d 100644 --- a/floris/tools/visualization.py +++ b/floris/tools/visualization.py @@ -125,7 +125,7 @@ def add_turbine_id_labels(fi: FlorisInterface, ax: plt.Axes, **kwargs): ax (plt.Axes): Axes object to add the labels. """ - # Rotate layout to intertial frame for plotting turbines relative to wind direction + # Rotate layout to inertial frame for plotting turbines relative to wind direction coordinates_array = np.array([ [x, y, 0.0] for x, y in list(zip(fi.layout_x, fi.layout_y)) diff --git a/floris/utilities.py b/floris/utilities.py index ce41648c5..583ff7659 100644 --- a/floris/utilities.py +++ b/floris/utilities.py @@ -229,9 +229,13 @@ def rotate_coordinates_rel_west( than North (0). The rotation happens about the centroid of the coordinates. Args: - wind_directions (NDArrayFloat): Series of wind directions to base the rotation + wind_directions (NDArrayFloat): Series of wind directions to base the rotation. coordinates (NDArrayFloat): Series of coordinates to rotate with shape (N coordinates, 3) - so that each element of the array coordinates[i] yields a three-component coordinate + so that each element of the array coordinates[i] yields a three-component coordinate. + x_center_of_rotation (float, optional): The x-coordinate for the rotation center of the + input coordinates. Defaults to None. + y_center_of_rotation (float, optional): The y-coordinate for the rotational center of the + input coordinates. Defaults to None. """ # Calculate the difference in given wind direction from 270 / West @@ -267,11 +271,28 @@ def rotate_coordinates_rel_west( def reverse_rotate_coordinates_rel_west( wind_directions, - grid_x, grid_y, + grid_x, + grid_y, grid_z, x_center_of_rotation, y_center_of_rotation ): + """ + This function reverses the rotation of the given grid so that the coordinates are aligned with + the original wind direction. The rotation happens about the centroid of the coordinates. + + Args: + wind_directions (NDArrayFloat): Series of wind directions to base the rotation. + coordinates (NDArrayFloat): Series of coordinates to rotate with shape (N coordinates, 3) + so that each element of the array coordinates[i] yields a three-component coordinate. + grid_x (NDArrayFloat): X-coordinates to be rotated. + grid_y (NDArrayFloat): Y-coordinates to be rotated. + grid_z (NDArrayFloat): Z-coordinates to be rotated. + x_center_of_rotation (float): The x-coordinate for the rotation center of the + input coordinates. + y_center_of_rotation (float): The y-coordinate for the rotational center of the + input coordinates. + """ # Calculate the difference in given wind direction from 270 / West # We are rotating in the other direction wind_deviation_from_west = -1.0 * wind_delta(wind_directions) From 331853f5daa36d5ab553e3cb669fe02afcb54830 Mon Sep 17 00:00:00 2001 From: bayc Date: Thu, 11 May 2023 19:00:28 -0600 Subject: [PATCH 31/36] fixing naming of x and y center of rotations --- examples/02_visualizations.py | 2 +- floris/simulation/grid.py | 22 +++++++++++----------- 2 files changed, 12 insertions(+), 12 deletions(-) diff --git a/examples/02_visualizations.py b/examples/02_visualizations.py index c812bd95d..669e91fa0 100644 --- a/examples/02_visualizations.py +++ b/examples/02_visualizations.py @@ -65,7 +65,7 @@ y_plane = fi.calculate_y_plane( x_resolution=200, z_resolution=100, - crossstream_dist=630.0, + crossstream_dist=0.0, yaw_angles=np.array([[[25.,0.,0.]]]), ) cross_plane = fi.calculate_cross_plane( diff --git a/floris/simulation/grid.py b/floris/simulation/grid.py index c4e33d832..a72b86db1 100644 --- a/floris/simulation/grid.py +++ b/floris/simulation/grid.py @@ -83,8 +83,8 @@ class Grid(ABC): x_sorted_inertial_frame: NDArrayFloat = field(init=False) y_sorted_inertial_frame: NDArrayFloat = field(init=False) z_sorted_inertial_frame: NDArrayFloat = field(init=False) - xc_rot: NDArrayFloat = field(init=False) - yc_rot: NDArrayFloat = field(init=False) + x_center_of_rotation: NDArrayFloat = field(init=False) + y_center_of_rotation: NDArrayFloat = field(init=False) def __attrs_post_init__(self) -> None: self.turbine_coordinates_array = np.array([c.elements for c in self.turbine_coordinates]) @@ -216,7 +216,7 @@ def set_grid(self) -> None: # the foot of the turbine where the tower meets the ground. # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z, self.xc_rot, self.yc_rot = rotate_coordinates_rel_west( + x, y, z, self.x_center_of_rotation, self.y_center_of_rotation = rotate_coordinates_rel_west( self.wind_directions, self.turbine_coordinates_array, ) @@ -287,8 +287,8 @@ def set_grid(self) -> None: grid_x=self.x_sorted, grid_y=self.y_sorted, grid_z=self.z_sorted, - x_center_of_rotation=self.xc_rot, - y_center_of_rotation=self.yc_rot, + x_center_of_rotation=self.x_center_of_rotation, + y_center_of_rotation=self.y_center_of_rotation, ) @define @@ -322,7 +322,7 @@ def set_grid(self) -> None: """ # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z, self.xc_rot, self.yc_rot = rotate_coordinates_rel_west( + x, y, z, self.x_center_of_rotation, self.y_center_of_rotation = rotate_coordinates_rel_west( self.wind_directions, self.turbine_coordinates_array, ) @@ -354,8 +354,8 @@ def set_grid(self) -> None: grid_x=self.x_sorted, grid_y=self.y_sorted, grid_z=self.z_sorted, - x_center_of_rotation=self.xc_rot, - y_center_of_rotation=self.yc_rot, + x_center_of_rotation=self.x_center_of_rotation, + y_center_of_rotation=self.y_center_of_rotation, ) @define @@ -395,7 +395,7 @@ def set_grid(self) -> None: Then, create the grid based on this wind-from-left orientation """ # These are the rotated coordinates of the wind turbines based on the wind direction - x, y, z, self.xc_rot, self.yc_rot = rotate_coordinates_rel_west( + x, y, z, self.x_center_of_rotation, self.y_center_of_rotation = rotate_coordinates_rel_west( self.wind_directions, self.turbine_coordinates_array, ) @@ -467,8 +467,8 @@ def set_grid(self) -> None: grid_x=self.x_sorted, grid_y=self.y_sorted, grid_z=self.z_sorted, - x_center_of_rotation=self.xc_rot, - y_center_of_rotation=self.yc_rot, + x_center_of_rotation=self.x_center_of_rotation, + y_center_of_rotation=self.y_center_of_rotation, ) # self.sorted_indices = self.x.argsort(axis=2) From ac54f090dab79863e479441e45d07705bbb6b35a Mon Sep 17 00:00:00 2001 From: Rafael M Mudafort Date: Fri, 12 May 2023 13:44:40 -0500 Subject: [PATCH 32/36] Resolve merge conflicts and revert extra changes Reducing the number of lines that were changed here to make reviewing the PR easier --- floris/simulation/flow_field.py | 9 +-------- floris/simulation/grid.py | 13 +++++++------ floris/tools/floris_interface.py | 7 +------ floris/tools/visualization.py | 2 +- tests/utilities_unit_test.py | 6 +++--- 5 files changed, 13 insertions(+), 24 deletions(-) diff --git a/floris/simulation/flow_field.py b/floris/simulation/flow_field.py index 19f60ee05..e69337083 100644 --- a/floris/simulation/flow_field.py +++ b/floris/simulation/flow_field.py @@ -136,14 +136,7 @@ def initialize_velocity_field(self, grid: Grid) -> None: # If heterogeneous flow data is given, the speed ups at the defined # grid locations are determined in either 2 or 3 dimensions. else: - bounds = np.array( - list( - zip( - self.het_config['x_locs'], - self.het_config['y_locs'], - ) - ) - ) + bounds = np.array(list(zip(self.het_config['x_locs'], self.het_config['y_locs']))) hull = ConvexHull(bounds) polygon = Polygon(bounds[hull.vertices]) path = mpltPath.Path(polygon.boundary.coords) diff --git a/floris/simulation/grid.py b/floris/simulation/grid.py index 76e4c1317..8e2325252 100644 --- a/floris/simulation/grid.py +++ b/floris/simulation/grid.py @@ -83,8 +83,6 @@ class Grid(ABC): x_sorted_inertial_frame: NDArrayFloat = field(init=False) y_sorted_inertial_frame: NDArrayFloat = field(init=False) z_sorted_inertial_frame: NDArrayFloat = field(init=False) - x_center_of_rotation: NDArrayFloat = field(init=False) - y_center_of_rotation: NDArrayFloat = field(init=False) cubature_weights: NDArrayFloat = field(init=False, default=None) def __attrs_post_init__(self) -> None: @@ -252,7 +250,7 @@ def set_grid(self) -> None: disc_area_radius, self.grid_resolution, dtype=floris_float_type, - axis=1, + axis=1 ) # Construct the turbine grids # Here, they are already rotated to the correct orientation for each wind direction @@ -472,6 +470,8 @@ class FlowFieldGrid(Grid): reference_turbine_diameter (:py:obj:`float`): The reference turbine's rotor diameter. grid_resolution (:py:obj:`int`): The number of points on each turbine """ + x_center_of_rotation: NDArrayFloat = field(init=False) + y_center_of_rotation: NDArrayFloat = field(init=False) def __attrs_post_init__(self) -> None: super().__attrs_post_init__() @@ -496,7 +496,7 @@ def set_grid(self) -> None: # These are the rotated coordinates of the wind turbines based on the wind direction x, y, z, self.x_center_of_rotation, self.y_center_of_rotation = rotate_coordinates_rel_west( self.wind_directions, - self.turbine_coordinates_array, + self.turbine_coordinates_array ) # Construct the arrays storing the grid points @@ -543,7 +543,8 @@ class FlowFieldPlanarGrid(Grid): planar_coordinate: float = field() x1_bounds: tuple = field(default=None) x2_bounds: tuple = field(default=None) - + x_center_of_rotation: NDArrayFloat = field(init=False) + y_center_of_rotation: NDArrayFloat = field(init=False) sorted_indices: NDArrayInt = field(init=False) unsorted_indices: NDArrayInt = field(init=False) @@ -569,7 +570,7 @@ def set_grid(self) -> None: # These are the rotated coordinates of the wind turbines based on the wind direction x, y, z, self.x_center_of_rotation, self.y_center_of_rotation = rotate_coordinates_rel_west( self.wind_directions, - self.turbine_coordinates_array, + self.turbine_coordinates_array ) max_diameter = np.max(self.reference_turbine_diameter) diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index 1b243dd45..c56832b33 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -182,7 +182,6 @@ def reinitialize( self, wind_speeds: list[float] | NDArrayFloat | None = None, wind_directions: list[float] | NDArrayFloat | None = None, - # wind_layout: list[float] | NDArrayFloat | None = None, wind_shear: float | None = None, wind_veer: float | None = None, reference_wind_height: float | None = None, @@ -194,9 +193,6 @@ def reinitialize( layout_y: list[float] | NDArrayFloat | None = None, turbine_type: list | None = None, turbine_library_path: str | Path | None = None, - # turbine_id: list[str] | None = None, - # wtg_id: list[str] | None = None, - # with_resolution: float | None = None, solver_settings: dict | None = None, time_series: bool = False, het_config=None, @@ -241,10 +237,9 @@ def reinitialize( ## Wake # if wake is not None: # self.floris.wake = wake - # if turbulence_intensity is not None: - # pass # TODO: this should be in the code, but maybe got skipped? # if turbulence_kinetic_energy is not None: # pass # TODO: not needed until GCH + if solver_settings is not None: floris_dict["solver"] = solver_settings diff --git a/floris/tools/visualization.py b/floris/tools/visualization.py index 47fbcf72d..c0d14c8d4 100644 --- a/floris/tools/visualization.py +++ b/floris/tools/visualization.py @@ -85,7 +85,7 @@ def plot_turbines_with_fi( from a FLORIS interface object Args: - fi (:py:class:`floris.tools.floris_interface.FlorisIntreface`): FlorisIntreface object. + fi (:py:class:`floris.tools.floris_interface.FlorisInterface`): FlorisInterface object. ax (:py:class:`matplotlib.pyplot.axes`): Figure axes. Defaults to None. color (str, optional): Color to plot turbines. Defaults to None. wd (list, optional): The wind direction to plot the turbines relative to. Defaults to None. diff --git a/tests/utilities_unit_test.py b/tests/utilities_unit_test.py index 40ea594a8..4ec7e9d3c 100644 --- a/tests/utilities_unit_test.py +++ b/tests/utilities_unit_test.py @@ -93,7 +93,7 @@ def test_rotate_coordinates_rel_west(): wind_directions = np.array([270.0]) x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west( wind_directions, - coordinates, + coordinates ) np.testing.assert_array_equal( X_COORDS, x_rotated[0,0] ) @@ -112,7 +112,7 @@ def test_rotate_coordinates_rel_west(): wind_directions = np.array([360.0]) x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west( wind_directions, - coordinates, + coordinates ) np.testing.assert_almost_equal( Y_COORDS, x_rotated[0,0] - np.min(x_rotated[0,0])) np.testing.assert_almost_equal( X_COORDS, y_rotated[0,0] - np.min(y_rotated[0,0])) @@ -124,7 +124,7 @@ def test_rotate_coordinates_rel_west(): wind_directions = np.array([90.0]) x_rotated, y_rotated, z_rotated, _, _ = rotate_coordinates_rel_west( wind_directions, - coordinates, + coordinates ) np.testing.assert_almost_equal( X_COORDS[-1:-4:-1], x_rotated[0,0] ) np.testing.assert_almost_equal( Y_COORDS, y_rotated[0,0] ) From fc47c2f12d7d13896a822d190514f09de62a5ade Mon Sep 17 00:00:00 2001 From: Rafael M Mudafort Date: Fri, 12 May 2023 14:13:51 -0500 Subject: [PATCH 33/36] Add docstring --- examples/16_heterogeneous_inflow.py | 2 +- floris/simulation/flow_field.py | 6 +++--- floris/utilities.py | 12 ++++++------ 3 files changed, 10 insertions(+), 10 deletions(-) diff --git a/examples/16_heterogeneous_inflow.py b/examples/16_heterogeneous_inflow.py index 3426028d7..121f2610b 100644 --- a/examples/16_heterogeneous_inflow.py +++ b/examples/16_heterogeneous_inflow.py @@ -37,7 +37,7 @@ # Initialize FLORIS with the given input file via FlorisInterface. -# Note that the hetergeneous flow is defined in the input file. The het_config +# Note that the heterogeneous flow is defined in the input file. The het_config # dictionary is defined as below. The speed ups are multipliers of the ambient wind speed, # and the x_locs and y_locs are the locations of the speed ups. # diff --git a/floris/simulation/flow_field.py b/floris/simulation/flow_field.py index e69337083..3e838c5ab 100644 --- a/floris/simulation/flow_field.py +++ b/floris/simulation/flow_field.py @@ -254,7 +254,7 @@ def calculate_speed_ups(self, het_map, x, y, z=None): def generate_heterogeneous_wind_map(self): """This function creates the heterogeneous interpolant used to calculate heterogeneous inflows. The interpolant is for computing wind speed based on an x and y location in the - flow field. This is computed using Scipy's LinearNDInterpolator and uses a fill value + flow field. This is computed using SciPy's LinearNDInterpolator and uses a fill value equal to the freestream for interpolated values outside of the user-defined heterogeneous map bounds. @@ -277,7 +277,7 @@ def generate_heterogeneous_wind_map(self): z = self.het_config['z_locs'] if z is not None: - # Compute the 3-dimensional interpolants for each wind diretion + # Compute the 3-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 = [ @@ -285,7 +285,7 @@ def generate_heterogeneous_wind_map(self): for speed_up in speed_ups ] else: - # Compute the 2-dimensional interpolants for each wind diretion + # 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 = [ diff --git a/floris/utilities.py b/floris/utilities.py index 1de20d88c..6e565d225 100644 --- a/floris/utilities.py +++ b/floris/utilities.py @@ -271,12 +271,12 @@ def rotate_coordinates_rel_west( def reverse_rotate_coordinates_rel_west( - wind_directions, - grid_x, - grid_y, - grid_z, - x_center_of_rotation, - y_center_of_rotation + wind_directions: NDArrayFloat, + grid_x: NDArrayFloat, + grid_y: NDArrayFloat, + grid_z: NDArrayFloat, + x_center_of_rotation: float, + y_center_of_rotation: float ): """ This function reverses the rotation of the given grid so that the coordinates are aligned with From 2f9ee1af56abc26d7b30dd9c85310dc64f95f89d Mon Sep 17 00:00:00 2001 From: Rafael M Mudafort Date: Tue, 16 May 2023 11:40:43 -0500 Subject: [PATCH 34/36] Reset example settings --- examples/03_making_adjustments.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/03_making_adjustments.py b/examples/03_making_adjustments.py index 2144262f9..750288d5a 100644 --- a/examples/03_making_adjustments.py +++ b/examples/03_making_adjustments.py @@ -76,7 +76,7 @@ 5.0 * fi.floris.farm.rotor_diameters[0][0][0] * np.arange(0, N, 1), 5.0 * fi.floris.farm.rotor_diameters[0][0][0] * np.arange(0, N, 1), ) -fi.reinitialize(layout_x=X.flatten(), layout_y=Y.flatten(), wind_directions=[360.0]) +fi.reinitialize(layout_x=X.flatten(), layout_y=Y.flatten(), wind_directions=[270.0]) horizontal_plane = fi.calculate_horizontal_plane(height=90.0) wakeviz.visualize_cut_plane( horizontal_plane, From cca4e62193cad52b90ea117973d8152cc4ca12b6 Mon Sep 17 00:00:00 2001 From: Rafael M Mudafort Date: Tue, 16 May 2023 11:41:50 -0500 Subject: [PATCH 35/36] Early return in error checking --- floris/simulation/flow_field.py | 36 ++++++++++++++++++--------------- 1 file changed, 20 insertions(+), 16 deletions(-) diff --git a/floris/simulation/flow_field.py b/floris/simulation/flow_field.py index 3e838c5ab..51dfc33de 100644 --- a/floris/simulation/flow_field.py +++ b/floris/simulation/flow_field.py @@ -81,28 +81,32 @@ def het_config_validator(self, instance: attrs.Attribute, value: dict | None) -> """Using the validator method to check that the het_config dictionary has the correct key-value pairs. """ - if value is not None: - # Check that the correct keys are supplied for the het_config dict - for k in ["speed_ups", "x_locs", "y_locs"]: - if k not in value.keys(): - raise ValueError( - "het_config must contain entries for 'speeds_ups', 'x_locs', and 'y_locs'," - f" with 'z_locs' as optional. Missing '{k}'." - ) - if "z_locs" not in value: - # If only a 2D case, add "None" for the z locations - value["z_locs"] = None + if value is None: + return + + # Check that the correct keys are supplied for the het_config dict + for k in ["speed_ups", "x_locs", "y_locs"]: + if k not in value.keys(): + raise ValueError( + "het_config must contain entries for 'speeds_ups', 'x_locs', and 'y_locs'," + f" with 'z_locs' as optional. Missing '{k}'." + ) + if "z_locs" not in value: + # If only a 2D case, add "None" for the z locations + value["z_locs"] = None @het_map.validator def het_map_validator(self, instance: attrs.Attribute, value: list | None) -> None: """Using this validator to make sure that the het_map has an interpolant defined for each wind direction. """ - if value is not None: - if self.n_wind_directions!= np.array(value).shape[0]: - raise ValueError( - "The het_map's wind direction dimension not equal to number of wind directions." - ) + if value is None: + return + + if self.n_wind_directions!= np.array(value).shape[0]: + raise ValueError( + "The het_map's wind direction dimension not equal to number of wind directions." + ) def __attrs_post_init__(self) -> None: From e620bd65bdea87ea1444954f5adc91aa1696c065 Mon Sep 17 00:00:00 2001 From: Rafael M Mudafort Date: Tue, 16 May 2023 11:59:39 -0500 Subject: [PATCH 36/36] Rename heterogeneous config dictionary MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit New name and key/value’s are more explicit --- examples/16_heterogeneous_inflow.py | 27 ++++---- examples/16b_heterogenaity_multiple_ws_wd.py | 23 +++---- examples/inputs/gch_heterogeneous_inflow.yaml | 8 +-- floris/simulation/flow_field.py | 62 +++++++++---------- floris/tools/floris_interface.py | 6 +- floris/tools/visualization.py | 4 +- 6 files changed, 64 insertions(+), 66 deletions(-) diff --git a/examples/16_heterogeneous_inflow.py b/examples/16_heterogeneous_inflow.py index 121f2610b..3f04d5bc4 100644 --- a/examples/16_heterogeneous_inflow.py +++ b/examples/16_heterogeneous_inflow.py @@ -37,16 +37,17 @@ # Initialize FLORIS with the given input file via FlorisInterface. -# Note that the heterogeneous flow is defined in the input file. The het_config +# Note that the heterogeneous flow is defined in the input file. The heterogenous_inflow_config # dictionary is defined as below. The speed ups are multipliers of the ambient wind speed, -# and the x_locs and y_locs are the locations of the speed ups. +# and the x and y are the locations of the speed ups. # -# het_config = { -# 'speed_ups': [[2.0, 1.0, 2.0, 1.0]], -# 'x_locs': [-300.0, -300.0, 2600.0, 2600.0], -# 'y_locs': [ -300.0, 300.0, -300.0, 300.0], +# heterogenous_inflow_config = { +# 'speed_multipliers': [[2.0, 1.0, 2.0, 1.0]], +# 'x': [-300.0, -300.0, 2600.0, 2600.0], +# 'y': [ -300.0, 300.0, -300.0, 300.0], # } + fi_2d = FlorisInterface("inputs/gch_heterogeneous_inflow.yaml") # Set shear to 0.0 to highlight the heterogeneous inflow @@ -88,24 +89,24 @@ # Define the speed ups of the heterogeneous inflow, and their locations. # For the 3-dimensional case, this requires x, y, and z locations. # The speed ups are multipliers of the ambient wind speed. -speed_ups = [[1.0, 1.0, 2.0, 2.0, 1.0, 1.0, 2.0, 2.0]] +speed_multipliers = [[1.0, 1.0, 2.0, 2.0, 1.0, 1.0, 2.0, 2.0]] x_locs = [-300.0, -300.0, -300.0, -300.0, 2600.0, 2600.0, 2600.0, 2600.0] y_locs = [-300.0, 300.0, -300.0, 300.0, -300.0, 300.0, -300.0, 300.0] z_locs = [540.0, 540.0, 0.0, 0.0, 540.0, 540.0, 0.0, 0.0] # Create the configuration dictionary to be used for the heterogeneous inflow. -het_config_3d = { - 'speed_ups': speed_ups, - 'x_locs': x_locs, - 'y_locs': y_locs, - 'z_locs': z_locs, +heterogenous_inflow_config = { + 'speed_multipliers': speed_multipliers, + 'x': x_locs, + 'y': y_locs, + 'z': z_locs, } # Initialize FLORIS with the given input file via FlorisInterface. # Note that we initialize FLORIS with a homogenous flow input file, but # then configure the heterogeneous inflow via the reinitialize method. fi_3d = FlorisInterface("inputs/gch.yaml") -fi_3d.reinitialize(het_config=het_config_3d) +fi_3d.reinitialize(heterogenous_inflow_config=heterogenous_inflow_config) # Set shear to 0.0 to highlight the heterogeneous inflow fi_3d.reinitialize(wind_shear=0.0) diff --git a/examples/16b_heterogenaity_multiple_ws_wd.py b/examples/16b_heterogenaity_multiple_ws_wd.py index c4822df72..43ac6f7eb 100644 --- a/examples/16b_heterogenaity_multiple_ws_wd.py +++ b/examples/16b_heterogenaity_multiple_ws_wd.py @@ -33,13 +33,6 @@ x_locs = [-300.0, -300.0, 2600.0, 2600.0] y_locs = [ -300.0, 300.0, -300.0, 300.0] -# Create the configuration dictionary to be used for the heterogeneous inflow. -het_config_2d = { - 'speed_ups': [[2.0, 1.0, 2.0, 1.0]], - 'x_locs': [-300.0, -300.0, 2600.0, 2600.0], - 'y_locs': [ -300.0, 300.0, -300.0, 300.0], -} - # Initialize FLORIS with the given input file via FlorisInterface. # Note the heterogeneous inflow is defined in the input file. fi = FlorisInterface("inputs/gch_heterogeneous_inflow.yaml") @@ -77,13 +70,17 @@ # To change the number of wind directions however it is necessary to make a matching # change to the dimensions of the het map -speed_ups = [[2.0, 1.0, 2.0, 1.0], [2.0, 1.0, 2.0, 1.0]] # Expand to two wind directions -het_config_2d = { - 'speed_ups': speed_ups, - 'x_locs': x_locs, - 'y_locs': y_locs, +speed_multipliers = [[2.0, 1.0, 2.0, 1.0], [2.0, 1.0, 2.0, 1.0]] # Expand to two wind directions +heterogenous_inflow_config = { + 'speed_multipliers': speed_multipliers, + 'x': x_locs, + 'y': y_locs, } -fi.reinitialize(wind_directions=[270., 275.], wind_speeds=[8.], het_config=het_config_2d) +fi.reinitialize( + wind_directions=[270.0, 275.0], + wind_speeds=[8.0], + heterogenous_inflow_config=heterogenous_inflow_config +) fi.calculate_wake() turbine_powers = np.round(fi.get_turbine_powers() / 1000.) print('With wind directions now set to 270 and 275 deg') diff --git a/examples/inputs/gch_heterogeneous_inflow.yaml b/examples/inputs/gch_heterogeneous_inflow.yaml index c0458c851..d7cffa0d5 100644 --- a/examples/inputs/gch_heterogeneous_inflow.yaml +++ b/examples/inputs/gch_heterogeneous_inflow.yaml @@ -27,18 +27,18 @@ farm: flow_field: air_density: 1.225 - het_config: - speed_ups: + heterogenous_inflow_config: + speed_multipliers: - - 2.0 - 1.0 - 2.0 - 1.0 - x_locs: + x: - -300. - -300. - 2600. - 2600. - y_locs: + y: - -300. - 300. - -300. diff --git a/floris/simulation/flow_field.py b/floris/simulation/flow_field.py index 51dfc33de..a8dbf7393 100644 --- a/floris/simulation/flow_field.py +++ b/floris/simulation/flow_field.py @@ -42,7 +42,7 @@ class FlowField(BaseClass): turbulence_intensity: float = field(converter=float) reference_wind_height: float = field(converter=float) time_series : bool = field(default=False) - het_config: dict = field(default=None) + heterogenous_inflow_config: dict = field(default=None) n_wind_speeds: int = field(init=False) n_wind_directions: int = field(init=False) @@ -76,24 +76,24 @@ def wind_directions_validator(self, instance: attrs.Attribute, value: NDArrayFlo """Using the validator method to keep the `n_wind_directions` attribute up to date.""" self.n_wind_directions = value.size - @het_config.validator - def het_config_validator(self, instance: attrs.Attribute, value: dict | None) -> None: - """Using the validator method to check that the het_config dictionary has the correct - key-value pairs. + @heterogenous_inflow_config.validator + def heterogenous_config_validator(self, instance: attrs.Attribute, value: dict | None) -> None: + """Using the validator method to check that the heterogenous_inflow_config dictionary has + the correct key-value pairs. """ if value is None: return - # Check that the correct keys are supplied for the het_config dict - for k in ["speed_ups", "x_locs", "y_locs"]: + # Check that the correct keys are supplied for the heterogenous_inflow_config dict + for k in ["speed_multipliers", "x", "y"]: if k not in value.keys(): raise ValueError( - "het_config must contain entries for 'speeds_ups', 'x_locs', and 'y_locs'," - f" with 'z_locs' as optional. Missing '{k}'." + "heterogenous_inflow_config must contain entries for 'speed_multipliers'," + f"'x', and 'y', with 'z' optional. Missing '{k}'." ) - if "z_locs" not in value: + if "z" not in value: # If only a 2D case, add "None" for the z locations - value["z_locs"] = None + value["z"] = None @het_map.validator def het_map_validator(self, instance: attrs.Attribute, value: list | None) -> None: @@ -110,7 +110,7 @@ def het_map_validator(self, instance: attrs.Attribute, value: list | None) -> No def __attrs_post_init__(self) -> None: - if self.het_config is not None: + if self.heterogenous_inflow_config is not None: self.generate_heterogeneous_wind_map() @@ -140,7 +140,10 @@ def initialize_velocity_field(self, grid: Grid) -> None: # If heterogeneous flow data is given, the speed ups at the defined # grid locations are determined in either 2 or 3 dimensions. else: - bounds = np.array(list(zip(self.het_config['x_locs'], self.het_config['y_locs']))) + bounds = np.array(list(zip( + self.heterogenous_inflow_config['x'], + self.heterogenous_inflow_config['y'] + ))) hull = ConvexHull(bounds) polygon = Polygon(bounds[hull.vertices]) path = mpltPath.Path(polygon.boundary.coords) @@ -263,38 +266,35 @@ def generate_heterogeneous_wind_map(self): map bounds. Args: - het_config (dict): The heterogeneous inflow configuration dictionary. + heterogenous_inflow_config (dict): The heterogeneous inflow configuration dictionary. The configuration should have the following inputs specified. - - **speed_ups** (list): A list of speed up factors that will multiply the specified - freestream wind speed. This 2-dimensional array should have an array of - multiplicative factors defined for each wind direction. - - **x_locs** (list): A list of x locations at which the speed up factors are - defined. - - **y_locs**: A list of y locations at which the speed up factors are - defined. - - **z_locs** (optional): A list of z locations at which the speed up factors are - defined. + - **speed_multipliers** (list): A list of speed up factors that will multiply + the specified freestream wind speed. This 2-dimensional array should have an + array of multiplicative factors defined for each wind direction. + - **x** (list): A list of x locations at which the speed up factors are defined. + - **y**: A list of y locations at which the speed up factors are defined. + - **z** (optional): A list of z locations at which the speed up factors are defined. """ - speed_ups = self.het_config['speed_ups'] - x = self.het_config['x_locs'] - y = self.het_config['y_locs'] - z = self.het_config['z_locs'] + speed_multipliers = self.heterogenous_inflow_config['speed_multipliers'] + x = self.heterogenous_inflow_config['x'] + y = self.heterogenous_inflow_config['y'] + z = self.heterogenous_inflow_config['z'] if z is not None: # Compute the 3-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, z)), speed_up, fill_value=1.0) - for speed_up in speed_ups + LinearNDInterpolator(list(zip(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)), speed_up, fill_value=1.0) - for speed_up in speed_ups + LinearNDInterpolator(list(zip(x, y)), multiplier, fill_value=1.0) + for multiplier in speed_multipliers ] self.het_map = in_region diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index c56832b33..4e4cc352e 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -195,7 +195,7 @@ def reinitialize( turbine_library_path: str | Path | None = None, solver_settings: dict | None = None, time_series: bool = False, - het_config=None, + heterogenous_inflow_config=None, ): # Export the floris object recursively as a dictionary floris_dict = self.floris.as_dict() @@ -219,8 +219,8 @@ def reinitialize( flow_field_dict["turbulence_intensity"] = turbulence_intensity if air_density is not None: flow_field_dict["air_density"] = air_density - if het_config is not None: - flow_field_dict["het_config"] = het_config + if heterogenous_inflow_config is not None: + flow_field_dict["heterogenous_inflow_config"] = heterogenous_inflow_config ## Farm if layout_x is not None: diff --git a/floris/tools/visualization.py b/floris/tools/visualization.py index c0d14c8d4..529e9f8da 100644 --- a/floris/tools/visualization.py +++ b/floris/tools/visualization.py @@ -391,8 +391,8 @@ def visualize_heterogeneous_cut_plane( points = np.array( list( zip( - fi.floris.flow_field.het_config['x_locs'], - fi.floris.flow_field.het_config['y_locs'], + fi.floris.flow_field.heterogenous_inflow_config['x'], + fi.floris.flow_field.heterogenous_inflow_config['y'], ) ) )