calculate_AIC_matrices!(body_aero::BodyAerodynamics, model::Model, 
                     core_radius_fraction,
                     va_norm_array, 
                     va_unit_array)

Calculate Aerodynamic Influence Coefficient matrices.

See also: BodyAerodynamics, Model

Returns: nothing

gamma_loop!(solver::Solver, AIC_x::Matrix{Float64},
          AIC_y::Matrix{Float64}, AIC_z::Matrix{Float64},
          panels::AbstractVector{<:Panel}, relaxation_factor::Float64; log=true)

Main iteration loop for calculating circulation distribution.

When solver.is_with_artificial_viscosity is set, the LOOP solver replaces the explicit target F(gamma) with the implicit Li/Gaunaa solution (I - diag(mu) L) gamma = F(gamma) before relaxation, stabilizing post-stall distributions. The Python reference gates this on precomputed per-panel stall angles; here the polar is an interpolation object, so the expensive linear solve is gated on any(mu > 0) instead, which is behaviourally identical.

build_spanwise_laplacian!(laplacian, n_panels)

Fill the n_panels × n_panels matrix laplacian with the discrete spanwise Laplacian used by the Li/Gaunaa post-stall artificial-viscosity regularization. Interior rows use the three-point stencil gamma[i-1] - 2 gamma[i] + gamma[i+1]; the tip rows use the second-order closures of Li, Gaunaa, Pirrung & Lønbæk (TORQUE 2026, Eq. 15) that enforce gamma -> 0 at the wing tips. Panels are assumed ordered consecutively along the span and approximately uniformly spaced. Returns the matrix unchanged (all zeros) when n_panels < 3.

local_lift_slope!(slopes, panels, alpha_dist, delta=deg2rad(0.5))

Per-panel local lift-curve slope dCl/dalpha, evaluated from each panel's own 2-D polar at its current effective angle of attack by central differences. The slope turns negative in post-stall, which is what activates the artificial-viscosity regularization in gamma_loop!.

apply_artificial_viscosity!(gamma, panels, alpha_dist, laplacian, viscosity_matrix,
                            lift_slope, mu_array, gamma_target, planform_area, factor)

