From 92ab0df38ba3a4fbf123b19b98caa57d2cffc873 Mon Sep 17 00:00:00 2001 From: Rafael M Mudafort Date: Thu, 11 May 2023 15:06:35 -0500 Subject: [PATCH] Save cubature weights on Grids Rather than the full cubature coefficient data --- floris/simulation/grid.py | 32 ++++++++++++------------ floris/simulation/solver.py | 14 +++++------ floris/simulation/turbine.py | 18 ++++++------- floris/tools/floris_interface.py | 2 +- tests/reg_tests/gauss_regression_test.py | 3 +-- 5 files changed, 34 insertions(+), 35 deletions(-) diff --git a/floris/simulation/grid.py b/floris/simulation/grid.py index 7a5b343d8..a79030280 100644 --- a/floris/simulation/grid.py +++ b/floris/simulation/grid.py @@ -79,7 +79,7 @@ class Grid(ABC): x_sorted: NDArrayFloat = field(init=False) y_sorted: NDArrayFloat = field(init=False) z_sorted: NDArrayFloat = field(init=False) - cubature_coefficients: dict = field(init=False, default=None) + cubature_weights: NDArrayFloat = field(init=False, default=None) def __attrs_post_init__(self) -> None: self.turbine_coordinates_array = np.array([c.elements for c in self.turbine_coordinates]) @@ -304,15 +304,17 @@ def set_grid(self) -> None: ) # Coefficients - self.cubature_coefficients = CubatureGrid.get_cubature_coefficients(self.grid_resolution) + cubature_coefficients = CubatureGrid.get_cubature_coefficients(self.grid_resolution) # Generate grid points - yv = np.kron(self.cubature_coefficients["r"], self.cubature_coefficients["q"]) - zv = np.kron(self.cubature_coefficients["r"], self.cubature_coefficients["t"]) + yv = np.kron(cubature_coefficients["r"], cubature_coefficients["q"]) + zv = np.kron(cubature_coefficients["r"], cubature_coefficients["t"]) - # Calculate weighing terms for the grid points - self.cubature_coefficients["CUBcoeff"] = \ - np.kron(self.cubature_coefficients["A"],np.ones((1, N))) * self.cubature_coefficients["B"] / np.pi + # Calculate weighting terms for the grid points + self.cubature_weights = ( + np.kron(cubature_coefficients["A"], np.ones((1, self.grid_resolution))) + * cubature_coefficients["B"] / np.pi + ) # Here, the coordinates are already rotated to the correct orientation for each # wind direction @@ -424,15 +426,13 @@ def get_cubature_coefficients(cls, N: int): q = [ 0.1564344650402308690101053, 0.4539904997395467915604084, 0.7071067811865475244008444, 0.8910065241883678623597096, 0.9876883405951377261900402, 0.9876883405951377261900402, 0.8910065241883678623597096, 0.7071067811865475244008444, 0.4539904997395467915604084, 0.1564344650402308690101053] # noqa: E501 A = [ 0.0592317212640472718785660, 0.1196571676248416170103229, 0.1422222222222222222222222, 0.1196571676248416170103229, 0.0592317212640472718785660, 0.0592317212640472718785660, 0.1196571676248416170103229, 0.1422222222222222222222222, 0.1196571676248416170103229, 0.0592317212640472718785660] # noqa: E501 - cubature_coefficients = { - "r": np.array(r, dtype=float), - "t": np.array(t, dtype=float), - "q": np.array(q, dtype=float), - "A": np.array(A, dtype=float), - "B": np.pi/N, - } - - return cubature_coefficients + return { + "r": np.array(r, dtype=float), + "t": np.array(t, dtype=float), + "q": np.array(q, dtype=float), + "A": np.array(A, dtype=float), + "B": np.pi/N, + } @define class FlowFieldGrid(Grid): diff --git a/floris/simulation/solver.py b/floris/simulation/solver.py index c3c680ca6..88f0fbdfd 100644 --- a/floris/simulation/solver.py +++ b/floris/simulation/solver.py @@ -99,7 +99,7 @@ def sequential_solver( average_velocities = average_velocity( flow_field.u_sorted, method=grid.average_method, - cubature_coefficients=grid.cubature_coefficients + cubature_weights=grid.cubature_weights ) ct_i = Ct( @@ -491,7 +491,7 @@ def cc_solver( turb_avg_vels = average_velocity( turb_inflow_field, method=grid.average_method, - cubature_coefficients=grid.cubature_coefficients + cubature_weights=grid.cubature_weights ) turb_Cts = Ct( turb_avg_vels, @@ -523,7 +523,7 @@ def cc_solver( different_average = average_velocity( flow_field.u_sorted, method=grid.average_method, - cubature_coefficients=grid.cubature_coefficients + cubature_weights=grid.cubature_weights ) axial_induction_i = axial_induction( average_velocities=different_average, @@ -739,7 +739,7 @@ def full_flow_cc_solver( turb_avg_vels = average_velocity( turbine_grid_flow_field.u_sorted, method=turbine_grid.average_method, - cubature_coefficients=turbine_grid.cubature_coefficients + cubature_weights=turbine_grid.cubature_weights ) turb_Cts = Ct( velocities=turb_avg_vels, @@ -889,7 +889,7 @@ def turbopark_solver( average_velocities = average_velocity( flow_field.u_sorted, method=grid.average_method, - cubature_coefficients=grid.cubature_coefficients + cubature_weights=grid.cubature_weights ) Cts = Ct( @@ -972,7 +972,7 @@ def turbopark_solver( average_velocities = average_velocity( flow_field.u_sorted, method=grid.average_method, - cubature_coefficients=grid.cubature_coefficients + cubature_weights=grid.cubature_weights ) yaw_ii = farm.yaw_angles_sorted[:, :, ii:ii+1, None, None] @@ -1211,7 +1211,7 @@ def empirical_gauss_solver( average_velocities = average_velocity( flow_field.u_sorted, method=grid.average_method, - cubature_coefficients=grid.cubature_coefficients + cubature_weights=grid.cubature_weights ) ct_i = Ct( diff --git a/floris/simulation/turbine.py b/floris/simulation/turbine.py index a924ae7b9..0d0fe0752 100644 --- a/floris/simulation/turbine.py +++ b/floris/simulation/turbine.py @@ -430,14 +430,14 @@ def simple_mean(array, axis=0): def cubic_mean(array, axis=0): return np.cbrt(np.mean(array ** 3.0, axis=axis)) -def simple_cubature(array, cubature_coefficients, axis=0): - weights = cubature_coefficients["CUBcoeff"].flatten() +def simple_cubature(array, cubature_weights, axis=0): + weights = cubature_weights.flatten() weights = weights * len(weights) / np.sum(weights) product = (array * weights[None, None, None, :, None]) return simple_mean(product, axis) -def cubic_cubature(array, cubature_coefficients, axis=0): - weights = cubature_coefficients["CUBcoeff"].flatten() +def cubic_cubature(array, cubature_weights, axis=0): + weights = cubature_weights.flatten() weights = weights * len(weights) / np.sum(weights) return np.cbrt(np.mean((array**3.0 * weights[None, None, None, :, None]), axis=axis)) @@ -445,7 +445,7 @@ def average_velocity( velocities: NDArrayFloat, method: str, ix_filter: NDArrayFilter | Iterable[int] | None = None, - cubature_coefficients: dict | None = None + cubature_weights: dict | None = None ) -> NDArrayFloat: """This property calculates and returns the cube root of the mean cubed velocity in the turbine's rotor swept area (m/s). @@ -471,16 +471,16 @@ def average_velocity( axis = tuple([3 + i for i in range(velocities.ndim - 3)]) if method == "simple-mean": - return simple_mean(velocities, axis=axis) + return simple_mean(velocities, axis) elif method == "cubic-mean": - return cubic_mean(velocities, axis=axis) + return cubic_mean(velocities, axis) elif method == "simple-cubature": - return simple_cubature(velocities, cubature_coefficients=cubature_coefficients, axis=axis) + return simple_cubature(velocities, cubature_weights, axis) elif method == "cubic-cubature": - return cubic_cubature(velocities, cubature_coefficients=cubature_coefficients, axis=axis) + return cubic_cubature(velocities, cubature_weights, axis) else: raise ValueError("Incorrect method given.") diff --git a/floris/tools/floris_interface.py b/floris/tools/floris_interface.py index 329b35c59..420bf1acd 100644 --- a/floris/tools/floris_interface.py +++ b/floris/tools/floris_interface.py @@ -655,7 +655,7 @@ def turbine_average_velocities(self) -> NDArrayFloat: return average_velocity( velocities=self.floris.flow_field.u, method=self.floris.grid.average_method, - cubature_coefficients=self.floris.grid.cubature_coefficients + cubature_weights=self.floris.grid.cubature_weights ) @property diff --git a/tests/reg_tests/gauss_regression_test.py b/tests/reg_tests/gauss_regression_test.py index e1b75c6b5..ca28bcc24 100644 --- a/tests/reg_tests/gauss_regression_test.py +++ b/tests/reg_tests/gauss_regression_test.py @@ -274,8 +274,7 @@ def test_regression_tandem(sample_inputs_fixture): farm_avg_velocities = average_velocity( floris.flow_field.u, - method="cubic-mean", - cubature_coefficients=floris.grid.cubature_coefficients + method="cubic-mean" ) farm_eff_velocities = rotor_effective_velocity( floris.flow_field.air_density,