Merge pull request #2375 from abhineet-gupta/f/tail_fin_unsteadymodel

Tail fin aerodynamics: address changes requested in pull request #1874

docs/source/user/aerodyn/input.rst

@@ -997,19 +997,19 @@ An example of tail fin input file is given below:
     Comment
     ======  General inputs =============================================================
     1         TFinMod     - Tail fin aerodynamics model {0: none, 1: polar-based, 2: USB-based} (switch)
-    0.5       TFinChord   - Tail fin chord (m) [used only when TFinMod=1]
-    0.3       TFinArea    - Tail fin planform area (m^2) [used only when TFinMod=1]
+    0.3       TFinArea    - Tail fin planform area (m^2)
     10.,0.,0. TFinRefP_n  - Undeflected position of the tail fin reference point wrt the tower top (m)
     0.,0.,0.  TFinAngles  - Tail fin chordline skew, tilt, and bank angles about the reference point (degrees)
     0         TFinIndMod  - Model for induced velocity calculation {0: none, 1:rotor-average} (switch)
     ====== Polar-based model ================================ [used only when TFinMod=1] 
-    1        TFinAFID - Index of Tail fin airfoil number [1 to NumAFfiles]
+    1         TFinAFID - Index of Tail fin airfoil number [1 to NumAFfiles]
+    0.5       TFinChord   - Tail fin chord (m)
     ====== Unsteady slender body model  ===================== [used only when TFinMod=2]
-    0.9           TFinKp        - Tail fin moment of area about reference point
-    0.3,0.1,0.1   TFinSigma     - Tail fin empirical constant for vortex separation functions
+    0.9           TFinKp        - Tail fin potential flow coefficient (-)
+    0.3,0.1,0.1   TFinSigma     - Tail fin empirical constant for vortex separation functions (1/deg)
     40,60,60      TFinAStar     - Tail fin initial angles for vortex separation functions (deg)
-    3.1416        TFinKv        - Tail fin vortex lift coefficient
-    1.3           TFinCDc       - Tail fin drag coefficient
+    3.1416        TFinKv        - Tail fin vortex lift coefficient (-)
+    1.3           TFinCDc       - Tail fin drag coefficient (-)
 
 General inputs
 ~~~~~~~~~~~~~~

docs/source/user/aerodyn/theory_tailfin.rst

@@ -34,7 +34,8 @@ The reference orientation (when the structure is un-deflected), the transformati
    \boldsymbol{R}_\text{tf,i} = \operatorname{EulerConstruct}(\theta_\text{bank}, \theta_\text{tilt}, \theta_\text{skew})
 
 For a common application with a vertical fin, the three angles are zero.
-:red:`TODO: The order of the angles might be different in the current implementation (3-2-1) instead of (1-2-3) above)`
+
+.. :red:`TODO: The order of the angles might be different in the current implementation (3-2-1) instead of (1-2-3) above)`
 
 
 **Velocities**