Apply one implicit Li/Gaunaa artificial-viscosity step to gamma in place and return true when it fired. The per-panel viscosity is mu_i = max(0, -factor * planform_area * Cl'_i / width_i^2), with the local lift slope Cl'_i from local_lift_slope! evaluated at alpha_dist. When any panel is post-stall (mu_i > 0), gamma is replaced by the solution of (I - diag(mu) L) gamma = gamma, with L the spanwise Laplacian in laplacian (see build_spanwise_laplacian!); otherwise gamma is left unchanged. The remaining arguments are preallocated work buffers reused across iterations.

frozen_wake(body_aero::BodyAerodynamics, va_distribution)

Update the filaments of the panels with frozen wake model. Uses one shared wake vector computed from area-weighted distributed inflow.

Replaces older filaments if present by checking length of filaments.

Arguments

• body_aero::BodyAerodynamics: see: BodyAerodynamics
• va_distribution::Matrix{Float64}: Array of velocity vectors at each panel

Returns

• nothing

calc_forces!(solver::Solver, body_aero::BodyAerodynamics;
             reference_point=solver.reference_point, moment_frac=0.1)

Assemble aerodynamic forces and moments from the circulation already converged by solve_base! and stored in solver. Split out of solve!.

calculate_cl(panel::Panel, alpha)
calculate_cl(panel::Panel, alpha, delta)

Calculate lift coefficient for given angle of attack alpha [rad]. The 3-arg form evaluates the (α, δ) polar (POLAR_MATRICES) at the passed flap deflection delta [rad] instead of the panel's stored delta; the 2-arg form forwards with panel.delta. Other aero models ignore delta.

Returns

• Float64: Lift coefficient (Cl)

calculate_cd(panel::Panel, alpha)
calculate_cd(panel::Panel, alpha, delta)

Calculate the drag coefficient for the given angle of attack. The 3-arg form evaluates the (α, δ) polar at the passed flap deflection delta; see calculate_cl.

calculate_cm(panel::Panel, alpha)
calculate_cm(panel::Panel, alpha, delta)

Calculate the pitching-moment coefficient for the given angle of attack. The 3-arg form evaluates the (α, δ) polar at the passed flap deflection delta; see calculate_cl.

calculate_cd_cm(panel::Panel, alpha)

Calculate drag and moment coefficients for the given angle of attack (calculate_cd and calculate_cm).

set_pitch_rate_dist!(body_aero, omega)

Fill body_aero.pitch_rate_dist from a rigid-body turn rate by projecting it onto each panel's own spanwise axis. Panels with different dihedral see different rates from the same omega.

calculate_relative_alpha_and_velocity(panel::Panel, induced_velocity)

Calculate relative angle of attack and relative velocity of the panel.

calculate_relative_alpha_and_relative_velocity(panel::Panel, induced_velocity::Vector{Float64})

Calculate the relative angle of attack and relative velocity of the panel.

Arguments

• panel::Panel: The panel object
• induced_velocity::Vector{Float64}: Induced velocity at the control point

Returns

• Tuple{Float64,Vector{Float64}}: Tuple containing:
  • alpha: Relative angle of attack of the panel (in radians)
  • relative_velocity: Relative velocity vector of the panel

update_effective_angle_of_attack_if_VSM(body_aero::BodyAerodynamics, gamma,
                                      core_radius_fraction,
                                      z_airf_array,
                                      x_airf_array,
                                      va_array,
                                      va_norm_array,
                                      va_unit_array)

Update angle of attack at aerodynamic center for VSM method.

Returns: nothing

calculate_stall_angle_list(panels::Vector{<:Panel};
                         begin_aoa=9.0,
                         end_aoa=22.0,
                         step_aoa=1.0,
                         stall_angle_if_none_detected=50.0,
                         cl_initial=-10.0)

Calculate stall angles for each panel.

Returns: Vector{Float64}: Stall angles in radians

wing_span_flip(wing) -> Int8

-1 when wing's sections run against its spanwise_direction, +1 otherwise: the flip every panel of the wing is reinitialized with (reinit!). One answer per wing, so neighbouring panels cannot disagree and invert a single normal by 180°.

calculate_circulation_distribution_elliptical_wing(body_aero::BodyAerodynamics, gamma_0=1.0)

Calculate circulation distribution for an elliptical wing.

Returns: nothing

_compute_reference_velocity_from_distribution(va_input, n_panels, panel_areas=nothing)

Return a single reference velocity vector from uniform or distributed inflow. For distributed inflow, the speed is area-weighted RMS and the direction is the area-weighted mean direction.

smooth_circulation!(damp, circulation, 
                  smoothness_factor::Float64, 
                  damping_factor::Float64)

Smooth circulation distribution if needed.

Returns:

• Tuple of smoothed circulation and boolean indicating if smoothing was applied

smooth_distribution!(dist, window_size)

Apply moving average smoothing to a distribution in-place. Uses a centered window of size window_size (must be odd).

Arguments

• dist::Vector{Float64}: Distribution to smooth (modified in-place)
• window_size::Int: Size of smoothing window (must be odd)

make_dual_shadow(solver, body_aero, ::Type{TD}) where TD

Build a BodyAerodynamics/Solver pair whose mutating buffers carry element type TD (typically a ForwardDiff.Dual). Wing geometry, polar coefficients and interpolators are reused from the Float64 originals; circulation/AIC/scratch buffers are freshly allocated as TD-typed.

Length smooth_norm adds in quadrature to keep a norm positive at zero.

smooth_norm(vec, floor=SMOOTH_FLOOR)

Norm floored in quadrature: positive and differentiable at the origin, and branch-free, which is what makes the functions below traceable.

panel_span_vector(le_1, te_1, le_2, te_2)

Quarter-chord vector from section 2 to section 1: the panel's spanwise axis, its length the span width.

panel_chord_weight(width_prev, width_own, width_next)

Section 1's share of the edge blend setting a panel's chord direction, weighted by panel spacing so a panel leans towards a narrower neighbour. nothing for a missing neighbour at a tip; a lone panel gets 0.5.

panel_axes(le_1, te_1, le_2, te_2, chord_weight=0.5, orient=1)

Airfoil frame and size of the panel between two sections, as (; x_airf, y_airf, z_airf, chord, width).

chord_weight (panel_chord_weight) enters as an offset from the midpoint rather than as w·p₁ + (1-w)·p₂: the two are equal, but the offset form leaves a constant term to fold, which a symbolic consumer builds several times faster. orient is ±1, flipping y_airf/z_airf so the frame does not depend on section ordering.

panel_axes(panel::Panel)

The panel's stored airfoil frame and size in the shape panel_axes returns, so a panel built by the geometry pass feeds the same functions as one built straight from section points.

effective_alpha(alpha, deficiency)

Angle the polars are read at: the geometric inflow angle less an unsteady lag. The geometric angle still turns the force, so a lag shifts the coefficients only.

panel_inflow(axes, va_1, va_2, v_ind, dva_1=nothing, dva_2=nothing,
             deficiency=0)

Flow a panel sees, as (; v_eff, alpha, alpha_eff, v_span, pitch_rate), where v_span is the effective velocity across the span. axes is a panel_axes result. dva_1/dva_2 are the sections' trailing minus leading edge apparent wind, giving the section_pitch_rate; nothing leaves it zero. deficiency feeds effective_alpha.

dynamic_pressure(rho_1, rho_2, v_span)

Panel dynamic pressure from its two section densities and the panel_inflow span-wise velocity, as that vector or its magnitude.

flow_curvature_cm(pitch_rate, chord, v_rel)

Thin-airfoil quarter-chord moment increment of a section pitching about its own spanwise axis: Δcm = -(π/4)·q̂ with q̂ = q·c/(2·v_rel), q positive nose-up. Pivot-independent, and lift needs no matching correction because the inflow is already sampled at three-quarter chord.

panel_force_directions(axes, alpha_dir, spanwise)

Lift and drag unit vectors, (; dir_lift, dir_drag). alpha_dir is the angle that turns the force, spanwise the wing's spanwise direction.

Panel pitching moment per unit span about y_airf, positive nose-up.

panel_couple_force(cm, q_dyn, chord, width, scale=1)

panel_moment as a force couple: the magnitude two opposed forces one chord apart need to produce it. A particle model places the moment this way, along the panel normal at its leading and trailing edge.

panel_loads(axes, dirs, q_dyn, cl, cd, cm, scale=1)

Panel load from its polar coefficients, panel_axes and panel_force_directions, as (; lift, drag, moment, force, pitching_moment). The first three are per unit span; force and pitching_moment are the whole panel's, scale included.

velocity_3D_bound_vortex(vel, filament::BoundFilament, XVP, 
                       gamma, core_radius_fraction, work_vectors)

Calculate induced velocity by a bound vortex filament at a point in space.

velocity_3D_trailing_vortex(vel, filament::BoundFilament, 
                          XVP, gamma, v_a, work_vectors)

Calculate induced velocity by a trailing vortex filament.

Arguments

• XVP: Control point coordinates
• gamma: Vortex strength
• v_a: Inflow velocity magnitude
• work_vectors: preallocated array of intermediate variables

Reference: Rick Damiani et al. "A vortex step method for nonlinear airfoil polar data as implemented in KiteAeroDyn".

velocity_3D_trailing_vortex_semiinfinite(filament::SemiInfiniteFilament, 
                                         Vf, XVP, GAMMA, v_a, work_vectors)

Calculate induced velocity by a semi-infinite trailing vortex filament.

calculate_velocity_induced_bound_2D!(U2D, panel::Panel, evaluation_point, work_vectors)

Calculate velocity induced by bound vortex filaments at the control point. Only needed for VSM, as LLT bound and filament align, thus no induced velocity.

Arguments

• U_2D: Resulting 2D velocity vector
• panel::Panel: Panel object
• evaluation_point: Point where induced velocity is evaluated
• work_vectors: Pre-allocated temporary variables

calculate_velocity_induced_single_ring_semiinfinite!(
    velind::MVec3,
    tempvel::MVec3,
    filaments,
    evaluation_point::MVec3,
    evaluation_point_on_bound::Bool,
    va_norm::Float64,
    va_unit::MVec3,
    gamma::Float64,
    core_radius_fraction::Float64,
    work_vectors::NTuple{10, MVec3}
)

Calculate the velocity induced by a vortex ring at a control point.

Arguments

• velind
• tempvel
• filaments
• evaluation_point::MVec3: Point where induced velocity is evaluated
• evaluation_point_on_bound::Bool: Whether evaluation point is on bound vortex
• va_norm::Float64: Norm of apparent velocity
• va_unit::MVec3: Unit vector of apparent velocity
• gamma::Float64: Circulation strength
• core_radius_fraction::Float64: Vortex core radius as fraction of panel width
• work_vectors::NTuple{10, MVec3} Pre-allocated temporary variables

Returns

• nothing

cross3!(result::AbstractVector{T}, a::AbstractVector{T}, b::AbstractVector{T}) where T

Compute cross product of 3D vectors in-place.

update_panel_properties!(panel_props::PanelProperties, section_list::Vector{Section}, n_panels::Int)

Update geometric properties for each panel.

Arguments

• section_list::Vector{Section}: List of Sections
• n_panels::Int: Number of Panels

Returns:

nothing, updates the PanelProperties in-place

build_interps(section_1, section_2, remove_nan) -> (cl, cd, cm, section_aero)

Build the averaged aerodynamic interpolations for the panel between two sections. Returns (cl_interp, cd_interp, cm_interp, section_aero), each nothing for models that do not use it (INVISCID, POLY). cl/cm clamp (Flat) past the alpha range but extrapolate linearly (Line) over delta; cd extrapolates linearly in both. The concrete return types parameterise Panel — see panel_interp_types.

panel_interp_types(section, remove_nan) -> (CL, CD, CM, CP)

Field types for a Panel whose aero model matches section. Since a wing has one aero model, every panel shares these, so Vector{Panel{...}} is concretely typed.

reinit!(wing::AbstractWing)

Reinitialize wing panel properties based on current refined_sections geometry.

This function only updates panel properties (chord, area, etc.) from the existing refinedsections. It does NOT refine the mesh - call refineaerodynamic_mesh!(wing) first if needed.

Note

After deformation via unrefined_deform!() or deform!(), call reinit! to update panel properties while preserving the deformed geometry.

reinit!(panel, section_1, section_2, aero_center, control_point, bound_point_1,
        bound_point_2, x_airf, y_airf, z_airf, delta, vec; kwargs...)

Reinitialize a panel's geometry, horseshoe filaments and aerodynamic interpolations.

flip reverses the section order so y_airf points along +spanwise_direction and z_airf to the airfoil upper surface, making the aero independent of section ordering. The caller owns it and must not derive it from the live geometry, which would invert normals mid-run; reinit!(::BodyAerodynamics) decides it once per wing.

rotated_te(le, te, y_hat, θ)

Compute the trailing-edge position after rotating the chord vector te - le around y_hat by θ radians (Rodrigues). Returns an SVector{3} (stack-allocated, zero heap allocs).

calculate_filaments_for_plotting(panel::Panel)

Calculate filaments for plotting with their positions and colors.

Returns

• Vector{Tuple{Vector{Float64}, Vector{Float64}, String}}: List of tuples containing:
  • First point (x1)
  • Second point (x2)
  • Color string

unrefined_deform!(wing::Wing, theta_angles=nothing, delta_angles=nothing)

Apply deformation angles defined per unrefined section.

Refined-section twist and TE-deflection values are computed by linear interpolation from the unrefined-section inputs using the precomputed refined_section_left_idx / refined_section_weight cache built at refinement time. Endpoint refined sections take the unrefined endpoint values exactly. The panel-level theta_dist / delta_dist arrays are then filled by averaging adjacent refined-section values, so downstream consumers (solver, body aerodynamics) see a per-panel value.

Arguments

• wing::Wing: Wing to deform (must have nondeformedsections, populated by refine! for manual/YAML wings or by OBJ refinement for OBJ-based wings).
• theta_angles::AbstractVector: Twist angles in radians, one per unrefined section. Pass nothing to leave twist unchanged.
• delta_angles::AbstractVector: TE deflection angles in radians, one per unrefined section. Pass nothing to leave deflection unchanged.

Keyword arguments

• smooth, smooth_window: accepted for backwards compatibility with callers of deform!, but ignored here — the linear interpolation between unrefined sections is already smooth, so no post-hoc smoothing is applied.

deform!(wing::Wing, theta_dist::AbstractVector, delta_dist::AbstractVector;
        smooth=false, smooth_window=nothing)

Deform wing by applying theta and delta distributions at the panel level.

Arguments

• wing::Wing: Wing to deform (must support deformation)
• theta_dist::AbstractVector: Twist angles for each panel (length = n_panels)
• delta_dist::AbstractVector: TE deflections for each panel (length = n_panels)
• smooth::Bool: Whether to apply smoothing (default: false)
• smooth_window::Union{Nothing, Int}: Smoothing window size (default: auto-calculated)

Effects

Updates wing.refinedsections with deformed geometry based on wing.nondeformed_sections

deform!(wing::Wing; smooth=false, smooth_window=nothing)

Apply stored thetadist and deltadist to deform the wing geometry. Converts panel angles (npanels) to section angles (npanels+1) by averaging adjacent panels.

Arguments

• wing::Wing: Wing to deform (must have nondeformedsections)
• smooth::Bool: Whether to apply smoothing to thetadist and deltadist (default: false)
• smooth_window::Union{Nothing, Int}: Smoothing window size (default: auto-calculated)

Effects

Updates wing.refinedsections based on wing.nondeformed_sections and stored distributions

compute_refined_panel_mapping!(wing::AbstractWing)

Compute the mapping from refined panels to unrefined sections by finding the closest unrefined section for each refined panel (based on section center distance). Maps each refined panel index to its corresponding unrefined section index (1 to nunrefinedsections). Works after refinement is complete.

compute_refined_section_interpolation!(wing::AbstractWing;
                                       reuse_aero_data=false)

Compute per-refined-section linear-interpolation weights from the unrefined sections. For refined section i, the interpolated value is:

out[i] = weight[i] * unrefined[left_idx[i]] +
         (1 - weight[i]) * unrefined[left_idx[i] + 1]

Positions are quarter-chord arc-length along the unrefined and refined sections. The first refined section is pinned to left_idx == 1, weight == 1 (returns unrefined[1] exactly) and the last refined section to left_idx == n_unref - 1, weight == 0 (returns unrefined[end] exactly).

reuse_aero_data keeps the surface tables the refined sections already hold instead of reblending them from the unrefined ones, the same preservation refine! applies to aero_data under use_prior_polar. A remesh that replaces the unrefined sections with a coarser set would otherwise resample the surface tables down to that set even while the polars stay at full resolution.

copy_sections(sections) -> Vector{Section}

Copy a vector of Sections into fresh objects with their own LE_point/TE_point storage; the read-only aero_data/section_aero tables are shared by reference.

copy_sections_to_refined!(wing; reuse_aero_data=false)

Copy unrefined sections to refined sections 1:1 (no interpolation). If refined_sections is empty, allocates via copy; otherwise reinitialises in-place. Warns if billowing was requested but there are no intermediate sections to billow.

_apply_refined_section_thetas!(wing, section_thetas)

Rotate each refined section's TE point around its LE point by the given per-section twist angle, using wing.non_deformed_sections as the reference geometry.

The rotation is a full Rodrigues rotation about the spanwise axis. A swept or dihedral section has a chord component along that axis, and dropping the axial term would scale it by cos(theta), shortening the chord and tilting it.

_panel_thetas_to_section_thetas!(section_thetas, panel_thetas)

Convert panel-level twist angles (length n) to refined-section angles (length n+1). Interior sections take the average of adjacent panels. Boundary sections use linear extrapolation from the two nearest panels, so a linear input produces a linear output across the full refined-section range.

_interpolate_unrefined_to_refined(wing, unrefined_values) -> Vector

Linearly interpolate unrefined_values (length n_unrefined_sections) to refined sections (length n_panels + 1) using the cached weights. Endpoints match exactly.

span_order_key(section) -> Float64

Spanwise coordinate that normalize_span_order! and refine!'s sort_sections order sections on.

normalize_span_order!(sections) -> sections

Reverse sections if they do not already run +y to -y, the order panel normals are built from. Use refine!'s sort_sections for a scrambled list.

can_reuse_prior_refined_surface_tables(wing) -> Bool

Whether every refined section already carries a SectionAero, so a remesh can keep them instead of reblending from the unrefined sections. False on a wing that has none, where the blend is what fills them in the first place.

refine_mesh_for_linear_cosine_distribution!(wing, idx, dist,
    n_sections, sections; endpoints, reuse_aero_data)

Refine wing mesh using linear or cosine spacing. Reads LE/TE directly from a Vector{Section} (zero matrix allocations).

refine_mesh_by_splitting_provided_sections!(wing; reuse_aero_data,
    billowing_percentage)

Refine mesh by splitting provided sections into desired number of panels.

When billowing_percentage > 0, rotates chord vectors around the leading edge with a sinusoidal profile to simulate fabric billowing between ribs.

refine_mesh_with_billowing!(wing; reuse_aero_data)

Refine wing mesh using SPLIT_PROVIDED spacing with TE billowing.

Between each pair of unrefined (rib) sections, chord vectors are rotated around the leading edge to simulate fabric billowing. The rotation amplitude is found iteratively so that the TE arc length matches the wing's billowing_percentage.

Delegates to refine_mesh_by_splitting_provided_sections!.

apply_billowing_to_pair!(sections, start_si, end_si, y_hat,
                         span_len, le_ref, te_left, te_right,
                         percentage)

Apply billowing to refined sections between two ribs by rotating each section's chord around y_hat. Uses Newton iteration to find the rotation amplitude that matches the target TE arc-length percentage, then applies the rotations in-place.

The angle profile is angle_max * sin(π t) (zero at ribs, maximum at centre). All angles are in radians.

Non-allocating: modifies sections[start_si:end_si].TE_point in-place using scalar arithmetic only.

billowing_arc_length(sections, start_si, end_si, y_hat,
                     span_len, le_ref, te_left, te_right,
                     angle_max)

Compute the TE arc length that would result from rotating each section's chord around y_hat by angle_max * sin(π t) radians, where t is the normalised spanwise position within the rib pair.

Non-allocating (uses scalar and SVector arithmetic only).

flip_created_coord_in_pairs_if_needed!(coord::Matrix{Float64})

Ensure coordinates are ordered from positive to negative along y-axis.

update_non_deformed_sections!(wing::AbstractWing)

Create nondeformedsections to match refinedsections. This enables deformation support for all wings (YAML and OBJ). Should be called after refinedsections are populated. Once populated, nondeformedsections serves as the undeformed reference geometry.

calculate_new_aero_data(sections, section_index, left_weight, right_weight)

Interpolate aerodynamic input between two adjacent sections (zero-copy variant).

calculate_new_aero_data(aero_model, aero_data, section_index,
                        left_weight, right_weight)

Interpolate aerodynamic input between two sections.

set_polar!(panel, alphas, cl, cd, cm; shape=nothing)

Rewrite a panel's polar table with values at alphas [rad], ascending, and rebuild its interpolations. The table keeps the shape the panel was built with: a POLAR_VECTORS panel takes the angles as they are, and a POLAR_MATRICES panel takes them at every delta it already spans, since a regenerated polar carries its deflection as shape rather than as a flap angle and so says the same thing at each. The panel's alpha_ref and alpha_window are taken from the angles, so a table covering only a window around one angle is held at its ends rather than extrapolated past them (see window_alpha).

Written in place: the knots and the three value vectors reuse the panel's own alpha_knots and cl_coeffs/cd_coeffs/cm_coeffs storage whenever the sample count is unchanged (polar_column!), and the interpolations read that storage rather than a copy of it, so writing it is the whole update and the interpolation objects stay put (refresh_polar!). Only a table that changes shape is rebuilt. Values may be views into a caller's own buffers; nothing here holds on to them.

shape is the KulfanParameters the values were generated from, stored on the panel as live_shape so a panel's polar and the shape behind it are set together and cannot drift apart. The panel holds the shape itself, not a copy, so a source that rewrites one in place has already updated every panel flying it.

refresh_polar!(old, knots, values) -> Extrapolation

An interpolation reading knots and values themselves, carrying old's extrapolation and knot container so it keeps the type the panel's field was parameterised with. interpolate! takes both arrays by reference, so once an interpolation reads a panel's own storage, writing that storage is the whole update and there is nothing to rebuild. A table that changes shape is rebuilt once, and reads in place from then on. The two dimensional form spreads the angles over every delta the panel spans.

polar_column!(stored, values) -> Vector{Float64}

One column of a rewritten polar table in the panel's own stored vector, which is returned, or in a fresh vector when the sample count has changed — the one case a rewritten table has to grow, and the one that costs the panel's interpolations a rebuild.

polar_model(interp) -> AeroModel

The model a rewritten table leaves the panel on: the one its interpolations were built with, since their shape is fixed by the panel's type.

Whether an interpolation's knot container is a view of knots rather than a copy.

The angles in the container the panel's interpolations were parameterised with.

polar_knots(alphas) -> AbstractVector

The angles of a 1D polar in the container whose knot search fits its length, see ScanKnots.

window_alpha(panel, alpha)

alpha [rad] clamped into the range the panel's polar was built over. A table generated over a window around one angle of attack says nothing past its ends, so it is held at its end values rather than extrapolated out of them; alpha_window of 0 means unbounded and leaves alpha alone, which is what a full-range polar wants.

load_polar_data(csv_file_path::String) -> Tuple{Union{Nothing, Tuple{Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}}}, Symbol}

