OpenFAST
Wind turbine multiphysics simulator
Functions/Subroutines | Variables
ifw_uniformwind Module Reference

This module contains all the data and procedures that define uniform wind files (formerly known as hub-height files). More...

Functions/Subroutines

subroutine, public ifw_uniformwind_init (InitData, ParamData, MiscVars, Interval, InitOutData, ErrStat, ErrMsg)
 A subroutine to initialize the UniformWind module. More...
 
subroutine, public ifw_uniformwind_calcoutput (Time, PositionXYZ, p, Velocity, DiskVel, m, ErrStat, ErrMsg)
 This routine and its subroutines calculate the wind velocity at a set of points given in PositionXYZ. More...
 
subroutine interpparams (Time, p, m, op)
 This subroutine linearly interpolates the parameters that are used to compute uniform wind. More...
 
subroutine getwindspeed (InputPosition, p, m, op, WindSpeed, ErrStat, ErrMsg)
 This subroutine linearly interpolates the columns in the uniform input file to get the values for the requested time, then uses the interpolated values to calclate the wind speed at a point in space represented by InputPosition. More...
 
real(reki) function, dimension(3) windinf_adhack_diskvel (t, p, m, ErrStat, ErrMsg)
 This function should be deleted ASAP. More...
 
subroutine, public ifw_uniformwind_end (ParamData, MiscVars, ErrStat, ErrMsg)
 This routine closes any open files and clears all data stored in UniformWind derived Types. More...
 
subroutine, public ifw_uniformwind_jacobianpinput (t, Position, CosPropDir, SinPropDir, p, m, dYdu)
 Routine to compute the Jacobians of the output (Y) function with respect to the inputs (u). More...
 
subroutine, public ifw_uniformwind_getop (t, p, m, OP_out)
 Routine to compute the Jacobians of the output (Y) function with respect to the inputs (u). More...
 

Variables

type(progdesc), parameter ifw_uniformwind_ver = ProgDesc( 'IfW_UniformWind', '', '' )
 

Detailed Description

This module contains all the data and procedures that define uniform wind files (formerly known as hub-height files).

This could more accurately be called a point wind file since the wind speed at any point is calculated by shear applied to the point where wind is defined. It is basically uniform wind over the rotor disk. The entire file is read on initialization, then the columns that make up the wind file are interpolated to the time requested, and wind is calculated based on the location in space.

the file contains header information (rows that contain "!"), followed by numeric data stored in 8 columns:

Column Description Variable Name Units
1 Time Time [s]
2 Horizontal wind speed V [m/s]
3 Wind direction Delta [deg]
4 Vertical wind speed VZ [m/s]
5 Horizontal linear shear HLinShr [-]
6 Vertical power-law shear VShr [-]
7 Vertical linear shear VLinShr [-]
8 Gust (horizontal) velocity VGust [m/s]

The horizontal wind speed at (X, Y, Z) is then calculated using the interpolated columns by

\begin{eqnarray} V_h & = & V \, \left( \frac{Z}{Z_{Ref}} \right) ^ {VShr} & \mbox{power-law wind shear} \\ & + & V \, \frac{H_{LinShr}}{RefWid} \, \left( Y \cos(Delta) + X \sin(Delta) \right) & \mbox{horizontal linear shear} \\ & + & V \, \frac{V_{LinShr}}{RefWid} \, \left( Z-Z_{Ref} \right) & \mbox{vertical linear shear} \\ & + & V_{Gust} & \mbox{gust speed} \end{eqnarray}

Function/Subroutine Documentation

◆ getwindspeed()

subroutine ifw_uniformwind::getwindspeed ( real(reki), dimension(3), intent(in)  InputPosition,
type(ifw_uniformwind_parametertype), intent(in)  p,
type(ifw_uniformwind_miscvartype), intent(inout)  m,
type(ifw_uniformwind_intrp), intent(in)  op,
real(reki), dimension(3), intent(out)  WindSpeed,
integer(intki), intent(out)  ErrStat,
character(*), intent(out)  ErrMsg 
)
private

