In this example we use the lifting line method to solve the aerodynamics of a $36^\circ$ swept-back wing with a cambered NACA 64(1)-612 airfoil from the NACA RM A50K27 Flying Wing report.
#=##############################################################################
# DESCRIPTION
36deg swept-back wing with a cambered NACA 64(1)-612 airfoil from the NACA
RM A50K27 Flying Wing report. This wing has an aspect ratio of 15 and tapper
ratio of 0.5. The experiment was run at Re_c = 2e6 and Mach 0.27.
This example uses the nonlinear lifting line method.
# AUTHORSHIP
* Author : Eduardo J. Alvarez
* Email : Edo.AlvarezR@gmail.com
* Created : Jan 2026
* License : MIT License
=###############################################################################
import FLOWPanel as pnl
import FLOWPanel: mean, norm, dot, cross
import PyPlot as plt
run_name = "ll-a50k27" # Name of this run
save_path = run_name # Where to save outputs
airfoil_path = joinpath(pnl.examples_path, "data") # Where to find 2D polars
paraview = false # Whether to visualize with Paraview
# ----------------- SIMULATION PARAMETERS --------------------------------------
magUinf = 91.3 # (m/s) freestream velocity
rho = 1.225 # (kg/m^3) air density
# ------------------ GEOMETRY PARAMETERS ---------------------------------------
# High-level parameters
b = 2*1.5490 # (m) wing span
# Discretization parameters
nelements = 40 # Number of stripwise elements per semi-span
discretization = [ # Multi-discretization of wing (seg length, ndivs, expansion, central)
(1.00, nelements, 1/10, false),
]
symmetric = true # Whether the wing is symmetric
# Chord distribution (nondim y-position, nondim chord)
chord_distribution = [
# 2*y/b c/b
0.0 0.132989
1.0 0.0664945
]
# Twist distribution (nondim y-position, twist)
twist_distribution = [
# 2*y/b twist (deg)
0.0 0.0
1.0 0.0
]
# Sweep distribution (nondim y-position, sweep)
sweep_distribution = [
# 2*y/b sweep (deg)
0.0 36.2603
1.0 36.2603
]
# Dihedral distribution (nondim y-position, dihedral)
dihedral_distribution = [
# 2*y/b dihedral (deg)
0.0 0.0
1.0 0.0
]
# Span-axis distribution: chordwise point about which the wing is twisted,
# swept, and dihedralized (nondim y-position, nondim chord-position)
spanaxis_distribution = [
# 2*y/b x/c
0.0 0.25
1.0 0.25
]
# Airfoil contour distribution (nondim y-position, polar, airfoil type)
airfoil_distribution = [
# 2*y/b polar file airfoil type
(0.00, "n64_1_A612-Re0p5e6-neuralfoil180-3.csv", pnl.SimpleAirfoil),
(1.00, "n64_1_A612-Re0p5e6-neuralfoil180-3.csv", pnl.SimpleAirfoil)
]
element_optargs = (; path = airfoil_path,
plot_polars = true,
extrapolatepolar = false, # Whether to extrapolate the 2D polars to ±180 deg
)
# ------------------ SOLVER PARAMETERS -----------------------------------------
deltasb = 1.0 # Blending distance, deltasb = 2*dy/b
deltajoint = 1.0 # Joint distance, deltajoint = dx/c
sigmafactor = 0.0 # Dragging line amplification factor (set to -1.0 for rebust post-stall method)
sigmaexponent = 1.0 # Dragging line amplification exponent (no effects if `sigmafactor==0.0`)
# Nonlinear solver
solver = pnl.SimpleNonlinearSolve.SimpleDFSane() # Indifferent to initial guess, but somewhat not robust
# solver = pnl.SimpleNonlinearSolve.SimpleTrustRegion() # Trust region needs a good initial guess, but it converges very reliably
solver_optargs = (;
abstol = 1e-9,
maxiters = 800,
)
align_joints_with_Uinfs = true # Whether to align joint bound vortices with the freestream
use_Uind_for_force = false # Whether to use Uind as opposed to selfUind for force postprocessing
# (`true` for more accurate spanwise cd distribution, but worse integrated CD)
X0 = [0.595, 0, 0] # Center about which to calculate moments
cref = 0.2472 # (m) reference chord
Aref = 0.957 # (m^2) reference area
nondim = 0.5*rho*magUinf^2*Aref # Normalization factor
# ------------------ GENERATE LIFTING LINE -------------------------------------
ll = pnl.LiftingLine{Float64}(
airfoil_distribution;
b, chord_distribution, twist_distribution,
sweep_distribution, dihedral_distribution,
spanaxis_distribution,
discretization,
symmetric,
deltasb, deltajoint, sigmafactor, sigmaexponent,
element_optargs,
plot_discretization = true,
)
display(ll)
# ------------------ OUTPUT GEOMETRY -------------------------------------------
if !isnothing(save_path)
str = pnl.save(ll, run_name; path=save_path, debug=true) # Use `debug=true` to output the effective horseshoes
if paraview
run(`paraview --data=$(joinpath(save_path, str))`)
end
end
# ----------------- AOA SWEEP --------------------------------------------------
# Sequence of sweeps to run
# NOTE: To help convergence and speed it up the sweep, we recommend starting
# each sweep at 0deg AOA since the sweep steps over points using the last
# converged solution as the initial guess
sweep_neg = range(0, -50, step=-0.5) # Sweep from 0 into deep negative stall (-50deg)
sweep_pos = range(0, 50, step=0.5) # Sweep from 0 into deep positive stall (50deg)
@time wingpolar = pnl.run_polarsweep(ll,
magUinf, rho, X0, cref, b;
Aref,
aoa_sweeps = (sweep_neg, sweep_pos),
solver,
solver_optargs,
align_joints_with_Uinfs,
use_Uind_for_force
)
(see the complete example under examples/liftingline_a50k27.jl to see how to postprocess the spanwise loading that is plotted below)
Wing Polar
Wing Polar Post-Stall
You can also automatically run this example and generate these plots with the following command:
import FLOWPanel as pnl
include(joinpath(pnl.examples_path, "liftingline_a50k27.jl"))