In this tutorial you will use PyWake, an open-source wind farm simulation tool developed at DTU (Pedersen et al. (2019)), to see how turbines affect each other through their wakes. A turbine extracts energy from the wind, so the air behind it is slower and more turbulent. Downstream turbines therefore produce less power than turbines in undisturbed flow.
After this tutorial you will be able to:
create a PyWake wind farm model from a site, a wind turbine and a wake model;
define a wind farm layout and the inflow conditions (wind speed, wind direction, turbulence intensity);
run the simulation and read the power of each turbine and of the whole farm;
visualise the wake flow field and explain how the wind direction changes the farm power.
1. Set up the model¶
A PyWake model (a wind farm model) combines three ingredients:
a site, describing the wind conditions.
UniformSitehas the same wind everywhere; here we set the turbulence intensity to 0.06;a wind turbine, with its rotor diameter, hub height and power curve. We use the IEA 10 MW reference turbine (198 m rotor diameter, 119 m hub height), with its power and thrust curves taken from the IEA Wind 740-10-MW reference offshore wind plants (Kainz et al. (2024));
a wake model, which computes the velocity deficit behind each turbine. We use the Gaussian wake model of Niayifar and Porté-Agel (Niayifar & Porté-Agel (2016)), in which the wake expands faster when the turbulence intensity is higher. Deficits from several turbines are combined with a sum of squares. These are the same choices as in the default wake model of the FLORIS tutorial, so both tutorials give similar wake losses.
The layout is given as lists of and coordinates in metres.
The turbine data are listed in the first code cell below. A turbine is described by its power curve and thrust coefficient curve, which give the power and the thrust of the rotor for each wind speed.
import numpy as np
# Power and thrust coefficient curves of the IEA 10 MW reference turbine, from the
# IEA Wind 740-10-MW Reference Offshore Wind Plants (file IEA37_10MW_turbine.yaml)
turbine_ws = np.array([
4.000, 4.515, 5.001, 5.457, 5.883, 6.278, 6.640, 6.968,
7.263, 7.523, 7.748, 7.938, 8.091, 8.208, 8.288, 8.331,
8.337, 8.368, 8.436, 8.540, 8.681, 8.859, 9.072, 9.320,
9.604, 9.921, 10.272, 10.656, 10.758, 11.518, 11.994, 12.499,
13.032, 13.592, 14.177, 14.786, 15.418, 16.070, 16.743, 17.434,
18.142, 18.865, 19.602, 20.351, 21.110, 21.877, 22.652, 23.432,
24.215, 25.000,
]) # Wind speed [m/s]
turbine_power = np.array([
387.5, 644.7, 938.0, 1257.9, 1604.1, 1969.5, 2340.1, 2708.6,
3067.1, 3408.7, 3723.7, 4003.4, 4239.6, 4425.8, 4556.6, 4628.1,
4638.4, 4690.0, 4804.9, 4985.7, 5236.9, 5564.4, 5975.8, 6480.5,
7089.7, 7816.4, 8690.6, 9716.7, 10000.0, 10000.0, 10000.0, 10000.0,
10000.0, 10000.0, 10000.0, 10000.0, 10000.0, 10000.0, 10000.0, 10000.0,
10000.0, 10000.0, 10000.0, 10000.0, 10000.0, 10000.0, 10000.0, 10000.0,
10000.0, 10000.0,
]) # Power [kW]
turbine_ct = np.array([
0.7701, 0.7763, 0.7824, 0.7820, 0.7802, 0.7772, 0.7719, 0.7768,
0.7768, 0.7768, 0.7768, 0.7768, 0.7768, 0.7768, 0.7768, 0.7768,
0.7768, 0.7768, 0.7768, 0.7768, 0.7768, 0.7768, 0.7768, 0.7768,
0.7768, 0.7768, 0.7675, 0.7651, 0.7587, 0.5056, 0.4310, 0.3708,
0.3209, 0.2788, 0.2432, 0.2128, 0.1868, 0.1645, 0.1454, 0.1289,
0.1147, 0.1024, 0.0918, 0.0825, 0.0745, 0.0675, 0.0613, 0.0559,
0.0512, 0.0470,
]) # Thrust coefficient [-]
rotor_diameter = 198.0 # [m]
hub_height = 119.0 # [m]
cut_in, cut_out = 4.0, 25.0 # Cut-in and cut-out wind speed [m/s]import matplotlib.pyplot as plt
import numpy as np
import ipywidgets as widgets
from ipywidgets import interact
from py_wake.flow_map import XYGrid
from py_wake.literature.gaussian_models import Niayifar_PorteAgel_2016
from py_wake.rotor_avg_models import RotorCenter
from py_wake.site import UniformSite
from py_wake.superposition_models import SquaredSum
from py_wake.wind_turbines import WindTurbine
from py_wake.wind_turbines.power_ct_functions import PowerCtTabular
site = UniformSite(ti=0.06)
wind_turbine = WindTurbine(
name="IEA 10 MW",
diameter=rotor_diameter,
hub_height=hub_height,
powerCtFunction=PowerCtTabular(
turbine_ws, turbine_power, power_unit="kW", ct=turbine_ct, ws_cutin=4.0, ws_cutout=25.0
),
)
wfm = Niayifar_PorteAgel_2016(
site,
wind_turbine,
a=[0.38, 0.004], # Wake expansion: k = 0.38 * TI + 0.004
superpositionModel=SquaredSum(),
rotorAvgModel=RotorCenter(),
)
# Turbine positions [m]
layout_x = [0, 500, 0, 128, 1000, 900, 1500, 1250]
layout_y = [0, 300, 750, 1400, 0, 567, 888, 1450]
print(f"Rotor diameter: {wind_turbine.diameter():.0f} m, hub height: {wind_turbine.hub_height():.0f} m")
# Plot style: accent colour palette, no grid
from cycler import cycler
plt.rcParams.update({
'axes.prop_cycle': cycler(color=['#017188', '#EC6842', '#00A6D6', '#6CC24A', '#A50034']),
'axes.grid': False,
'axes.spines.top': False,
'axes.spines.right': False,
'lines.linewidth': 2,
})---------------------------------------------------------------------------
ModuleNotFoundError Traceback (most recent call last)
Cell In[2], line 6
2 import numpy as np
3 import ipywidgets as widgets
4 from ipywidgets import interact
5
----> 6 from py_wake.flow_map import XYGrid
7 from py_wake.literature.gaussian_models import Niayifar_PorteAgel_2016
8 from py_wake.rotor_avg_models import RotorCenter
9 from py_wake.site import UniformSite
ModuleNotFoundError: No module named 'py_wake'Let’s check the layout. The numbers are the turbine indices, which we will use later to read the power per turbine. The circles show the rotor size.
def plot_layout(ax):
D = wind_turbine.diameter()
for i, (xi, yi) in enumerate(zip(layout_x, layout_y)):
ax.add_patch(plt.Circle((xi, yi), D / 2, fill=False, color="k"))
ax.text(xi, yi + D / 2 + 30, str(i), ha="center", fontsize=9, color="C1")
fig, ax = plt.subplots()
plot_layout(ax)
ax.set_title("Wind farm layout")
ax.set_xlabel("x [m]")
ax.set_ylabel("y [m]")
ax.set_xlim([-500, 2000])
ax.set_ylim([-500, 2000])
ax.set_aspect("equal", adjustable="box")
plt.show()The power curve gives the power of one turbine as a function of the wind speed. Below the cut-in wind speed (4 m/s) the turbine does not produce power. Above the rated wind speed the power stays at the rated power of 10 MW. The thrust coefficient describes how strongly the rotor slows down the wind, and so how strong the wake is.
fig, axs = plt.subplots(1, 2, figsize=(9, 3.5))
axs[0].plot(turbine_ws, turbine_power / 1e3)
axs[0].set_ylabel("Power [MW]")
axs[1].plot(turbine_ws, turbine_ct, color="C1")
axs[1].set_ylabel("Thrust coefficient [-]")
for ax in axs:
ax.axvline(cut_in, color="k", linestyle=":")
ax.axvline(cut_out, color="k", linestyle=":")
ax.set_xlabel("Wind speed [m/s]")
ax.set_xlim([0, 26])
fig.suptitle("IEA 10 MW reference turbine (dotted: cut-in and cut-out)")
plt.tight_layout()
plt.show()2. Run a single wind condition¶
The inflow is defined by:
the wind speed
wsat hub height [m/s];the wind direction
wd[deg], defined as the direction the wind comes from (meteorological convention): is a westerly wind, blowing towards the east (+), and is a northerly wind, blowing towards the south.
The turbulence intensity is set in the site (0.06). We call the wind farm model like a function with the layout and the conditions. PyWake can evaluate many wind directions and speeds at once; the result sim_res stores the power in W with dimensions (turbine, wind direction, wind speed).
sim_res = wfm(layout_x, layout_y, wd=[270], ws=[8])
turbine_powers = sim_res.Power.sel(wd=270, ws=8) / 1e6
farm_power = sim_res.Power.sel(wd=270, ws=8).sum() / 1e6
for i, p in enumerate(turbine_powers.values):
print(f"Turbine {i}: {p:.2f} MW")
print(f"Farm power: {farm_power.item():.2f} MW")3. Visualise the flow field¶
To visualise the wakes we ask PyWake for a flow map: the wind speed on a horizontal grid. By default the map is computed at hub height. The resolution is kept low so that the calculation is fast.
grid = XYGrid(x=np.linspace(-500, 2000, 100), y=np.linspace(-500, 2000, 100))
flow_map = sim_res.flow_map(grid, wd=270, ws=8)
fig, ax = plt.subplots()
flow_map.plot_wake_map(ax=ax, cmap="bone")
ax.set_title("Wind speed at hub height (wind from the west)")
ax.set_xlabel("x [m]")
ax.set_ylabel("y [m]")
plt.show()4. Interactive: effect of the wind direction¶
Use the sliders to change the wind direction and the wind speed. The farm power is shown in the plot title. The calculation is only updated when you release the slider.
grid_interactive = XYGrid(x=np.linspace(-500, 2000, 50), y=np.linspace(-500, 2000, 50))
def update(wind_direction, wind_speed):
sim_res = wfm(layout_x, layout_y, wd=[wind_direction], ws=[wind_speed])
farm_power = sim_res.Power.sum().item() / 1e6
flow_map = sim_res.flow_map(grid_interactive, wd=wind_direction, ws=wind_speed)
fig, ax = plt.subplots()
flow_map.plot_wake_map(ax=ax, cmap="bone")
ax.set_title(f"Wind from {wind_direction:.0f}$^\\circ$, {wind_speed:.1f} m/s: farm power = {farm_power:.2f} MW")
ax.set_xlabel("x [m]")
ax.set_ylabel("y [m]")
ax.set_xlim([-500, 2000])
ax.set_ylim([-500, 2000])
ax.set_aspect("equal", adjustable="box")
plt.show()
interact(
update,
wind_direction=widgets.FloatSlider(
min=0, max=360, step=30, value=270, description="Direction [deg]",
continuous_update=False, style={"description_width": "initial"},
),
wind_speed=widgets.FloatSlider(
min=4, max=14, step=1, value=8, description="Speed [m/s]",
continuous_update=False, style={"description_width": "initial"},
),
);5. Farm power versus wind direction¶
Instead of moving the slider one value at a time, we can pass many wind directions to PyWake in a single call. Here we sweep the wind direction from to at constant wind speed.
wind_directions = np.arange(0, 360, 1.0)
sim_sweep = wfm(layout_x, layout_y, wd=wind_directions, ws=[8])
farm_power = sim_sweep.Power.sum("wt").sel(ws=8) / 1e6
fig, ax = plt.subplots()
ax.plot(wind_directions, farm_power)
ax.set_xlabel("Wind direction [deg]")
ax.set_ylabel("Farm power [MW]")
ax.set_title("Farm power at 8 m/s")
plt.show()- Pedersen, M. M., van der Laan, P., Friis-Møller, M., Rinker, J., & Réthoré, P.-E. (2019). PyWake: An open-source wind farm simulation tool. Technical University of Denmark. https://gitlab.windenergy.dtu.dk/TOPFARM/PyWake
- Kainz, S., Quick, J., Souza de Alencar, M., Sanchez Perez Moreno, S., Dykes, K., Bay, C. J., Zaaijer, M. B., & Bortolotti, P. (2024). The IEA Wind 740-10-MW Reference Offshore Wind Plants [Techreport]. IEA Wind Task 55, National Laboratory of the Rockies. https://www.nlr.gov/docs/fy24osti/87923.pdf
- Niayifar, A., & Porté-Agel, F. (2016). Analytical Modeling of Wind Farms: A New Approach for Power Prediction. Energies, 9(9), 741. 10.3390/en9090741