diff --git a/CHANGELOG.md b/CHANGELOG.md index 5c872947..7f004675 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -15,6 +15,12 @@ fit: on the SK100's 45 mesh sections nine came out with `.dat` coordinates at 1e2 to 1e4 instead of 0..1, and any consumer mapping the mesh onto a structure rejected the geometry. +- A deflection the 2D solver returned no contour for is no longer written as an + all-`NaN` airfoil `.dat`: `write_section_aero` skips that column, so `isfile` stays + the honest test for a generated deflection. `read_section_aero` warns and returns + `nothing` for a contour that does not hold its `Cp` table's nodes, where it used to + read the `NaN`s back and interpolate a `NaN` airfoil shape at *every* deflection, + 0° included. ## VortexStepMethod v5.1.0 2026-09-11 diff --git a/docs/src/private_functions.md b/docs/src/private_functions.md index cd8f3db7..d7807c42 100644 --- a/docs/src/private_functions.md +++ b/docs/src/private_functions.md @@ -116,7 +116,7 @@ window_alpha assemble_polar_matrix load_matrix_polar_data read_aero_matrix -read_dat +read_dat_coordinates read_node_table write_node_rows convert_node_table @@ -203,7 +203,6 @@ sigmoid create_2d_polars lei_poly_coeffs resolve_airfoil -read_dat_coordinates write_dat write_polar_csv write_polar_matrix_csv diff --git a/src/airfoil_aero/AirfoilAero.jl b/src/airfoil_aero/AirfoilAero.jl index 16f10813..2e7eea0f 100644 --- a/src/airfoil_aero/AirfoilAero.jl +++ b/src/airfoil_aero/AirfoilAero.jl @@ -9,8 +9,8 @@ using Xfoil using Printf: @sprintf using ..VortexStepMethod: SectionAero, interpolate_matrix_nans!, delta_suffix, write_node_rows, section_surface, set_polar!, - KulfanParameters, calculate_cl, calculate_cd, - calculate_cm + read_dat_coordinates, KulfanParameters, calculate_cl, + calculate_cd, calculate_cm include("kulfan.jl") include("deform.jl") diff --git a/src/airfoil_aero/airfoil_io.jl b/src/airfoil_aero/airfoil_io.jl index 0918a2b7..10cefff5 100644 --- a/src/airfoil_aero/airfoil_io.jl +++ b/src/airfoil_aero/airfoil_io.jl @@ -1,26 +1,3 @@ -""" - read_dat_coordinates(path) -> (x, y) - -Read airfoil coordinates from a Selig-format `.dat` file, skipping header and -comment lines. -""" -function read_dat_coordinates(path::String) - x = Float64[] - y = Float64[] - for line in eachline(path) - s = strip(line) - (isempty(s) || !(isdigit(s[1]) || s[1] == '-' || s[1] == '.')) && continue - parts = split(s) - length(parts) >= 2 || continue - xp = tryparse(Float64, parts[1]) - yp = tryparse(Float64, parts[2]) - (xp === nothing || yp === nothing) && continue - push!(x, xp) - push!(y, yp) - end - return x, y -end - """ write_dat(filepath, name, x, y) -> filepath diff --git a/src/airfoil_aero/section_aero_gen.jl b/src/airfoil_aero/section_aero_gen.jl index d8dfb161..9c515f88 100644 --- a/src/airfoil_aero/section_aero_gen.jl +++ b/src/airfoil_aero/section_aero_gen.jl @@ -119,7 +119,9 @@ end Write a [`SectionAero`](@ref) as human-readable files. The airfoil contours share `dat_prefix`: `{dat_prefix}.dat` (contour at `delta=0`) plus -`{dat_prefix}_{delta_suffix(δ)}.dat` per non-zero deflection. The per-node `Cp`/`cf` +`{dat_prefix}_{delta_suffix(δ)}.dat` per non-zero deflection. A deflection whose +contour is all `NaN` is skipped, and `{dat_prefix}.dat` falls back to the first +deflection that has one. The per-node `Cp`/`cf` tables share `table_prefix` (defaults to `dat_prefix`): `{table_prefix}_cp.{ext}` and `{table_prefix}_cf.{ext}`, where `ext` is `table_format`, either `:csv` (default, readable) or `:arrow` (binary, an order of magnitude faster to load). Pass a separate @@ -134,15 +136,15 @@ function write_section_aero(dat_prefix::AbstractString, aero::SectionAero; throw(ArgumentError("table_format must be :csv or :arrow, got :$table_format")) mkpath(dirname(dat_prefix)) mkpath(dirname(table_prefix)) + finite(jd) = any(isfinite, view(aero.x, :, jd)) for (jd, d) in enumerate(aero.delta_range) - iszero(d) && continue + (iszero(d) || !finite(jd)) && continue write_dat("$(dat_prefix)_$(delta_suffix(d)).dat", "section", aero.x[:, jd], aero.y[:, jd]) end - # never write `$dat_prefix.dat` all-NaN (read_dat_coordinates would drop it to empty) - finite(jd) = any(isfinite, view(aero.x, :, jd)) jz = findfirst(iszero, aero.delta_range) - jdat = jz !== nothing && finite(jz) ? jz : findfirst(finite, eachindex(aero.delta_range)) + jdat = jz !== nothing && finite(jz) ? jz : + findfirst(finite, eachindex(aero.delta_range)) write_dat("$dat_prefix.dat", "section", aero.x[:, jdat], aero.y[:, jdat]) cp_path = "$(table_prefix)_cp.$table_format" cf_path = "$(table_prefix)_cf.$table_format" diff --git a/src/section_aero.jl b/src/section_aero.jl index eb84ebb0..af52413a 100644 --- a/src/section_aero.jl +++ b/src/section_aero.jl @@ -118,21 +118,24 @@ function delta_suffix(delta) end """ - read_dat(path) -> (x, y) + read_dat_coordinates(path) -> (x, y) Read Selig `.dat` airfoil coordinates (two whitespace-separated columns), skipping the -name/header line and any non-numeric lines. +name/header line, comments, and any row that is not a finite pair. The single `.dat` +reader; [`write_dat`](@ref VortexStepMethod.AirfoilAero.write_dat) is its writer. """ -function read_dat(path::AbstractString) +function read_dat_coordinates(path::AbstractString) x = Float64[] y = Float64[] - for ln in eachline(String(path)) - p = split(strip(ln)) - length(p) >= 2 || continue - xv, yv = tryparse(Float64, p[1]), tryparse(Float64, p[2]) - (xv === nothing || yv === nothing) && continue - push!(x, xv) - push!(y, yv) + for line in eachline(String(path)) + fields = split(strip(line)) + length(fields) >= 2 || continue + xp = tryparse(Float64, fields[1]) + yp = tryparse(Float64, fields[2]) + (xp === nothing || yp === nothing) && continue + (isfinite(xp) && isfinite(yp)) || continue + push!(x, xp) + push!(y, yp) end return x, y end @@ -212,10 +215,11 @@ convert_node_table(src::AbstractString, dst::AbstractString) = read_section_aero(dat_file, cp_file, cf_file) -> Union{Nothing, SectionAero} Assemble a [`SectionAero`](@ref) from the human-readable files: the airfoil contour -(`dat_file`, plus `{stem}_{delta_suffix(δ)}.dat` per non-zero deflection) and the per-node `Cp` -and `cf` tables in the `.dat` node order, CSV or Arrow as their suffix says (see -[`read_node_table`](@ref)). Returns `nothing` if any file is missing. This is the single -loader for both provided and generated aero. +(`dat_file`, plus `{stem}_{delta_suffix(δ)}.dat` per non-zero deflection) and the +per-node `Cp` and `cf` tables in the `.dat` node order, CSV or Arrow as their suffix +says (see [`read_node_table`](@ref)). Returns `nothing` if one of the three named files +is missing, and `nothing` with a warning if a contour does not hold exactly the tables' +nodes. This is the single loader for both provided and generated aero. """ function read_section_aero(dat_file::AbstractString, cp_file::AbstractString, cf_file::AbstractString) @@ -229,8 +233,13 @@ function read_section_aero(dat_file::AbstractString, cp_file::AbstractString, x = fill(NaN, n_node, length(delta_range)) y = fill(NaN, n_node, length(delta_range)) for (jd, d) in enumerate(delta_range) - xd, yd = read_dat(iszero(d) ? String(dat_file) : - "$(stem)_$(delta_suffix(d)).dat") + contour = iszero(d) ? String(dat_file) : "$(stem)_$(delta_suffix(d)).dat" + xd, yd = isfile(contour) ? read_dat_coordinates(contour) : (Float64[], Float64[]) + if length(xd) != n_node + @warn "$contour holds $(length(xd)) finite coordinates, not the $n_node " * + "nodes of its Cp table; this airfoil gets no surface aero." + return nothing + end x[:, jd] .= xd y[:, jd] .= yd end diff --git a/test/airfoil_aero/test_airfoil_aero.jl b/test/airfoil_aero/test_airfoil_aero.jl index cd09c3f5..ac4b9363 100644 --- a/test/airfoil_aero/test_airfoil_aero.jl +++ b/test/airfoil_aero/test_airfoil_aero.jl @@ -14,17 +14,6 @@ seg_dist(px, py, ax, ay, bx, by) = begin hypot(px - (ax + t * vx), py - (ay + t * vy)) end -read_dat_coords(path) = begin - x = Float64[]; y = Float64[] - for ln in eachline(path) - s = strip(ln) - (isempty(s) || !(isdigit(s[1]) || s[1] == '-' || s[1] == '.')) && continue - p = split(s); length(p) >= 2 || continue - push!(x, parse(Float64, p[1])); push!(y, parse(Float64, p[2])) - end - x, y -end - @testset "Kulfan fit and NeuralFoil" begin @testset "Round-trip recovers known parameters" begin truth = KulfanParameters(fill(0.2, 8), fill(-0.2, 8), 0.0, 0.0) @@ -37,7 +26,7 @@ end end dat = joinpath(@__DIR__, "data", "test_airfoil.dat") - xr, yr = read_dat_coords(dat) + xr, yr = read_dat_coordinates(dat) params = fit_kulfan_parameters(xr, yr) @testset "Fit matches aerosandbox get_kulfan_parameters" begin @@ -165,6 +154,29 @@ end @test_throws ArgumentError write_section_aero(prefix, aero; table_format=:parquet) end +@testset "a deflection with no contour is neither written nor loaded as NaN" begin + alpha_range = deg2rad.([-5.0, 0.0, 5.0]) + delta_range = deg2rad.([0.0, 1.0]) + xc = [1.0, 0.5, 0.0, 0.5, 1.0] + yc = [0.0, 0.06, 0.0, -0.04, 0.0] + n_node = length(xc) + cp = [float(100i + 10ia + jd) for i in 1:n_node, + ia in eachindex(alpha_range), jd in eachindex(delta_range)] + aero = SectionAero(alpha_range, delta_range, hcat(xc, fill(NaN, n_node)), + hcat(yc, fill(NaN, n_node)), cp, cp ./ 1000) + + prefix = joinpath(mktempdir(), "af") + dat, cp_csv, cf_csv = write_section_aero(prefix, aero) + deflected = "$(prefix)_d1.dat" + @test !isfile(deflected) + @test isapprox(first(read_dat_coordinates(dat)), xc; atol=1e-6) + @test (@test_logs (:warn,) read_section_aero(dat, cp_csv, cf_csv)) === nothing + + write_dat(deflected, "section", fill(NaN, n_node), fill(NaN, n_node)) + @test isempty(first(read_dat_coordinates(deflected))) + @test (@test_logs (:warn,) read_section_aero(dat, cp_csv, cf_csv)) === nothing +end + @testset "generate_section_aero builds a surface table" begin truth = KulfanParameters(fill(0.15, 8), fill(-0.15, 8), 0.1, 0.0) alpha_range = deg2rad.(-4.0:2.0:4.0) @@ -181,7 +193,7 @@ end _, _, cp_hi_a, _ = section_surface(aero, alpha_range[end], delta_range[1]) @test maximum(abs.(cp_lo_a .- cp_hi_a)) > 0.1 - xdat, ydat = read_dat_coords(joinpath(@__DIR__, "data", "test_airfoil.dat")) + xdat, ydat = read_dat_coordinates(joinpath(@__DIR__, "data", "test_airfoil.dat")) aero2 = generate_section_aero(NeuralFoilSolver(model_size="medium"), xdat, ydat; alpha_range, delta_range=deg2rad.([0.0, 5.0]), reynolds_number=5e5) @test aero2 isa SectionAero @@ -237,7 +249,7 @@ end end @testset "generate_polar_from_coordinates POLAR_VECTORS sweep" begin - x, y = read_dat_coords(joinpath(@__DIR__, "data", "test_airfoil.dat")) + x, y = read_dat_coordinates(joinpath(@__DIR__, "data", "test_airfoil.dat")) csv = joinpath(mktempdir(), "polar.csv") sols = generate_polar_from_coordinates(x, y, csv; Re=5e5, alpha_range=-4:2:4, solver=NeuralFoilSolver(model_size="medium")) @@ -249,7 +261,7 @@ end end @testset "turn_trailing_edge! legacy crease cleanup" begin - x, y = read_dat_coords(joinpath(@__DIR__, "data", "test_airfoil.dat")) + x, y = read_dat_coordinates(joinpath(@__DIR__, "data", "test_airfoil.dat")) crease_frac = 0.7 for angle in (deg2rad(10.0), deg2rad(-10.0)) xd, yd = collect(float.(x)), collect(float.(y)) @@ -273,7 +285,7 @@ end end @testset "generate_airfoils fits the wrapped contour it is handed" begin - x_raw, y_raw = read_dat_coords(joinpath(@__DIR__, "data", "test_airfoil.dat")) + x_raw, y_raw = read_dat_coordinates(joinpath(@__DIR__, "data", "test_airfoil.dat")) x_fit, y_fit = shrink_wrap(x_raw, y_raw, ShrinkWrap(clearance=0.0)) _, fitted_y = kulfan_to_coordinates( fit_kulfan_parameters(x_fit, y_fit, LeastSquaresFit())) @@ -282,7 +294,7 @@ end Re=5e5, alpha_range=-2:2:2, aero_solver=NeuralFoilSolver(model_size="medium"), verbose=false) @test ok == [1] - _, written_y = read_dat_coords(joinpath(out, "airfoils", "1.dat")) + _, written_y = read_dat_coordinates(joinpath(out, "airfoils", "1.dat")) # A second shrink wrap inflates the section by its clearance, 0.006 — 60x this bound. @test maximum(abs, collect(extrema(written_y)) .- collect(extrema(fitted_y))) < 1e-4