GitHub

Tutorial

This tutorial covers the basics of FLOWFarm. For more specifics refer to the How-to guide.

This tutorial discusses how to do the following with FLOWFarm:

  • (1) setting up a problem description
  • (2) setting up an analysis model set
  • (3) running analyses
  • (4) calculating and visualizing a flow field

Details for setting up an optimization will depend heavily on the optimization package you are using, your objective, and your design variables. Optimization examples using various packages are provided in the example scripts located in the test directory.

(1) Problem description

The problem description involves a definition of the wind farm, the turbines, and the wind resource. While this tutorial uses the same design across all the wind turbines and mostly equal properties across all wind flow states, all turbines and flow states can be unique.

For API demonstration purposes, we have directly assigned all values. However, values may be loaded from .csv and/or .yaml files.

Set up the running environment

using FLOWFarm

Initialize the wind farm design

# set initial turbine x and y locations
turbinex = [-240.0, -240.0, -240.0, 0.0, 0.0, 0.0, 240.0, 240.0, 240.0]
turbiney = [-240.0, 0.0, 240.0, -240.0, 0.0, 240.0, -240.0, 0.0, 240.0]

# get the number of turbines
nturbines = length(turbinex)

# set turbine base heights
turbinez = zeros(nturbines)

# set turbine yaw values
turbineyaw = zeros(nturbines)

# set wind farm boundary parameters in meters (we won't really need this until we optimize)
boundarycenter = [0.0,0.0]
boundaryradius = hypot(300, 300)
nothing

Initialize wind turbine design

# set turbine design parameters (these values correspond to the Vestas V80 turbine)
rotordiameter = zeros(nturbines) .+ 80.0   # m
hubheight = zeros(nturbines) .+ 70.0           # m
cutinspeed = zeros(nturbines) .+ 4.0           # m/s
cutoutspeed = zeros(nturbines) .+ 25.0         # m/s
ratedspeed = zeros(nturbines) .+ 16.0          # m/s
ratedpower = zeros(nturbines) .+ 2.0E6         # W
generatorefficiency = ones(nturbines)
nothing

Determine effective inflow speeds sampling method

Rotor swept area sample points are normalized by the rotor radius. These arrays define which points on the rotor swept area should be used to estimate the effective inflow wind speed for each wind turbine. Values of 0.0 are at the rotor hub, 1.0 is at the blade tip, z is vertical, and y is horizontal. These points track the rotor when yawed. A single sample point will always be placed at the hub. More points can be arranged in either a grid pattern or a sunflower packing pattern with various options. See doc strings for more information.

# get the sample points
nsamplepoints = 50
rotorsamplepointsy, rotorsamplepointsz = FLOWFarm.rotor_sample_points(nsamplepoints, method="sunflower")
([-0.07410836062639709, 0.015218865228956489, 0.13673646546308013, -0.2618432050129595, 0.25440179343349323, -0.08653476830049676, -0.16701958726816815, 0.3656303187083795, -0.38303743548855274, 0.1856809233100207  …  -0.4817598863166705, 0.8830705168069845, -0.8243825104580116, 0.32586882131681083, 0.3574932865011348, -0.8663755084100319, 0.9262022154807192, -0.49530244582727084, -0.20782082756188439, 0.8152470554454784], [0.06788932895734043, -0.1734111197246969, 0.17834850579052913, -0.046316376106210956, -0.16182959681996725, 0.321905024791926, -0.3215860519081978, 0.13352760597105143, 0.15811228520113338, -0.3967893479390924  …  -0.765564647902427, 0.24201301768723907, 0.42305949349566085, -0.8789751930978891, 0.8781733593639224, -0.41059152161141643, -0.28555804214300595, 0.8451458139004672, -0.9677853498733453, 0.5791133210240265])

The sunflower sampling method:

Setting up the wind resource

The wind resource determines the properties of the flowfield at all wind states. A wind state is any combination of wind speed, wind direction, turbulence intensity, etc...

# set flow parameters
windspeed = 8.0        # m/2
airdensity = 1.1716    # kg/m^3
ambientti = 0.1      # %
shearexponent = 0.15
ndirections = 5
winddirections = collect(range(0, 2*pi*(1-1/ndirections), length=ndirections))   # radians
windspeeds = ones(ndirections).*windspeed   # m/s
windprobabilities = ones(ndirections).*(1.0/ndirections)       # %
ambienttis = ones(ndirections).*ambientti  # %
measurementheight = ones(ndirections).*hubheight[1] # m

