Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Update in Computing Default Unsteady Airfoil Coefficients #920

Merged
merged 6 commits into from
Mar 1, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
19 changes: 13 additions & 6 deletions docs/source/user/aerodyn/input.rst
Original file line number Diff line number Diff line change
Expand Up @@ -238,7 +238,8 @@ value, *f’*, will be calculated. When ``FLookup`` is set to TRUE, *f’*
is determined via a lookup into the static lift-force coefficient and
drag-force coefficient data. **Using best-fit exponential equations
(``FLookup = FALSE``) is not yet available, so ``FLookup`` must be
``TRUE`` in this version of AeroDyn.** Note, ``FLookup`` is not used when ``UAMod=5``.
``TRUE`` in this version of AeroDyn.** Note, ``FLookup`` is not used
when ``UAMod=4`` or ``UAMod=5``.

``UAStartRad`` is the starting rotor radius where dynamic stall
will be turned on. Enter a number between 0 and 1, representing a fraction of rotor radius,
Expand Down Expand Up @@ -434,6 +435,7 @@ The default value of 0.2 if provided for convenience.

``NonDimArea`` is the nondimensional airfoil area (normalized by the
local ``BlChord`` squared), but is currently unused by AeroDyn.

``NumCoords`` is the number of points to define the exterior shape of
the airfoil, plus one point to define the aerodynamic center, and
determines the number of rows in the subsequent table; ``NumCoords``
Expand All @@ -449,6 +451,7 @@ used by OpenFAST for blade surface visualization when enabled.

``BL_file`` is the name of the file containing boundary-layer characteristics
of the profile. It is ignored if the aeroacoustic module is not used.
This parameter may also be omitted from the file if the aeroacoustic module is not used.

Specify the number of Reynolds number- or aerodynamic-control
setting-dependent tables of data for the given airfoil via the
Expand All @@ -472,18 +475,22 @@ or calculating it based on the polar coefficient data in the airfoil table:

- ``alpha1`` specifies the AoA (in degrees) larger than ``alpha0``
for which *f* equals 0.7; approximately the positive stall angle;
This parameter is used when ``flookup=false`` and when determining
a default value of ``Cn1``.

- ``alpha2`` specifies the AoA (in degrees) less than ``alpha0``
for which *f* equals 0.7; approximately the negative stall angle;
This parameter is used when ``flookup=false`` and when determining
a default value of ``Cn2``.

- ``alphaUpper`` specifies the AoA (in degrees) of the upper boundary of
fully-attached region of the cn or cl curve. It is used to
compute the separation function when ``UAMod=5``.
fully-attached region of the cn or cl curve. It is used to
compute the separation function when ``UAMod=5``.

- ``alphaLower`` specifies the AoA (in degrees) of the lower boundary of
fully-attached region of the cn or cl curve. It is used to
compute the separation function when ``UAMod=5``. (The separation function
will have a value of 1 between ``alphaUpper`` and ``alphaLower``.)
fully-attached region of the cn or cl curve. It is used to
compute the separation function when ``UAMod=5``. (The separation function
will have a value of 1 between ``alphaUpper`` and ``alphaLower``.)

- ``eta_e`` is the recovery factor and typically has a value in the
range [0.85 to 0.95] for ``UAMod = 1``; if the keyword ``DEFAULT`` is
Expand Down
89 changes: 89 additions & 0 deletions docs/source/user/aerodyn/theory_ua.rst
Original file line number Diff line number Diff line change
Expand Up @@ -139,7 +139,29 @@ as the function reaches :math:`0` on both sides of :math:`\alpha_0`,
then :math:`f_s^{st}` is kept at the constant value :math:`0`.


**Note that for UAMod=5, a different separation function is formed.**
We define an offset for the :math:`C_n` function, ``cn_offset``, where
:math:`C_{n,offset}=\frac{C_n\left(\alpha^{Lower}\right)+C_n\left(\alpha^{Upper}\right)}{2}`. Then, the separation function
is a value between 0 and 1, given by the following equation:

.. math::

f_s^{st}(\alpha) = \left[ 2 \max\left\{\frac{1}{4} , \sqrt{\frac{C_n^{st}(\alpha) - C_{n,offset}}{C_n^{fullyAttached}(\alpha)-C_{n,offset}}} \right\} -1 \right]^2

with the fully-attached :math:`C_n` curve defined as :math:`C_n` between :math:`alpha^{Lower}` and :math:`alpha^{Upper}` and linear functions outside of that range:

.. math::

C_n^{fullyAttached}(\alpha) = \begin{cases} C_n\left(\alpha^{Upper}\right) + C_n^{slope}\left(\alpha^{Upper}\right) \cdot \left(\alpha-\alpha^{Upper}\right) & \alpha>\alpha^{Upper} \\
C_n(\alpha) & \alpha^{Lower}<=\alpha<=\alpha^{Upper} \\
C_n\left(\alpha^{Lower}\right) + C_n^{slope}\left(\alpha^{Lower}\right) \cdot \left(\alpha-\alpha^{Lower}\right) & \alpha<\alpha^{Lower} \end{cases}

