OWENSAero

Environment(rho::TF,mu::TF,Vx::TAF #Vinf is Vx,Vy::TAF,Vz::TAF,Vtwist::TAF,windangle::TF #radians,DSModel::TS,AModel::TS,awwarm::TVF,steplast::TAI,idxRPI::TAI,Vwakeold::TVF2,BVDynamicFlagL::TAI,BVDynamicFlagD::TAI,alphalast::TAF2,suction::TB) Environment(rho,mu,Vx,Vy,Vz,Vtwist,windangle,DSModel,AModel,awwarm) = Environment(rho,mu,Vx,Vy,Vz,Vtwist,windangle,DSModel,AModel,awwarm,zeros(Int,1),zeros(Int,length(Vx)),deepcopy(Vx),zeros(Int,1),zeros(Int,1),zeros(Real,1),false) Environment(rho,mu,Vx,DSModel,AModel,awwarm) = Environment(rho,mu,Vx,zeros(Real,size(Vx)),zeros(Real,size(Vx)),zeros(Real,size(Vx)),0.0,DSModel,AModel,awwarm,zeros(Int,1),zeros(Int,length(Vx)),deepcopy(Vx),zeros(Int,1),zeros(Int,1),zeros(Real,1),false)

Contains specications for turbine slice environment/operating conditions as well as some backend memory for dynamic stall and unsteady calculations

Inputs

• rho::TF: Working fluid density (kg/m^3)
• mu::TF: Working fluid viscosity (standard SI units)
• V_x::TAF Vinf is Vx for simple simulations (m/s), array corresponding to each azimuthal position
• V_y::TAF: y input velocity (m/s), array corresponding to each azimuthal position
• V_z::TAF: z input velocity (m/s), array corresponding to each azimuthal position
• V_twist::TAF: rotational velocity from active twist (rad/s), array corresponding to each azimuthal position
• windangle::TF: angle of mean oncoming wind (rad)
• DSModel::TS: dynamic stall model ("BV" or "none" or "LB" - once it is finished)
• AModel::TS: aero model used ("DMS" or "AC")
• aw_warm::TVF: warm start induction factor array, first half corresponding to u, second half to v
• steplast::TAI: prior simulation step index, used for unsteady wake propogation
• idx_RPI::TAI: used to specify the azimuthal indices needed for a partial solve (i.e. not every azimuthal index), such as is used in the RPI method
• V_wake_old::TVF2: Prior step's mean wake velocity (m/s)
• BV_DynamicFlagL::TAI: Boeing-vertol dynamic stall lift flag
• BV_DynamicFlagD::TAI: Boeing-vertol dynamic stall drag flag
• alpha_last::TAF2: Boeing-vertol dynamic stall prior step's angle of attack
• suction::TB: DMS flag for alternate induction model

Outputs:

• none:

Turbine(R::TF,r::TAF,z::TF,chord::TAF3,twist::TAF5,delta::TAF,omega::TAF4,B::TI,af::TFN,ntheta::TI,r_delta_influence::TB,centerX::TAF2,centerY::TAF2)
Turbine(R,r,z,chord,twist,delta,omega,B,af,ntheta,r_delta_infl) = Turbine(R,r,z,chord,twist,delta,omega,B,af,ntheta,r_delta_infl,zeros(Real,size(R)),zeros(Real,size(R)))
Turbine(R,r,chord,twist,delta,omega,B,af,ntheta,r_delta_infl) = Turbine(R,r,1.0,chord,twist,delta,omega,B,af,ntheta,r_delta_infl,zeros(Real,size(R)),zeros(Real,size(R)))

Contains specications for turbine slice (geometry, location, airfoil)

Inputs

• R::TF: Nominal turbine radius (m)
• r::TAF: Array of local radaii corresponding to each azimuthal position for the slice, allows for active blade deformation (m)
• z::TF: Vertical location of slice (only used when calling inflow-wind turbulent input)(m)
• chord::TAF3: Array of chord corresponding to each azimuthal position, allows for active blade deformation (m)
• twist::TAF5: Array of blade twist corresponding to each azimuthal position, allows for active blade deformation (rad)
• delta::TAF: Array of blade slope corresponding to each azimuthal position, allows for active blade deformation (rad)
• omega::TAF4: Array of rotational rate corresponding to each azimuthal position, allows for active blade deformation (rad/s)
• B::TI: Number of blades
• af::TFN: Airfoil function - see tests for example of how to create
• ntheta::TI: Number of azimuthal discretizations
• r_delta_influence::TB: Specification of whether local radius and blade slope are used in the influence coefficients for the actuator cylinder method
• centerX::TAF2: Turbine center x location (only used if multiple turbines are modeled)
• centerY::TAF2: Turbine center y location (only used if multiple turbines are modeled)

Outputs:

• none:

UnsteadyParams(RPI::TB,tau::TAF,ifw::TB,IECgust::TB,nominalVinf::TF,G_amp::TF,gustX0::TF,gustT::TF) UnsteadyParams(RPI,tau,ifw) = UnsteadyParams(RPI,tau,ifw,false,1.0,0.0,1.0,1.0)

Contains specications for turbine slice unsteady inputs

Inputs

• RPI::TB: Flag to specify if RPI is being used
• tau::TAF: Unsteady method wake propogation weighting [3.0,0.3]
• ifw::TB: Flag to specify if inflow-wind is being used
• IECgust::TB: Flag to specify if the simple sin-cos gust profile in the x-direction will be used
• nominalVinf::TF: Nominal velocity used to calculate the IEC gust size (m/s)
• G_amp::TF: IEC gust amplitude (m/s)
• gustX0::TF: IEC gust normalized starting point (x-location divided by reference radius)
• gustT::TF: IEC gust duration (s)

Outputs:

• none:

AC(turbines, env; w=zeros(Real,2turbines[1].ntheta), idx_RPI=1:2turbine.ntheta, solve=true, ifw=false)

see ?steady for detailed i/o description

Double multiple streamtube model

integrand used for computing AIJ

integrand used for computing Ay

Boeing_Vertol(af,alpha,adotnorm,umach,Re,aoaStallPos,aoaStallNeg,AOA0,tc,BV_DynamicFlagL,BV_DynamicFlagD; family_factor = 0.0)

Boeing-Vertol Dynamic Stall Model. All angles are in rad unless explicitely stated otherwise (e.g. alpha_d) Arguments

• af::airfoil_data4D: airfoil function callable by: CL, CD, CM = af(aoa,Re,mach,family_factor)
• alpha::Float64: Static Angle of Attack (at 0.75 chord)
• adotnorm::Float64: Normalized Change in Angle of Attack adotc/(2U)
• umach::Float64: Blade mach number
• Re::Float64: Blade Reynolds number
• aoaStallPos::Float64: Positive Stall Angle (onset)
• aoaStallNeg::Float64: Negative Stall Angle (onset)
• AOA0::Float64: Zero Lift AOA
• tc::Float64: Thickness to chord ratio
• BV_DynamicFlagL::Int: lagged dynamic stall state for lift
• BV_DynamicFlagD::Int: lagged dynamic stall state for drag
• family_factor::float64: factor indexing airfoil family, if used

DMS(turbine, env; w=0, idx_RPI=1:turbine.ntheta, solve=true)

see ?steady for detailed i/o description

Double multiple streamtube model

integrand used for computing DxII

integrand used for computing DxIJ

integrand used for computing Dx

Unsteady_Step(turbine,env,us_param,mystep)

calls inflow wind init

Inputs

• turbine::Turbine: turbine input for slice see ?Turbine
• env::Env: environment input for slice see ?Env
• us_param::UnsteadyParams: unsteady inputs for slice see ?UnsteadyParams
• mystep::int: continuous index cooresponding to the azimuthal discretation - i.e. for ntheta of 30 step 1 is the first step of rev 1, sep 31 is the first step of rev 2, etc. Keeps track of temporal locaion

Outputs:

• CP: This slice's coefficient of performance at this step
• Th: This slice's thrust coefficient at this step
• Q: Torque (N0m) at this step
• Rp: Radial force per height (N) at this step
• Tp: Tangential force per height (N) at this step
• Zp: Vertical force per height (N) at this step
• Vloc: Local velocity array for each azimuthal position (includes induction) (m/s) at this step
• CD: This slice's drag coefficient at this step
• CT: This slice's thrust coefficient (should equal drag, but may no depending on usage or solver status) at this step
• amean: Mean turbine induction in the streamwise direction at this step
• astar: Solved induction factors for each azimuthal location. First half are streamwise (u), second are cross-steam (v) at this step
• alpha: Local angle of attack array for each azimuthal position (includes induction) (rad) at this step
• cl: Local lift coefficient used for each azimuthal position at this step
• cd_af: Local drag coefficient used for each azimuthal position at this step
• thetavec: Azimuthal location of each discretization (rad)
• Re: Reynolds number for each azimuthal position at this step

integrand used for computing WxII

integrand used for computing WxIJ

advanceTurb(tnew;ts=2*pi/(turbslices[1].omega[1]*turbslices[1].ntheta))

Runs a previously initialized aero model (see ?setupTurb) in the unsteady mode (can be repeateadly called, or called for a specific time, or repeatedly called for sections of time)

Inputs

• tnew::float: new time (s); will run from last time specified from the last call, to the current time specified, or from t=ts if the first time called
• ts::float: optional, desired timestep. Will run at finer timesteps than the azimuthal discretization without interfering with wake propogation. While possible, it is not recommended to run with timesteps larger than the azimuthal discretization (hence the optional nature and automatic calculation)

Outputs:

• CP: Turbine coefficient of performance
• Rp: Array(B,Nslices,nsteps) of radial force (N) where nsteps = max(1,round(Int,(tnew-timelast)/ts))
• Tp: Array(B,Nslices,n_steps) of tangential force (N)
• Zp: Array(B,Nslices,n_steps) of vertical force (N)
• alpha: Array(B,Nslices,n_steps) of angle of attack (rad)
• cl: Array(B,Nslices,n_steps) of airfoil cl used
• cd_af: Array(B,Nslices,n_steps) of airfoil cd used
• Vloc: Array(B,Nslices,n_steps) of airfoil local velocity used
• Re: Array(B,Nslices,n_steps) of airfoil Reynolds number used
• thetavec: Azimuthal discretization location (rad)
• ntheta: number of azimuthal discretizations used
• Fx_base: Array(ntheta)Turbine base Fx (N)
• Fy_base: Array(ntheta)Turbine base Fy (N)
• Fz_base: Array(ntheta)Turbine base Fz (N)
• Mx_base: Array(ntheta)Turbine base Mx (N-m)
• My_base: Array(ntheta)Turbine base My (N-m)
• Mz_base: Array(ntheta)Turbine base Mz (N-m)
• power: Array(ntheta)Turbine power (watts)
• power2: Turbine average power for the revolution (watts)
• torque: Array(ntheta)Turbine torque (N-m) (alternative calculation method from Mz-base)

deformTurb(azi;newOmega=-1,newVinf=-1,bldx=-1, bldz=-1, bld_twist=-1, steady=false)

Equivalent to an update states call, mutating the internal aerodynamic inputs within the unsteady model.

Inputs

• azi: Current azimuth position of the turbine in radians (continuously growing with numbers of revolutions)
• bld_x: Blade structural x shape, size(NBlade,any), any as it is splined against bld_z and the aero discretization
• bld_z: Blade structural z shape, size(NBlade,any), any as it is splined against bld_x and the aero discretization
• bld_twist: Blade structural twist, size(NBlade,any), any as it is splined against bld_z and the aero discretization. Note that in the calcs, this will be in addition to the aero twist offset already applied in initialization.
• accel_flap_in: Blade structural acceleration in the flap direction, size(NBlade,any), any as it is splined against bld_z and the aero discretization
• accel_edge_in: Blade structural acceleration in the edge direction, size(NBlade,any), any as it is splined against bld_z and the aero discretization
• steady::bool: if steady is true, it just updates a single step. TODO: verify this is correct

Outputs:

• none:

Internal, assembles the matrices of multiple turbine systems into a combined system centerX, centerY: array of x,y coordinates for centers of the VAWTs in the farm radii: corresponding array of their radii

Internal, integration for a periodic function where end points don't reach ends (uses trapezoidal method)

applies for both Ay and Rx depending on which function ifunc(x, y, phi) is passed in

Internal, precomputes influence coefficient matricies and saves them as HDF5 files

Internal, calculates the radial force used in the residual function as well as the turbine performance when converged

readaerodyn(filename)

create airfoil lookup for a file with only one reynolds number

Inputs

• filename::string: file path/name to airfoil file formatted like in the test folder

Outputs:

• af::function: cl, cd = af(alpha,re,mach) with alpha in rad

readaerodyn_BV(filename)

create airfoil lookup function with boeing vertol dynamic stall model for a file with only one reynolds number

Inputs

• filename::string: file path/name to airfoil file formatted like in the test folder

Outputs:

• af::function: cl, cd = afBV(alpha,Re,M,env,Vtwist,c,dt,U;solvestep=false) with alpha in rad, OWENSAero.Env, V_twist in rad/s, c chord in m, dt in sec, U Vloc in m/s, solvestep true during solve loop

readaerodyn_BV_NEW(filename;DSModel="BV")

for a file with multiple reynolds numbers create airfoil lookup function with boeing vertol dynamic stall model and wrap interpolation

Inputs

• filename::string: file path/name to airfoil file formatted like in the test folder
• DSModel::string: "BV" or "none"

Outputs:

• af::function: cl, cd = afBV(alpha,Re,M,env,Vtwist,c,dt,U;solvestep=false) with alpha in rad, OWENSAero.Env, V_twist in rad/s, c chord in m, dt in sec, U Vloc in m/s, solvestep true during solve loop
• af::function: cl, cd = af(alpha,re,mach) with alpha in rad

Internal, sets up the residual function

setupTurb(bldx,bldz,B,chord,omega,Vinf; Height = maximum(bldz), Radius = maximum(bldx), eta = 0.25, twist = 0.0, #or array{Float,Nslices} rho = 1.225, mu = 1.7894e-5, RPI = true, tau = [0.3,3.0], ntheta = 30, Nslices = 30, #TODO: make this different from ntheta ifw = false, DSModel = "BV", AModel = "DMS", windangleD = 0.0, afname = "(path)/airfoils/NACA0015RE3E5.dat", #TODO: analytical airfoil as default turbsimfilename = "(path)/data/ifw/turbDLC1p313mps330mseed1.bts", ifwlibfile = joinpath(dirname(@FILE), "../bin/libifwcbinding"), AMflag = false, buoyflag = false, rotAccelflag = false, AMCoeffCa = 1.0)

Initializes aerodynamic models and sets up backend persistent memory to simplify intermittent calling within coupled solver loops

Inputs

• bld_x: Blade x shape
• bld_z: Blade z shape
• B: Number of blades
• chord: chord length (m)
• omega: rotation rate in rad/s. size(1) or size(ntheta), pass in an array(Real,ntheta) when propogating automatic gradients
• Vinf: Inflow velocity
• Height: turbine total height (m) typically maximum(bldz) unless only the shape and not size of bldz is being used
• Radius: turbine nominal radius (m) typically maximum(bldx) unless only shape and not size of bldx is used
• eta: blade mount point ratio, i.e. 0.25 would be at the quarter chord
• twist: 0.0, #or array{Float,Nslices}
• rho: working fluid density (kg/m^3)
• mu: working fluid dynamic viscosity (Pa*s)
• RPI: RPI method flag
• tau: Unsteady wake propogation time constants [0.3,3.0],
• ntheta: Number of azimuthal discretizations
• Nslices: Number of vertical slices of the turbine
• ifw: flag for inflow wind
• DSModel: Dynamic stall model "BV" or "none" or "LB" when we get it working
• AModel: Aerodynamic model "DMS" or "AC"
• windangle_D: Inflow wind angle (degrees)
• afname: airfoil path and name e.g. "(path)/airfoils/NACA0015RE3E5.dat"
• turbsim_filename: turbsim path and name e.g. "(path)/data/ifw/turbDLC1p313mps330mseed1.bts",
• ifw_libfile: inflow wind dynamic library location e.g. joinpath(dirname(@FILE), "../../../openfast/build/modules/inflowwind/libifwcbinding"))
• AM_flag::bool: flag to turn on added mass effects
• buoy_flag::bool: flag to turn on buoyancy forces
• rotAccel_flag::bool: flag to turn on the rotational acceleration portion of added mass for a crossflow turbine
• AM_Coeff_Ca::float: added mass coefficient, typically 1.0

Outputs:

• none:

steady(turbine::Turbine, env::Env; w=zeros(Real,2*turbine.ntheta), idx_RPI=1:2*turbine.ntheta,solve=true,ifw=false)

Calculates steady state aerodynamics for a single VAWT slice

Inputs

• turbine::Turbine: Turbine struct, see ?Turbine for details
• env::Env: Env struct, see ?Env for details
• w::Array(<:Real): Optional, used if solve=false, induction factor array, first half corresponding to u, second half to v
• idx_RPI::Array(<:Int): Optional, used to specify the azimuthal indices needed for a partial solve (i.e. not every azimuthal index), such as is used in the RPI method
• solve::Bool: Optional, False is used when you want the model outputs for a given set of induction factors without resolving them.
• ifw::Bool: Optional, used to tell the Vinf lookup to attempt to use the dynamic inflow wind library, requires preprocessing as is shown in the test cases.

Outputs:

• CP: This slice's coefficient of performance
• Th: This slice's thrust coefficient
• Q: Torque (N0m)
• Rp: Radial force per height (N)
• Tp: Tangential force per height (N)
• Zp: Vertical force per height (N)
• Vloc: Local velocity array for each azimuthal position (includes induction) (m/s)
• CD: This slice's drag coefficient
• CT: This slice's thrust coefficient (should equal drag, but may no depending on usage or solver status)
• amean: Mean turbine induction in the streamwise direction
• astar: Solved induction factors for each azimuthal location. First half are streamwise (u), second are cross-steam (v)
• alpha: Local angle of attack array for each azimuthal position (includes induction) (rad)
• cl: Local lift coefficient used for each azimuthal position
• cd_af: Local drag coefficient used for each azimuthal position
• thetavec: Azimuthal location of each discretization (rad)
• Re: Reynolds number for each azimuthal position

steadyTurb(omega,Vinf)

Runs a previously initialized aero model (see ?setupTurb) in the steady state mode

Inputs

• omega::float: turbine rotation rate (rad/s)
• Vinf::float: turbine steady inflow velocity (m/s)

Outputs:

• CP: Turbine coefficient of performance
• Rp: Array(B,Nslices,ntheta) of radial force (N)
• Tp: Array(B,Nslices,ntheta) of tangential force (N)
• Zp: Array(B,Nslices,ntheta) of vertical force (N)
• alpha: Array(B,Nslices,ntheta) of angle of attack (rad)
• cl: Array(B,Nslices,ntheta) of airfoil cl used
• cd_af: Array(B,Nslices,ntheta) of airfoil cd used
• Vloc: Array(B,Nslices,ntheta) of airfoil local velocity used
• Re: Array(B,Nslices,ntheta) of airfoil Reynolds number used
• thetavec: Azimuthal discretization location (rad)
• ntheta: number of azimuthal discretizations used
• Fx_base: Array(ntheta)Turbine base Fx (N)
• Fy_base: Array(ntheta)Turbine base Fy (N)
• Fz_base: Array(ntheta)Turbine base Fz (N)
• Mx_base: Array(ntheta)Turbine base Mx (N-m)
• My_base: Array(ntheta)Turbine base My (N-m)
• Mz_base: Array(ntheta)Turbine base Mz (N-m)
• power: Array(ntheta)Turbine power (watts)
• power2: Turbine average power for the revolution (watts)
• torque: Array(ntheta)Turbine torque (N-m) (alternative calculation method from Mz-base)

INTERNAL streamtube(a,theta,turbine,env;output_all=false,Vxwake=nothing,solvestep=false)

Double multiple streamtube individual streamtube calculation

Output:

if outputall return Th, Q, Rp, Tp, Zp, Vloc, CD, CT, alpha, cl, cdaf, Re else return CD-CT # Residual, section 2.4 end

Internal, trapezoidal integration of y w.r.t. x