# initialize the wind shear model
windshearmodel = FLOWFarm.PowerLawWindShear(shearexponent)

# initialize the wind resource definition
windresource = FLOWFarm.DiscretizedWindResource(winddirections, windspeeds, windprobabilities,
measurementheight, airdensity, ambienttis, windshearmodel)
DiscretizedWindResource{Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}, Float64, PowerLawWindShear{Float64, Float64, String}}([0.0, 1.2566370614359172, 2.5132741228718345, 3.7699111843077517, 5.026548245743669], [8.0, 8.0, 8.0, 8.0, 8.0], [0.2, 0.2, 0.2, 0.2, 0.2], [70.0, 70.0, 70.0, 70.0, 70.0], 1.1716, [0.1, 0.1, 0.1, 0.1, 0.1], PowerLawWindShear{Float64, Float64, String}(0.15, 0.0, "first"))

(2) Analysis models

A model set requires a Wake Deficit Model, Wake Deflection Model, Wake Combination Model, and a Local Turbulence Intensity Model. There are several options for each model type. To facilitate research studies, any of the models in each type can be used with any of the models in any other type. However, behavior is not guaranteed. It is recommended that common, validated, model combinations be used in most cases.

Model types and options are:

  • Deficit Models: JensenTopHat, JensenCosine, MultiZone, GaussOriginal, GaussYaw, GaussYawVariableSpread, GaussSimple.
  • Deflection Models: GaussYawDeflection, GaussYawVariableSpreadDeflection, JiminezYawDeflection, MultizoneDeflection.
  • Combination Models: LinearFreestreamSuperposition, SumOfSquaresFreestreamSuperposition SumOfSquaresLocalVelocitySuperposition, LinearLocalVelocitySuperposition.
  • Turbulence Models: LocalTIModelNoLocalTI, LocalTIModelMaxTI.

The model set can be set up as follows:

Initialize the power model (other options can be found in the Reference section).

powermodel = FLOWFarm.PowerModelPowerCurveCubic()
PowerModelPowerCurveCubic{Int64}(2)

The user can define different power models for different wind turbines, but here we use the same power model for every turbine. 

powermodels = Vector{typeof(powermodel)}(undef, nturbines)
for i = 1:nturbines
    powermodels[i] = powermodel
end

Initialize thrust model(s) (other options can be found in the Reference section). 

ctmodel = FLOWFarm.ThrustModelConstantCt(0.65)
ctmodels = Vector{typeof(ctmodel)}(undef, nturbines)
for i = 1:nturbines
    ctmodels[i] = ctmodel
end

Next, we can set up wake and related models. Here we will use the default values provided in FLOWFarm. However, it is important to use the correct model parameters. More information and references are provided in the doc strings attached to each model.

  • The wake deficit model estimates the reduction in wind speed caused by the wake of a wind turbine.
  • The wake deflection model predicts the lateral displacement of the wake center due to yaw misalignment or ambient wind shear.
  • The wake combination model defines how the deficits from multiple wakes are aggregated to compute the total velocity deficit at a given point.
  • Local turbulence intensity models estimate the turbulence intensity at each turbine or evaluation point, providing improved inputs to the wake deficit and deflection models when applicable.
wakedeficitmodel = FLOWFarm.GaussYawVariableSpread()
wakedeflectionmodel = FLOWFarm.GaussYawDeflection()
wakecombinationmodel = FLOWFarm.LinearLocalVelocitySuperposition()
localtimodel = FLOWFarm.LocalTIModelMaxTI()
LocalTIModelMaxTI{Float64, Float64, Float64, Float64}(2.32, 0.154, 0.3837, 0.003678)

Initialize the model set. 

modelset = FLOWFarm.WindFarmModelSet(wakedeficitmodel, wakedeflectionmodel, wakecombinationmodel, localtimodel)
nothing

(3) Running the analysis

Now that the wind farm and analysis models have been defined, we can calculate AEP. The output is in Watt-hours.

aep = FLOWFarm.calculate_aep(turbinex, turbiney, turbinez, rotordiameter,
    hubheight, turbineyaw, ctmodels, generatorefficiency, cutinspeed,
    cutoutspeed, ratedspeed, ratedpower, windresource, powermodels, modelset,
    rotor_sample_points_y=rotorsamplepointsy, rotor_sample_points_z=rotorsamplepointsz)