Note that to avoid numerical issues at the :math:`\pm180` degree boundary, this function changes slope when the separation function is 0 above :math:`alpha^{Upper}` and below :math:`alpha^{Lower}`.
This allow the fully-attached linear sections to be periodic and avoid numerical issues with large magnitudes of angle of attack.


**Inviscid and fully separated lift coefficient:**

The inviscid lift coefficient is
:math:`C_{l,\text{inv}}= C_{l,\alpha} (\alpha-\alpha_0)`.
The fully separated lift coefficient may
Expand Down Expand Up @@ -465,13 +487,79 @@ Inputs
See :numref:`ad_ua_inputs` for a description of the inputs necessary in the AeroDyn primary file (e.g. ``UAMod``)

See :numref:`airfoil_data_input_file` for a more comprehensive description of all the inputs present in the profile input file.
Their default values are described in :numref:`UA_AFI_defaults`

See :numref:`ua_notations` for a list of notations and definitions specific to unsteady aerodynamic inputs.

An example of profile data (containing some of the unsteady aerodynamic parameters) is available here
:download:`(here) <examples/ad_polar_example.dat>`.


.. _UA_AFI_defaults:

Calculating Default Airfoil Coefficients
----------------------------------------

The default value for ``cd0`` is the minimum value of the :math:`C_d` curve between :math:`\pm20` degrees angle of attack.
:math:`\alpha_{c_{d0}}` is defined to be the angle of attack where ``cd0`` occurs.

After computing ``cd0``, the :math:`C_n` curve is computed by

.. math::
C_{n}(\alpha) = C_l(\alpha) \cos\alpha + \left(C_d(\alpha) - c_{d0}\right) \sin\alpha

The slope of the :math:`C_n` curve is computed as follows:

.. math::
C_{n}^{Slope}\left(\frac{\alpha_{i+1} + \alpha_i}{2}\right) = \frac{C_n(\alpha_{i+1}) - C_n(\alpha_i)}{\alpha_{i+1} - \alpha_i}
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is this computed for a range of angle of attack and then few indices? How is this range determined? I'm guessing it's +/- 20 deg?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The slope is computed for NumAlf-1 points (i.e., between all consecutive AoA specified in the polar). The smoothed values of this curve can have some weird behavior at the ends, but alphaUpper and alphaLower are far enough from the ends that it shouldn't matter.


:math:`C_{n,smooth}^{Slope}` is a smoothed version of :math:`C_{n}^{Slope}`, calculated using a triweight kernel with a window of 2 degrees.


.. math::
C_{l}^{Slope}\left(\frac{\alpha_{i+1} + \alpha_i}{2}\right) = \frac{C_l(\alpha_{i+1}) - C_l(\alpha_i)}{\alpha_{i+1} - \alpha_i}


Using :math:`C_{n,smooth}^{Slope}`, ``alphaUpper`` and ``alphaLower`` are computed:

``alphaUpper`` is the smallest angle of attack value between :math:`\alpha_{c_{d0}}` and 20 degrees where the :math:`C_{n,smooth}^{Slope}` curve has started to decrease to 90% of its maximum slope.

.. math::
C_{n,smooth}^{Slope}\left(\alpha^{Upper}\right) < 0.9 \max_{\alpha \in \left[\alpha_{c_{d0}}, \alpha^{Upper}\right]} C_{n,smooth}^{Slope}\left( \alpha \right)


``alphaLower`` is the largest angle of attack value between -20 degrees and :math:`\alpha_{c_{d0}}` where the :math:`C_{n,smooth}^{Slope}` curve has started to decrease to 90% of its maximum slope.

.. math::
C_{n,smooth}^{Slope}\left(\alpha^{Lower}\right) < 0.9 \max_{\alpha \in \left[\alpha^{Lower}, \alpha_{c_{d0}}\right]} C_{n,smooth}^{Slope}\left( \alpha \right)

``Cn1`` is the value of :math:`C_n(\alpha)` at the smallest value of :math:`\alpha` where :math:`\alpha >= \alpha^{Upper}` and the separation function, :math:`f_{st}(\alpha)` = 0.7.

``Cn2`` is the value of :math:`C_n(\alpha)` at the largest value of :math:`\alpha` where :math:`\alpha <= \alpha^{Lower}` and the separation function, :math:`f_{st}(\alpha)` = 0.7.

``Cn_offset`` is the average value of the :math:`C_n` curve at ``alphaUpper`` and ``alphaLower``:

.. math::
C_{n}^{offset} = \frac{C_n\left(\alpha^{Lower}\right) + C_n\left(\alpha^{Upper}\right)}{2}

``C_nalpha`` is defined as the maximum slope of the smoothed :math:`C_n` curve, :math:`C_{n,smooth}^{Slope}` between :math:`\pm20` degrees angle of attack.
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

