Skip to content

Truncate the SVD in the least-squares Kulfan fit so crowded stations give bounded weights - #342

Open
1-Bort-1 wants to merge 2 commits into
mainfrom
agent/295-the-kulfan-least-squares-fit-has-no-rank
Open

1-Bort-1 wants to merge 2 commits into
mainfrom
agent/295-the-kulfan-least-squares-fit-has-no-rank

Conversation

@1-Bort-1

@1-Bort-1 1-Bort-1 commented Sep 16, 2026

Copy link
Copy Markdown
Contributor

TL;DR

fit_kulfan_parameters(x, y, LeastSquaresFit()) now solves with pinv(A; rtol=1e-4) instead of A \ y, dropping singular values below 1e-4 of the largest. A contour whose stations crowd into a narrow band of the chord made the solve return weights that kulfan_to_coordinates blew up to 1e4 scale, and nothing downstream checks magnitude.

What was wrong

When most stations on one surface sit inside a narrow band of x, the upper Bernstein columns become numerically dependent (cond 1e12–1e15). The QR solve then fits the non-representable part of the contour with huge, cancelling weights, which only cancel on the stations they were fitted against. I reproduced this without the SK100 mesh: 118 upper stations in [0.30, 0.31] on a known CST shape with a 0.005 ripple on top gives a resampled max|y| of 13946 against a true 0.077. That is the same size as the 13013 in the issue.

Why truncation and not the KulfanBasis ridge

The obvious fix was the Tikhonov ridge KulfanBasis uses. It does not work here because "Fit matches aerosandbox get_kulfan_parameters" pins the fitted weights to AeroSandbox's at 1e-9. A ridge λ changes each weight by about λ·cond². At a normal cond of about 1e3, that means either λ ≤ 1e-16, which does nothing for a crowded fit, or every fit drifts off AeroSandbox. A hard SVD cutoff returns exactly the same answer whenever cond(A) < 1e4, and the minimum-norm answer when it is higher.

Measured by the smallest kept singular value relative to the largest (cond of A):

contour cond \ max |y| resampled pinv(rtol=1e-4) change in weights vs \
test_airfoil.dat 916 0.0964 0.0964 6.8e-15
spread stations, ripple 1120 0.0768 0.0768 6.7e-15
118 stations in 0.01 band, ripple 1.2e15 13900 0.119
118 stations in 0.03 band, ripple 1.5e12 629000 1.12

The cutoff allows a cond of 1e4, about 9× the worst normal cond I found (1120). At 1e-3 it would already start cutting test_airfoil.dat, whose cond is 916.

Where I would push back

  • Truncation bounds the weights, but it cannot recover shape the stations never sampled. The 0.03-band row still comes out 15× too thick, only no longer 1e4×. The real cause is the resampler: smoothed_curvature has no clamp on its density, and resample_arc crowds the stations. That fix, and the thickness check in write_section_aero, are options 2 and 3 of the issue. Neither is in this branch. With both, a degenerate section would error instead of being written.
  • The fit gives no warning when it truncates. A @warn would bring crowding to the surface at the point it happens. I left it out to keep this to one idea.
  • KulfanBasis keeps its ridge. Its job is to bound a deflection projection, not to match a reference fit, so the two remedies do different things. I searched src/ for ridge|svdvals|pinv|Tikhonov and found no third regularised solve.
  • Cleanup around the change: the trailing-edge refit reuses A[:, 1:end-1] instead of rebuilding the hcat. The result is the same.

Verification

  • Reproduced first: the new testset on unchanged src/, Evaluated: 13946.420825128731 < 0.1631808615770601
  • test/airfoil_aero/test_airfoil_aero.jl "Fit to stations crowded into a narrow band stays airfoil-sized": red before, green after (juliaserver)
  • Full test/airfoil_aero/test_airfoil_aero.jl passes, including the AeroSandbox reference and round-trip testsets. test/airfoil_aero/test_live_polar.jl and test/solver/test_backend_comparison.jl pass too (juliaserver, 0 fail / 0 error)
  • Local full suite (agent ci-local, Julia 1.12.7): PASS (15 min)
  • GitHub CI: Documentation, macOS 1.12, ubuntu 1.12 (on re-run) and ubuntu 1.13 pass. windows 1.12 fails on one assertion, test/solver/test_forwarddiff.jl:87 (POLAR_MATRICES), on both the first run and the re-run. That is test_forwarddiff.jl's POLAR_MATRICES check flakes at ~4%, and the oper #287, not this branch. Both windows runs give 0.04206403446037785 with norm_fd = 3.2851809771070055, the exact value test_forwarddiff.jl's POLAR_MATRICES check flakes at ~4%, and the oper #287 records for windows when the runner builds polar table A without this change. The first ubuntu 1.12 run showed the matching ubuntu value, 0.0420681631196805.
  • REUSE: not used in this repo. Up to date with main.
  • Benchmark: n/a
  • Risk: a real contour with cond between 1e3 and 1e4 would now fit slightly differently from AeroSandbox. Of the contours I checked, none came near that range.

Scope

+22 / −4 across 3 files: a two-line solve change plus a constant in src/airfoil_aero/kulfan.jl, one 12-line testset, and a changelog entry.

