"""Canonical CST airfoil runtime family."""
from __future__ import annotations
from typing import override
import numpy as np
from buffalo_core.numeric import as_float_array, as_float_scalar
from buffalo_core.typing import FloatArray, FloatInput, FloatScalar
from buffalo_wings.airfoil.internal.schema.analytic import (
CstAirfoilSpec,
CstSurfaceSpec,
)
from .cst_airfoil_side import CstAirfoilSide
from .cst_geometry import CstGeometry
[docs]
class CstAirfoil(CstGeometry):
"""
Canonical CST airfoil with fixed airfoil class exponents.
Notes
-----
This runtime represents the airfoil-specific CST case with ``n1 = 0.5``
and ``n2 = 1.0`` on both sides.
Unlike :class:`CstGeometry`, the curve parameterization uses the
canonical side parameter ``x = s**2`` with full-airfoil curve parameter
``s = |u|``.
This keeps canonical curve-parameter derivatives finite at the leading
edge while preserving the same geometric shape as the corresponding
general CST
geometry under the mapped parameter ``tau = sign(u) * u**2``.
"""
_upper: CstAirfoilSide
_lower: CstAirfoilSide
def __init__(
self,
*,
upper: CstAirfoilSide,
lower: CstAirfoilSide,
trailing_edge_thickness: FloatScalar = 0.0,
) -> None:
"""
Initialize a canonical CST airfoil from two canonical sides.
Parameters
----------
upper : CstAirfoilSide
Upper-side canonical CST definition.
lower : CstAirfoilSide
Lower-side canonical CST definition.
trailing_edge_thickness : FloatScalar, default=0.0
Explicit trailing-edge gap as a fraction of chord.
"""
super().__init__(
upper=upper,
lower=lower,
trailing_edge_thickness=trailing_edge_thickness,
)
self._upper = CstAirfoilSide(
shape=upper.shape,
delta_te=self._upper.delta_te,
)
self._lower = CstAirfoilSide(
shape=lower.shape,
delta_te=self._lower.delta_te,
)
@property
def upper(self) -> CstAirfoilSide:
"""
Return the upper-side canonical CST geometry.
This property exposes the upper-side canonical CST definition.
"""
return self._upper
@property
def lower(self) -> CstAirfoilSide:
"""
Return the lower-side canonical CST geometry.
This property exposes the lower-side canonical CST definition.
"""
return self._lower
@property
def spec(self) -> CstAirfoilSpec:
"""
Return the schema definition used to create this airfoil.
The returned schema reproduces the current runtime coefficients
and trailing-edge thickness.
"""
return CstAirfoilSpec(
trailing_edge_thickness=float(self.trailing_edge_thickness),
upper=CstSurfaceSpec(
a=[float(value) for value in self.upper.coefficients]
),
lower=CstSurfaceSpec(
a=[float(value) for value in self.lower.coefficients]
),
)
[docs]
@override
def xy_from_u(self, u: FloatInput) -> tuple[FloatArray, FloatArray]:
"""
Calculate canonical CST airfoil coordinates.
Parameters
----------
u : FloatInput
Signed airfoil parameter values in ``[-1, 1]``.
Returns
-------
tuple[FloatArray, FloatArray]
Tuple ``(x, y)`` of ``float64`` arrays.
Notes
-----
This uses ``x = u**2`` on each surface branch rather than
the general CST geometry parameterization ``x = |u|``.
"""
u_array = as_float_array(u)
u_abs_array = as_float_array(np.abs(u_array))
x = as_float_array(u_abs_array**2)
y = np.empty_like(x)
upper_mask = u_array >= 0.0
lower_mask = ~upper_mask
if np.any(upper_mask):
_, upper_y = self.upper.xy_canonical(u_abs_array[upper_mask])
y[upper_mask] = upper_y
if np.any(lower_mask):
_, lower_y = self.lower.xy_canonical(u_abs_array[lower_mask])
y[lower_mask] = lower_y
return x, as_float_array(y)
[docs]
@override
def xy_u(self, u: FloatInput) -> tuple[FloatArray, FloatArray]:
"""
Calculate first derivatives with respect to the canonical parameter.
Parameters
----------
u : FloatInput
Signed airfoil parameter values in ``[-1, 1]``.
Returns
-------
tuple[FloatArray, FloatArray]
Tuple ``(dx/du, dy/du)`` of ``float64`` arrays.
Notes
-----
At listed breakpoints, this method returns the minus-side
derivative so array-valued evaluations remain single-valued.
"""
u_array = as_float_array(u)
x_u = np.empty_like(u_array)
y_u = np.empty_like(u_array)
analytic_mask = np.ones_like(u_array, dtype=bool)
breakpoints = self.breakpoints()
flat_u = u_array.ravel()
flat_x_u = x_u.ravel()
flat_y_u = y_u.ravel()
flat_mask = analytic_mask.ravel()
for index, value in enumerate(flat_u):
breakpoint_index = self._breakpoint_index(
as_float_scalar(value),
breakpoints=breakpoints,
)
if breakpoint_index is not None:
minus, _ = self.xy_u_breakpoint(index=breakpoint_index)
flat_x_u[index] = minus[0]
flat_y_u[index] = minus[1]
flat_mask[index] = False
if np.any(analytic_mask):
analytic_u = u_array[analytic_mask]
s_array = as_float_array(np.abs(analytic_u))
upper_mask = analytic_u >= 0.0
lower_mask = ~upper_mask
analytic_x_u = np.empty_like(analytic_u)
analytic_y_u = np.empty_like(analytic_u)
if np.any(upper_mask):
upper_x_u, upper_y_u = self.upper.xy_canonical_t(
s_array[upper_mask]
)
analytic_x_u[upper_mask] = upper_x_u
analytic_y_u[upper_mask] = upper_y_u
if np.any(lower_mask):
lower_x_u, lower_y_u = self.lower.xy_canonical_t(
s_array[lower_mask]
)
analytic_x_u[lower_mask] = -lower_x_u
analytic_y_u[lower_mask] = -lower_y_u
x_u[analytic_mask] = analytic_x_u
y_u[analytic_mask] = analytic_y_u
return as_float_array(x_u), as_float_array(y_u)
[docs]
@override
def xy_uu(self, u: FloatInput) -> tuple[FloatArray, FloatArray]:
"""
Calculate second derivatives with respect to the canonical parameter.
Parameters
----------
u : FloatInput
Signed airfoil parameter values in ``[-1, 1]``.
Returns
-------
tuple[FloatArray, FloatArray]
Tuple ``(d2x/du2, d2y/du2)`` of ``float64`` arrays.
Notes
-----
At listed breakpoints, this method returns the minus-side second
derivative so array-valued evaluations remain single-valued.
"""
u_array = as_float_array(u)
x_uu = np.empty_like(u_array)
y_uu = np.empty_like(u_array)
analytic_mask = np.ones_like(u_array, dtype=bool)
breakpoints = self.breakpoints()
flat_u = u_array.ravel()
flat_x_uu = x_uu.ravel()
flat_y_uu = y_uu.ravel()
flat_mask = analytic_mask.ravel()
for index, value in enumerate(flat_u):
breakpoint_index = self._breakpoint_index(
as_float_scalar(value),
breakpoints=breakpoints,
)
if breakpoint_index is not None:
minus, _ = self.xy_uu_breakpoint(index=breakpoint_index)
flat_x_uu[index] = minus[0]
flat_y_uu[index] = minus[1]
flat_mask[index] = False
if np.any(analytic_mask):
analytic_u = u_array[analytic_mask]
s_array = as_float_array(np.abs(analytic_u))
upper_mask = analytic_u >= 0.0
lower_mask = ~upper_mask
analytic_x_uu = np.empty_like(analytic_u)
analytic_y_uu = np.empty_like(analytic_u)
if np.any(upper_mask):
upper_x_uu, upper_y_uu = self.upper.xy_canonical_tt(
s_array[upper_mask]
)
analytic_x_uu[upper_mask] = upper_x_uu
analytic_y_uu[upper_mask] = upper_y_uu
if np.any(lower_mask):
lower_x_uu, lower_y_uu = self.lower.xy_canonical_tt(
s_array[lower_mask]
)
analytic_x_uu[lower_mask] = lower_x_uu
analytic_y_uu[lower_mask] = lower_y_uu
x_uu[analytic_mask] = analytic_x_uu
y_uu[analytic_mask] = analytic_y_uu
return as_float_array(x_uu), as_float_array(y_uu)
[docs]
def xy_u_breakpoint(
self,
*,
index: int,
) -> tuple[
tuple[FloatScalar, FloatScalar], tuple[FloatScalar, FloatScalar]
]:
"""
Return one-sided first derivatives at one breakpoint index.
Notes
-----
Endpoint breakpoints return the same exact boundary value for both
entries.
The interior leading-edge breakpoint returns finite canonical
one-sided derivatives in the ``x = s**2`` parameterization.
"""
u_breakpoint = self.breakpoints()[index]
if u_breakpoint <= -1.0:
boundary_values = self._xy_u_boundary(u=-1.0)
return boundary_values, boundary_values
if u_breakpoint >= 1.0:
boundary_values = self._xy_u_boundary(u=1.0)
return boundary_values, boundary_values
lower_x_u, lower_y_u = self.lower.xy_canonical_t(0.0)
upper_x_u, upper_y_u = self.upper.xy_canonical_t(0.0)
minus_values = (
-as_float_scalar(lower_x_u),
-as_float_scalar(lower_y_u),
)
plus_values = (
as_float_scalar(upper_x_u),
as_float_scalar(upper_y_u),
)
return minus_values, plus_values
[docs]
def xy_uu_breakpoint(
self,
*,
index: int,
) -> tuple[
tuple[FloatScalar, FloatScalar], tuple[FloatScalar, FloatScalar]
]:
"""
Return one-sided second derivatives at one breakpoint index.
Notes
-----
Endpoint breakpoints return the same exact boundary value for both
entries.
The interior leading-edge breakpoint returns finite canonical
one-sided second derivatives in the ``x = s**2`` parameterization.
"""
u_breakpoint = self.breakpoints()[index]
if u_breakpoint <= -1.0:
boundary_values = self._xy_uu_boundary(u=-1.0)
return boundary_values, boundary_values
if u_breakpoint >= 1.0:
boundary_values = self._xy_uu_boundary(u=1.0)
return boundary_values, boundary_values
lower_x_uu, lower_y_uu = self.lower.xy_canonical_tt(0.0)
upper_x_uu, upper_y_uu = self.upper.xy_canonical_tt(0.0)
minus_values = (
as_float_scalar(lower_x_uu),
as_float_scalar(lower_y_uu),
)
plus_values = (
as_float_scalar(upper_x_uu),
as_float_scalar(upper_y_uu),
)
return minus_values, plus_values
def _xy_u_boundary(
self, *, u: FloatScalar
) -> tuple[FloatScalar, FloatScalar]:
"""
Return one exact boundary first derivative value.
Parameters
----------
u : FloatScalar
Endpoint airfoil parameter, either negative on the lower trailing
edge or positive on the upper trailing edge.
Returns
-------
tuple[FloatScalar, FloatScalar]
Exact one-sided boundary derivative ``(x_u, y_u)``.
"""
if u < 0.0:
x_u, y_u = self.lower.xy_canonical_t(1.0)
return -as_float_scalar(x_u), -as_float_scalar(y_u)
x_u, y_u = self.upper.xy_canonical_t(1.0)
return as_float_scalar(x_u), as_float_scalar(y_u)
def _xy_uu_boundary(
self, *, u: FloatScalar
) -> tuple[FloatScalar, FloatScalar]:
"""
Return one exact boundary second derivative value.
Parameters
----------
u : FloatScalar
Endpoint airfoil parameter, either negative on the lower trailing
edge or positive on the upper trailing edge.
Returns
-------
tuple[FloatScalar, FloatScalar]
Exact one-sided boundary second derivative ``(x_uu, y_uu)``.
"""
if u < 0.0:
x_uu, y_uu = self.lower.xy_canonical_tt(1.0)
return as_float_scalar(x_uu), as_float_scalar(y_uu)
x_uu, y_uu = self.upper.xy_canonical_tt(1.0)
return as_float_scalar(x_uu), as_float_scalar(y_uu)