Optimization of a wind farm layout
In the following code, I have optimized aep of a power plant. Now I want to change the objective function and, instead of aep, maximize the Power output of each turbine (sim_res.Power
).
To make the code a bit clear, I have to say that: A variable (c) has been defined in which we have two values, 0 and 1, and its length is equal to the number of wind turbines (wt=80). If c=1, we say that that wind turbine is On, otherwise it is off, so power production of that turbine will be zero. In the current code.
Now instead, I want to change "c" ; For each wind turbine, I wanna define some bins for wind speed (WS), and wind direction (wd) for instance for wd the criteria is:[0,360] and we define 12 bins and also for ws:[3,25] is criteria and we divide to 5 bins. Now I wanna change the C by using these bins and we just investigate the mean value of each bin and according to that, we can analyze the values on that range. For instance, for the first bin of ws(=[0,30]), mean=15 so if C for mean is zero we say that wt is off and for each wd in this range, the turbine is off. Also,I want to optimize Power output instead of aep, and I have 3 arguments for it: wt, wd, ws (sim_res.Power
is a function of wt, wd, ws)
C has 3 arguments while we have 80 turbines, so for using it in funC
it makes a problem since x_selected
and y_selected
is for 80 points while c is defined for 80 * 360(wd) * 23 (ws) points. So I do not know how to match these numbers together to maximize sim_res.Power
.
"""
References:
https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.minimize.html
https://github.com/DTUWindEnergy/PyWake
"""
import time
from py_wake.examples.data.hornsrev1 import V80
from py_wake.examples.data.hornsrev1 import Hornsrev1Site # We work with the Horns Rev 1 site, which comes already set up with PyWake.
from py_wake import BastankhahGaussian
from scipy.optimize import minimize
import numpy as np
def funC(x, y, c):
"""
Turns on/off the use of wind turbine depending on the value of c.
scipy generates c real values in the range [0, 1] as specified by the bounds including 0.2 etc.
If c is 0.5 or more turbine will be used otherwise turbine will not be used.
"""
x_selected = x[c >= 0.5]
y_selected = y[c >= 0.5]
return (x_selected, y_selected)
def wt_simulation(c):
"""
This is our objective function. It will return the aep=annual energy production in GWh.
We will maximize aep.
"""
site = Hornsrev1Site()
x, y = site.initial_position.T
windTurbines = V80()
wf_model = BastankhahGaussian(site, windTurbines)
x_new, y_new = funC(x, y, c)
# run wind farm simulation
sim_res = wf_model(
x_new, y_new, # wind turbine positions
h=None, # wind turbine heights (defaults to the heights defined in windTurbines)
type=0, # Wind turbine types
wd=None, # Wind direction (defaults to site.default_wd (0,1,...,360 if not overriden))
ws=None, # Wind speed (defaults to site.default_ws (3,4,...,25m/s if not overriden))
)
aep_output = sim_res.aep().sum() # we maximize aep
return -float(aep_output) # negate because of scipy minimize
def solve():
t0 = time.perf_counter()
wt = 80 # for V80
x0 = np.ones(wt) # initial value
bounds = [(0, 1) for _ in range(wt)]
res = minimize(wt_simulation, x0=x0, bounds=bounds)
print(f'success status: {res.success}')
print(f'aep: {-res.fun}') # negate to get the true maximum aep
print(f'c values: {res.x}\n')
print(f'elapse: {round(time.perf_counter() - t0)}s')
# start
solve()