Closes #295 · task VortexStepMethod.jl-295

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 <noreply@anthropic.com>
@1-Bort-1 1-Bort-1 added agent:running Agent task state agent:ci Agent task state and removed agent:running Agent task state labels Sep 16, 2026

@1-Bort-1 1-Bort-1 left a comment

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Independent review (advisory)

Verdict: APPROVE WITH COMMENTS · 2 inline, 0 off the diff

Good

  • The diff matches the card: git diff --stat HEAD~1 shows +22/−4 across CHANGELOG, kulfan.jl and one testset, with nothing unplanned
  • Switching A \ y_norm to pinv(A; rtol) changes no result where cond(A) < 1e4, since both return the least-squares solution; the AeroSandbox 1e-9 testset still covers this path
  • The refit A[:, 1:end-1] is the same matrix as the old hcat without te_col (the last column), so the named cleanup is behaviour-preserving and saves a second allocation
  • LinearAlgebra is already loaded in AirfoilAero.jl, and the test file already imports class_function and bernstein_basis, so no new imports are needed
  • No AD path goes through the fit: grepping src/airfoil_aero finds no ForwardDiff or Dual, so the SVD in pinv cannot break a Dual-number caller
  • Every caller (geometry_gen, live_polar, section_aero_gen, common.jl) runs at build or generation time, not in the solve loop, so the SVD's extra cost does not matter
  • The card explains why a ridge was rejected (it would break the AeroSandbox pin), and a grep for a third regularised solve found none, so there is no duplicate path
  • The new testset is named for the behaviour it protects, and the card shows it failing before the fix (13946 < 0.163)

Not good

  • CHANGELOG.md:18 — The entry says a crowded contour now fits to an airfoil-sized shape, but the card's own table shows the 0.03-band case still comes out 15× too thick. The entry should say the weights are bounded, not that the shape is correct, or users will trust sections that are still wrong.
  • src/airfoil_aero/kulfan.jl:150 — Truncation happens silently, so a degenerate contour now gives a plausible-looking but wrong section instead of an obviously broken one, and nothing downstream checks it. Without the resampler or thickness fixes in this branch or a tracked issue for them, the failure gets harder to spot.
  • The cond survey covers test_airfoil.dat plus synthetic contours, but the real inputs are shrink-wrapped kite sections (deform_section, panel_kulfan_parameters); the TUDELFT_V3 and ram-air sections should be checked against the 1e4 cutoff
  • A hard cutoff is discontinuous: a section whose cond crosses 1e4 during a deform_section delta sweep would jump in weights between neighbouring deflections, and nothing warns
  • Closes #295 shuts the issue while options 2 and 3 (clamping resampler density, a thickness check in write_section_aero) stay undone and are not linked to a follow-up issue
  • The docstring and CHANGELOG repeat the literal 1e-4 instead of naming KULFAN_FIT_RTOL, so they go stale if the constant changes
  • The testset defines a local surface(weights, xs) closure with logic in it; §6 prefers a hoisted helper, though this one is a one-liner
  • The local CI mirror result was not read and GitHub CI is still pending, so the only green runs are on juliaserver

claude, rubric CLEAN_CODE.md. A different lab from the implementer
on purpose: a reviewer sharing its blind spots would not flag its mistakes.

Comment thread CHANGELOG.md
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.

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

MINOR: The entry says a crowded contour now fits to an airfoil-sized shape, but the card's own table shows the 0.03-band case still comes out 15× too thick. The entry should say the weights are bounded, not that the shape is correct, or users will trust sections that are still wrong.


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

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

MINOR: Truncation happens silently, so a degenerate contour now gives a plausible-looking but wrong section instead of an obviously broken one, and nothing downstream checks it. Without the resampler or thickness fixes in this branch or a tracked issue for them, the failure gets harder to spot.

@1-Bort-1

Copy link
Copy Markdown
Contributor Author

CI: Julia 1.12 - windows-latest - x64 - pull_request failed on 245ab9620.

https://github.com/OpenSourceAWE/VortexStepMethod.jl/actions/runs/35158113156/job/105002277668

@1-Bort-1 1-Bort-1 added agent:queued Agent task state agent:running Agent task state and removed agent:ci Agent task state agent:queued Agent task state labels Sep 16, 2026
@codecov

codecov Bot commented Sep 16, 2026

Copy link
Copy Markdown

Codecov Report

✅ All modified and coverable lines are covered by tests.

📢 Thoughts on this report? Let us know!

@1-Bort-1 1-Bort-1 added agent:queued Agent task state agent:running Agent task state agent:ci Agent task state and removed agent:running Agent task state agent:queued Agent task state labels Sep 16, 2026
@1-Bort-1

Copy link
Copy Markdown
Contributor Author

Local full suite: PASS (7 min, Julia 1.12.7, one cell of the matrix)

@1-Bart-1 1-Bart-1 left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Add the warning

@1-Bort-1 1-Bort-1 added agent:review Agent task state and removed agent:ci Agent task state labels Sep 17, 2026
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

agent:review Agent task state

Projects

None yet

Development

Successfully merging this pull request may close these issues.

The Kulfan least-squares fit has no rank guard, so a crowded resample writes an airfoil at 1e4 scale with a normal residual

2 participants