From the definition of Cnalpha or Clalpha it seems the local slope is used (between i and i+1). I wonder if it would be better to use the slope between i and i0 where i0 if the index of alpha0. That the definition used in the Riso-R-1354 report equation 16.

Copy link
Contributor Author

@bjonkman bjonkman Nov 30, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We could change that. I was trying to make it somewhat consistent with the method for computing C_lalpha (which UAMod4 says has to be the max slope), but it's used in different UA models, so could certainly be done differently.


``C_lalpha`` is defined as the maximum slope of the (un-smoothed) :math:`C_l` curve, :math:`C_{l}^{Slope}` between :math:`\pm20` degrees angle of attack.

The default ``alpha0`` is computed as the zero-crossing of a line with a slope equal to ``C_lalpha`` that goes through the :math:`C_l` curve at :math:`\alpha = \frac{\alpha^{Upper} + \alpha^{Lower}}{2}`

.. math::
\alpha_0 = \frac{\alpha^{Upper} + \alpha^{Lower}}{2} - \frac{C_l\left(\frac{\alpha^{Upper} + \alpha^{Lower}}{2}\right) }{C_{l,\alpha}}

``Cm0`` is the value of the :math:`C_m` curve at ``alpha0``: :math:`C_{m,0} = C_m\left(\alpha_0\right)`. If the :math:`C_m` polar values have not been included, :math:`C_{m,0} =0`.

``alpha1`` is the angle of attack above ``alphaUpper`` where the separation function, :math:`f_s^{st}` is 0.7.

``alpha2`` is the angle of attack below ``alphaLower`` where the separation function, :math:`f_s^{st}` is 0.7.

``Cn1`` is the value of the :math:`C_n` curve at ``alpha1``.

``Cn2`` is the value of the :math:`C_n` curve at ``alpha2``.



Expand All @@ -485,6 +573,7 @@ To activate these outputs with `cmake`, compile using ``-DCMAKE_Fortran_FLAGS="-




Driver
------

Expand Down
2 changes: 1 addition & 1 deletion glue-codes/fast-farm/src/FASTWrapper.f90
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,7 @@ MODULE FASTWrapper

PRIVATE

TYPE(ProgDesc), PARAMETER :: FWrap_Ver = ProgDesc( 'FASTWrapper', 'v1.00.00', '7-Feb-2017' ) !< module date/version information
TYPE(ProgDesc), PARAMETER :: FWrap_Ver = ProgDesc( 'FASTWrapper', '', '' ) !< module date/version information

REAL(DbKi), PARAMETER :: t_initial = 0.0_DbKi ! Initial time

Expand Down
4 changes: 2 additions & 2 deletions glue-codes/fast-farm/src/FAST_Farm_IO.f90
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@ module FAST_Farm_IO
IMPLICIT NONE
! The maximum number of output channels which can be output by the code.
INTEGER(IntKi), PARAMETER :: Farm_MaxOutPts = 9423
TYPE(ProgDesc), PARAMETER :: Farm_Ver = ProgDesc( 'FAST.Farm', 'v0.01.02a-jmj', '04-Jul-2017' ) !< module date/version information
TYPE(ProgDesc), PARAMETER :: Farm_Ver = ProgDesc( 'FAST.Farm', '', '' ) !< module date/version information
! ===================================================================================================
! NOTE: The following lines of code were generated by a Matlab script called "Write_ChckOutLst.m"
! using the parameters listed in the "OutListParameters3.xlsx" Excel file. Any changes to these
Expand Down Expand Up @@ -11064,7 +11064,7 @@ SUBROUTINE WriteFarmOutputToFile( t_global, farm, ErrStat, ErrMsg )
REAL(DbKi) :: OutTime ! Used to determine if output should be generated at this simulation time
INTEGER(IntKi) :: ErrStat2
CHARACTER(ErrMsgLen) :: ErrMSg2
CHARACTER(*), PARAMETER :: RoutineName = 'WriteOutputToFile'
CHARACTER(*), PARAMETER :: RoutineName = 'WriteFarmOutputToFile'
CHARACTER(200) :: Frmt ! A string to hold a format specifier
CHARACTER(farm%p%TChanLen) :: TmpStr ! temporary string to print the time output as text
CHARACTER(ChanLen) :: TmpStr2 ! temporary string to print the output as text
Expand Down
2 changes: 1 addition & 1 deletion modules/aerodyn/src/AeroAcoustics_IO.f90
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@ MODULE AeroAcoustics_IO

implicit none

type(ProgDesc), parameter :: AA_Ver = ProgDesc( 'AeroAcoustics', 'v1.00.00', '18-Aug-2016' )
type(ProgDesc), parameter :: AA_Ver = ProgDesc( 'AeroAcoustics', '', '' )
character(*), parameter :: AA_Nickname = 'AA'


Expand Down
Loading