4.1752381500085173e9

We can also get the AEP in each wind-direction using the following.

state_aeps = FLOWFarm.calculate_state_aeps(turbinex, turbiney, turbinez, rotordiameter,
        hubheight, turbineyaw, ctmodels, generatorefficiency, cutinspeed,
        cutoutspeed, ratedspeed, ratedpower, windresource, powermodels, modelset,
        rotor_sample_points_y=rotorsamplepointsy, rotor_sample_points_z=rotorsamplepointsz,
        hours_per_year=365.25*24.0, weighted=true)
5-element Vector{Float64}:
 4.238504382525765e8
 9.75999028763892e8
 8.989462395040804e8
 8.934110908062623e8
 9.830313526817064e8

If we instead set weighted=false then we would get the power in each direction in Watts.

If we want to get the individual turbine powers in each directions, we use the following.

turbine_powers_by_direction = FLOWFarm.calculate_state_turbine_powers(turbinex, turbiney, turbinez, rotordiameter,
    hubheight, turbineyaw, ctmodels, generatorefficiency, cutinspeed,
    cutoutspeed, ratedspeed, ratedpower, windresource, powermodels, modelset,
    rotor_sample_points_y=rotorsamplepointsy, rotor_sample_points_z=rotorsamplepointsz)
5×9 Matrix{Float64}:
  2544.31   6417.76  71624.0   2544.31  …   2544.31   6417.76  71624.0
 50455.7   50455.7   60389.0  60173.9      71624.0   71624.0   71624.0
 71624.0   41504.8   35763.4  71624.0      71624.0   71616.9   71616.9
 71624.0   71618.3   71618.3  71624.0      71624.0   40492.8   35288.1
 71624.0   71624.0   71624.0  61031.2      50765.9   50765.9   61207.3

The output is a 2D array where each entry represents the power of a turbine, with dimensions corresponding to the number of wind directions and turbines (ndirections x nturbines).

(4) Calculating a flow field

It is helpful to visualize the whole flow-field, not just the turbine powers. The input wind_farm_state_id determines which state from the wind resource is being used to calculate the flow-field. 

# define how many points should be in the flow field
xres = 1000
yres = 1000
zres = 1

# define flow field domain
maxy = boundaryradius*1.5
miny = -boundaryradius*1.5
maxx = boundaryradius*1.5
minx = -boundaryradius*1.5

# set up point grid for flow field
xrange = minx:(maxx-minx)/xres:maxx
yrange = miny:(maxy-miny)/yres:maxy
zrange = hubheight[1]

# run flowfarm
ffvelocities = FLOWFarm.calculate_flow_field(xrange, yrange, zrange,
    modelset, turbinex, turbiney, turbinez, turbineyaw,
    rotordiameter, hubheight, ctmodels, rotorsamplepointsy, rotorsamplepointsz,
    windresource, wind_farm_state_id=5)
1×1001×1001 Array{Float64, 3}:
[:, :, 1] =
 8.0  8.0  8.0  8.0  8.0  8.0  8.0  8.0  …  8.0  8.0  8.0  8.0  8.0  8.0  8.0

[:, :, 2] =
 8.0  8.0  8.0  8.0  8.0  8.0  8.0  8.0  …  8.0  8.0  8.0  8.0  8.0  8.0  8.0

[:, :, 3] =
 8.0  8.0  8.0  8.0  8.0  8.0  8.0  8.0  …  8.0  8.0  8.0  8.0  8.0  8.0  8.0

;;; … 

[:, :, 999] =
 7.88965  7.88497  7.88014  7.87516  …  8.0  8.0  8.0  8.0  8.0  8.0  8.0

[:, :, 1000] =
 7.88796  7.88323  7.87836  7.87333  …  8.0  8.0  8.0  8.0  8.0  8.0  8.0

[:, :, 1001] =
 7.88626  7.88148  7.87656  7.87148  …  8.0  8.0  8.0  8.0  8.0  8.0  8.0

ffvelocities is a three-dimensional array containing flow velocities over the specified domain. The dimensions correspond to the spatial ranges as follows: the first dimension maps to zrange, the second to xrange, and the third to yrange.