diff --git a/py_wake/superposition_models.py b/py_wake/superposition_models.py index 06406738bc8de8ac9ac4a7567a4636dec58510ae..86159e805e0227b596695e1cc48adb09f0a499ff 100644 --- a/py_wake/superposition_models.py +++ b/py_wake/superposition_models.py @@ -52,6 +52,7 @@ class AddedTurbulenceSuperpositionModel(): class SquaredSum(SuperpositionModel, AddedTurbulenceSuperpositionModel): def __call__(self, value_jxxx): + assert not np.any(value_jxxx < 0), "SquaredSum only works for deficit - not speedups" return np.sqrt(np.sum(np.power(value_jxxx, 2), 0)) diff --git a/py_wake/tests/test_deficit_models/test_deficit_models.py b/py_wake/tests/test_deficit_models/test_deficit_models.py index eb3192e2669b0b9247d643636e478e8098fefab0..d34d1db10b752d53a20dd224097da91f35fbf1c7 100644 --- a/py_wake/tests/test_deficit_models/test_deficit_models.py +++ b/py_wake/tests/test_deficit_models/test_deficit_models.py @@ -16,7 +16,7 @@ from py_wake.examples.data.iea37 import iea37_path from py_wake.examples.data.iea37._iea37 import IEA37_WindTurbines, IEA37Site from py_wake.examples.data.iea37.iea37_reader import read_iea37_windfarm from py_wake.flow_map import HorizontalGrid, XYGrid, YZGrid -from py_wake.superposition_models import SquaredSum, WeightedSum +from py_wake.superposition_models import SquaredSum, WeightedSum, LinearSum from py_wake.tests import npt from py_wake.tests.test_files import tfp from py_wake.turbulence_models.gcl_turb import GCLTurbulence @@ -62,10 +62,10 @@ class GCLLocalDeficit(GCLDeficit): (IEA37SimpleBastankhahGaussianDeficit(), read_iea37_windfarm(iea37_path + 'iea37-ex16.yaml')[2]), (FugaDeficit(LUT_path=tfp + 'fuga/2MW/Z0=0.00408599Zi=00400Zeta0=0.00E+00.nc', smooth2zero_x=250, smooth2zero_y=50), - (404441.6306021485, [9912.33731, 9762.05717, 12510.14066, 15396.76584, 23017.66483, - 27799.7161, 43138.41606, 49623.79059, 24979.09001, 15460.45923, - 16723.02619, 35694.35526, 77969.14805, 19782.41376, 13721.45739, - 8950.79218])), + (404307.913694007, [9962.5736, 9801.78637, 12608.917, 15452.89633, 22577.55462, + 27901.06282, 43479.02414, 49825.74736, 25105.68548, 15542.67624, + 16918.79822, 35640.98079, 76856.89565, 19752.83273, 13882.09085, + 8998.39151])), (GCLDeficit(), (370863.6246093183, [9385.75387, 8768.52105, 11450.13309, 14262.42186, 21178.74926, 25751.59502, 39483.21753, 44573.31533, 23652.09976, 13924.58752, @@ -111,8 +111,9 @@ def test_IEA37_ex16(deficitModel, aep_ref): site = IEA37Site(16) x, y = site.initial_position.T windTurbines = IEA37_WindTurbines() + superpositionModel = [SquaredSum, LinearSum][isinstance(deficitModel, FugaDeficit)] wf_model = PropagateDownwind(site, windTurbines, wake_deficitModel=deficitModel, - superpositionModel=SquaredSum(), turbulenceModel=GCLTurbulence()) + superpositionModel=superpositionModel(), turbulenceModel=GCLTurbulence()) aep_ilk = wf_model(x, y, wd=np.arange(0, 360, 22.5), ws=[9.8]).aep_ilk(normalize_probabilities=True) aep_MW_l = aep_ilk.sum((0, 2)) * 1000 @@ -196,7 +197,7 @@ def test_huge_distance_blockage(deficitModel): (IEA37SimpleBastankhahGaussianDeficit(), [3.32, 4.86, 7.0, 8.1, 7.8, 7.23, 6.86, 6.9, 7.3, 7.82, 8.11, 8.04, 7.87, 7.79, 7.85, 8.04, 8.28]), (FugaDeficit(LUT_path=tfp + 'fuga/2MW/Z0=0.00408599Zi=00400Zeta0=0.00E+00.nc'), - [7.06, 7.87, 8.77, 8.85, 8.52, 7.96, 7.49, 7.55, 8.06, 8.58, 8.69, 8.45, 8.18, 8.05, 8.15, 8.41, 8.68]), + [7.08, 7.89, 8.78, 8.88, 8.55, 7.99, 7.51, 7.55, 8.05, 8.5, 8.62, 8.46, 8.21, 8.08, 8.18, 8.45, 8.72]), (GCLDeficit(), [2.39, 5.01, 7.74, 8.34, 7.95, 7.58, 7.29, 7.32, 7.61, 7.92, 8.11, 8.09, 7.95, 7.83, 7.92, 8.1, 8.3]), (GCLLocalDeficit(), @@ -214,8 +215,8 @@ def test_deficitModel_wake_map(deficitModel, ref): site = IEA37Site(16) x, y = site.initial_position.T windTurbines = IEA37_WindTurbines() - - wf_model = PropagateDownwind(site, windTurbines, wake_deficitModel=deficitModel, superpositionModel=SquaredSum(), + superpositionModel = [SquaredSum, LinearSum][isinstance(deficitModel, FugaDeficit)] + wf_model = PropagateDownwind(site, windTurbines, wake_deficitModel=deficitModel, superpositionModel=superpositionModel(), turbulenceModel=GCLTurbulence()) x_j = np.linspace(-1500, 1500, 200) @@ -473,10 +474,10 @@ def test_deficitModel_wake_map_convection_all2all(deficitModel, ref): (IEA37SimpleBastankhahGaussian, read_iea37_windfarm(iea37_path + 'iea37-ex16.yaml')[2]), (lambda *args, **kwargs: Fuga(tfp + 'fuga/2MW/Z0=0.00408599Zi=00400Zeta0=0.00E+00.nc', *args, **kwargs), - (404422.9302289747, [9911.82745, 9761.75297, 12509.30601, 15396.26107, 23016.74091, - 27798.80471, 43135.53798, 49622.24427, 24977.80518, 15459.8399, - 16721.53971, 35692.36221, 77966.92736, 19781.30917, 13720.23771, - 8950.43363])), + (404411.56948244764, [9963.95691, 9804.76072, 12613.39798, 15457.51387, 22582.08811, + 27909.40005, 43494.4758, 49840.86701, 25109.17142, 15546.80155, + 16923.11586, 35647.55255, 76875.57937, 19756.4749, 13885.63352, + 9000.77984])), (GCL, (370863.6246093183, [9385.75387, 8768.52105, 11450.13309, 14262.42186, 21178.74926, 25751.59502, 39483.21753, 44573.31533, 23652.09976, 13924.58752, @@ -505,7 +506,7 @@ def test_IEA37_ex16_windFarmModel(windFarmModel, aep_ref): with warnings.catch_warnings(): warnings.filterwarnings('ignore', category=DeprecationWarning) wf_model = windFarmModel(site, windTurbines, turbulenceModel=GCLTurbulence()) - wf_model.superpositionModel = SquaredSum() + wf_model.superpositionModel = [SquaredSum(), LinearSum()][isinstance(wf_model.wake_deficitModel, FugaDeficit)] aep_ilk = wf_model(x, y, wd=np.arange(0, 360, 22.5), ws=[9.8]).aep_ilk(normalize_probabilities=True) aep_MW_l = aep_ilk.sum((0, 2)) * 1000 diff --git a/py_wake/tests/test_superposition_models.py b/py_wake/tests/test_superposition_models.py index 56ea50e65d1a008af67eeea684a99d35070a4e0c..0ad664e8f4b065ea0a20524aa1c61bed7a0df9d7 100644 --- a/py_wake/tests/test_superposition_models.py +++ b/py_wake/tests/test_superposition_models.py @@ -172,3 +172,10 @@ def test_superpositionModels_rotor_avg_model(superpositionModel): sim_res_a2a = wfm_a2a(wt9_x, wt9_y, operating=operating) sim_res_pdw = wfm_pdw(wt9_x, wt9_y, operating=operating) npt.assert_array_almost_equal(sim_res_a2a.WS_eff.sel(wt=8), sim_res_pdw.WS_eff.sel(wt=8)) + + +def test_squaredSum_speedups(): + wfm = All2AllIterative(site=UniformSite(), windTurbines=V80(), wake_deficitModel=NOJDeficit(), + blockage_deficitModel=SelfSimilarityDeficit(exclude_wake=False), superpositionModel=SquaredSum()) + with pytest.raises(AssertionError, match='SquaredSum only works for deficit - not speedups'): + wfm([0, 200, 400], [0, 0, 0], wd=270) diff --git a/py_wake/tests/test_wind_farm_models/test_enginering_wind_farm_model.py b/py_wake/tests/test_wind_farm_models/test_enginering_wind_farm_model.py index 767b146e7f3c8cc117f8e9a03eaf15939e9aabec..bd4eec66e4639479c5a687205255beef55ff9744 100644 --- a/py_wake/tests/test_wind_farm_models/test_enginering_wind_farm_model.py +++ b/py_wake/tests/test_wind_farm_models/test_enginering_wind_farm_model.py @@ -614,3 +614,16 @@ def test_PropagateUpDownIterative_blockage_deficitModels(blockage_deficitModel): RankineHalfBody, Rathmann, RathmannScaled, VortexCylinder]: atol = 0.01 npt.assert_allclose(pudi(x, y, wd=270).WS_eff, all2all(x, y, wd=270).WS_eff, atol=atol) + + +def test_PropagateUpDownIterative_SquaredSum(): + class MyBlockage(SelfSimilarityDeficit): + def calc_blockage_deficit(self, dw_ijlk, **kwargs): + return np.abs(SelfSimilarityDeficit.calc_blockage_deficit(self, dw_ijlk, **kwargs)) + wfm = PropagateUpDownIterative(site=UniformSite(), + windTurbines=V80(), + wake_deficitModel=BastankhahGaussianDeficit(use_effective_ws=True), + blockage_deficitModel=MyBlockage(exclude_wake=False, use_effective_ws=True), + superpositionModel=SquaredSum()) + with pytest.raises(AssertionError, match='PropagateDownwind does not work with SquaredSum when the blockage model has both up and downstream effects'): + wfm([0, 200, 400], [0, 0, 0], wd=270) diff --git a/py_wake/wind_farm_models/engineering_models.py b/py_wake/wind_farm_models/engineering_models.py index 36e5623eb8b62843658d34b9b235013d60e68aab..84d53dffd310c531b7fc3ea2f38c751c833bb7b5 100644 --- a/py_wake/wind_farm_models/engineering_models.py +++ b/py_wake/wind_farm_models/engineering_models.py @@ -1,7 +1,7 @@ from abc import abstractmethod from numpy import newaxis as na from py_wake import np -from py_wake.superposition_models import SuperpositionModel, LinearSum, WeightedSum, CumulativeWakeSum +from py_wake.superposition_models import SuperpositionModel, LinearSum, WeightedSum, CumulativeWakeSum, SquaredSum from py_wake.wind_farm_models.wind_farm_model import WindFarmModel from py_wake.deflection_models.deflection_model import DeflectionModel from py_wake.utils.gradients import cabs @@ -483,6 +483,9 @@ class PropagateUpDownIterative(EngineeringWindFarmModel): alt_model = [self.superpositionModel, LinearSum()][isinstance( self.superpositionModel, (WeightedSum, CumulativeWakeSum))] self.blockage_superpositionModel = self.blockage_deficitModel.superpositionModel or alt_model + msg = "PropagateDownwind does not work with SquaredSum when the blockage model has both up and downstream effects" + assert self.blockage_deficitModel.upstream_only or not isinstance( + self.blockage_superpositionModel, SquaredSum), msg WS_eff_ilk_last = WS_ilk for j in tqdm(range(I), disable=I <= 1 or not self.verbose, desc="Calculate flow interaction", unit="wt"): @@ -1069,7 +1072,7 @@ def main(): # NOJ wake and selfsimilarity blockage noj_ss = All2AllIterative(site, windTurbines, wake_deficitModel=NOJDeficit(), superpositionModel=SquaredSum(), - blockage_deficitModel=SelfSimilarityDeficit()) + blockage_deficitModel=SelfSimilarityDeficit(upstream_only=True)) # Zong convection superposition zongp = PropagateDownwind(site, windTurbines, wake_deficitModel=ZongGaussianDeficit(), superpositionModel=WeightedSum(),