This subroutine linearly interpolates the columns in the uniform input file to get the values for the requested time, then uses the interpolated values to calclate the wind speed at a point in space represented by InputPosition.

16-Apr-2013 - A. Platt, NREL. Converted to modular framework. Modified for NWTC_Library 2.0

Parameters
[in]inputpositioninput information: positions X,Y,Z
[in]pParameters
[in,out]mMisc variables (stores last index into array time for efficiency)
[in]opoperating point values; interpolated UniformWind parameters for this time (for glue-code linearization operating point)
[out]errstaterror status
[out]errmsgThe error message
[out]windspeedreturn velocities (U,V,W)
  1. Calculate the wind speed at this time (if z<0, return an error):

Let

\begin{eqnarray} V_h & = & V \, \left( \frac{Z}{Z_{ref}} \right) ^ {V_{shr}} & \mbox{power-law wind shear} \\ & + & V \, \frac{H_{LinShr}}{RefWid} \, \left( Y \cos(Delta) + X \sin(Delta) \right) & \mbox{horizontal linear shear} \\ & + & V \, \frac{V_{LinShr}}{RefWid} \, \left( Z - Z_{ref} \right) & \mbox{vertical linear shear} \\ & + & V_{Gust} & \mbox{gust speed} \end{eqnarray}

Then the returned wind speed, \(Vt\), is
\(Vt_u = V_h \, \cos(Delta) \)
\(Vt_v = -V_h \, \sin(Delta) \)
\(Vt_w = V_z \)
using input positions \(X,Y,Z\) and interpolated values for time-dependent input-file parameters \(V, Delta, V_z, H_{LinShr}, V_{Shr}, V_{LinShr}, V_{Gust}\).

◆ ifw_uniformwind_calcoutput()

subroutine, public ifw_uniformwind::ifw_uniformwind_calcoutput ( real(dbki), intent(in)  Time,
real(reki), dimension(:,:), intent(in)  PositionXYZ,
type(ifw_uniformwind_parametertype), intent(in)  p,
real(reki), dimension(:,:), intent(inout)  Velocity,
real(reki), dimension(3), intent(out)  DiskVel,
type(ifw_uniformwind_miscvartype), intent(inout)  m,
integer(intki), intent(out)  ErrStat,
character(*), intent(out)  ErrMsg 
)

This routine and its subroutines calculate the wind velocity at a set of points given in PositionXYZ.

The UVW velocities are returned in Velocity

Note
This routine does not satisfy the Modular framework.
Date
16-Apr-2013 - A. Platt, NREL. Converted to modular framework. Modified for NWTC_Library 2.0
Parameters
[in]timetime from the start of the simulation
[in]positionxyzArray of XYZ coordinates, 3xN
[in]pParameters
[in,out]velocityVelocity output at Time (Set to INOUT so that array does not get deallocated)
[out]diskvelHACK for AD14: disk velocity output at Time
[in,out]mMisc variables for optimization (not copied in glue code)
[out]errstaterror status
[out]errmsgThe error message
  1. Linearly interpolate parameters in time (or use nearest-neighbor to extrapolate) (compare with nwtc_num::interpstpreal)

◆ ifw_uniformwind_end()

subroutine, public ifw_uniformwind::ifw_uniformwind_end ( type(ifw_uniformwind_parametertype), intent(inout)  ParamData,
type(ifw_uniformwind_miscvartype), intent(inout)  MiscVars,
integer(intki), intent(out)  ErrStat,
character(*), intent(out)  ErrMsg 
)

This routine closes any open files and clears all data stored in UniformWind derived Types.

Note
This routine does not satisfy the Modular framework. The InputType is not used, rather an array of points is passed in.
Date
: 16-Apr-2013 - A. Platt, NREL. Converted to modular framework. Modified for NWTC_Library 2.0
Parameters
[in,out]paramdataParameters
[in,out]miscvarsMisc variables for optimization (not copied in glue code)
[out]errstatdetermines if an error has been encountered
[out]errmsgMessage about errors

