Calculate AEP
PyWAsP can estimate the Annual Energy Production (AEP) of wind turbines.
For now, this includes estimating the gross and potential AEP (as per the workflow below), with only limited support, for now, for estimating other losses and P50 and P90.
We will start by importing needed libraries and some test data from the Serra Santa Luzia site.
In [1]: import numpy as np
In [2]: import pandas as pd
In [3]: import xarray as xr
In [4]: import matplotlib.pyplot as plt
In [5]: import windkit as wk
In [6]: import pywasp as pw
In [7]: bwc = wk.read_bwc("./source/tutorials/data/SerraSantaLuzia.omwc", "EPSG:4326")
In [8]: bwc = wk.spatial.reproject(bwc, "EPSG:32629")
In [9]: elev_map = wk.read_vector_map(
...: "./source/tutorials/data/SerraSantaLuzia.map",
...: crs="EPSG:32629",
...: map_type="elevation"
...: )
...:
In [10]: lc_map, lc_tbl = wk.read_vector_map(
....: "./source/tutorials/data/SerraSantaLuzia.map",
....: crs="EPSG:32629",
....: map_type="roughness"
....: )
....:
In [11]: topo_map = pw.wasp.TopographyMap(elev_map, lc_map, lc_tbl)
We have used the site in the WAsP Flow Model page, but here we will
go further and include some turbines. Their locations are stored in a .csv,
so we will create a xarray.Dataset to hold their spatial information: x, y, and hub height.
We can plot their locations together with the elevation map to get an idea of where they are located:
In [12]: wtg_locs = pd.read_csv('./source/tutorials/data/turbine_positions.csv')
In [13]: output_locs = wk.spatial.create_dataset(
....: wtg_locs.Easting.values,
....: wtg_locs.Northing.values,
....: wtg_locs['Hub height'].values,
....: crs="EPSG:32629"
....: )
....:
In [14]: fig, ax = plt.subplots()
In [15]: ax.set(xlim=(510000, 521500), ylim=(4616000, 4627500))
Out[15]: [(510000.0, 521500.0), (4616000.0, 4627500.0)]
In [16]: elev_map.plot("elev", ax=ax);
In [17]: ax.scatter(
....: output_locs.west_east.values,
....: output_locs.south_north.values,
....: c="k",
....: marker="x",
....: s=15,
....: zorder=2,
....: );
....:
To estimate the AEP, we will fist estimate wind climate at the hub height of each turbine:
In [18]: pwc = pw.wasp.generalize_and_downscale(output_locs, bwc, topo_map)
In [19]: print(pwc)
<xarray.Dataset>
Dimensions: (sector: 12, point: 15)
Coordinates:
* sector (sector) float64 0.0 30.0 60.0 90.0 ... 270.0 300.0 330.0
height (point) int64 50 50 50 50 50 50 50 50 50 50 50 50 50 50 50
south_north (point) int64 4622313 4622199 4622336 ... 4624252 4624142
west_east (point) int64 513914 514161 514425 ... 515808 516060 516295
crs int8 0
sector_ceil (sector) float64 15.0 45.0 75.0 105.0 ... 285.0 315.0 345.0
sector_floor (sector) float64 345.0 15.0 45.0 75.0 ... 255.0 285.0 315.0
Dimensions without coordinates: point
Data variables:
A (sector, point) float32 7.428 7.399 7.054 ... 7.52 7.75 8.305
k (sector, point) float32 2.029 2.037 2.041 ... 2.127 2.123
wdfreq (sector, point) float32 0.06445 0.06485 ... 0.07584 0.0792
site_elev (point) float32 459.7 460.0 460.0 427.1 ... 520.0 520.0 520.0
air_density (point) float32 1.163 1.163 1.163 1.167 ... 1.156 1.156 1.156
wspd (point) float32 7.434 7.257 6.978 7.08 ... 7.548 7.818 8.102
power_density (point) float32 414.6 385.3 340.9 368.9 ... 426.1 481.9 532.1
Attributes:
Conventions: CF-1.8
history: 2024-02-16T10:49:06:\twindkit==0.7.1.dev17+ga2ef68d\t w...
title: WAsP site effects
Package name: windkit
Package version: 0.7.1.dev17+ga2ef68d
Creation date: 2024-02-16T10:49:09
Object type: Met fields
author: DTU CI Config
author_email: pywasp@dtu.dk
institution: DTU Wind
Gross AEP
Before we estiamte the AEP we have to get a Wind Turbine Generator (WTG) corresponding to the turbine models used.
We use windkit.read_wtg() to read the .wtg file.
In [20]: wtg = wk.read_wtg("./source/tutorials/data/Bonus_1_MW.wtg")
In [21]: print(wtg)
<xarray.Dataset>
Dimensions: (mode: 1, wind_speed: 22)
Coordinates:
* wind_speed (wind_speed) float64 4.0 5.0 ... 24.0 25.0
* mode (mode) int64 0
Data variables:
power_output (mode, wind_speed) float64 2.41e+04 ... 1e+06
thrust_coefficient (mode, wind_speed) float64 0.915 ... 0.161
air_density (mode) float64 1.225
stationary_thrust_coefficient (mode) float64 0.161
wind_speed_cutin (mode) float64 4.0
wind_speed_cutout (mode) float64 25.0
rated_power (mode) float64 1e+06
name <U10 'Bonus 1 MW'
rotor_diameter float64 54.2
hub_height float64 50.0
regulation_type int64 2
Attributes:
Conventions: CF-1.8
Package name: windkit
Package version: 0.7.1.dev17+ga2ef68d
Creation date: 2024-02-16T10:49:10
author: DTU CI Config
author_email: pywasp@dtu.dk
institution: DTU Wind
The wtg holds power and thrust curves for a number of wind speeds and operational modes (in this case just 1).
Now we can estimate the gross AEP using gross_aep():
In [22]: aep = pw.wasp.gross_aep(pwc, wtg)
In [23]: print(aep)
<xarray.Dataset>
Dimensions: (point: 15)
Coordinates:
height (point) int64 50 50 50 50 50 50 50 50 50 50 50 50 50 50 50
south_north (point) int64 4622313 4622199 4622336 ... 4624252 4624142
west_east (point) int64 513914 514161 514425 ... 515808 516060 516295
crs int8 0
Dimensions without coordinates: point
Data variables:
gross_AEP (point) float32 2.905 2.767 2.552 2.629 ... 2.996 3.189 3.407
Attributes:
Conventions: CF-1.8
history: 2024-02-16T10:49:06:\twindkit==0.7.1.dev17+ga2ef68d\t w...
title: WAsP site effects
Package name: windkit
Package version: 0.7.1.dev17+ga2ef68d
Creation date: 2024-02-16T10:49:10
Object type: Anual Energy Production
author: DTU CI Config
author_email: pywasp@dtu.dk
institution: DTU Wind
Instead of a WTG, a more flexible option is to use the windkit.WindTurbines object instead. It holds
a collection of WTGs and associated wind turbine locations. This is useful when estimating the AEP of a wind farm
with multiple turbine types:
In [24]: wind_turbines = wk.WindTurbines([(wtg, output_locs)])
In [25]: aep = pw.wasp.gross_aep(pwc, wind_turbines)
In [26]: print(aep)
<xarray.Dataset>
Dimensions: (point: 15)
Coordinates:
height (point) int64 50 50 50 50 50 50 50 50 50 50 50 50 50 50 50
south_north (point) int64 4622313 4622199 4622336 ... 4624252 4624142
west_east (point) int64 513914 514161 514425 ... 515808 516060 516295
crs int8 0
Dimensions without coordinates: point
Data variables:
gross_AEP (point) float32 2.905 2.767 2.552 2.629 ... 2.996 3.189 3.407
Attributes:
Conventions: CF-1.8
history: 2024-02-16T10:49:06:\twindkit==0.7.1.dev17+ga2ef68d\t w...
title: WAsP site effects
Package name: windkit
Package version: 0.7.1.dev17+ga2ef68d
Creation date: 2024-02-16T10:49:10
Object type: Anual Energy Production
author: DTU CI Config
author_email: pywasp@dtu.dk
institution: DTU Wind
Here, a single type of turbine is used, but it is possible to add more turbines to the collection by passing
several tuples of WTGs and locations to the windkit.WindTurbines constructor.
Potential AEP
To add wind farm effects (such as wakes and blockage) to the AEP estimation and obtain the
“potential AEP”, PyWAsP integrates with PyWake. PyWAsP defines a number
of named wind farm models (see potential_aep()), including “PARK2_onshore” and “PARK2_offshore”,
that are indentical to the wake model in the WAsP GUI. Here we will use the onshore version:
In [27]: aep_wakes = pw.wasp.potential_aep(pwc, wind_turbines, wind_farm_model="PARK2_onshore")
In [28]: print(aep_wakes)
<xarray.Dataset>
Dimensions: (sector: 12, point: 15)
Coordinates:
* sector (sector) float64 0.0 30.0 ... 300.0 330.0
height (point) int64 50 50 50 50 ... 50 50 50 50
south_north (point) int64 4622313 4622199 ... 4624142
west_east (point) int64 513914 514161 ... 516295
crs int8 0
* point (point) int64 0 1 2 3 4 ... 10 11 12 13 14
Data variables: (12/15)
potential_AEP_sector (point, sector) float64 0.149 ... 0.2285
gross_AEP_sector (point, sector) float64 0.1498 ... 0.2289
AEP_deficit_sector (point, sector) float64 0.004942 ... 0.0...
wspd_sector (point, sector) float64 6.582 ... 7.356
wspd_eff_sector (point, sector) float64 6.582 ... 7.356
wspd_deficit_sector (point, sector) float64 0.0 ... 0.0
... ...
gross_AEP (point) float64 2.917 2.778 ... 3.196 3.414
AEP_deficit (point) float64 0.04692 0.04875 ... 0.03774
wspd (point) float64 7.434 7.257 ... 7.818 8.102
wspd_eff (point) float64 7.309 7.14 ... 7.611 7.958
wspd_deficit (point) float64 0.01856 0.01824 ... 0.01793
turbulence_intensity_eff (point) float64 0.0 0.0 0.0 ... 0.0 0.0 0.0
Attributes:
Conventions: CF-1.8
history: 2024-02-16T10:49:06:\twindkit==0.7.1.dev17+ga2ef68d\t w...
title: WAsP site effects
Package name: windkit
Package version: 0.7.1.dev17+ga2ef68d
Creation date: 2024-02-16T10:49:10
Object type: Anual Energy Production
author: DTU CI Config
author_email: pywasp@dtu.dk
institution: DTU Wind
Wind farm flow map
To map an area around a wind farm PyWAsP has wind_farm_flow_map(), which
takes a resource grid as input and enriches it with wind farm variables including
effects of nearby turbines:
In [29]: output_grid = wk.spatial.create_dataset(
....: np.linspace(511700.0, 517700.0, 41),
....: np.linspace(4619500.0, 4625500.0, 41),
....: [50.0],
....: crs="EPSG:32629",
....: struct="cuboid",
....: )
....:
In [30]: pwc_grid = pw.wasp.downscale(gwc, topo_map, output_grid, genwc_interp="nearest")
In [31]: flow_map = pw.wasp.wind_farm_flow_map(
....: pwc_grid,
....: wind_turbines,
....: output_grid,
....: wind_farm_model="PARK2_onshore",
....: )
....:
In [32]: print(flow_map)
<xarray.Dataset>
Dimensions: (sector: 12, height: 1, south_north: 41,
west_east: 41)
Coordinates:
* sector (sector) float64 0.0 30.0 ... 300.0 330.0
* height (height) float64 50.0
* south_north (south_north) float64 4.62e+06 ... 4.626...
* west_east (west_east) float64 5.117e+05 ... 5.177e+05
crs int8 0
Data variables: (12/15)
potential_AEP_sector (height, south_north, west_east, sector) float64 ...
gross_AEP_sector (height, south_north, west_east, sector) float64 ...
AEP_deficit_sector (height, south_north, west_east, sector) float64 ...
wspd_sector (height, south_north, west_east, sector) float64 ...
wspd_eff_sector (height, south_north, west_east, sector) float64 ...
wspd_deficit_sector (height, south_north, west_east, sector) float64 ...
... ...
gross_AEP (height, south_north, west_east) float64 ...
AEP_deficit (height, south_north, west_east) float64 ...
wspd (height, south_north, west_east) float64 ...
wspd_eff (height, south_north, west_east) float64 ...
wspd_deficit (height, south_north, west_east) float64 ...
turbulence_intensity_eff (height, south_north, west_east) float64 ...
Attributes:
Conventions: CF-1.8
history: 2024-02-16T10:49:10:\twindkit==0.7.1.dev17+ga2ef68d\twk...
Package name: windkit
Package version: 0.7.1.dev17+ga2ef68d
Creation date: 2024-02-16T10:49:44
Object type: Wind Farm Flow Map
author: DTU CI Config
author_email: pywasp@dtu.dk
institution: DTU Wind
title: WAsP site effects
We can plot the wind speed deficit for one direction to show the wake effects:
In [33]: flow_map["wspd_deficit_sector"].sel(sector=210).plot(vmax=0.1)
Using custom PyWake wind farm models
Instead of using a named wind farm model, it is possible to construct a
PyWake wind farm model manually, using functools.partial(), and use it for adding wind farm effects in PyWAsP:
In [34]: from functools import partial
In [35]: import py_wake
In [36]: wind_farm_model = partial(
....: py_wake.wind_farm_models.engineering_models.All2AllIterative,
....: wake_deficitModel=py_wake.deficit_models.noj.NOJLocalDeficit(
....: a=[0, 0.06],
....: ct2a=py_wake.deficit_models.utils.ct2a_mom1d,
....: use_effective_ws=True,
....: use_effective_ti=False,
....: rotorAvgModel=py_wake.rotor_avg_models.AreaOverlapAvgModel(),
....: ),
....: superpositionModel=py_wake.superposition_models.LinearSum(),
....: blockage_deficitModel=py_wake.deficit_models.Rathmann(),
....: deflectionModel=None,
....: turbulenceModel=None,
....: )
....:
In [37]: aep_wakes = pw.wasp.potential_aep(pwc, wtg, wind_farm_model=wind_farm_model)
In [38]: print(aep_wakes)
<xarray.Dataset>
Dimensions: (sector: 12, point: 15)
Coordinates:
* sector (sector) float64 0.0 30.0 ... 300.0 330.0
height (point) int64 50 50 50 50 ... 50 50 50 50
south_north (point) int64 4622313 4622199 ... 4624142
west_east (point) int64 513914 514161 ... 516295
crs int8 0
* point (point) int64 0 1 2 3 4 ... 10 11 12 13 14
Data variables: (12/15)
potential_AEP_sector (point, sector) float64 0.149 ... 0.229
gross_AEP_sector (point, sector) float64 0.1498 ... 0.2289
AEP_deficit_sector (point, sector) float64 0.005177 ... -0....
wspd_sector (point, sector) float64 6.582 ... 7.356
wspd_eff_sector (point, sector) float64 6.581 ... 7.364
wspd_deficit_sector (point, sector) float64 9.693e-05 ... -0...
... ...
gross_AEP (point) float64 2.917 2.778 ... 3.196 3.414
AEP_deficit (point) float64 0.05398 0.05819 ... 0.04531
wspd (point) float64 7.434 7.257 ... 7.818 8.102
wspd_eff (point) float64 7.287 7.112 ... 7.575 7.925
wspd_deficit (point) float64 0.02182 0.02262 ... 0.02194
turbulence_intensity_eff (point) float64 0.0 0.0 0.0 ... 0.0 0.0 0.0
Attributes:
Conventions: CF-1.8
history: 2024-02-16T10:49:06:\twindkit==0.7.1.dev17+ga2ef68d\t w...
title: WAsP site effects
Package name: windkit
Package version: 0.7.1.dev17+ga2ef68d
Creation date: 2024-02-16T10:49:44
Object type: Anual Energy Production
author: DTU CI Config
author_email: pywasp@dtu.dk
institution: DTU Wind