From 245ab962004b02b826782f1400ae57c2c247afb0 Mon Sep 17 00:00:00 2001 From: 1-Bort-1 <323661610+1-Bort-1@users.noreply.github.com> Date: Thu, 17 Sep 2026 00:30:47 +0200 Subject: [PATCH 1/2] Truncate the SVD in the least-squares Kulfan fit Stations crowded into a narrow band of the chord make the upper Bernstein columns numerically dependent, and A \ y returned weights that resample the contour at 1e4 scale. Singular values below 1e-4 of the largest are now dropped; fits of well-spread stations (cond ~1e3) are unchanged to 1e-14. Refs #295 Co-Authored-By: Claude Opus 5 --- CHANGELOG.md | 4 ++++ src/airfoil_aero/kulfan.jl | 10 ++++++---- test/airfoil_aero/test_airfoil_aero.jl | 12 ++++++++++++ 3 files changed, 22 insertions(+), 4 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 4c00c754..7669769f 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -13,6 +13,10 @@ - The `VSMSolution` docstring gives `lift_dist`, `drag_dist` and `panel_moment_dist` in the per-unit-span units they hold, [N/m] and [Nm/m], instead of [N] and [Nm]. +- `fit_kulfan_parameters` with `LeastSquaresFit` drops singular values below `1e-4` + times the largest, so a contour whose stations crowd into a narrow band of the chord + fits to an airfoil-sized shape instead of weights that resample it to 1e4 scale. + Fits of well-spread stations are unchanged. ## VortexStepMethod v5.1.1 2026-09-12 diff --git a/src/airfoil_aero/kulfan.jl b/src/airfoil_aero/kulfan.jl index dce90044..c7488b0b 100644 --- a/src/airfoil_aero/kulfan.jl +++ b/src/airfoil_aero/kulfan.jl @@ -106,6 +106,8 @@ function leading_edge_basis(x::AbstractVector{T}, n_weights::Int) where T return x .* max.(1 .- x, zero(T)).^(n_weights + 0.5) end +const KULFAN_FIT_RTOL = 1e-4 # [-] singular values dropped below this times the largest + """ fit_kulfan_parameters(x::Vector, y::Vector, method::KulfanFitMethod) fit_kulfan_parameters(x::Vector, y::Vector; n_weights=8) @@ -128,7 +130,8 @@ end Least-squares fit matching AeroSandbox's `get_kulfan_parameters`: both surfaces share a single least-squares system with a shared leading-edge weight and a -trailing-edge thickness. +trailing-edge thickness. Singular values below `1e-4` times the largest are dropped, so +stations crowded into part of the chord give bounded weights. """ function fit_kulfan_parameters(x::Vector{T}, y::Vector{T}, method::LeastSquaresFit) where T @@ -144,12 +147,11 @@ function fit_kulfan_parameters(x::Vector{T}, y::Vector{T}, te_col = ifelse.(is_upper, xv ./ 2, .-xv ./ 2) A = hcat((.!is_upper) .* CS, is_upper .* CS, le_col, te_col) - coeffs = A \ y_norm + coeffs = pinv(A; rtol=KULFAN_FIT_RTOL) * y_norm TE_thickness = coeffs[end] if TE_thickness < 0 - A = hcat((.!is_upper) .* CS, is_upper .* CS, le_col) - coeffs = A \ y_norm + coeffs = pinv(A[:, 1:end-1]; rtol=KULFAN_FIT_RTOL) * y_norm TE_thickness = zero(T) end diff --git a/test/airfoil_aero/test_airfoil_aero.jl b/test/airfoil_aero/test_airfoil_aero.jl index ac4b9363..e35a18a5 100644 --- a/test/airfoil_aero/test_airfoil_aero.jl +++ b/test/airfoil_aero/test_airfoil_aero.jl @@ -43,6 +43,18 @@ end @test params.TE_thickness ≈ 0.0 atol = 1e-12 end + @testset "Fit to stations crowded into a narrow band stays airfoil-sized" begin + spread = (1 .- cos.(range(0, pi, 121))) ./ 2 + crowded = vcat(0.0, range(0.30, 0.31, 118), 1.0) + surface(weights, xs) = class_function(xs) .* (bernstein_basis(xs, 7) * weights) + y_upper = surface(fill(0.2, 8), crowded) .+ 0.005 .* sin.(40pi .* crowded) + y_lower = surface(fill(-0.1, 8), spread) + x = vcat(reverse(crowded), spread[2:end]) + y = vcat(reverse(y_upper), y_lower[2:end]) + _, y_fit = kulfan_to_coordinates(fit_kulfan_parameters(x, y)) + @test maximum(abs, y_fit) < 2 * maximum(abs, y) + end + @testset "Shrink-wrap encloses points with clearance" begin xn, yn, _ = normalize_airfoil(collect(float.(xr)), collect(float.(yr))) cloud_to_wrap(xw, yw) = minimum( From b725371e008293fb13aecf988b75ac2442e833cf Mon Sep 17 00:00:00 2001 From: 1-Bort-1 <323661610+1-Bort-1@users.noreply.github.com> Date: Fri, 18 Sep 2026 15:15:54 +0200 Subject: [PATCH 2/2] Warn when the least-squares Kulfan fit drops singular values The truncated solve moves into truncated_least_squares, which returns how many singular values it dropped; fit_kulfan_parameters warns once per fit when any were. Co-Authored-By: Claude Opus 5 --- CHANGELOG.md | 6 +++--- docs/src/private_functions.md | 1 + src/airfoil_aero/kulfan.jl | 22 +++++++++++++++++++--- test/airfoil_aero/test_airfoil_aero.jl | 7 ++++++- 4 files changed, 29 insertions(+), 7 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index caae262a..13e1e426 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -23,9 +23,9 @@ - The `VSMSolution` docstring gives `lift_dist`, `drag_dist` and `panel_moment_dist` in the per-unit-span units they hold, [N/m] and [Nm/m], instead of [N] and [Nm]. - `fit_kulfan_parameters` with `LeastSquaresFit` drops singular values below `1e-4` - times the largest, so a contour whose stations crowd into a narrow band of the chord - fits to an airfoil-sized shape instead of weights that resample it to 1e4 scale. - Fits of well-spread stations are unchanged. + times the largest and warns when it does, so a contour whose stations crowd into a + narrow band of the chord gets bounded weights instead of ones that resample it to 1e4 + scale. Fits of well-spread stations are unchanged. - Inside its vortex core, `velocity_3D_trailing_vortex!` induces an azimuthal velocity instead of a radial one. Only points within the millimetre-scale Oseen core of a panel's chordwise trailing segment were affected. diff --git a/docs/src/private_functions.md b/docs/src/private_functions.md index 5f7eec56..5f954cb1 100644 --- a/docs/src/private_functions.md +++ b/docs/src/private_functions.md @@ -150,6 +150,7 @@ chord_residual! bernstein_basis class_function leading_edge_basis +truncated_least_squares compute_optimal_x_points normalize_airfoil get_lower_upper diff --git a/src/airfoil_aero/kulfan.jl b/src/airfoil_aero/kulfan.jl index c7488b0b..7d7f2c30 100644 --- a/src/airfoil_aero/kulfan.jl +++ b/src/airfoil_aero/kulfan.jl @@ -108,6 +108,19 @@ end const KULFAN_FIT_RTOL = 1e-4 # [-] singular values dropped below this times the largest +""" + truncated_least_squares(A, b) + +Minimum-norm least-squares solution of `A * coeffs = b` with the singular values of `A` +below `KULFAN_FIT_RTOL` times the largest dropped. Returns `(coeffs, n_dropped)`. +""" +function truncated_least_squares(A::AbstractMatrix, b::AbstractVector) + F = svd(A) + kept = F.S .> KULFAN_FIT_RTOL * F.S[1] + coeffs = F.V[:, kept] * ((F.U[:, kept]' * b) ./ F.S[kept]) + return coeffs, count(!, kept) +end + """ fit_kulfan_parameters(x::Vector, y::Vector, method::KulfanFitMethod) fit_kulfan_parameters(x::Vector, y::Vector; n_weights=8) @@ -131,7 +144,7 @@ end Least-squares fit matching AeroSandbox's `get_kulfan_parameters`: both surfaces share a single least-squares system with a shared leading-edge weight and a trailing-edge thickness. Singular values below `1e-4` times the largest are dropped, so -stations crowded into part of the chord give bounded weights. +stations crowded into part of the chord give bounded weights, and a warning says so. """ function fit_kulfan_parameters(x::Vector{T}, y::Vector{T}, method::LeastSquaresFit) where T @@ -147,13 +160,16 @@ function fit_kulfan_parameters(x::Vector{T}, y::Vector{T}, te_col = ifelse.(is_upper, xv ./ 2, .-xv ./ 2) A = hcat((.!is_upper) .* CS, is_upper .* CS, le_col, te_col) - coeffs = pinv(A; rtol=KULFAN_FIT_RTOL) * y_norm + coeffs, n_dropped = truncated_least_squares(A, y_norm) TE_thickness = coeffs[end] if TE_thickness < 0 - coeffs = pinv(A[:, 1:end-1]; rtol=KULFAN_FIT_RTOL) * y_norm + coeffs, n_dropped = truncated_least_squares(A[:, 1:end-1], y_norm) TE_thickness = zero(T) end + n_dropped > 0 && @warn "Kulfan fit dropped $n_dropped of $(length(coeffs)) singular \ + values: the stations leave part of the shape unconstrained, as when they crowd \ + into a narrow band of the chord." lower_weights = coeffs[1:n_weights] upper_weights = coeffs[n_weights+1:2n_weights] diff --git a/test/airfoil_aero/test_airfoil_aero.jl b/test/airfoil_aero/test_airfoil_aero.jl index cb6da398..65329629 100644 --- a/test/airfoil_aero/test_airfoil_aero.jl +++ b/test/airfoil_aero/test_airfoil_aero.jl @@ -51,8 +51,13 @@ end y_lower = surface(fill(-0.1, 8), spread) x = vcat(reverse(crowded), spread[2:end]) y = vcat(reverse(y_upper), y_lower[2:end]) - _, y_fit = kulfan_to_coordinates(fit_kulfan_parameters(x, y)) + crowded_fit = @test_logs (:warn, r"dropped \d+ of \d+ singular values") ( + fit_kulfan_parameters(x, y)) + _, y_fit = kulfan_to_coordinates(crowded_fit) @test maximum(abs, y_fit) < 2 * maximum(abs, y) + + y_spread = vcat(reverse(surface(fill(0.2, 8), spread)), y_lower[2:end]) + @test_logs fit_kulfan_parameters(vcat(reverse(spread), spread[2:end]), y_spread) end @testset "Shrink-wrap encloses points with clearance" begin