.. _DynamicDerivatives:
.. |deg|    unicode:: U+000B0 .. DEGREE SIGN
   :trim:

Calculating Dynamic Derivatives using Sliding Interfaces
=========================================================

This tutorial demonstrates how to obtain dynamic stability derivatives by oscillating a wing enclosed in a sliding interface.

Geometry and Mesh
------------------
In this example, the geometry is a wing with NACA0012 airfoil. 
The chord :math:`c` = 1 m and the span :math:`b` = 2 m. 
So the :code:`momentLength` = [:math:`b`/2, :math:`c`, :math:`b`/2] = [1, 1, 1]. 
The :code:`momentCenter` is set at [0, 0, 0], which is the airfoil quarter-chord point and midspan of the wing. 
As shown in the following figure, there is a cylindrical sliding interface enclosing the wing. 
The radius of the cylinder is 1 m and the length is 2.5 m. 
The rotational center is also located at [0, 0, 0] and the axis of rotation is along the y-axis.

.. _wing_in_sliding_interface:
.. figure:: Figures/wing_in_sliding_interface.png
   :align: center

   Wing enclosed in sliding interface

The mesh was generated using the automated meshing workflow. The `geometry file <https://simcloud-public-1.s3.amazonaws.com/tutorials/calculating_dynamic_derivatives_using_sliding_interface/geometry.csm>`_, `surfaceMesh.json <https://simcloud-public-1.s3.amazonaws.com/tutorials/calculating_dynamic_derivatives_using_sliding_interface/surfaceMesh.json>`_ and `volumeMesh.json <https://simcloud-public-1.s3.amazonaws.com/tutorials/calculating_dynamic_derivatives_using_sliding_interface/volumeMesh_singleSlidingInterface.json>`_ are provided so the readers can easily generate the mesh in their account. For more information about this automated meshing workflow please see :ref:`Automated Meshing <automatedMeshing>`.

Case Setup
-----------
The entire simulation contains two cases:

#. Steady case with a fixed sliding interface for initializing the flow field.
#. Unsteady case with an oscillating sliding interface for collecting the data.

The steady case is quite trivial since the sliding interface is stationary.
For readers' convenience, the `Flow360Steady.json <https://simcloud-public-1.s3.amazonaws.com/tutorials/calculating_dynamic_derivatives_using_sliding_interface/flow360_singleSlidingInterface-steady.json>`_ is provided.
As for the unsteady case, when preparing the `Flow360Unsteady.json <https://simcloud-public-1.s3.amazonaws.com/tutorials/calculating_dynamic_derivatives_using_sliding_interface/flow360_singleSlidingInterface-mixed-2.5flowthruPerRad-2deg-400steps.json>`_, the key is to correctly set up the :code:`thetaDegrees` input expression:

.. math::
    \text{thetaDegrees} = A * \sin (\text{omega} * \text{t})

Where :math:`A` is the amplitude of oscillation in degrees, 
:code:`omega` is the nondimensional angular frequency of the oscillation,
and :code:`t` is the nondimensional time. Note that :code:`omega` and :code:`t` are nondimensional variables of the Flow360 solver.

In this case we set :math:`A=2` |deg|, so the pitching angle of the wing is oscillating between :math:`\pm 2` |deg|.
:code:`omega` is typically computed from the number of flowthroughs per radian, denoted by :math:`n`.
In this example, we have :math:`n = 2.5` flowthru/rad, which means when :code:`omega*t` changes 1 radian, the air will flow through 2.5 chord lengths.
In this tutorial, the freestream velocity :math:`U_\infty=50` m/s.
The physical angular frequency is given by,

.. math::
    \omega = \frac{U_\infty}{n \cdot c}
           = \frac{50 \ \text{m/s}}{2.5 \times 1 \ \text{m}}
           = 20 \ \text{rad/s}

Since we have the freestream speed of sound :math:`C_\infty = 340.294` m/s and gridUnit = 1 m.
The nondimensional angular frequency can be written as,

