Airfoil Generation
Parameter Types
Most of the parameterizations below have associated composite types whose fields are the parameters used in defining the airfoil geometries. Each of these composite types is defined using the @kwdef
macro such that the user does not need to remember the order of the fields, but can simply define the fields as though they were keyword arguments. In general, few of the fields are given defaults with the exception of things like trailing edge gap or y-positions, which are always defaulted to zero. In addition, some parameterization methods have specific values inherent to their methods. These are exposed to the user for convenience, but are also defaulted to the values inherent in the parameterization method.
Coordinates
All coordinates are given clockwise from the trailing edge, even if the coordinates are split between upper and lower sides. In other words, the coordinates are given from the lower trailing edge to the leading edge, then the leading edge back to the upper side trailing edge. In general, if the coordinates are given in an upper and lower split, the leading edge point is repeated. Note that all coordinates are given with x being in the chord-wise direction and y being orthogonal to x (just like a standard Cartesian x,y plot). Airfoil geometries are all given normalized to the chord length, so the x values will generally have the leading edge at x=0 and the trailing edge at x=1, though there may be some slight variation depending on camber and trailing edge gap.
NACA Parameterizations
NACA 4-series
We begin with the standard NACA 4-series airfoil defined by maximum camber, position of maximum camber, and maximum thickness.
using FLOWFoil.AirfoilTools
using Plots
using LaTeXStrings
parameters = NACA4(;
maximum_camber=2.0,
maximum_camber_position=4.0,
maximum_thickness=12.0,
blunt_trailing_edge=false,
)
FLOWFoil.AirfoilTools.NACA4{Float64, Float64, Float64, Bool}(2.0, 4.0, 12.0, false)
The NACA4
parameter type is the only one with fully defined defaults, which happen to default to the NACA 2412 airfoil with a sharp trailing edge due to its ubiquity.
We can then determine the x,y coordinates from the parameters.
x, y = naca4(parameters)
plot(x, y; aspectratio=1, xlabel=L"x", ylabel=L"y", label="")
If you are already familiar with a specific parameterization method, you can also forego defining the parameter object and call the method directly.
x, y = naca4(2.0, 4.0, 12.0)
plot(
x,
y;
aspectratio=1,
xlabel=L"x",
ylabel=L"y",
label="Direct",
linewidth=4,
linestyle=:dash,
)
If you have a set of coordinates, and would like to fit NACA 4-series parameters to them, you can use the determine_naca4
function.
fit_parameters = determine_naca4(x,y)
x, y = naca4(fit_parameters)
plot!(x, y; label="Fit")
FLOWFoil.AirfoilTools.NACA4
— TypeNACA4
Fields
maximum_camber::Float=2.0
: maximum camber in % chordmaximum_camber_position::Float=4.0
: x-position of maximum camber point in 1/10 chordmaximum_thickness::Float=12.0
: maximum thickness in % chordblunt_trailing_edge::Bool=false
: flag for whether to use the blunt trailing edge NACA 4-series definition or not.
FLOWFoil.AirfoilTools.naca4
— Functionnaca4(parameters::NACA4; N=161, x=nothing, split=false)
Compute x, y airfoil coordinates for N nodes, based on NACA 4-Series Parameterization.
Arguments
parameters::NACA4
: NACA 4-series parameters
Keyword Arguments
N::Int=161
: Total number of coordinates to use. This values should be odd, but if not, the number of points returned will be N-1.x::AbstractArray{Float}
: x coordinates (cosine spaced coordinates used by default)split::Bool=false
: Flag wheter to split into upper and lower halves.
Returns
If split
== false:
x::AbstractArray{Float}
: Vector of x coordinates, clockwise from trailing edge.y::AbstractArray{Float}
: Vector of y coordinates, clockwise from trailing edge.
If split
== true:
xl::AbstractArray{Float}
: Vector of lower half of x coordinates from trailing edge to leading edge.xu::AbstractArray{Float}
: Vector of upper half of x coordinates from leading edge to trailing edge.yl::AbstractArray{Float}
: Vector of lower half of y coordinates from trailing edge to leading edge.yu::AbstractArray{Float}
: Vector of upper half of y coordinates from leading edge to trailing edge.
naca4(c=2.0, p=4.0, t=12.0; N=161, x=nothing, blunt_trailing_edge=false, split=false)
Compute x, y airfoil coordinates for N nodes, based on NACA 4-Series Parameterization.
Arguments
c::Float
: Maximum camber value (percent of chord)p::Float
: Position along chord (in 10ths of chord) where maximum naca4_camber liest::Float
: Maximum thickness of airfoil in percent chord
Keyword Arguments
N::Int
: Total number of coordinates to use. This values should be odd, but if not, the number of points returned will be N-1.x::AbstractArray{Float}
: x-coordinates (cosine spaced coordinates used by default)blunt_trailing_edge::Bool
: Flag whether trailing edge is blunt or notsplit::Bool
: Flag wheter to split into upper and lower halves.
Returns
If split
== false:
x::AbstractArray{Float}
: Vector of x coordinates, clockwise from trailing edge.y::AbstractArray{Float}
: Vector of y coordinates, clockwise from trailing edge.
If split
== true:
xl::AbstractArray{Float}
: Vector of lower half of x coordinates from trailing edge to leading edge.xu::AbstractArray{Float}
: Vector of upper half of x coordinates from leading edge to trailing edge.yl::AbstractArray{Float}
: Vector of lower half of y coordinates from trailing edge to leading edge.yu::AbstractArray{Float}
: Vector of upper half of y coordinates from leading edge to trailing edge.
FLOWFoil.AirfoilTools.determine_naca4
— Functiondetermine_naca4(x,y)
Calculate NACA 4-series parameters based on input x,y coordinates.
Arguments
x::AbstractArray{Float}
: vector of x coordinates for airfoily::AbstractArray{Float}
: vector of y coordinates for airfoil
Keyword Arguments
blunt_trailing_edge::Bool=false
: Flag whether trailing edge is blunt or not
Returns
parameters::NACA4
: a parameter object of type NACA4.
NACA 65-series
The NACA 65-series of airfoils was historically used for various linear cascade experiments done by NACA (see "NACA Report No 824 Summary of Airfoil Data" by Ira H. Abbott, Albert E. Von Doenhoff, and Louis S. Stivers, Jr. ). We have implemented some NACA 65-series parameterizations associated with various NACA experimental studies. In general, the 65-series requires a target lift coefficient, a mean line designation, and a series number.
using FLOWFoil.AirfoilTools
using Plots
using LaTeXStrings
#define design lift coefficient
clo = 1.0
#define mean line designation
a = 1.0
#define series number
series_number = "3-018"
x, y = AirfoilTools.naca65(clo, a, series_number)
plot(x, y; aspectratio=1, xlabel=L"x", ylabel=L"y", label="")
There exists a special scaled case for the NACA 65-010 family of airfoils. This special case only requires the user to input a design lift coefficient. using
#define design lift coefficient
clo = 1.0
x, y = AirfoilTools.naca65_scaled(clo)
plot(x, y; aspectratio=1, xlabel=L"x", ylabel=L"y", label="")
FLOWFoil.AirfoilTools.naca65
— FunctionAssumes x in non-dimensional range [0.0,1.0]
Description from NACA Research Memorandum L51G31: "Systematic Two-dimensional Cascade Tests of NACA 65-Series Compressor Blades at Low Speeds:"
The 65-series compressor blade family is formed by combining a basic thickness form with cambered mean lines. The basic thickness form used is the NACA 65(216)-010 thickness form with the ordinates increased by 0.0015 times the chordwise stations to provide slightly, increased thickness toward the trailing edge. In the scaled case, it was not derived for 10-percent thickness but was scaled down from the NACA 65,2-016 airfoil. The scaling procedure gives the best results whep it is restricted to maximum thickness changes of a few percent. The NACA 65-010 basic thickness has also been derived. These thickness forms differ slightly but are considered to be interchangeable.
The basic mean line used is the a=1.0 mean line. The amount of camber is for the design lift coefficient for the isolated airfoil with cl_o
of 1.0. Both ordinates and slopes are scaled directly to obtain other cambers. Cambered blade sections are obtained by applying the thickness perpendicular to the mean line at stations laid out along the chord line. In the designation the camber is given by the first number after the dash in tenths of cl_o
. For example, the NACA 65-810 and NACA 65-(12)10 blade sections are cambered for cl_o
= 0.8 and cl_o
= 1.2, respectively.
Arguments
clo::TF
: Design lift coefficient in tenths of chord. Usually first number after the 2nd dash (ie NACA 65-3-818 would input 0.8 for clo)a::TF
: Mean-line designation, fraction of chord from leading edge over which design load is uniform.series_number::String
: digits of the airfoil series family with clo of 0 (ie NACA 65-3-818 would be "3-010")
Keyword Arguments
x::Vector{TF} = nothing
: input x values if specificing x values manuallysplit::Boolean = false
: if true, then the output will be split between top and bottom coordinatesextra_blending::Boolean = false
: If desired number of points is large (> 300ish) then set to true and it will add some extra blending if desired. Note: This is generally not needed!
Returns
x::Vector{TF}
: x coordinatesy::Vector{TF}
: y coordinates
FLOWFoil.AirfoilTools.naca65_scaled
— Functionnaca65_scaled(clo; N=161, x=nothing, split=false, extra_blending = false)
Specific version of naca65 for the NACA 65-010 series which uses special mean line and thickness form coordinates. It scales those coordinates based on the cl_o
value.
Arguments
clo::TF
: Design lift coefficient in tenths of chord. Usually first number after the 2nd dash (ie NACA 65-3-818 would input 0.8 for clo)
Keyword Arguments
x::Vector{TF} = nothing
: input x values if specificing x values manuallysplit::Boolean = false
: if true, then the output will be split between top and bottom coordinatesextra_blending::Boolean = false
: If desired number of points is large (> 300ish) then set to true and it will add some extra blending if desired. Note: This is generally not needed!
ouptuts
x::Vector{TF}
: x coordinatesy::Vector{TF}
: y coordinates
Conformal Mapping Parameterizations
Conformally mapped airfoils are defined from analytical solutions for flow about a cylinder. Currently we have the Joukowsky conformal mapping implemented, which features a cusped trailing edge. In development is the Karman-Trefftz conformal mapping, which allows for a non-zero trailing edge angle.
Joukowsky
The Joukowsky map is based on the center and radius of a circle
using FLOWFoil.AirfoilTools
using Plots
using LaTeXStrings
center = [-0.1; 0.1]
radius = 1.0
pc = plot(
cosd.(0:360) * radius .+ center[1],
sind.(0:360) * radius .+ center[2];
aspectratio=1,
label="",
xlabel=L"\xi",
ylabel=L"\eta",
framestyle=:origin
)
scatter!(pc, [center[1]], [center[2]]; color=1, label="")
parameters = Joukowsky(; center, radius)
x, y = joukowsky(parameters; N=161)
paf = plot(x, y; aspectratio=1, xlabel=L"x", ylabel=L"y", label="")
plot(pc, paf; size=(900, 300))
Since conformally mapped airfoils are defined from analytic solutions, we can also use those solutions for validation as done in the Additional Examples.
angle_of_attack = 5.0
surface_velocity, surface_pressure_coefficient, cl = joukowsky_flow(
center, radius, angle_of_attack
)
plot(
x[2:(end - 1)],
surface_pressure_coefficient[2:(end - 1)];
xlabel=L"x",
ylabel=L"c_p",
yflip=true,
label="",
)
FLOWFoil.AirfoilTools.Joukowsky
— TypeJoukowsky
Fields
center::AbstractArray{Float}
: [x y] location of center of circle relative to originradius::Float
: radius of circle
FLOWFoil.AirfoilTools.joukowsky
— Functionjoukowsky(parameters::Joukowsky; N=361, fortest=false, normalize=true, split=false)
Joukowsky airfoil parameterization.
Arguments
parameters::Joukowsky
: Joukowsky parameters
Keyword Arguments
N::Int=361
: Total number of coordinates to use. Can be even or odd, but it is recommended to be odd for a clear leading edge point.fortest::Bool=false
: Flag to output non-coordinate paramters used in 'joukowsky_flow()'normalize::Bool=true
: Flag whether to normalize to unit chord and translate the leading edge to zero.split::Bool=false
: Flag wheter to split output into upper and lower surfaces.
Returns
IF split == false
x::AbstractArray{Float}
: Array of x coordinatesy::AbstractArray{Float}
: Array of y coordinates
IF split == true
xu::AbstractArray{Float}
: Array of upper half of x coordinatesxl::AbstractArray{Float}
: Array of lower half of x coordinatesyu::AbstractArray{Float}
: Array of upper half of y coordinatesyl::AbstractArray{Float}
: Array of lower half of y coordinates
joukowsky(center, radius; N=361, fortest=false, normalize=true, split=false)
Joukowsky airfoil parameterization.
Arguments
center::AbstractArray{Float}
: [x y] location of center of circle relative to originradius::Float
: radius of circle
Keyword Arguments
N::Int=361
: Total number of coordinates to use. Can be even or odd, but it is recommended to be odd for a clear leading edge point.fortest::Bool=false
: Flag to output non-coordinate paramters used in 'joukowsky_flow()'normalize::Bool=true
: Flag whether to normalize to unit chord and translate the leading edge to zero.split::Bool=false
: Flag wheter to split output into upper and lower surfaces.
Returns
IF split == false
x::AbstractArray{Float}
: Array of x coordinatesy::AbstractArray{Float}
: Array of y coordinates
IF split == true
xu::AbstractArray{Float}
: Array of upper half of x coordinatesxl::AbstractArray{Float}
: Array of lower half of x coordinatesyu::AbstractArray{Float}
: Array of upper half of y coordinatesyl::AbstractArray{Float}
: Array of lower half of y coordinates
FLOWFoil.AirfoilTools.joukowsky_flow
— Functionjoukowsky_flow(center, radius, alpha; N=361)
Calculate the analytic surface velocities and pressures as well as lift coefficient for a joukowsky airfoil.
Arguments
center::AbstractArray{Float}
: [x y] location of circle center relative to originradius::Float
: Radius of circlealpha::Float
: Angle of attack in degrees
Keyword Arguments
N::Int=361
: Total number of coordinates to use. Can be even or odd, but it is recommended to be odd for a clear leading edge point.
Returns
vsurf::AbstractArray{Float}
: Magnitude of surface velocities at the nodescpsurf::AbstractArray{Float}
: Surface pressures at the nodescl::Float
: Lift coefficient
Class Shape Transformation (CST) Paramterizations
Standard Kulfan CST
We implement the CST parameterization presented by Kulfan, where airfoil shapes are defined based on upper and lower coefficients.
using FLOWFoil.AirfoilTools
using Plots
using LaTeXStrings
parameters = CST(;
upper_coefficients=[0.2; 0.3; 0.2; 0.2], lower_coefficients=[-0.1; 0.1; 0.0; 0.0]
)
x, y = cst(parameters)
plot(x, y; aspectratio=1, xlabel=L"x", ylabel=L"y", label="")
This CST method also has a fitting function.
x, y = naca4()
plot(
x,
y;
aspectratio=1,
xlabel=L"x",
ylabel=L"y",
label="Coordinates",
linewidth=4,
linestyle=:dash,
)
fit_parameters = determine_cst(x, y; n_upper_coefficients=4, n_lower_coefficients=4)
x, y = cst(fit_parameters)
plot!(x, y; label="Fit")
FLOWFoil.AirfoilTools.CST
— TypeCST
Fields
upper_coefficients::AbstractArray{Float}
: Vector of coefficients defining the upper side.lower_coefficients::AbstractArray{Float}
: Vector of coefficients defining the lower side.trailing_edge_yu::Float=0.0
: y-position of the upper side trailing edge.trailing_edge_yl::Float=0.0
: y-position of the lower side trailing edge.N1::Float=0.5
: inherent parameter for round-nosed airfoils.N2::Float=1.0
: inherent parameter for sharp trailing edge (with optional blunt trailing edge) airfoils.
FLOWFoil.AirfoilTools.cst
— Functioncst(
parameters::CST;
N::Integer=80,
x=split_cosine_spacing(N),
split=false,
)
Obtain airfoil coordiantes (clockwise from trailing edge) from the class shape transformation (CST) parameterization.
Arguments
parameters::CST
: CST parameters for airfoil.
Keyword Arguments
N::Integer=80
: number of points to use for each sidex::AbstractArray{Float}=split_cosine_spacing(N)
: x-coordinates to use.trailing_edge_yu::Float=0.0
: upper side trailing edge gaptrailing_edge_yl::Float=0.0
: lower side trailing edge gapN1::Float=0.5
: Class shape parameter 1N2::Float=1.0
: Class shape parameter 2split::Bool=false
: if true, returns upper and lower coordinates separately as xl, xu, yl, yu rather than just x, y.
Returns
If split
== false
x::AbstractArray{Float}
: Vector of x coordinates, clockwise from trailing edge.y::AbstractArray{Float}
: Vector of y coordinates, clockwise from trailing edge.
If split
== true
xl::AbstractArray{Float}
: Vector of lower half of x coordinates from trailing edge to leading edge.xu::AbstractArray{Float}
: Vector of upper half of x coordinates from leading edge to trailing edge.yl::AbstractArray{Float}
: Vector of lower half of y coordinates from trailing edge to leading edge.yu::AbstractArray{Float}
: Vector of upper half of y coordinates from leading edge to trailing edge.
cst(
upper_coefficients,
lower_coefficients;
N::Integer=80,
x=split_cosine_spacing(N),
trailing_edge_yu=0.0,
trailing_edge_yl=0.0,
N1=0.5,
N2=1.0,
split=false,
)
Obtain airfoil coordiantes (clockwise from trailing edge) from the class shape transformation (CST) parameterization.
Arguments
upper_coefficients::AbstractArray{Float}
: Vector of CST coefficients for upper side of airfoil.lower_coefficients::AbstractArray{Float}
: Vector of CST coefficients for lower side of airfoil.
Keyword Arguments
N::Integer=80
: number of points to use for each sidex::AbstractArray{Float}=split_cosine_spacing(N)
: x-coordinates to use.trailing_edge_yu::Float=0.0
: upper side trailing edge gaptrailing_edge_yl::Float=0.0
: lower side trailing edge gapN1::Float=0.5
: Class shape parameter 1N2::Float=1.0
: Class shape parameter 2split::Bool=false
: if true, returns upper and lower coordinates separately as xl, xu, yl, yu rather than just x, y.
Returns
If split
== false
x::AbstractArray{Float}
: Vector of x coordinates, clockwise from trailing edge.y::AbstractArray{Float}
: Vector of y coordinates, clockwise from trailing edge.
If split
== true
xl::AbstractArray{Float}
: Vector of lower half of x coordinates from trailing edge to leading edge.xu::AbstractArray{Float}
: Vector of upper half of x coordinates from leading edge to trailing edge.yl::AbstractArray{Float}
: Vector of lower half of y coordinates from trailing edge to leading edge.yu::AbstractArray{Float}
: Vector of upper half of y coordinates from leading edge to trailing edge.
FLOWFoil.AirfoilTools.determine_cst
— Functiondetermine_cst(
x,
y;
n_upper_coefficients::Integer=8,
n_lower_coefficients::Integer=8,
trailing_edge_yu=0.0,
trailing_edge_yl=0.0,
N1=0.5,
N2=1.0,
)
Determine best-fit CST parameters using a least squares solve.
Arguments
x::AbstractArray{Float}
: vector of x-coordinates.y::AbstractArray{Float}
: vector of y-coordinates.
Keyword Arguments
n_upper_coefficients::Integer=8
: number of upper side coefficients to fitn_lower_coefficients::Integer=8
: number of lower side coefficients to fittrailing_edge_yu::Float=0.0
: y coordiante of the trailing edge upper surfacetrailing_edge_yl::Float=0.0
: y coordinate of the trailing edge lower surfaceN1::Float=0.5
: Class shape parameter 1N2::Float=1.0
: Class shape parameter 2
Returns
paramters::CST
: CST paramters for airfoil.
determine_cst(
xl,
xu,
yl,
yu;
n_upper_coefficients::Integer=8,
n_lower_coefficients::Integer=8,
N1=0.5,
N2=1.0,
)
Determine best-fit CST parameters for upper and lower sides of airfoil using a least squares solve.
Arguments
xl::AbstractArray{Float}
: vector of lower side x-coordinates.xu::AbstractArray{Float}
: vector of upper side x-coordinates.yl::AbstractArray{Float}
: vector of lower side y-coordinates.yu::AbstractArray{Float}
: vector of upper side y-coordinates.
Keyword Arguments
n_upper_coefficients::Integer=8
: number of upper side coefficients to fitn_lower_coefficients::Integer=8
: number of lower side coefficients to fittrailing_edge_yu::Float=0.0
: y coordiante of the trailing edge upper surfacetrailing_edge_yl::Float=0.0
: y coordinate of the trailing edge lower surfaceN1::Float=0.5
: Class shape parameter 1N2::Float=1.0
: Class shape parameter 2
Returns
parameters::CST
: CST paramters for airfoil.
Circular Arc Camber CST
We also implement a CST-based circular arc camber line airfoil parameterization that may be useful for axial cascade geometry. This parameterization comes from "Aerodynamics of Low Reynolds Number Axial Compressor Sections" by Maffioli et al., and is defined from a camber, maximum thickness, and position of maximum thickness.
using FLOWFoil.AirfoilTools
using Plots
using LaTeXStrings
parameters = CircularArcCST(;
maximum_camber=0.05, maximum_thickness_postition = 0.3, maximum_thickness = 0.12
)
x, y = circular_arc_cst(parameters)
plot(x, y; aspectratio=1, xlabel=L"x", ylabel=L"y", label="")
Spline Parameterizations
Basic B-Spline
The basic B-Spline parameterization comes from "Universal Airfoil Parametrization Using B-Splines" by Rajnarayan, Ning, and Mehr. It is a cubic B-Spline parameterization based on leading edge radius, trailing edge camber and wedge angle, and optional trailing edge gap distance.
using FLOWFoil.AirfoilTools
using Plots
using LaTeXStrings
parameters = BasicBSpline(;
leading_edge_radius=0.015, trailing_edge_camber_angle=12.0, wedge_angle=12.0
)
x, y = basic_bspline(parameters)
plot(x, y; aspectratio=1, xlabel=L"x", ylabel=L"y", label="")
FLOWFoil.AirfoilTools.BasicBSpline
— TypeBasicBSpline
Fields
leading_edge_radius::Float
: leading edge radiustrailing_edge_camber_angle::Float
: trailing edge camber angle (angle of chordline from horizontal at trailing edge).wedge_angle::Float
: Wedge angle (angle between upper and lower surfaces at trailing edge).trailing_edge_gap::Float=0.0
: distance between upper and lower surfaces at trailing edge. A value of zero indicates a sharp trailing edge.third_ctrlpt_position::Float=1.0/3.0
: the position of the third control point. This is an inherent value in the parameterization and if changed, the other parameters will not behave as they are defined here.
FLOWFoil.AirfoilTools.basic_bspline
— Functionbasic_bspline(parameters::BasicBSpline; N=160, split=false, return_nurbs=false)
Obtain airfoil coordinates from a B-Spline parameterization method.
Arguments
parameters::BasicBSpline
: BasicBSpline parameters.
Keyword Arguments
N::Integer=160
: number of points to use when defining the airfoilsplit::Bool
: flag whether to output upper and lower coordinates separatelyreturn_nurbs::Bool
: flag whether to output spline knots and control points as well
Returns
if split=false
x::AbstractArray{Float}
: x-coordinates from lower TE clockwise to upper TEy::AbstractArray{Float}
: y-coordinates from lower TE clockwise to upper TE
if split=true
xu::AbstractArray{Float}
: array of x-coordinates for the upper half of the airfoil (LE to TE)yu::AbstractArray{Float}
: array of y-coordinates for the upper half of the airfoil (LE to TE)xl::AbstractArray{Float}
: array of x-coordinates for the lower half of the airfoil (LE to TE)yl::AbstractArray{Float}
: array of y-coordinates for the lower half of the airfoil (LE to TE)
if return_nurbs=true, also return:
NURBSu::NURBS.NURBScurve
: upper spline objectNURBSl::NURBS.NURBScurve
: lower spline object
gbs(leading_edge_radius, trailing_edge_camber_angle, wedge_angle; perturbations=nothing, trailing_edge_gap=0, degree=3, third_ctrlpt_position=1/3, weights=nothing, split=false, return_nurbs=false)
Obtain airfoil coordinates from a B-Spline parameterization method.
Arguments
leading_edge_radius::Float
: Leading Edge Radiustrailing_edge_camber_angle::Float
: Trailing Edge Camber Angle (degrees)wedge_angle::Float
: Wedge Angle (degrees)
Keyword Arguments
trailing_edge_gap::Float=0
: Trailing Edge Gapthird_ctrlpt_position::Float=1/3
: The x postion of the third control point.split::Bool=false
: flag whether to output upper and lower coordinates separatelyreturn_nurbs::Bool=false
: flag whether to output spline object
Returns
if split=false
x::AbstractArray{Float}
: x-coordinates from lower TE clockwise to upper TEy::AbstractArray{Float}
: y-coordinates from lower TE clockwise to upper TE
if split=true
xu::AbstractArray{Float}
: array of x-coordinates for the upper half of the airfoil (LE to TE)yu::AbstractArray{Float}
: array of y-coordinates for the upper half of the airfoil (LE to TE)xl::AbstractArray{Float}
: array of x-coordinates for the lower half of the airfoil (LE to TE)yl::AbstractArray{Float}
: array of y-coordinates for the lower half of the airfoil (LE to TE)
if return_nurbs=true, also return
NURBSu::NURBS.NURBScurve
: upper spline objectNURBSl::NURBS.NURBScurve
: lower spline object
PARametric SECtion (PARSEC) Parameterizations
Standard PARSEC
We implement the PARSEC method developed by Sobieczky (using this reference). There are 11 parameters in the traditional PARSEC parameterization method. They are as follows:
Parameter | Definition |
---|---|
$r_{LE}$ | Leading Edge Radius |
$x_u$ | Chordwise position of maximum thickness of upper surface |
$x_l$ | Chordwise position of maximum thickness of lower surface |
$y_u$ | y-coordinate at maximum thickness of upper surface |
$y_l$ | y-coordinate at maximum thickness of lower surface |
$y_{xx_u}$ | Second derivative of upper surface curve at point of maximum thickness |
$y_{xx_l}$ | Second derivative of lower surface curve at point of maximum thickness |
$\alpha_{TE}$ | Trailing edge angle |
$\beta_{TE}$ | Boat-tail angle |
$y_{TE}$ | y-position of center of trailing edge |
$\Delta y_{TE}$ | y-distance between upper and lower surface trailing edge points |
using FLOWFoil.AirfoilTools
using Plots
using LaTeXStrings
parameters = PARSEC(;
leading_edge_radius=0.015,
maximum_thickness_xu=0.33,
maximum_thickness_xl=0.20,
maximum_thickness_yu=0.08,
maximum_thickness_yl=-0.04,
curvature_u=-0.63,
curvature_l=0.30,
trailing_edge_angle=-0.05,
boattail_angle=-0.15,
trailing_edge_y=0.0,
trailing_edge_gap=0.0,
)
x, y = parsec(parameters)
plot(x, y; aspectratio=1, xlabel=L"x", ylabel=L"y", label="")
The PARSEC method also has a fit functionality implemented.
x, y = naca4()
plot(
x,
y;
aspectratio=1,
xlabel=L"x",
ylabel=L"y",
label="Coordinates",
linewidth=4,
linestyle=:dash,
)
fit_parameters = determine_parsec(x, y)
x, y = parsec(fit_parameters)
plot!(x, y; label="Fit")
FLOWFoil.AirfoilTools.PARSEC
— TypePARSEC
Fields
leading_edge_radius::Float
: leading edge radiusmaximum_thickness_xu::Float
: x-position of maximum thickness for upper side.maximum_thickness_xl::Float
: x-position of maximum thickness for lower side.maximum_thickness_yu::Float
: value of maximum thickness (from zero) for upper side.maximum_thickness_yl::Float
: value of maximum thickness (from zero) for lower side.curvature_u::Float
: curvature at point of maximum thickness on upper side.curvature_l::Float
: curvature at point of maximum thickness on lower side.trailing_edge_angle::Float
: angle from chordline to horizontal at trailing edge.boattail_angle::Float
: angle from chordline to upper/lower surfaces (half of wedge angle).trailing_edge_gap::Float=0.0
: total gap distance between upper and lower surfaces at treailing edge.trailing_edge_y::Float=0.0
: y-position of midpoint between upper and lower surfaces at trailing edge.
FLOWFoil.AirfoilTools.parsec
— Functionparsec(p::PARSEC; N::Integer=80, split=false)
Calculate the x,y airfoil coordinates for both top and bottom surfaces using standard PARSEC Parameterization method.
Use parsec() for modified PARSEC implementation.
Arguments
p::PARSEC
: PARSEC paramters including:leading_edge_radius
: Leading edge radiusmaximum_thickness_xu
: chordwise position of maximum thickness of upper sidemaximum_thickness_xl
: chordwise position of maximum thickness of lower sidemaximum_thickness_yu
: y-coordinate at maximum thickness of upper sidemaximum_thickness_yl
: y-coordinate at maximum thickness of lower sidecurvature_u
: second derivative of surface geometry at maximum thickness of upper sidecurvature_l
: second derivative of surface geometry at maximum thickness of lower sidetrailing_edge_angle
: trailing edge angleboattail_angle
: boat-tail angletrailing_edge_gap
: y-position of center of trailing edgetrailing_edge_y
: y-distance between upper and lower surface trailing edge points
Keyword Arguments
N::Integer=80
: Number of x stations along chordsplit::Bool
: Flag wheter to split into upper and lower halves.
Returns
If split
== false:
x::AbstractArray{Float}
: Vector of x coordinates, clockwise from trailing edge.y::AbstractArray{Float}
: Vector of y coordinates, clockwise from trailing edge.
If split
== true:
xl::AbstractArray{Float}
: Vector of lower half of x coordinates from trailing edge to leading edge.xu::AbstractArray{Float}
: Vector of upper half of x coordinates from leading edge to trailing edge.yl::AbstractArray{Float}
: Vector of lower half of y coordinates from trailing edge to leading edge.yu::AbstractArray{Float}
: Vector of upper half of y coordinates from leading edge to trailing edge.
parsec(p::AbstractArray{Float}; N::Integer=80, split=false)
Calculate the x,y airfoil coordinates for both top and bottom surfaces using standard PARSEC Parameterization method.
Use parsec() for modified PARSEC implementation.
Arguments
p::AbstractArray{Float}
: PARSEC paramters including:leading_edge_radius
: Leading edge radiusmaximum_thickness_xu
: chordwise position of maximum thickness of upper sidemaximum_thickness_xl
: chordwise position of maximum thickness of lower sidemaximum_thickness_yu
: y-coordinate at maximum thickness of upper sidemaximum_thickness_yl
: y-coordinate at maximum thickness of lower sidecurvature_u
: second derivative of surface geometry at maximum thickness of upper sidecurvature_l
: second derivative of surface geometry at maximum thickness of lower sidetrailing_edge_angle
: trailing edge angleboattail_angle
: boat-tail angletrailing_edge_gap
: y-position of center of trailing edgetrailing_edge_y
: y-distance between upper and lower surface trailing edge points
Keyword Arguments
N::Integer=80
: Number of x stations along chordsplit::Bool
: Flag wheter to split into upper and lower halves.
Returns
If split
== false:
x::AbstractArray{Float}
: Vector of x coordinates, clockwise from trailing edge.y::AbstractArray{Float}
: Vector of y coordinates, clockwise from trailing edge.
If split
== true:
xl::AbstractArray{Float}
: Vector of lower half of x coordinates from trailing edge to leading edge.xu::AbstractArray{Float}
: Vector of upper half of x coordinates from leading edge to trailing edge.yl::AbstractArray{Float}
: Vector of lower half of y coordinates from trailing edge to leading edge.yu::AbstractArray{Float}
: Vector of upper half of y coordinates from leading edge to trailing edge.
FLOWFoil.AirfoilTools.determine_parsec
— Functiondetermine_parsec(x,y)
Uses LsqFit to go from x-y coordinates to standard PARSEC parameters.
Arguments
x::AbstractArray{Float}
: vector of x coordinatesy::AbstractArray{Float}
: vector of y coordinates
Returns
parameters::PARSEC
: an parameter object of type PARSEC
Modified PARSEC
Also implemented in AirfoilTools is a modified PARSEC parameterization that gives direct control to the trailing edge surfaces of the upper and lower sides. In order to create a more intuitive parameterization, the trailing edge parameters, both position and angles are replaced to be directly the trailing edge positions of the upper and lower surface and trailing edge surface angles.
Original Parameter | Modified Parameter | New Definition |
---|---|---|
$\alpha_{TE}$ | $\theta_{TE_u}$ | Upper surface trailing edge tangent angle |
$\beta_{TE}$ | $\theta_{TE_l}$ | Upper surface trailing edge tangent angle |
$y_{TE}$ | $y_{TE_u}$ | y-coordinate of upper surface trailing edge |
$\Delta y_{TE}$ | $y_{TE_l}$ | y-coordinate of lower surface trailing edge |
using FLOWFoil.AirfoilTools
using Plots
using LaTeXStrings
parameters = ModifiedPARSEC(;
leading_edge_radius=0.015,
maximum_thickness_xu=0.33,
maximum_thickness_xl=0.20,
maximum_thickness_yu=0.08,
maximum_thickness_yl=-0.04,
curvature_u=-0.63,
curvature_l=0.30,
trailing_edge_tangent_u=-0.2,
trailing_edge_tangent_l=0.1,
trailing_edge_yu=0.0,
trailing_edge_yl=0.0,
)
x, y = modified_parsec(parameters)
plot(x, y; aspectratio=1, xlabel=L"x", ylabel=L"y", label="")
The Modified PARSEC method also has a fit functionality implemented.
x, y = naca4()
plot(
x,
y;
aspectratio=1,
xlabel=L"x",
ylabel=L"y",
label="Coordinates",
linewidth=4,
linestyle=:dash,
)
fit_parameters = determine_modified_parsec(x, y)
x, y = modified_parsec(fit_parameters)
plot!(x, y; label="Fit")
FLOWFoil.AirfoilTools.ModifiedPARSEC
— TypeModifiedPARSEC
Fields
leading_edge_radius::Float
: leading edge radiusmaximum_thickness_xu::Float
: x-position of maximum thickness for upper side.maximum_thickness_xl::Float
: x-position of maximum thickness for lower side.maximum_thickness_yu::Float
: value of maximum thickness (from zero) for upper side.maximum_thickness_yl::Float
: value of maximum thickness (from zero) for lower side.curvature_u::Float
: curvature at point of maximum thickness on upper side.curvature_l::Float
: curvature at point of maximum thickness on lower side.trailing_edge_tangent_u::Float
: angle of surface at upper side trailing edge.trailing_edge_tangent_l::Float
: angle of surface at lower side trailing edge.trailing_edge_yu::Float=0.0
: y-position of upper side trailing edge.trailing_edge_yl::Float=0.0
: y-position of lower side trailing edge.
FLOWFoil.AirfoilTools.modified_parsec
— Functionparsec(p::ModifiedPARSEC; N::Int=80, split=false)
Calculate the x,y airfoil coordinates for both top and bottom surfaces using modified PARSEC Parameterization method.
Use parsec() for standard PARSEC implementation. This modified version employs direct values for trailing edge position and angles for each surface.
Arguments
p::ModifiedPARSEC
: ModifiedPARSEC paramters including:leading_edge_radius
: Leading edge radiusmaximum_thickness_xu
: chordwise position of maximum thickness of upper sidemaximum_thickness_xl
: chordwise position of maximum thickness of lower sidemaximum_thickness_yu
: y-coordinate at maximum thickness of upper sidemaximum_thickness_yl
: y-coordinate at maximum thickness of lower sidecurvature_u
: second derivative of surface geometry at maximum thickness of upper sidecurvature_l
: second derivative of surface geometry at maximum thickness of lower sidetrailing_edge_tangent_u
: trailing edge tangent angle of upper sidetrailing_edge_tangent_l
: trailing edge tangent angle of lower sidetrailing_edge_yu
: y-position of trailing edge of upper sidetrailing_edge_yl
: y-position of trailing edge of lower side
Keyword Arguments
N::Integer=80
: Number of x stations along chordsplit::Bool
: Flag wheter to split into upper and lower halves.
Returns
If split
== false:
x::AbstractArray{Float}
: Vector of x coordinates, clockwise from trailing edge.y::AbstractArray{Float}
: Vector of y coordinates, clockwise from trailing edge.
If split
== true:
xl::AbstractArray{Float}
: Vector of lower half of x coordinates from trailing edge to leading edge.xu::AbstractArray{Float}
: Vector of upper half of x coordinates from leading edge to trailing edge.yl::AbstractArray{Float}
: Vector of lower half of y coordinates from trailing edge to leading edge.yu::AbstractArray{Float}
: Vector of upper half of y coordinates from leading edge to trailing edge.
parsec(p::AbstractArray{Float}; N::Int=80, split=false)
Calculate the x,y airfoil coordinates for both top and bottom surfaces using modified PARSEC Parameterization method.
Use parsec() for standard PARSEC implementation. This modified version employs direct values for trailing edge position and angles for each surface.
Arguments
p::AbstractArray{Float}
: ModifiedPARSEC paramters including:leading_edge_radius
: Leading edge radiusmaximum_thickness_xu
: chordwise position of maximum thickness of upper sidemaximum_thickness_xl
: chordwise position of maximum thickness of lower sidemaximum_thickness_yu
: y-coordinate at maximum thickness of upper sidemaximum_thickness_yl
: y-coordinate at maximum thickness of lower sidecurvature_u
: second derivative of surface geometry at maximum thickness of upper sidecurvature_l
: second derivative of surface geometry at maximum thickness of lower sidetrailing_edge_tangent_u
: trailing edge tangent angle of upper sidetrailing_edge_tangent_l
: trailing edge tangent angle of lower sidetrailing_edge_yu
: y-position of trailing edge of upper sidetrailing_edge_yl
: y-position of trailing edge of lower side
Keyword Arguments
N::Integer=80
: Number of x stations along chordsplit::Bool
: Flag wheter to split into upper and lower halves.
Returns
If split
== false:
x::AbstractArray{Float}
: Vector of x coordinates, clockwise from trailing edge.y::AbstractArray{Float}
: Vector of y coordinates, clockwise from trailing edge.
If split
== true:
xl::AbstractArray{Float}
: Vector of lower half of x coordinates from trailing edge to leading edge.xu::AbstractArray{Float}
: Vector of upper half of x coordinates from leading edge to trailing edge.yl::AbstractArray{Float}
: Vector of lower half of y coordinates from trailing edge to leading edge.yu::AbstractArray{Float}
: Vector of upper half of y coordinates from leading edge to trailing edge.
FLOWFoil.AirfoilTools.determine_modified_parsec
— Functiondetermine_modified_parsec(x,y)
Uses LsqFit to go from x-y coordinates to modified ModifiedPARSEC parameters.
Arguments
x::AbstractArray{Float}
: vector of x coordinatesy::AbstractArray{Float}
: vector of y coordinates
Returns
parameters::ModifiedPARSEC
: an parameter object of type ModifiedPARSEC
Contributing
We welcome the addition of other common parameterizations, as well as improvements/additions to currently implemented options. Additions should have outputs consistent with current parameterizations.