From ebb6080045f69de4791a5f0805bcff6f6c034aa8 Mon Sep 17 00:00:00 2001
From: "Mads M. Pedersen" <mmpe@dtu.dk>
Date: Tue, 28 Jan 2025 14:03:21 +0000
Subject: [PATCH] Jimenez local wake deflection
---
py_wake/deflection_models/deflection_model.py | 3 +-
py_wake/deflection_models/jimenez.py | 24 +++++++++
py_wake/flow_map.py | 2 +-
.../test_deflection_models.py | 5 +-
.../test_deflection_models/test_jimenez.py | 49 ++++++++++++++++++-
py_wake/utils/model_utils.py | 4 +-
6 files changed, 81 insertions(+), 6 deletions(-)
diff --git a/py_wake/deflection_models/deflection_model.py b/py_wake/deflection_models/deflection_model.py
index 3a5e8ddb1..c6ed7dbcf 100644
--- a/py_wake/deflection_models/deflection_model.py
+++ b/py_wake/deflection_models/deflection_model.py
@@ -6,6 +6,7 @@ from numpy import newaxis as na
class DeflectionModel(ABC):
+ additional_args = set()
@property
def args4model(self):
@@ -13,7 +14,7 @@ class DeflectionModel(ABC):
@property
def args4deflection(self):
- return method_args(self.calc_deflection)
+ return method_args(self.calc_deflection) | self.additional_args
@abstractmethod
def calc_deflection(self, dw_ijl, hcw_ijl, dh_ijl, **kwargs):
diff --git a/py_wake/deflection_models/jimenez.py b/py_wake/deflection_models/jimenez.py
index bcbd38354..fbbbbeacb 100644
--- a/py_wake/deflection_models/jimenez.py
+++ b/py_wake/deflection_models/jimenez.py
@@ -23,6 +23,30 @@ class JimenezWakeDeflection(DeflectionIntegrator):
return np.tan(alpha_ijlkx)
+class JimenezLocalWakeDeflection(DeflectionIntegrator):
+ """Implemented according to
+ Jiménez, Ã., Crespo, A. and Migoya, E. (2010), Application of a LES technique to characterize
+ the wake deflection of a wind turbine in yaw. Wind Energ., 13: 559-572. doi:10.1002/we.380
+ """
+
+ def __init__(self, N=20, a=[0.38, 4e-3], use_effective_ti=False):
+ DeflectionIntegrator.__init__(self, N)
+ self.a = a
+ self.TI_key = ['TI_ilk', 'TI_eff_ilk'][use_effective_ti]
+
+ @property
+ def additional_args(self):
+ return {self.TI_key}
+
+ def get_deflection_rate(self, theta_ilk, ct_ilk, D_src_il, dw_ijlkx, **kwargs):
+ denominator_ilk = np.cos(theta_ilk)**2 * np.sin(theta_ilk) * (ct_ilk / 2)
+ TI_ref_ilk = kwargs[self.TI_key]
+ beta_ilk = TI_ref_ilk * self.a[0] + self.a[1]
+ nominator_ijlkx = (1 + (beta_ilk / D_src_il[:, :, na])[:, na, :, :, na] * np.maximum(dw_ijlkx, 0))**2
+ alpha_ijlkx = denominator_ilk[:, na, :, :, na] / nominator_ijlkx
+ return np.tan(alpha_ijlkx)
+
+
def main():
if __name__ == '__main__':
from py_wake import Fuga
diff --git a/py_wake/flow_map.py b/py_wake/flow_map.py
index 30e25e6db..0f228e854 100644
--- a/py_wake/flow_map.py
+++ b/py_wake/flow_map.py
@@ -198,7 +198,7 @@ class FlowMap(FlowBox):
if plot_windturbines:
self.plot_windturbines(normalize_with=normalize_with, ax=ax)
- plt.axis('equal')
+ ax.axis('equal')
return c
def plot_windturbines(self, normalize_with=1, ax=None):
diff --git a/py_wake/tests/test_deflection_models/test_deflection_models.py b/py_wake/tests/test_deflection_models/test_deflection_models.py
index 6e36c0e3d..4c688a02d 100644
--- a/py_wake/tests/test_deflection_models/test_deflection_models.py
+++ b/py_wake/tests/test_deflection_models/test_deflection_models.py
@@ -175,9 +175,10 @@ def test_upstream_deflection():
ref = {'None': 0.0,
'FugaDeflection': 1.6,
'GCLHillDeflection': -1.6,
- 'JimenezWakeDeflection': -18}[l]
+ 'JimenezWakeDeflection': -18.4,
+ 'JimenezLocalWakeDeflection': -18.4}[l]
- # npt.assert_almost_equal(ref, blockage_center, err_msg=l)
+ npt.assert_almost_equal(ref, blockage_center, err_msg=l)
if plot:
fm.WS_eff.squeeze().plot(ax=ax, label=l)
plt.figure()
diff --git a/py_wake/tests/test_deflection_models/test_jimenez.py b/py_wake/tests/test_deflection_models/test_jimenez.py
index 1e38c5238..899152dd2 100644
--- a/py_wake/tests/test_deflection_models/test_jimenez.py
+++ b/py_wake/tests/test_deflection_models/test_jimenez.py
@@ -1,7 +1,7 @@
from py_wake.deficit_models.gaussian import IEA37SimpleBastankhahGaussianDeficit
from py_wake.examples.data.hornsrev1 import V80
from py_wake.examples.data.iea37._iea37 import IEA37Site
-from py_wake.deflection_models.jimenez import JimenezWakeDeflection
+from py_wake.deflection_models.jimenez import JimenezWakeDeflection, JimenezLocalWakeDeflection
from py_wake.flow_map import YZGrid, XYGrid, XZGrid
import matplotlib.pyplot as plt
from py_wake import np
@@ -9,6 +9,7 @@ import pytest
from py_wake.tests import npt
from py_wake.wind_farm_models.engineering_models import PropagateDownwind
from py_wake.superposition_models import SquaredSum
+from py_wake.turbulence_models.crespo import CrespoHernandez
@pytest.mark.parametrize('yaw,tilt,cx,cz', [
@@ -45,3 +46,49 @@ def test_yaw_tilt(yaw, tilt, cx, cz):
npt.assert_allclose(wake_center.y, cx, atol=.1)
npt.assert_allclose(wake_center.h, cz, atol=.24)
+
+
+def test_JimenezLocal():
+ x, y = [0, 500], [0, 0]
+ windTurbines = V80()
+ yaw = [0, 30]
+ tilt = 0
+ D = windTurbines.diameter()
+ plot = 0
+ if plot:
+ ax = plt.gca()
+ for ti, ref_amb, ref_eff in [(.06, -35, -29), (.1, -32, -28),
+ ((.1 - 4e-3) / .38, -26, -24) # beta=.1
+ ]:
+ for ref, use_effective_ti in [(ref_amb, False), (ref_eff, True)]:
+ site = IEA37Site(16, ti=ti)
+ wfm = PropagateDownwind(site, windTurbines, IEA37SimpleBastankhahGaussianDeficit(),
+ superpositionModel=SquaredSum(),
+ deflectionModel=JimenezLocalWakeDeflection(use_effective_ti=use_effective_ti),
+ turbulenceModel=CrespoHernandez())
+
+ x_lst = np.arange(-200, 1010, 10)
+ y_lst = np.linspace(-100, 100, 201)
+ z_lst = np.arange(10, 200, 1)
+
+ fm_xy = wfm(x, y, wd=270, ws=10, yaw=yaw, tilt=tilt).flow_map(XYGrid(x=x_lst))
+ fm_yz = wfm(x, y, wd=270, ws=10, yaw=yaw, tilt=tilt).flow_map(YZGrid(x=500 + 5 * D, y=y_lst, z=z_lst))
+ wake_center = fm_yz.WS_eff.where(fm_yz.WS_eff == fm_yz.WS_eff.min(), drop=True).squeeze()
+
+ if plot:
+ wc = fm_xy.min_WS_eff(x=x_lst)
+ c = ax.plot(wc.x, wc, label=f'{ti}, {("amb", "eff")[use_effective_ti]}')[0].get_color()
+ ax.plot(wake_center.x, wake_center.y, '.', color=c)
+
+ axes = plt.subplots(2, 1)[1].flatten()
+ fm_xy.plot_wake_map(ax=axes[0])
+ axes[0].axvline(5 * D)
+ axes[0].axvline(5 * D + 500)
+ axes[0].plot(wake_center.x, wake_center.y, '.')
+ fm_yz.plot_wake_map(ax=axes[1])
+ axes[1].plot(wake_center.y, wake_center.h, '.')
+ # print(wake_center.y.item())
+ npt.assert_allclose(wake_center.y, ref, atol=.1)
+ if plot:
+ ax.legend()
+ plt.show()
diff --git a/py_wake/utils/model_utils.py b/py_wake/utils/model_utils.py
index 58cda481b..6ecf0abf4 100644
--- a/py_wake/utils/model_utils.py
+++ b/py_wake/utils/model_utils.py
@@ -263,7 +263,9 @@ def get_model_input(wfm, x, y, ws=10, wd=270, **kwargs):
('tilt_ilk', 'tilt'),
('WS_ilk', 'WS'),
('WS_eff_ilk', 'WS_eff'),
- ('ct_ilk', 'CT')]
+ ('ct_ilk', 'CT'),
+ ('TI_ilk', 'TI'),
+ ('TI_eff_ilk', 'TI_eff')]
if n in sim_res})
args['IJLK'] = (1, len(x), len(wd), len(ws))
return args
--
GitLab