Load aerodynamic polar data from a CSV file using only readlines.

The CSV file must contain a header row with columns for alpha, cl, cd, and cm (case-insensitive, order arbitrary). Each subsequent row should contain numeric values for these columns.

Arguments

• csv_file_path::String: Path to the CSV file containing polar data.

Returns

• A tuple (aero_data, model_type) where:
  • aero_data: A tuple of vectors (alpha, cl, cd, cm) if the file is valid, or nothing if invalid or missing.
    • alpha: Angle of attack in degrees (converted to radians internally).
    • cl: Lift coefficient.
    • cd: Drag coefficient.
    • cm: Moment coefficient.
  • model_type: POLAR_VECTORS if data is loaded, or INVISCID if not.

Behavior

• If the file is missing, empty, or invalid, a warning is issued and (nothing, INVISCID) is returned.

Example

# Create a YAML-based wing from configuration file
wing = Wing(
    "path/to/wing_config.yaml";
    n_panels=40,
    n_groups=4
)

load_matrix_polar_data(cl_path, cd_path, cm_path) -> (aero_data, POLAR_MATRICES)

Read the three (alpha × delta) coefficient matrices (see read_aero_matrix) and assemble the POLAR_MATRICES aero_data = (alpha, delta, cl, cd, cm).

read_aero_matrix(filepath::AbstractString) -> (Matrix{Float64}, Vector{Float64}, Vector{Float64})

