Allow weightedSum and cumulative Sum in combination with blockage and PropagateUpDownIterative

py_wake/tests/test_deficit_models/test_blockage_models.py

@@ -4,6 +4,7 @@ import matplotlib.pyplot as plt
 from py_wake import np
 from py_wake.deficit_models import VortexCylinder
 from py_wake.deficit_models.deficit_model import BlockageDeficitModel
+from py_wake.deficit_models.gaussian import BastankhahGaussianDeficit
 from py_wake.deficit_models.hybridinduction import HybridInduction
 from py_wake.deficit_models.no_wake import NoWakeDeficit
 from py_wake.deficit_models.noj import NOJDeficit
@@ -16,11 +17,11 @@ from py_wake.examples.data.hornsrev1 import Hornsrev1Site
 from py_wake.examples.data.iea37._iea37 import IEA37Site, IEA37_WindTurbines
 from py_wake.flow_map import XYGrid
 from py_wake.rotor_avg_models.rotor_avg_model import CGIRotorAvg, RotorCenter
-from py_wake.superposition_models import LinearSum
+from py_wake.superposition_models import LinearSum, WeightedSum, CumulativeWakeSum
 from py_wake.tests import npt
 from py_wake.turbulence_models.stf import STF2017TurbulenceModel
 from py_wake.utils.model_utils import get_models
-from py_wake.wind_farm_models.engineering_models import All2AllIterative
+from py_wake.wind_farm_models.engineering_models import All2AllIterative, PropagateUpDownIterative
 from py_wake.tests.test_wind_farm_models.test_enginering_wind_farm_model import OperatableV80
 from py_wake.deficit_models.utils import ct2a_mom1d
 
@@ -216,3 +217,21 @@ def test_All2AllIterative_Blockage_Deficit_RotorAvg():
             plt.show()
         npt.assert_almost_equal(sim_res.WS_eff[2, 0, 0].item(), ref)
     plt.close('all')
+
+
+@pytest.mark.parametrize('superpositionModel', [LinearSum(), WeightedSum(), CumulativeWakeSum()])
+def test_PropagateUpDownIterative_Blockage(setup, superpositionModel):
+    site, windTurbines = setup
+
+    def get_aep(wfm_cls):
+        wfm = wfm_cls(
+            site=site,
+            windTurbines=windTurbines,
+            wake_deficitModel=BastankhahGaussianDeficit(use_effective_ws=True),
+            superpositionModel=superpositionModel,
+            blockage_deficitModel=SelfSimilarityDeficit2020(use_effective_ws=True),
+        )
+
+        sim_res = wfm(x=[0, 500, 1000, 1500], y=[0, 0, 0, 0], wd=270, WS_eff=0)
+        return sim_res.aep().sum().item()
+    npt.assert_allclose(get_aep(All2AllIterative), get_aep(PropagateUpDownIterative), atol=0.03)

py_wake/wind_farm_models/engineering_models.py

@@ -606,10 +606,11 @@ class PropagateUpDownIterative(EngineeringWindFarmModel):
                                           'WS_eff_xxx': np.array(WS_eff_mk),
                                           'ct_xxx': np.array(ct_jlk),
                                           'D_xx': np.array(D_mk)})
-
-                WS_eff_lk = WS_mk[m] - self.superpositionModel.superpose_deficit(**sp_kwargs)
+                WS_eff_lk = WS_mk[m]
+                if self.direction == 'down':
+                    WS_eff_lk = WS_eff_lk - self.superpositionModel.superpose_deficit(**sp_kwargs)
                 if self.blockage_deficitModel:
-                    WS_eff_lk -= self.blockage_superpositionModel(get_value2WT(blockage_nk))
+                    WS_eff_lk = WS_eff_lk - self.blockage_superpositionModel(get_value2WT(blockage_nk))
                 WS_eff_mk.append(WS_eff_lk)
 
                 if self.turbulenceModel:
@@ -706,7 +707,7 @@ class PropagateUpDownIterative(EngineeringWindFarmModel):
                 # ======================================================================================================
                 # Calculate deficit
                 # ======================================================================================================
-                if isinstance(self.superpositionModel, (WeightedSum, CumulativeWakeSum)):
+                if isinstance(self.superpositionModel, (WeightedSum, CumulativeWakeSum)) and self.direction == 'down':
                     # only cw needs to be rotor averaged as remaining super position input is
                     # the same all over the rotor
                     if self.wake_deficitModel.rotorAvgModel:
@@ -724,9 +725,9 @@ class PropagateUpDownIterative(EngineeringWindFarmModel):
                         deficit = np.zeros_like(sigma_sqr)
                 else:
                     deficit, blockage = self._calc_deficit(**model_kwargs)
+                    if self.blockage_deficitModel:
+                        blockage_nk.append(blockage[0])
                 deficit_nk.append(deficit[0])
-                if self.blockage_deficitModel:
-                    blockage_nk.append(blockage[0])
 
                 if self.turbulenceModel:
 
@@ -784,6 +785,7 @@ class PropagateDownwind(PropagateUpDownIterative):
         EngineeringWindFarmModel.__init__(self, site, windTurbines, wake_deficitModel, superpositionModel, rotorAvgModel,
                                           blockage_deficitModel=None, deflectionModel=deflectionModel,
                                           turbulenceModel=turbulenceModel, inputModifierModels=inputModifierModels)
+        self.direction = 'down'
 
     def _calc_deficit(self, dw_ijlk, **kwargs):
         return EngineeringWindFarmModel._calc_deficit(self, dw_ijlk, **kwargs)
@@ -853,6 +855,7 @@ class All2AllIterative(EngineeringWindFarmModel):
             blockage_deficitModel = self.blockage_deficitModel
             self.blockage_deficitModel = None
             dw_order_indices_ld = self.site.distance.dw_order_indices(wd)[:, 0]
+            self.direction = 'down'
             WS_eff_ilk = PropagateUpDownIterative._propagate_deficit(
                 self, wd, dw_order_indices_ld, WD_ilk=WD_ilk, WS_ilk=WS_ilk, TI_ilk=TI_ilk,
                 WS_eff_ilk=WS_eff_ilk, TI_eff_ilk=TI_eff_ilk, D_i=D_i, I=I, L=L, K=K, **kwargs)[0]