@@ -45,14 +46,15 @@ The following velocity vectors (3D vectors in global coordinates) are defined (s
   Undisturbed Wind speed vector at the reference point 
 - :math:`\boldsymbol{V}_\text{dist}`:
   Disturbed wind speed vector at the reference point (the disturbed wind contains the influence of the tower on the flow). AeroDyn has internal methods to compute :math:`\boldsymbol{V}_\text{dist}` from :math:`\boldsymbol{V}_\text{wind}`.
-  :red:`For now, we use "wind" but in the future we might use "dist". In the theory below we would simply replace all the "wind" by "dist"`.
 - :math:`\boldsymbol{V}_\text{elast}`:
   Structural translational velocity vector at the reference point
 - :math:`\boldsymbol{V}_\text{ind}`:
   Induced velocity from the wake at the reference point (assumed to be zero for now)
 - :math:`\boldsymbol{\omega}`:
   Structural rotational velocity of the fin
 
+..  :red:`For now, we use "wind" but in the future we might use "dist". In the theory below we would simply replace all the "wind" by "dist"`.
+
 All velocities (except for :math:`\boldsymbol{V}_\text{ind}` and :math:`\boldsymbol{V}_\text{dist}` which are computed internally by AeroDyn) are provided as input to the aerodynamic solver.
 The relative wind experienced by the airfoil is given by:
 
@@ -145,7 +147,8 @@ The rotor-averaged induced velocity can also be used as an estimate (`TFinIndMod
     \boldsymbol{V}_\text{ind}=\frac{1}{n_B n_r}\sum_{i_b=1..n_B} \sum_{i_r=1..n_r}  \boldsymbol{V}_{\text{ind},\text{blade}}[i_b, i_r] 
 
 Where :math:`\boldsymbol{V}_{\text{ind},\text{blade}}[i_b, i_r]` is the induced velocity vector for blade :math:`i_b` and at the radial node :math:`i_r`.
-:red:`NOTE: This averaging corresponds to what is done for the disk-average of the inflow in AeroDyn. In the future, we can use something weighted by the radius, or using precomputed coefficients, as done by Envision`.
+
+.. :red:`NOTE: This averaging corresponds to what is done for the disk-average of the inflow in AeroDyn. In the future, we can use something weighted by the radius, or using precomputed coefficients, as done by Envision`.
 
 More advanced models could set the induced velocity to zero when outside of the wake boundary, or include a tower-shadow-like wake model. Such option is not yet available.
 
@@ -163,15 +166,20 @@ Unsteady slender body model
 
 The unsteady aerodynamics of the tail fin is modeled based on Unsteady Slender Body Theory.
 The theory is extended to include the effect of high yaw angle :cite:`ad-hammam_NREL:2023`.
+To simplify the implementation, it is assumed that that arm length of the tail fin is much greater than the chord and the characteristic time (chord/wind speed) is small.
 
 The normal force on the tail fin can be described as the sum of three contributions (potential lift, vortex lift, and drag), weighted by separation functions :math:`x_i` as:
 
-.. math::  :label: TFUSBForce
-
+.. math::  
+    :label: tfusbforce
+    
     N = \frac{\rho}{2} A_{tf} \bigg(  K_p x_1 V_{\text{rel},x} V_{\text{rel},y} +  \Big[x_2 K_v+(1- x_3)C_{Dc} \Big] V_{\text{rel},y}\big|V_{\text{rel},y}\big|\bigg)
 
-
 where :math:`\rho` is the density of air, :math:`A_{tf}` is the tail fin area, :math:`K_p` is the potential lift coefficient and :math:`K_v` is the vortex lift coefficient, and :math:`C_{Dc}` is the drag coefficient.
+Note that the sign convnetion of OpenFAST is slightly different than used in :cite:`ad-hammam_NREL:2023`.
+This is reflected in Equation :eq:`tfusbforce`.
+
+
 :math:`x_i` are the separation functions calculated using a quasi-steady approximation as:
 
 .. math::  :label: TFUSBxiEquation
@@ -180,6 +188,6 @@ where :math:`\rho` is the density of air, :math:`A_{tf}` is the tail fin area, :
 
 
 where :math:`\sigma_i` are empirical constants characterizing the decay of separation functions, :math:`\gamma_{tf}` is the yaw angle of the tail fin with respect to the free-stream wind (:math:`V_{\text{wind}}`), :math:`\alpha^*_i` are the characteristics angles for separation functions.
-:math:`x_i` takes on a value between 0 and 1, and are used to activate or deactivate a the contribution of potential lift, vortex lift and draft to the normal force on the tail fin.
+:math:`x_i` takes on a value between 0 and 1, and are used to activate or deactivate a the contribution of potential lift, vortex lift and drag to the normal force on the tail fin.
 
 The normal force is assumed to act at the user defined reference point on the tail fin and the moment of the normal force is calculated accordingly.

modules/aerodyn/src/AeroDyn.f90

@@ -1,4 +1,4 @@
-!**********************************************************************************************************************************
+!*********************************************************************************************************************************
 ! LICENSING
 ! Copyright (C) 2015-2016  National Renewable Energy Laboratory
 ! Copyright (C) 2016-2021  Envision Energy USA, LTD
@@ -393,10 +393,10 @@ subroutine AD_Init( InitInp, u, p, x, xd, z, OtherState, y, m, Interval, InitOut
    do iR = 1, nRotors
       p%rotors(iR)%TFinAero         = InputFileData%rotors(iR)%TFinAero
       p%rotors(iR)%TFin%TFinMod     = InputFileData%rotors(iR)%TFin%TFinMod
-      p%rotors(iR)%TFin%TFinChord   = InputFileData%rotors(iR)%TFin%TFinChord
       p%rotors(iR)%TFin%TFinArea    = InputFileData%rotors(iR)%TFin%TFinArea
       p%rotors(iR)%TFin%TFinIndMod  = InputFileData%rotors(iR)%TFin%TFinIndMod
       p%rotors(iR)%TFin%TFinAFID    = InputFileData%rotors(iR)%TFin%TFinAFID
+      p%rotors(iR)%TFin%TFinChord   = InputFileData%rotors(iR)%TFin%TFinChord
       p%rotors(iR)%TFin%TFinKp      = InputFileData%rotors(iR)%TFin%TFinKp
       p%rotors(iR)%TFin%TFinSigma   = InputFileData%rotors(iR)%TFin%TFinSigma
       p%rotors(iR)%TFin%TFinAStar   = InputFileData%rotors(iR)%TFin%TFinAStar
@@ -4962,9 +4962,9 @@ SUBROUTINE TFin_CalcOutput(p, p_AD, u, RotInflow, m, y, ErrStat, ErrMsg )
       ! Unsteady aerodynamic model
 
       ! Calculate separation function (quasi-steady)
-      x1 = 1.0_Reki/(1.0_Reki+exp(p%TFin%TFinSigma(1)*((ABS(gamma_tf)*180.0_ReKi/pi)-p%TFin%TFinAStar(1)))) 
-      x2 = 1.0_Reki/(1.0_Reki+exp(p%TFin%TFinSigma(2)*((ABS(gamma_tf)*180.0_ReKi/pi)-p%TFin%TFinAStar(2)))) 
-      x3 = 1.0_Reki/(1.0_Reki+exp(p%TFin%TFinSigma(3)*((ABS(gamma_tf)*180.0_ReKi/pi)-p%TFin%TFinAStar(3))))
+      x1 = 1.0_Reki/(1.0_Reki+exp(p%TFin%TFinSigma(1)*((ABS(gamma_tf)*R2D)-p%TFin%TFinAStar(1)))) 
+      x2 = 1.0_Reki/(1.0_Reki+exp(p%TFin%TFinSigma(2)*((ABS(gamma_tf)*R2D)-p%TFin%TFinAStar(2)))) 
+      x3 = 1.0_Reki/(1.0_Reki+exp(p%TFin%TFinSigma(3)*((ABS(gamma_tf)*R2D)-p%TFin%TFinAStar(3))))
 
       ! Calculate unsteady force on tail fin
       force_tf(2) = 0.5_ReKi * p%AirDens * p%TFin%TFinArea * &

modules/aerodyn/src/AeroDyn_IO.f90

@@ -1677,14 +1677,14 @@ SUBROUTINE ReadTailFinInputs(FileName, TFData, UnEc, ErrStat, ErrMsg)
    !====== General inputs ============================================================
    call ParseCom(FileInfo_in, iLine, DummyLine                          , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;
    call ParseVar(FileInfo_In, iLine, 'TFinMod'   , TFData%TFinMod       , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;
-   call ParseVar(FileInfo_In, iLine, 'TFinChord' , TFData%TFinChord     , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;
    call ParseVar(FileInfo_In, iLine, 'TFinArea'  , TFData%TFinArea      , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;
    call ParseAry(FileInfo_In, iLine, 'TFinRefP_n', TFData%TFinRefP_n, 3 , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;
    call ParseAry(FileInfo_In, iLine, 'TFinAngles', TFData%TFinAngles, 3 , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;
    call ParseVar(FileInfo_In, iLine, 'TFinIndMod', TFData%TFinIndMod    , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;
    !====== Polar-based model ================================ [used only when TFinMod=1]
    call ParseCom(FileInfo_in, iLine, DummyLine                          , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;
    call ParseVar(FileInfo_In, iLine, 'TFinAFID'  , TFData%TFinAFID      , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;
+   call ParseVar(FileInfo_In, iLine, 'TFinChord' , TFData%TFinChord     , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;
    !====== Unsteady slender body model ===================== [used only when TFinMod=2]
    call ParseCom(FileInfo_in, iLine, DummyLine                          , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;
    call ParseVar(FileInfo_In, iLine, 'TFinKp'  , TFData%TFinKp      , ErrStat2, ErrMsg2, UnEc); if (Failed()) return;

modules/aerodyn/src/AeroDyn_Registry.txt

@@ -55,29 +55,29 @@ param    ^              -  IntKi    AD_MaxBl_Out            -  3  -  "Maximum nu
 
 # Tail Fin parameters 
 typedef  ^  TFinParameterType   IntKi   TFinMod    -    -   0 "Tail fin aerodynamics model {0=none, 1=polar-based, 2=USB-based}" (switch)
-typedef  ^  TFinParameterType   ReKi    TFinChord  -    -   - "Tail fin chord [used only when TFinMod=1]" m
-typedef  ^  TFinParameterType   ReKi    TFinArea   -    -   - "Tail fin planform area [used only when TFinMod=1]" m^2
+typedef  ^  TFinParameterType   ReKi    TFinArea   -    -   - "Tail fin planform area" m^2
 typedef  ^  TFinParameterType   IntKi   TFinIndMod -    -   - "Model for induced velocity calculation {0=none, 1=rotor-average}" (switch)
 typedef  ^  TFinParameterType   IntKi   TFinAFID   -    -   - "Index of Tail fin airfoil number [1 to NumAFfiles]" -
-typedef  ^  TFinParameterType   ReKi   TFinKp   -    -   - "Tail fin potential lift coefficient for unsteady aerodynamics [used only when TFMod=2]" -
-typedef  ^  TFinParameterType   ReKi   TFinSigma   3    -   - "Tail fin empirical constants characterizing the decay of separation functions [used only when TFMod=2]" -
-typedef  ^  TFinParameterType   ReKi   TFinAStar   3    -   - "Tail fin characteristics angles for separation functions [used only when TFMod=2]" deg
-typedef  ^  TFinParameterType   ReKi   TFinKv   -    -   - "Tail fin vortex lift coefficient for unsteady aerodynamics [used only when TFMod=2]" -
-typedef  ^  TFinParameterType   ReKi   TFinCDc   -    -   - "Tail fin drag coefficient for unsteady aerodynamics [used only when TFMod=2]" -
+typedef  ^  TFinParameterType   ReKi    TFinChord  -    -   - "Tail fin chord" m
+typedef  ^  TFinParameterType   ReKi   TFinKp   -    -   - "Tail fin potential flow coefficient for unsteady aerodynamics" -
+typedef  ^  TFinParameterType   ReKi   TFinSigma   3    -   - "Tail fin empirical constants characterizing the decay of separation functions" 1/deg
+typedef  ^  TFinParameterType   ReKi   TFinAStar   3    -   - "Tail fin characteristics angles for separation functions" deg
+typedef  ^  TFinParameterType   ReKi   TFinKv   -    -   - "Tail fin vortex lift coefficient for unsteady aerodynamics" -
+typedef  ^  TFinParameterType   ReKi   TFinCDc   -    -   - "Tail fin drag coefficient for unsteady aerodynamics" -
 
 # Tail Fin input file 
 typedef  ^  TFinInputFileType   IntKi   TFinMod    -    -   0 "Tail fin aerodynamics model {0=none, 1=polar-based, 2=USB-based}" (switch)
-typedef  ^  TFinInputFileType   ReKi    TFinChord  -    -   - "Tail fin chord [used only when TFinMod=1]" m
-typedef  ^  TFinInputFileType   ReKi    TFinArea   -    -   - "Tail fin planform area [used only when TFinMod=1]" m^2
+typedef  ^  TFinInputFileType   ReKi    TFinArea   -    -   - "Tail fin planform area" m^2
 typedef  ^  TFinInputFileType   ReKi    TFinRefP_n 3    -   - "Undeflected position of the tail fin reference point wrt the tower top" m 
 typedef  ^  TFinInputFileType   ReKi    TFinAngles 3    -   - "Tail fin chordline skew, tilt, and bank angles about the reference point" (deg)
 typedef  ^  TFinInputFileType   IntKi   TFinIndMod -    -   - "Model for induced velocity calculation {0=none, 1=rotor-average}" (switch)
 typedef  ^  TFinInputFileType   IntKi   TFinAFID   -    -   - "Index of Tail fin airfoil number [1 to NumAFfiles]" -
-typedef  ^  TFinInputFileType   ReKi    TFinKp   -    -   - "Tail fin potential lift coefficient for unsteady aerodynamics [used only when TFMod=2]" -
-typedef  ^  TFinInputFileType   ReKi    TFinSigma   3    -   - "Tail fin empirical constants characterizing the decay of separation functions [used only when TFMod=2]" -
-typedef  ^  TFinInputFileType   ReKi    TFinAStar 3    -   - "Tail fin characteristics angles for separation functions [used only when TFMod=2]" deg
-typedef  ^  TFinInputFileType   ReKi   TFinKv   -    -   - "Tail fin vortex lift coefficient for unsteady aerodynamics [used only when TFMod=2]" -
-typedef  ^  TFinInputFileType   ReKi   TFinCDc   -    -   - "Tail fin drag coefficient for unsteady aerodynamics [used only when TFMod=2]" -
+typedef  ^  TFinInputFileType   ReKi    TFinChord  -    -   - "Tail fin chord" m
+typedef  ^  TFinInputFileType   ReKi    TFinKp   -    -   - "Tail fin potential flow coefficient for unsteady aerodynamics" -
+typedef  ^  TFinInputFileType   ReKi    TFinSigma   3    -   - "Tail fin empirical constants characterizing the decay of separation functions" 1/deg
+typedef  ^  TFinInputFileType   ReKi    TFinAStar 3    -   - "Tail fin characteristics angles for separation functions" deg
+typedef  ^  TFinInputFileType   ReKi   TFinKv   -    -   - "Tail fin vortex lift coefficient for unsteady aerodynamics" -
+typedef  ^  TFinInputFileType   ReKi   TFinCDc   -    -   - "Tail fin drag coefficient for unsteady aerodynamics" -