.. math::
    \text{omega} = \frac{\omega}{C_\infty/L_\text{gridUnit}}
                 \approx 0.05877271

Therefore, the :code:`slidingInterfaces` section is:

.. code-block:: json
    
    "slidingInterfaces": [
        {
            "stationaryPatches": ["stationaryBlock/slidingInterface-0"],
            "rotatingPatches": ["rotatingBlock-0/slidingInterface-0"],
            "axisOfRotation": [0,1,0],
            "centerOfRotation": [0,0,0],
            "volumeName": "rotatingBlock-0",
            "thetaDegrees": "2.00 * sin(0.05877271 * t)"
        }
    ]

Another critical value we need to set up carefully is the nondimensional :code:`timeStepSize`.
In this example, each period is spilt into 100 steps, therefore,

.. math::
    \text{timeStepSize} = 0.01 \cdot \frac{2\pi}{\text{omega}}

Note that the :code:`omega` in the above equation is the nondimensional angular frequency.

Submission
-----------
The surface/volume meshes and steady/unsteady cases can be easily submitted via our `PythonAPI <https://docs.flexcompute.com/projects/flow360/en/latest/pythonAPIreference/pythonAPIreference.html>`_:

.. code-block:: python

    surfaceMeshId  = flow360client.NewSurfaceMeshFromGeometry('path/to/geometry.csm', 'path/to/surfaceMesh.json')
    volumeMeshId   = flow360client.NewMeshFromSurface(surfaceMeshId=surfaceMeshId, config='path/to/volumeMesh.json')
    steadyCaseId   = flow360client.NewCase(meshId=volumeMeshId, config='path/to/SteadyFlow360.json', caseName='dynamic-derivatives-steady')
    unsteadyCaseId = flow360client.NewCase(meshId=volumeMeshId, config='path/to/UnsteadyFlow360.json', caseName='dynamic-derivatives-unsteady', parentId=steadyCaseId)

Postprocessing
---------------
Once the case is completed, the user can download the "total_forces_v2.csv", where the first column is :code:`physical_step`.
Note that the flowfield was initialized from a steady case, so the nondimensional time can be written as,

.. math::
    \text{t} = \text{physical_step} * \text{timeStepSize}


In this example, the angle of attack :math:`\alpha` is always identical to the pitch angle :math:`\theta`,

.. math::
    \alpha = \theta = A * \sin (\text{omega} * \text{t})

Therefore, :math:`\dot{\alpha}` is always identical to the pitch rate :math:`q`.
The :math:`\dot{\alpha}` at time step :code:`i` can be obtained by finite differencing :math:`\alpha` at the :code:`i-1` and :code:`i+1` time step

.. math::
    \dot{\alpha}_i = \frac{\alpha_{i+1} - \alpha_{i-1}}{2 * \text{timeStepSize}} \frac{C_\infty}{L_\text{gridUnit}} \frac{c}{U_\infty}

Note that the unit of :math:`\dot{\alpha}` in the equation above is rad/flowthru.

As shown in the following plot, when :math:`\alpha=0`, :math:`\dot{\alpha}` will reach the extrema.

.. _time_history_CMy_alpha_alpha_dot:
.. figure:: Figures/timeHistory.png
   :width: 60%
   :align: center 

   Time History of CMy, :math:`\alpha` and :math:`\dot{\alpha}`. The blue and red dots show when :math:`\dot{\alpha}` reaches the extrema.

If we plot CMy versus :math:`\dot{\alpha}` and linearly fit the points with :math:`\alpha=0` (i.e. min/max :math:`\dot{\alpha}`) the slope of the fitted line is the dynamic derivative:

.. math::
    \frac{dCMy}{d\dot{\alpha}} = \frac{dCMy}{dq} = -0.511

.. _dCMy_dq:
.. figure:: Figures/CMy.png
   :width: 60%
   :align: center

   Dynamic Derivative :math:`dCMy/d\dot{\alpha}` i.e. :math:`dCMy/dq`.