Read an aerodynamic coefficient matrix from CSV with angle labels. Returns the coefficient matrix and corresponding angle ranges.

Returns

• matrix: Matrix of aerodynamic coefficients
• alpha_range: Vector of angle of attack values in radians
• delta_range: Vector of flap deflection angles in radians

read_dat(path) -> (x, y)

Read Selig .dat airfoil coordinates (two whitespace-separated columns), skipping the name/header line and any non-numeric lines.

read_node_table(path) -> (alpha, delta, values)

Read a per-node aero table into radian alpha/delta vectors (one entry per row) and a nrow × n_node value matrix. The file suffix picks the format: .arrow (columns alpha, delta in degrees and a per-row list column values), anything else CSV (header alpha, delta, n0, n1, …, angles in degrees).

write_node_rows(path, alpha, delta, values) -> path

Write a per-node aero table from radian alpha/delta vectors (one entry per row) and a nrow × n_node value matrix. The file suffix picks the format, matching read_node_table: .arrow for the binary form, anything else CSV. Angles are written in degrees either way.

convert_node_table(src, dst) -> dst

Rewrite a per-node aero table in the format dst's suffix names. Use to move an existing dataset between CSV and Arrow without re-running the (slow) airfoil solver that produced it.

delta_suffix(delta) -> String

Filename tag for a non-zero trailing-edge deflection delta [rad], e.g. d5, dm3 (−3°), d2p5 (2.5°). Uses millidegree precision (negatives as m, the decimal point as p) so sub-degree deflections get distinct .dat files instead of colliding. Shared by write_section_aero and read_section_aero.

interpolate_matrix_nans!(matrix::Matrix{Float64}; prn=true)

Replace NaN values in a matrix by interpolating from nearest non-NaN neighbors. Uses an expanding search radius until valid neighbors are found.

Arguments

• matrix: Matrix containing NaN values to be interpolated

remove_vector_nans(aero_data)

Remove the indices from aero_data where a NaN is found.

generate_polar_data(solver, body_aero, angle_range;
    angle_type="angle_of_attack", angle_of_attack=0.0,
    side_slip=0.0, v_a=10.0)

Sweep over angle_range (degrees), solving at each angle. Returns a named tuple (polar_data, cmx, cmy, cmz, rey) where polar_data is [angle, cl, cd, cs, gamma_dist, cl_dist, cd_dist, cs_dist, reynolds].

extract_literature_polar_data(raw_data, path;
    angle_type="angle_of_attack")

Extract polar data from literature CSV data (as returned by readdlm). Returns a named tuple with polar_data (array of [angle, cl, cd, cs] vectors) and cmx, cmy, cmz vectors.

Supports both tuple format (table, header) from readdlm(...; header=true) and matrix format where the first row contains headers.

parse_literature_column(col)

Parse a column of mixed-type data (Real or String) into Float64, returning NaN for unparseable values.

interpolate_section_aero_to_refined!(wing)

Set each refined section's SectionAero by spanwise-interpolating the unrefined sections' surface tables (contour, Cp, cf) using the refined→unrefined mapping (refined_section_left_idx, refined_section_weight). No-op when the unrefined sections carry no surface aero.

validate_section_aero(sections)

Enforce the all-or-none surface-aero rule and a uniform node resolution: either every section in sections carries SectionAero or none does, and all present tables share the same node count. Throws ArgumentError on a violation.

copy_examples()