◆ ifw_uniformwind_getop()

subroutine, public ifw_uniformwind::ifw_uniformwind_getop ( real(dbki), intent(in)  t,
type(ifw_uniformwind_parametertype), intent(in)  p,
type(ifw_uniformwind_miscvartype), intent(inout)  m,
real(reki), dimension(2), intent(out)  OP_out 
)

Routine to compute the Jacobians of the output (Y) function with respect to the inputs (u).

The partial derivative dY/du is returned. This submodule does not follow the modularization framework.

Parameters
[in]tCurrent simulation time in seconds
[out]op_outoperating point (HWindSpeed and PLexp
[in]pParameters
[in,out]mMisc/optimization variables
  1. Linearly interpolate parameters in time at operating point (or use nearest-neighbor to extrapolate) (compare with nwtc_num::interpstpreal)

◆ ifw_uniformwind_init()

subroutine, public ifw_uniformwind::ifw_uniformwind_init ( type(ifw_uniformwind_initinputtype), intent(in)  InitData,
type(ifw_uniformwind_parametertype), intent(out)  ParamData,
type(ifw_uniformwind_miscvartype), intent(out)  MiscVars,
real(dbki), intent(in)  Interval,
type(ifw_uniformwind_initoutputtype), intent(out)  InitOutData,
integer(intki), intent(out)  ErrStat,
character(*), intent(out)  ErrMsg 
)

A subroutine to initialize the UniformWind module.

It reads the uniform wind file and stores the data in an array to use later. It requires an initial reference height (hub height) and width (rotor diameter), both in meters, which are used to define the volume where wind velocities will be calculated. This information is necessary because of the way the shears are defined.

Note
This routine does not conform to the framework. The InputType has been replaced with just the PositionXYZ array.
Date
16-Apr-2013 - A. Platt, NREL. Converted to modular framework. Modified for NWTC_Library 2.0
Parameters
[in]initdataInput data for initialization
[out]paramdataParameters
[out]miscvarsMisc variables for optimization (not copied in glue code)
[out]initoutdataInitial output
[in]intervalWe don't change this.
[out]errstatdetermines if an error has been encountered
[out]errmsgA message about the error

◆ ifw_uniformwind_jacobianpinput()

subroutine, public ifw_uniformwind::ifw_uniformwind_jacobianpinput ( real(dbki), intent(in)  t,
real(reki), dimension(3), intent(in)  Position,
real(reki), intent(in)  CosPropDir,
real(reki), intent(in)  SinPropDir,
type(ifw_uniformwind_parametertype), intent(in)  p,
type(ifw_uniformwind_miscvartype), intent(inout)  m,
real(r8ki), dimension(3,6), intent(inout)  dYdu 
)

Routine to compute the Jacobians of the output (Y) function with respect to the inputs (u).

The partial derivative dY/du is returned. This submodule does not follow the modularization framework.

Parameters
[in]tCurrent simulation time in seconds
[in]positionXYZ Position at which to find velocity (operating point)
[in]cospropdircosine of InflowWind propagation direction
[in]sinpropdirsine of InflowWind propagation direction
[in]pParameters
[in,out]mMisc/optimization variables
[in,out]dyduPartial derivatives of output functions (Y) with respect to the inputs (u)
  1. Linearly interpolate parameters in time at operating point (or use nearest-neighbor to extrapolate) (compare with nwtc_num::interpstpreal)
  2. Calculate \( \frac{\partial Y_{Output \, Equations}}{\partial u_{inputs}} = \begin{bmatrix} \frac{\partial Vt_u}{\partial X} & \frac{\partial Vt_u}{\partial Y} & \frac{\partial Vt_u}{\partial Z} \\ \frac{\partial Vt_v}{\partial X} & \frac{\partial Vt_v}{\partial Y} & \frac{\partial Vt_v}{\partial Z} \\ \frac{\partial Vt_w}{\partial X} & \frac{\partial Vt_w}{\partial Y} & \frac{\partial Vt_w}{\partial Z} \\ \end{bmatrix} \)

\( \frac{\partial Vt_u}{\partial X} = \left[\cos(PropagationDir)\cos(Delta) - \sin(PropagationDir)\sin(Delta) \right] V \, \frac{H_{LinShr}}{RefWid} \, \sin(Delta) \cos(PropagationDir) \)

\( \frac{\partial Vt_v}{\partial X} = \left[-\sin(PropagationDir)\cos(Delta) - \cos(PropagationDir)\sin(Delta) \right] V \, \frac{H_{LinShr}}{RefWid} \, \sin(Delta) \cos(PropagationDir) \)

\( \frac{\partial Vt_w}{\partial X} = 0 \)

\( \frac{\partial Vt_u}{\partial Y} = \left[\cos(PropagationDir)\cos(Delta) - \sin(PropagationDir)\sin(Delta) \right] V \, \frac{H_{LinShr}}{RefWid} \, \cos(Delta) \cos(PropagationDir) \)

\( \frac{\partial Vt_v}{\partial Y} = \left[-\sin(PropagationDir)\cos(Delta) - \cos(PropagationDir)\sin(Delta) \right] V \, \frac{H_{LinShr}}{RefWid} \, \cos(Delta) \cos(PropagationDir) \)

\( \frac{\partial Vt_w}{\partial Y} = 0 \)

\( \frac{\partial Vt_u}{\partial Z} = \left[\cos(PropagationDir)\cos(Delta) - \sin(PropagationDir)\sin(Delta) \right] V \, \left[ \frac{V_{shr}}{Z_{ref}} \left( \frac{Z}{Z_{ref}} \right) ^ {V_{shr}-1} + \frac{V_{LinShr}}{RefWid} \right] \)

\( \frac{\partial Vt_v}{\partial Z} = \left[-\sin(PropagationDir)\cos(Delta) - \cos(PropagationDir)\sin(Delta) \right] V \, \left[ \frac{V_{shr}}{Z_{ref}} \left( \frac{Z}{Z_{ref}} \right) ^ {V_{shr}-1} + \frac{V_{LinShr}}{RefWid} \right] \)

\( \frac{\partial Vt_w}{\partial Z} = 0 \)

\( \frac{\partial Vt_w}{\partial V} = 0 \)

\( \frac{\partial Vt_w}{\partial VShr} = 0 \)

\( \frac{\partial Vt_w}{\partial PropDir} = 0 \)

◆ interpparams()

subroutine ifw_uniformwind::interpparams ( real(dbki), intent(in)  Time,
type(ifw_uniformwind_parametertype), intent(in)  p,
type(ifw_uniformwind_miscvartype), intent(inout)  m,
type(ifw_uniformwind_intrp), intent(out)  op 
)
private

This subroutine linearly interpolates the parameters that are used to compute uniform wind.

Parameters
[in]timetime from the start of the simulation
[in]pParameters
[in,out]mMisc variables (index)
[out]opinterpolated V values at input TIME

◆ windinf_adhack_diskvel()

real(reki) function, dimension(3) ifw_uniformwind::windinf_adhack_diskvel ( real(dbki), intent(in)  t,
type(ifw_uniformwind_parametertype), intent(in)  p,
type(ifw_uniformwind_miscvartype), intent(inout)  m,
integer(intki), intent(out)  ErrStat,
character(*), intent(out)  ErrMsg 
)
private

This function should be deleted ASAP.

Its purpose is to reproduce results of AeroDyn 12.57; when a consensus on the definition of "average velocity" is determined, this function will be removed.

Parameters
[in]tTime
[in]pParameters
[in,out]mmisc/optimization data (storage for efficiency index)
[out]errstaterror status from this function
[out]errmsgerror message from this function