Copy all example scripts to the folder "examples" (it will be created if it doesn't exist).

deform_kulfan!(out, basis, base, upper_deflection, lower_deflection, residual)
    -> out

deform_kulfan written into the shape out already is, with residual as the scratch chord_residual! fills — the form a live polar refreshes a panel's shape with, since it neither allocates nor replaces the object the panel points at.

chord_residual!(residual, basis, deflection) -> residual

chord_residual written into storage the caller owns.

bernstein_basis(x::AbstractVector, n::Int)

Compute Bernstein polynomial basis matrix.

Returns matrix of shape (length(x), n+1) where each column is a Bernstein polynomial B_i(x) = binomial(n,i) * x^i * (1-x)^(n-i)

class_function(x; N1=0.5, N2=1.0)

CST class function for airfoils: C(x) = x^N1 * (1-x)^N2

Standard values N1=0.5, N2=1.0 give round leading edge and pointed trailing edge.

leading_edge_basis(x::AbstractVector, n_weights::Int)

Leading-edge-modification basis used by AeroSandbox/NeuralFoil: x * (1 - x)^(n_weights + 0.5). Added identically to both surfaces.

compute_optimal_x_points(n) -> Vector{Float64}

NeuralFoil's boundary-layer station x/c: midpoints of a uniform [0,1] grid.

normalize_airfoil(x::Vector, y::Vector)

Normalize airfoil coordinates to unit chord with LE at origin.

Returns normalized (x, y) and transformation parameters.

get_lower_upper(x, y, crease_frac) -> (lower, upper)

Find y-coordinates where upper/lower surfaces intersect the hinge line.

turn_trailing_edge!(angle, x, y, lower_turn, upper_turn, crease_frac; thickness_frac=nothing)

Deflect airfoil trailing edge by rotating coordinates behind the crease line. Positive angle deflects downward.

The pivot's through-thickness location is set by thickness_frac (0 = bottom surface, 0.5 = mid, 1 = top): y_pivot = lower_turn + thickness_frac*(upper_turn - lower_turn). When thickness_frac is given, points are only rotated (no crease cleanup) — intended for a following Kulfan refit that wraps the crease robustly. When thickness_frac === nothing (legacy), the pivot is the top surface for downward deflection and the bottom for upward, and folded-over crease points are removed and the band averaged so the coordinates stay XFoil-paneable.

distance_parabolas!(d, f, n, v, z)

One pass of the Felzenszwalb–Huttenlocher distance transform: writes into d the lower envelope of the parabolas (q - p)^2 + f[p] over p, for q in 1:n. v and z are scratch (length n and n + 1).

squared_distance_transform!(field) -> field

In-place 2D squared Euclidean distance transform in grid-index units. On input field holds 0 (or a local squared offset) at seeds and a large finite value elsewhere; on output every node holds its squared distance to the nearest seed.

grid_sampler(field, x0, y0, cell) -> f(px, py)

Bilinear interpolant of field at physical points; node (i, j) sits at (x0 + (i-1)*cell, y0 + (j-1)*cell). Queries are clamped to the grid.

flood_outside(blocked) -> BitMatrix

Mark the nodes reachable from the grid border without entering blocked (4-connected flood fill); unreached nodes are the solid plus its enclosed holes.

trace_level_set(field, level) -> Vector{Vector{NTuple{2,Float64}}}

Marching-squares contours of field .== level, each returned as a closed loop of grid-frame vertices (unit = one cell, node 1 at 1.0). Saddle cells are resolved by the cell-center average. Loops touching the grid border are dropped.

largest_linking_gap(x, y) -> Float64

Largest edge of the Euclidean minimum spanning tree of the points (Prim, O(n²)) — the longest hop needed to keep the cloud connected. The rolling ball must not fall through it.

resample_arc(ax, ay, n, curvature_weight) -> (x, y)

Resample the polyline (ax, ay) at n stations cosine-clustered in a measure that blends arclength with curvature_weight extra length per radian of turning, so both ends attract points and sharp features are resolved by several panels regardless of their size (like XFoil's curvature-attracted paneling). Endpoints are preserved.

smooth_turning!(turn) -> turn

Diffuse the per-node turning-angle density along the contour (in place). A sharp corner — e.g. a deflected section's hinge, rounded to clearance — otherwise dumps its whole turn into one node, and the curvature-weighted resampling then collapses a few panels to near-zero length, which XFoil's viscous solver cannot handle. Spreading the turn over a short band refines a group of panels gradually instead, the way XFoil's PANGEN bunches panels on a smoothed curvature. The total turn is conserved, so leading-edge clustering is preserved.

load_neuralfoil_model(model_size::String="xlarge"; weights_dir=nothing)

Load NeuralFoil neural network weights from .npz files.

Arguments

• model_size: One of "xxsmall", "xsmall", "small", "medium", "large", "xlarge", "xxlarge", "xxxlarge"
• weights_dir: Directory containing the .npz files (defaults to package data dir)

Returns

• NeuralFoilModel: Loaded model ready for evaluation

neuralfoil_section(params, alpha, Re; kwargs...) -> NamedTuple

Full NeuralFoil evaluation returning integrated coefficients and the surface pressure distribution reconstructed from the predicted edge-velocity ratios (Cp = 1 - (ue/vinf)^2) at NeuralFoil's N fixed station x/c.

Returns (; alpha, cl, cd, cm, confidence, x, cp_upper, cp_lower, ue_upper, ue_lower), with x of length N and the four matrices sized N × n_alpha. A reconstruction should interpolate ue and square afterwards: it is linear in arc length through a stagnation point, where Cp is quadratic.

neuralfoil_fused_output(params, alpha, Re; model_size, weights_dir,
                        n_crit, xtr_upper, xtr_lower) -> Matrix

Forward pass with the symmetry embedding: evaluate the network, evaluate the top/bottom-flipped case, flip its outputs back, and average. Returns the fused output matrix (n_outputs × n_cases), with the Mahalanobis penalty already applied to the confidence logit (row 1).

fused_output(x, model) -> Matrix

Symmetry-fused forward pass over a prepared input matrix, see neuralfoil_fused_output. Takes the inputs already built so a caller that assembles its own batch does not go back through prepare_inputs.

fused_output!(work, x, model) -> AbstractMatrix

fused_output run entirely inside work, allocating nothing for a batch the workspace has room for. The result is a view of the workspace's own storage and stays valid until the next pass through it.

decode_surface_velocity(y) -> (station_x, ue_upper, ue_lower)

The edge-velocity ratios a fused output matrix carries, at NeuralFoil's own N fixed stations. Each ue matrix is N × n_cases. Interpolate these rather than Cp = 1 - ue²: ue is linear in arc length through a stagnation point, where Cp is quadratic.

decode_surface_velocity!(ue_upper, ue_lower, y) -> (ue_upper, ue_lower)

decode_surface_velocity into matrices the caller owns. The station axis is fixed by the network and is carried by NeuralFoilWorkspace instead, so a live pressure refresh reads a pass out without allocating anything for it.

surface_velocity_rows(y) -> (n_stations, upper, lower)

The number of boundary-layer stations a fused output matrix carries and the row ranges its upper and lower edge-velocity ratios sit in. Viewing those rows is how a caller that already holds the output reads the velocities without copying them out (decode_surface_velocity is the copying form).

decode_coefficients(y) -> (cl, cd, cm, confidence)

Turn a fused network output matrix into the integrated coefficients per case, undoing NeuralFoil's output scaling. The single place that scaling is written down.

decode_coefficients!(cl, cd, cm, confidence, y) -> (cl, cd, cm, confidence)

decode_coefficients into vectors the caller owns, so a live polar reads a refresh out of the network without allocating four vectors for it.

nn_forward(x::AbstractMatrix, model::NeuralFoilModel)

Forward pass through the neural network.

Arguments

• x: Input matrix of shape (ninputs, ncases)
• model: Loaded NeuralFoil model

Returns

• Output matrix of shape (noutputs, ncases)

nn_forward!(layers, x, model, n_cases=size(x, 2)) -> AbstractMatrix

Forward pass writing each layer's activations into layers, returning a view of the last one over the n_cases columns used. The matrix products go through BLAS in place, so a pass over a batch the buffers were sized for allocates nothing.

add_bias!(out, bias, activate) -> out

Add a layer's bias to its activations in place, passing them through swish unless this is the output layer, which carries no activation.

layer_buffers(model, n_cases) -> Vector{Matrix{Float32}}

One activation matrix per network layer, n_cases wide: the storage a forward pass writes through. NeuralFoilWorkspace holds two sets of them, one per symmetry.

How many cases a workspace was sized for.

prepare_inputs(params::KulfanParameters, alpha, Re;
               n_crit=9.0, xtr_upper=1.0, xtr_lower=1.0)
prepare_inputs(params::AbstractVector{KulfanParameters}, alpha, Re; kwargs...)

Prepare neural network inputs from Kulfan parameters and flow conditions. The vector form takes one shape per case, so a whole wing's panels go through the network in a single forward pass.

Arguments

• params: Kulfan CST parameters
• alpha: Angle of attack in degrees (scalar or vector)
• Re: Reynolds number (scalar or vector)
• n_crit: Critical amplification factor (default 9.0)
• xtr_upper: Forced transition location on upper surface (0-1)
• xtr_lower: Forced transition location on lower surface (0-1)

Returns

• Input matrix of shape (25, n_cases)

fill_case_input!(x, case, params::KulfanParameters, alpha_deg, Re, n_crit,
                 xtr_upper, xtr_lower)

Write one NeuralFoil input column: rows 1–18 the Kulfan shape, rows 19–25 the flow condition (alpha_deg in degrees). The single point where the network's input layout is defined, shared by every prepare_inputs method.

flip_inputs(x::AbstractMatrix)

Flip inputs for symmetry embedding (swap upper/lower, negate alpha).

flip_inputs!(x_flip, x) -> x_flip

flip_inputs written into storage the caller owns: the upper and lower weights swap with a sign change, the leading-edge weight and sin(2·alpha) negate, and the two transition locations swap. Everything else carries over.

flip_outputs(y::AbstractMatrix)

Flip outputs back after evaluating with flipped inputs.

flipped_row(row, n_stations) -> (source, sign)

Which row of a network output a top/bottom-flipped case reads row from, and with which sign. The single place the mirror symmetry of the output layout is written down, shared by flip_outputs and fuse_flipped!.

CL and CM change sign, the two transition locations swap, and of the six n_stations-long boundary-layer blocks the upper and lower halves swap — the momentum thickness and shape factor as they are, the edge velocity mirrored, so it changes sign.

fuse_flipped!(y, y_flip) -> y

Average a direct output with its flipped counterpart in place, undoing the flip on the way (flipped_row). Reading y_flip while writing y is what lets the fusion land in storage that is already there instead of a third matrix.

squared_mahalanobis_distance(x::AbstractMatrix, model::NeuralFoilModel)

Compute squared Mahalanobis distance from training distribution, one per case.

This is used to penalize predictions far from the training data.

mahalanobis_case(x, case, model, centered) -> Float64

One case's squared Mahalanobis distance, taking centered as the scratch the centred input column is written into rather than allocating one per case.

penalize_confidence!(y, x, model, centered) -> y

Subtract every case's Mahalanobis penalty from the confidence logit, row 1 of the network output. This is what makes a shape far from the training distribution report a low confidence, and it is applied to each symmetry before the two are fused.

swish(x)

SiLU/Swish activation function: x * sigmoid(x) = x / (1 + exp(-x))

sigmoid(x)

Sigmoid activation function.

create_2d_polars(; dat_path, cl_polar_path, cd_polar_path, cm_polar_path, wind_vel,
                 area, width, crease_frac, alpha_range, delta_range,
                 solver=XFoilSolver(), remove_nan=true)

Generate 2D airfoil-section (alpha, delta) POLAR_MATRICES CSVs for an airfoil .dat, deflecting the trailing edge over delta_range and sweeping alpha_range (both radians) with solver. Pass solver=NeuralFoilSolver() to use NeuralFoil instead of XFoil — both emit the same three CSV files, so VSM loads them identically.

The Reynolds number is wind_vel * (area / width) / ν (kinematic viscosity of air). crease_frac is the chordwise hinge location (0–1). Writes lift, drag, and moment coefficient matrices to the three output paths.

lei_poly_coeffs(tube_diameter, camber) -> (cl_coeffs, cd_coeffs, cm_coeffs)

Breukels leading-edge-inflatable α-polynomial coefficients (α in degrees) for a section with normalized tube_diameter and camber. cl_coeffs is a cubic (4 coefficients), cd_coeffs/cm_coeffs are quadratic (3), each in ascending order (constant term first), ready to pass straight to evalpoly. Feed the result to a core POLY section as aero_data = (cl_coeffs, cd_coeffs, cm_coeffs).

resolve_airfoil(type, info, out_dir, id; Re, alpha_range) -> (new_type, new_info)

Resolve one awesIO wing_airfoils entry to a core-loadable form. breukels_regression (t, kappa) → poly coeffs (via lei_poly_coeffs); neuralfoil (dat_file_path, …) → a polars CSV (via generate_polar_from_dat) written under out_dir; polars/poly/inviscid pass through. masure_regression is not yet supported. info file paths should already be absolute.

read_dat_coordinates(path) -> (x, y)

Read airfoil coordinates from a Selig-format .dat file, skipping header and comment lines.

write_dat(filepath, name, x, y) -> filepath

Write airfoil coordinates to a Selig-format .dat file.

write_polar_csv(filepath, sols::Vector{SectionSolution})

Write a solver sweep (from any AbstractAirfoilSolver) to a POLAR_VECTORS CSV (alpha, Cd, Cs, Cl, Cm; alpha in degrees). Non-converged angles (NaN) are skipped, so this works for both NeuralFoil and XFoil sweeps.

write_polar_csv(filepath, result::NeuralFoilResult)

Write a NeuralFoil result to a POLAR_VECTORS CSV (alpha, Cd, Cs, Cl, Cm).

write_polar_matrix_csv(filepath, alpha_range, delta_range, cl, cd, cm)

Write an (alpha × delta) sweep to a long-format POLAR_MATRICES CSV with columns alpha, delta, Cl, Cd, Cm (both angles in degrees), one row per grid point. The delta column is what marks the file as a matrix polar to the loader. alpha_range and delta_range are in radians; cl/cd/cm are length(alpha) × length(delta).

write_aero_matrix(filepath, matrix, alpha_range, delta_range, label) -> filepath

Write an (alpha × delta) coefficient matrix to a labelled CSV: the header row holds the flap deflections (δ=…°), the first column the angles of attack (α=…°), both in degrees. alpha_range/delta_range are radians; label (e.g. "C_l") names the coefficient. The submodule-side writer; read_aero_matrix (main package) reads it back.

write_node_table(path, aero, values) -> path

Write a per-node aero table (Cp or cf, shaped n_node × n_alpha × n_delta), one row per (alpha, delta) with angles in degrees. The file suffix picks the format: .arrow (columns alpha, delta and a per-row list column values, an order of magnitude faster to load), anything else a human-readable CSV with header alpha, delta, n0, n1, … (node columns in the contour node order). read_node_table reads both.

flat_plate_cf(xc, Re) -> Float64

Approximate local skin-friction coefficient at chord fraction xc (0..1) for Reynolds number Re (Re_x = Re·xc), taking the larger of the two standard flat-plate correlations: the Blasius laminar solution cf = 0.664·Re_x^(-1/2) and Prandtl's one-seventh-power-law turbulent estimate cf = 0.027·Re_x^(-1/7) (see e.g. White, Viscous Fluid Flow; Schlichting & Gersten, Boundary-Layer Theory). A stand-in per-node cf for backends that do not expose one (NeuralFoil); XFoil returns the exact distribution via bldump.

neuralfoil_contour_solution(alpha, res, i, x, y, le, arc, cf) -> SectionSolution

Assemble a full-contour SectionSolution for case i of a NeuralFoil sweep res: reconstruct Cp on the contour nodes from the predicted edge velocities (contour_pressure), carrying the precomputed cf.

contour_arc(x, y, le) -> Vector

Signed arc length of every contour node from the leading edge node le: positive along the upper surface toward the trailing edge, negative along the lower. The natural coordinate across a blunt nose, where a small chordwise step is a long one along the skin.

contour_arc!(arc, x, y, le) -> arc

contour_arc written into arc, which is resized to the contour. A live pressure refresh walks panel after panel through one buffer this way.

arc_at_chord(x, arc, indices, fractions) -> Vector

Signed arc length at each chord fraction, read off the contour nodes indices — the one surface the fractions belong to.

arc_at_chord!(out, scratch, x, arc, indices, fractions) -> out

arc_at_chord written into out, taking the sorted copy of the surface from scratch rather than allocating one per call.

sorted_surface!(scratch, x, arc, indices) -> (surface_x, surface_arc)

One surface's nodes as chord fraction and signed arc length, ascending in chord, in the scratch's own storage. Sorted by insertion, which a surface listed from one end to the other is either already in or exactly reversed from, and which needs no scratch of its own.

sort_pairs!(sorted, carried)

Insertion-sort sorted ascending, carrying carried along with it. Stable, in place, and linear on input that is already ordered, which both of its callers hand it.

trailing_edge_speed(ue_upper, ue_lower, upper_arc, lower_arc, upper_te, lower_te)
    -> Float64

One edge speed for both surfaces at the trailing edge: each surface's last two stations extrapolated to its own trailing-edge arc length, the two magnitudes averaged, so a sharp edge carries a single pressure (Kutta).

velocity_knots(station_x, ue_upper, ue_lower, x, arc, le) -> (knots, values)

The whole contour's edge velocity as one curve of (signed arc length, signed ue/u∞), strictly increasing in arc length.

NeuralFoil reports nothing over the first and last 1/2N of chord. Signing the lower surface negative joins both surfaces into one curve through the stagnation point, so the unsampled nose is interpolated between the innermost station on each side rather than extrapolated off the end of one. The stagnation point enters as its own knot, where ue interpolated linearly between those two stations reaches zero — the stagnation-point-flow result, ue being linear in arc length there.

velocity_knots!(scratch, station_x, ue_upper, ue_lower, x, arc, le)
    -> (knots, values)

velocity_knots assembled in scratch. The two vectors returned are the scratch's own and are overwritten by the next call.

contour_pressure(station_x, ue_upper, ue_lower, x, arc, le) -> Vector{Float64}

Cp at every node of a closed contour x with signed arc length arc and nose at le, from one case's edge-velocity ratios at NeuralFoil's station_x. Builds the one-curve edge velocity (velocity_knots), interpolates it with a shape-preserving monotone cubic in arc length, and squares it into Cp = 1 - ue².

contour_pressure!(cp, scratch, station_x, ue_upper, ue_lower, x, arc, le) -> cp

contour_pressure written into cp, with scratch carrying the edge-velocity curve. Only the monotone interpolation itself still allocates, once per contour.

fill_node_nans!(grid, i)

Fill NaNs in node i's (alpha, delta) matrix grid[i, :, :] with interpolate_matrix_nans!; an all-NaN node is left untouched.

live_xfoil_solver(live::LivePolars) -> XFoilSolver

XFoil settings matching the ones the live polars were evaluated at: the same critical amplification factor, and free transition on both surfaces, which is what refresh_live_polars! feeds the network. A comparison run at other settings measures the settings rather than the network.

xfoil_ramp(alpha, step) -> Vector{Float64}

Angles [rad] from zero out to alpha in step increments, alpha itself last. The viscous march refuses a blunt inflated section jumped straight to its angle, and analyze_sweep reinitialises at zero and walks outward, so handing it the whole ramp is what gets the reference solved at all.

read_faces(filename)

Read vertices and faces from an OBJ file.

Arguments

• filename::String: Path to .obj file

Returns

• Tuple of (vertices, faces) where:
  • vertices: Vector of 3D coordinates [x,y,z]
  • faces: Vector of triangle vertex indices

slice_mesh_at_plane(vertices, faces, point, normal; tol=1e-6)

Slice a mesh with the plane through point with unit normal, returning the 3D segments (p1, p2) where the surface crosses it. The plane may have any orientation, so the cut can follow a curved or swept span.

order_segments_to_contour(segments; tol=1e-4)

Order line segments into a continuous contour.

Arguments

• segments: Vector of (p1, p2) tuples representing line segments

Returns

• Vector of [x, z] points forming the ordered contour

station_indices(march, n; wingtip_distance=0.0, min_chord_frac=0.01) -> Vector{Int}

Indices of the march_edges stations nearest n targets spread over the leading-edge arc length. Stations whose chord has closed to less than min_chord_frac of the longest one are left out of that range first, so a wing tapering to a point puts its outermost sections on the last stations that still have an airfoil to slice rather than on the point itself. The remaining first and last targets sit a further wingtip_distance (arc length) inboard.

airfoil_frame(LE_point, TE_point, span_tangent) -> (x_af, y_af, z_af)

Local airfoil axes for a slice, built so the slice plane contains the chord. With x̂ = [1,0,0] and the LE span_tangent: ẑ = x̂ × span (up), ŷ = ẑ × x̂ (spanwise slice normal — the span projected into the y-z plane), x̂_c = ŷ × ẑ (chord, oriented LE → TE). Slice the mesh with y_af. Returns nothing if the span is parallel to global-x.

build_section(vertices, faces, le, te, point, tangent) -> section or nothing

Slice the mesh at one marched station (airfoil_frame drops the chordwise tilt) and project it, producing the full (; LE_point, TE_point, span_dir, contour3d, x_airfoil, y_airfoil), or nothing.

contour_to_airfoil(contour::Vector{Vector{Float64}})

Convert a 2D contour to normalized airfoil coordinates.

Assumes contour is in [x, z] format where x is chordwise, z is thickness direction.

Returns

• (x, y): Normalized airfoil coordinates with chord = 1, LE at origin

plane_contour_to_airfoil(contour3d, LE_point, x_af, z_af)

Project a 3D slice contour into the local airfoil frame (chord axis x_af, up axis z_af) and return normalized Selig coordinates (x, y), or nothing if degenerate.

reorder_airfoil_selig(x::Vector, y::Vector)

Reorder airfoil coordinates to Selig format (TE upper -> LE -> TE lower).

This is the standard format expected by most airfoil tools.

densify_contour(contour, max_edge) -> Vector{Vector{Float64}}

Insert evenly spaced points along contour edges longer than max_edge. Coarse mesh triangles otherwise leave large hops in the slice cloud, which force the shrink wrap's auto-raised rolling ball far up and over-smooth the wrapped airfoil.

create_interpolations(vertices, circle_center_z, radius, gamma_tip)

Create interpolation functions for leading/trailing edges and area.

Arguments

• vertices: Vector of 3D point coordinates
• circle_center_z: Z-coordinate of circle center
• radius: Circle radius
• gamma_tip: Maximum angular extent

Returns

• Tuple of (leinterp, teinterp, area_interp) interpolation functions
• Where leinterp and teinterp are tuples themselves, containing the x, y and z interpolations

find_circle_center_and_radius(vertices)

Find the center and radius of the kite's curvature circle.

Arguments

• vertices: Vector of 3D point coordinates

Returns

• Tuple of (zcenter, radius, gammatip) where:
  • z_center: Z-coordinate of circle center
  • radius: Circle radius
  • gamma_tip: Angle of the kite tip from z-axis

march_edges(vertices, faces; step) -> (; le, te, point, tangent, arclen)

March the leading edge outward from mid-span in both directions in steps of arc length step. Each cut is a vertical spanwise plane (both the chordwise and vertical components of the running LE tangent dropped from the normal) so a tip that curls downward can't tilt the plane toward horizontal, where its min-chord "LE" pick would jump across the wing. Marching stops when the leading edge stops advancing spanwise. Cuts sample mesh edges, so the picks are robust to vertex density. Returns, ordered along the span, the LE/TE points, each cut's plane origin and tangent, and the cumulative LE arc length. Build the airfoil for a chosen station with build_section.

calculate_inertia_tensor(vertices, faces, mass, com)

Calculate the inertia tensor for a triangulated surface mesh, assuming a thin shell with uniform surface density.

Arguments

• vertices: Vector of 3D point coordinates representing mesh vertices
• faces: Vector of triangle indices, each defining a face of the mesh
• mass: Total mass of the shell in kg
• com: Center of mass coordinates [x,y,z]

Method

Uses the thin shell approximation where:

1. Mass is distributed uniformly over the surface area
2. Each triangle contributes to the inertia based on its area and position
3. For each triangle vertex p, contribution to diagonal terms is: area * (sum(p²) - p_i²)
4. For off-diagonal terms: area * (-p_i * p_j)
5. Final tensor is scaled by mass/(3*total_area) to get correct units

Returns

• 3×3 matrix representing the inertia tensor in kg⋅m²

center_to_com!(vertices, faces)

Calculate center of mass of a mesh and translate vertices so that COM is at origin.

Arguments

• vertices: Vector of 3D point coordinates
• faces: Vector of vertex indices for each face (can be triangular or non-triangular)

Returns

• Vector representing the original center of mass before translation

Notes

• Non-triangular faces are automatically triangulated into triangles
• Assumes uniform surface density

airfoils_from_yaml(geometry_file) -> Vector of (; id, x, y, x_raw, y_raw)

Read each airfoil's coordinates from the .dat files referenced by a YAML geometry's wing_airfoils. x, y come from dat_file (the fitted airfoil); x_raw, y_raw come from raw_dat_file if present (the raw sliced points), and are empty otherwise. Relative paths resolve against the geometry file's directory.

write_geometry_yaml(path, section_rows, airfoil_rows)

Write a geometry YAML via write_yaml, one line per section/airfoil row. section_rows are [airfoil_id, LE_x, LE_y, LE_z, TE_x, TE_y, TE_z]; airfoil_rows are [airfoil_id, type, info_dict] where info_dict holds dat_file, csv_file_path, and optionally raw_dat_file, cp_file, cf_file.

resolve_aero_geometry(yaml_in, out_dir; verbose=true) -> yaml_out

Read an awesIO-style geometry YAML and resolve every wing_airfoils entry to a core-loadable form via resolve_airfoil (breukels_regression → poly, neuralfoil → polars CSV, others pass through), writing generated CSVs and a resolved geometry.yaml under out_dir. wing_sections (incl. any VUP up-vectors) pass through unchanged. Load the result with Wing(yaml_out).

table_path_prefix(geometry_path, output_dir) -> String

Path from the geometry YAML's directory to the table directory, empty when they are the same. Table references resolve against the YAML's own directory, so a YAML written outside output_dir has to carry this hop.

prefix_table_paths!(airfoil_rows, prefix) -> airfoil_rows

Prepend prefix to every relative table reference in the info_dict of each row, so a geometry YAML written outside the table directory still resolves them. A no-op on an empty prefix.

plot_airfoil_fit(x, y; kwargs...)

Plot a single airfoil with its Kulfan CST fit. Implemented in ObjAdapterMakieExt; load Makie (or GLMakie) to use it.

migrate_node_tables(yaml_path, output_dir, table_format; verbose=true)

Rewrite a generated dataset's per-node Cp/cf tables in table_format and point geometry.yaml at them, when they are not in that format already. This is what lets an existing directory change format without re-running the airfoil solver that produced it — the polars are the slow part and they are untouched. The source tables are left in place.

output_dir is what the YAML's relative table references resolve against, which is its own directory — the rule the geometry loader follows.

display_named(fig, name) -> fig

Show fig in the window registered under name, opening a new titled window for a name that has none yet: plotting the same title again replaces its predecessor, while a new title gets its own window. Backends without titled screens, such as CairoMakie, fall back to a plain display.

span_axis(position, title, ylabel) -> Axis

Axis for a spanwise distribution, +y on the left, matching the kite seen from the front and the +y to -y order its sections are stored in.

create_geometry_plot_makie(body_aero::BodyAerodynamics, title,
                           view_elevation, view_azimuth; zoom=1.8)

Create a 3D Makie plot of wing geometry including panels and filaments.

Arguments

• body_aero: struct of type BodyAerodynamics
• title: plot title
• view_elevation: initial view elevation angle [°]
• view_azimuth: initial view azimuth angle [°]

Keyword arguments

• zoom: zoom factor (default: 1.8)

plot_line_segment_makie!(ax, segment, color, label; width=3)

Plot a line segment in 3D with arrow using Makie.

Arguments

• ax: Makie Axis3
• segment: Array of two points defining the segment
• color: Color of the segment
• label: Label for the legend

Keyword Arguments

• width: Line width (default: 3)

set_axes_equal_makie!(ax, panels; zoom=1.8)

Set 3D Makie axis to equal scale based on panel data.

Arguments

• ax: Makie Axis3
• panels: Array of panels
• zoom: zoom factor (default: 1.8)

map_airfoil_3d(le, te, tangent, x, y) -> 3×N or nothing

Map normalized airfoil coordinates into 3D through the airfoil frame of the station given by its LE/TE points and leading-edge tangent (chord scale = |TE - LE|).

fitted_airfoil_3d(section, wrap_method; delta=0.0, crease_frac=0.75) -> 3×N or nothing

Shrink-wrap a section's sliced contour with wrap_method and map it back into 3D through the section's local airfoil frame, for overlaying on the 3D slice diagnostic. A nonzero delta (degrees) deflects the trailing edge and re-wraps (deform_section), showing the geometry the solvers consume. Returns nothing for a degenerate slice.

generated_slices(out_dir, delta, fit_pts) -> (slices, le, te)

Read the stations of a generated obj_to_yaml output directory and their written .dat airfoils — raw slice, wrap, and the delta-degree deformed wrap when it was generated — assembled for plot_slices_3d. Nothing is re-sliced or re-wrapped; only the Kulfan fits of the stored coordinates are recomputed (via fit_pts), exactly as the polar pipeline fits them.

airfoil_skin_geometry(body; R_b_w=nothing, T_b_w=nothing) -> (vertices, faces, ribs)

Lofted airfoil skin of a BodyAerodynamics: each section's contour is fitted between the panel's corner_points by a 2D similarity, TE pinned so a deflection bulges the fore body up.

The contour is the panel's own live_shape when it has one — the very KulfanParameters object the live polar source deformed and handed to the airfoil solver, not a re-derivation of it, so a deformation bug shows up in the picture instead of being papered over. Otherwise it is section_surface at the panel's delta, which reflects delta only when the geometry carries per-delta slices (obj_to_yaml with a delta_range); with δ=0-only data it renders undeflected, a deliberate cue that the deflected slices are missing. Transformed to world by R_b_w/T_b_w. vertices/faces triangulate the skin between consecutive equal-node sections; ribs is one closed contour polyline per section. Sections without contour data are skipped.

panel_contour(panel, section) -> (x, y)

The airfoil contour to draw a panel with: coordinates of its live_shape when a live polar source has put one there, else the section's tabulated contour at the panel's delta. The live branch reads the stored shape object itself, so what is drawn is what was flown.

panel_normal(panel) -> Point3f

Body-frame airfoil-upper-surface unit normal of a panel, from its corner_points and sign-aligned to z_airf (both body-frame, so the sign relation survives the kite's attitude changes). Falls back to +z for a degenerate quad.

plate_hinge_local(crease_frac, delta) -> (lx, ly)

Chordwise/thickness fractions (chord = 1, LE at 0) of the flap hinge for a plate whose trailing edge is deflected by delta [rad] about crease_frac and then re-pinned so the TE stays at the panel TE corner. A downward deflection therefore lifts the hinge (ly > 0, the "bulge up"). crease_frac outside (0, 1) disables the kink (hinge stays on the chord line).

panel_plate_geometry(panel; R_b_w=nothing, T_b_w=nothing) -> Vector{Point3f}

The 6 vertices [LE_1, hinge_1, TE_1, TE_2, hinge_2, LE_2] of a panel's flat-plate skin, kinked at plate_hinge_local by the panel's delta/crease_frac (TE pinned to the corners, hinge bulging up). Triangulated by PLATE_FACES; with delta == 0 the two quads are coplanar (the original flat quad). Transformed to world by R_b_w/T_b_w.

PLATE_FACES

The 4 triangles that mesh the 6 panel_plate_geometry vertices of a panel's flat-plate skin, two per side of the hinge.

plot!(ax, panel::VortexStepMethod.Panel; use_observables=false, kwargs...)

Plot a single Panel as a flat-plate mesh, kinked at the flap hinge by the panel's delta/crease_frac (see panel_plate_geometry); with delta == 0 this is the flat quad LE1-TE1-TE2-LE2.

If use_observables=true, creates observables for dynamic updates.

plot!(ax, body::VortexStepMethod.BodyAerodynamics; use_observables=false,
      airfoils=false, kwargs...)

Plot a BodyAerodynamics object. By default draws each panel as a flat quad; with airfoils=true instead draws the lofted airfoil skin (one see-through wing-shaped mesh with a contour rib line per section, see airfoil_skin_geometry).

If use_observables=true, creates observables for dynamic updates keyed by (bodyid, panelindex). Otherwise, creates static plots (original behavior).

plot!(body::VortexStepMethod.BodyAerodynamics; R_b_w=nothing, T_b_w=nothing)

Update existing body aerodynamics plot observables with current geometry. This updates all panels in the body using their current corner_points.

Requires that plot(body; use_observables=true) or plot!(ax, body; use_observables=true) was called first to create the observables.