From 316c187467bd8e0cb666b22ee79423e1190466a6 Mon Sep 17 00:00:00 2001 From: aeronauty Date: Wed, 8 Apr 2026 21:43:43 +0100 Subject: [PATCH 1/2] feat(flexfoil): add forced transition (xstrip) support Add XFOIL-equivalent forced transition via xstrip_upper/xstrip_lower parameters (x/c chord fraction, 0-1; default 1.0 = free transition). Thread through all layers: Rust core -> XFOIL-faithful path -> coupling solver -> Python bindings -> WASM bindings -> Web UI. The low-level trchek2_full() already handled x_forced; this change plumbs the user-specified xstrip value down to it via xstrip_to_xiforc conversion (XFOIL's XIFSET equivalent). Natural e^N transition takes precedence when it occurs before the forced location. Includes 12 new Python tests (binding + API level), Trip Upper/Lower inputs in the Solve panel, and updated docs/README. Co-Authored-By: Claude Opus 4.6 (1M context) --- crates/rustfoil-coupling/src/global_newton.rs | 11 +- crates/rustfoil-coupling/src/march.rs | 343 +++++++++++------- crates/rustfoil-python/src/lib.rs | 16 +- crates/rustfoil-solver/src/viscous/config.rs | 10 + crates/rustfoil-solver/src/viscous/viscal.rs | 28 +- crates/rustfoil-wasm/src/lib.rs | 22 +- crates/rustfoil-xfoil/src/assembly.rs | 17 +- crates/rustfoil-xfoil/src/config.rs | 36 ++ crates/rustfoil-xfoil/src/march.rs | 59 ++- crates/rustfoil-xfoil/src/oper.rs | 2 +- crates/rustfoil-xfoil/tests/march_parity.rs | 4 +- crates/rustfoil-xfoil/tests/support.rs | 6 +- docs-site/docs/python-api.mdx | 25 +- .../src/components/panels/SolvePanel.tsx | 34 +- flexfoil-ui/src/lib/searchIndex.ts | 4 +- flexfoil-ui/src/lib/wasm.ts | 4 + flexfoil-ui/src/stores/airfoilStore.ts | 6 + packages/flexfoil-python/README.md | 18 +- .../src/flexfoil/_rustfoil.pyi | 6 + .../flexfoil-python/src/flexfoil/airfoil.py | 20 +- packages/flexfoil-python/tests/test_api.py | 65 ++++ packages/flexfoil-python/tests/test_solver.py | 48 +++ 22 files changed, 614 insertions(+), 170 deletions(-) diff --git a/crates/rustfoil-coupling/src/global_newton.rs b/crates/rustfoil-coupling/src/global_newton.rs index 2ff166c0..2a6f8c8c 100644 --- a/crates/rustfoil-coupling/src/global_newton.rs +++ b/crates/rustfoil-coupling/src/global_newton.rs @@ -138,6 +138,11 @@ pub struct GlobalNewtonSystem { pub current_iteration: usize, /// ANTE: base thickness contribution at blunt trailing edge (XFOIL WGAP(1)) pub ante: f64, + /// Forced transition arc-length on upper surface (from XSTRIP). + /// `None` falls back to TE arc length (default XFOIL behaviour). + pub xiforc_upper: Option, + /// Forced transition arc-length on lower surface (from XSTRIP). + pub xiforc_lower: Option, } impl GlobalNewtonSystem { @@ -189,6 +194,8 @@ impl GlobalNewtonSystem { // Current iteration - set in build_global_system current_iteration: 0, ante: 0.0, + xiforc_upper: None, + xiforc_lower: None, } } @@ -790,9 +797,9 @@ impl GlobalNewtonSystem { let is_transition_interval = !s1.is_turbulent && !s1.is_wake && s2.is_turbulent && !s2.is_wake; let x_forced = if surface == 0 { - stations.get(self.iblte_upper).map(|s| s.x) + self.xiforc_upper.or_else(|| stations.get(self.iblte_upper).map(|s| s.x)) } else { - stations.get(self.iblte_lower).map(|s| s.x) + self.xiforc_lower.or_else(|| stations.get(self.iblte_lower).map(|s| s.x)) }; let live_transition = if is_transition_interval { Some(trchek2_full( diff --git a/crates/rustfoil-coupling/src/march.rs b/crates/rustfoil-coupling/src/march.rs index c221325d..a2fb9b3f 100644 --- a/crates/rustfoil-coupling/src/march.rs +++ b/crates/rustfoil-coupling/src/march.rs @@ -387,6 +387,10 @@ pub struct MarchConfig { pub deps: f64, /// Enable debug tracing output pub debug_trace: bool, + /// Forced transition arc-length location (XFOIL's XIFORC). + /// `None` means default behaviour: force transition at trailing edge only. + /// `Some(xi)` forces transition at the specified arc-length coordinate. + pub xiforc: Option, } impl Default for MarchConfig { @@ -400,6 +404,7 @@ impl Default for MarchConfig { senswt: 1000.0, deps: 5.0e-6, debug_trace: false, + xiforc: None, } } } @@ -422,6 +427,39 @@ impl MarchConfig { } } +// ============================================================================ +// XSTRIP → XIFORC Conversion +// ============================================================================ + +/// Convert an XSTRIP value (x/c chord fraction, 0–1) to an arc-length +/// coordinate (XIFORC) by interpolating the station arrays. +/// +/// If `xstrip >= 1.0`, returns the arc-length of the last station (TE), +/// matching XFOIL's default behaviour. Stations must have `x` (arc length) +/// and `x_coord` (physical x) populated. +pub fn xstrip_to_xiforc(stations: &[BlStation], xstrip: f64) -> f64 { + let te_x = stations.last().map(|s| s.x).unwrap_or(1.0); + if xstrip >= 1.0 - 1e-9 || stations.len() < 2 { + return te_x; + } + // Walk downstream from stagnation (stations are ordered by arc length). + // On each surface x_coord increases monotonically from LE towards TE. + for i in 1..stations.len() { + let xc1 = stations[i - 1].x_coord; + let xc2 = stations[i].x_coord; + if (xc1 <= xstrip && xstrip <= xc2) || (xc2 <= xstrip && xstrip <= xc1) { + let dx = xc2 - xc1; + if dx.abs() < 1e-15 { + return stations[i - 1].x; + } + let frac = (xstrip - xc1) / dx; + return stations[i - 1].x + frac * (stations[i].x - stations[i - 1].x); + } + } + // Fallback: xstrip beyond the station range → TE + te_x +} + // ============================================================================ // Stagnation Point Initialization // ============================================================================ @@ -1787,11 +1825,10 @@ pub fn march_fixed_ue( config.tolerance, config.hlmax, config.htmax, - None, + config.xiforc, Some(0), Some(i), ); - // If Newton didn't converge, use best estimate (station still has result) // Note: With tolerance=1e-5 (matching XFOIL), Newton should converge for normal cases @@ -2139,118 +2176,144 @@ pub fn march_surface( // Do not clamp Hk here; inverse-mode handles near-separation behavior. // Check for transition (laminar only) using TRCHEK2 implicit integration + // Also check for user-specified forced transition (XSTRIP/XIFORC). + let x_forced = config.xiforc.filter(|&xi| xi < surface_te_x - 1e-9); if is_laminar && result.x_transition.is_none() { - let trchek_result = trchek2_stations( - prev_station.x, - station.x, - prev_station.hk, - prev_station.theta, - prev_station.r_theta, - prev_station.delta_star, - prev_station.u, - prev_station.ampl, - station.hk, - station.theta, - station.r_theta, - station.delta_star, - station.u, - config.ncrit, - msq, - re, - ); - station.ampl = trchek_result.ampl2; - emit_trchek2_station_iter( - side, - station_idx + 2, - prev_station, - &station, - &trchek_result, - config.ncrit, - ); + // --- Forced transition check (XSTRIP) --- + // If user specified a forced transition location and this station + // has reached or passed it, force transition here. + let forced = x_forced.is_some_and(|xf| station.x >= xf && prev_station.x < xf); + if forced { + let xf = x_forced.unwrap(); + if config.debug_trace { + println!( + "[march_surface] FORCED TRANSITION (xstrip) at station {}, x={:.4}, xiforc={:.4}", + station_idx, station.x, xf + ); + } + is_laminar = false; + station.is_laminar = false; + station.is_turbulent = true; + result.x_transition = Some(xf); + result.transition_index = Some(station_idx + 2); - if config.debug_trace && station.ampl > 0.1 && station.ampl < 12.0 { - println!( - "[march_surface] i={} x={:.4} Hk={:.3} Rθ={:.0} N={:.3} dN={:.4}", - i, + let mut initial_guess = prev_station.clone(); + initial_guess.x = station.x; + initial_guess.u = station.u; + initial_guess.ampl = station.ampl; + initial_guess.ctau = 0.05; + + let (turb_station, _turb_converged, _turb_dmax) = newton_solve_station_transition( + prev_station, + Some(&initial_guess), + x[i], + station.u, + re, + msq, + false, + false, + config.ncrit, + config.max_iter, + config.tolerance, + config.hlmax, + config.htmax, + x_forced, + Some(side), + Some(station_idx + 2), + ); + + station = turb_station; + } else { + // --- Natural transition via e^N method --- + let trchek_result = trchek2_stations( + prev_station.x, station.x, + prev_station.hk, + prev_station.theta, + prev_station.r_theta, + prev_station.delta_star, + prev_station.u, + prev_station.ampl, station.hk, + station.theta, station.r_theta, - station.ampl, - station.ampl - prev_station.ampl + station.delta_star, + station.u, + config.ncrit, + msq, + re, + ); + station.ampl = trchek_result.ampl2; + emit_trchek2_station_iter( + side, + station_idx + 2, + prev_station, + &station, + &trchek_result, + config.ncrit, ); - } - - // Check for transition via e^N method ONLY - // - // IMPORTANT: XFOIL does NOT force transition when Hk exceeds a threshold. - // It only uses the e^N method (amplification exceeds Ncrit) for transition. - // When Hk exceeds hlmax, XFOIL switches to inverse mode to constrain H, - // but keeps the flow laminar unless N exceeds Ncrit. - // - // Previously, RustFoil had "separation-induced transition" that forced - // transition when Hk > 4.3. This was WRONG and caused premature transition, - // preventing proper stall prediction. High Hk values (4-10) are normal - // during laminar separation bubbles and should be handled by inverse mode, - // not by forcing transition. - - if trchek_result.transition { - // e^N transition detected - result.x_transition = trchek_result.xt; - result.transition_index = Some(station_idx + 2); - - if station.is_turbulent { - is_laminar = false; - } else { - is_laminar = false; - station.is_laminar = false; - station.is_turbulent = true; - - let mut initial_guess = prev_station.clone(); - initial_guess.x = station.x; - initial_guess.u = station.u; - initial_guess.ampl = station.ampl; - // XFOIL xbl.f line 895: IF(IBL.EQ.ITRAN(IS)) CTI = 0.05 - initial_guess.ctau = 0.05; - let (turb_station, _turb_converged, _turb_dmax) = newton_solve_station_transition( - prev_station, - Some(&initial_guess), - x[i], - station.u, - re, - msq, - false, - false, - config.ncrit, - config.max_iter, - config.tolerance, - config.hlmax, - config.htmax, - None, - Some(side), - Some(station_idx + 2), + if config.debug_trace && station.ampl > 0.1 && station.ampl < 12.0 { + println!( + "[march_surface] i={} x={:.4} Hk={:.3} Rθ={:.0} N={:.3} dN={:.4}", + i, + station.x, + station.hk, + station.r_theta, + station.ampl, + station.ampl - prev_station.ampl ); + } - station = turb_station; - station.ampl = trchek_result.ampl2; + // e^N transition: only triggers when amplification exceeds Ncrit. + // XFOIL does NOT force transition based on Hk threshold. + if trchek_result.transition { + // If forced transition is set but further downstream, e^N wins + result.x_transition = trchek_result.xt; + result.transition_index = Some(station_idx + 2); + + if station.is_turbulent { + is_laminar = false; + } else { + is_laminar = false; + station.is_laminar = false; + station.is_turbulent = true; + + let mut initial_guess = prev_station.clone(); + initial_guess.x = station.x; + initial_guess.u = station.u; + initial_guess.ampl = station.ampl; + initial_guess.ctau = 0.05; + + let (turb_station, _turb_converged, _turb_dmax) = newton_solve_station_transition( + prev_station, + Some(&initial_guess), + x[i], + station.u, + re, + msq, + false, + false, + config.ncrit, + config.max_iter, + config.tolerance, + config.hlmax, + config.htmax, + x_forced, + Some(side), + Some(station_idx + 2), + ); + + station = turb_station; + station.ampl = trchek_result.ampl2; + } } } } - // NOTE: XFOIL does NOT force transition when Hk > hlmax. It only uses hlmax - // to switch between direct and inverse mode. When Hk exceeds hlmax, the solver - // switches to inverse mode where Hk is prescribed, allowing laminar flow to - // continue through near-separation conditions. Transition is only triggered by: - // 1. e^N method (amplification exceeds Ncrit) - // 2. Forced transition at a user-specified location (XIFORC) - // - // In XFOIL, when no forced transition strip is set (XSTRIP >= 1.0), XIFORC is - // set to the TE arc length. This means laminar flow is forced to transition - // at the trailing edge if no natural transition has occurred. - // - // Implement forced TE transition: if this is the last airfoil station (or very - // close to TE, arc length > 0.98) and still laminar, force transition. - // This matches XFOIL's behavior of forcing transition at TE. + // Forced TE transition: if no user-specified strip or e^N has triggered, + // and we've reached the trailing edge, force transition here. + // This matches XFOIL's default XIFORC = TE arc length. let at_surface_te = station.x >= surface_te_x - 1e-9; if is_laminar && at_surface_te && result.x_transition.is_none() { if config.debug_trace { @@ -2565,28 +2628,12 @@ pub fn march_mixed_du( let mut transition_detected_this_iter = false; // TRCHEK inside Newton loop (XFOIL xbl.f:1081-1088) + // Also check for user-specified forced transition (XSTRIP/XIFORC). if !simi && !turb && !wake && !tran { - let trchek = trchek2_stations( - prev.x, - stations[i].x, - prev.hk, - prev.theta, - prev.r_theta, - prev.delta_star, - prev.u, - prev.ampl, - stations[i].hk, - stations[i].theta, - stations[i].r_theta, - stations[i].delta_star, - stations[i].u, - config.ncrit, - msq, - re, - ); - ami = trchek.ampl2; - stations[i].ampl = ami; - if trchek.transition { + // Forced transition check: if user set a strip and this station + // has reached or passed it, force transition immediately. + let forced = config.xiforc.is_some_and(|xf| stations[i].x >= xf && prev.x < xf); + if forced { let tr_full = trchek2_full( prev.x, stations[i].x, @@ -2602,7 +2649,7 @@ pub fn march_mixed_du( stations[i].r_theta, prev.ampl, config.ncrit, - None, + config.xiforc, msq, re, ); @@ -2614,6 +2661,56 @@ pub fn march_mixed_du( stations[i].ctau = 0.03; cti = 0.03; } + } else { + let trchek = trchek2_stations( + prev.x, + stations[i].x, + prev.hk, + prev.theta, + prev.r_theta, + prev.delta_star, + prev.u, + prev.ampl, + stations[i].hk, + stations[i].theta, + stations[i].r_theta, + stations[i].delta_star, + stations[i].u, + config.ncrit, + msq, + re, + ); + ami = trchek.ampl2; + stations[i].ampl = ami; + if trchek.transition { + let tr_full = trchek2_full( + prev.x, + stations[i].x, + prev.theta, + stations[i].theta, + prev.delta_star, + stations[i].delta_star, + prev.u, + stations[i].u, + prev.hk, + stations[i].hk, + prev.r_theta, + stations[i].r_theta, + prev.ampl, + config.ncrit, + config.xiforc, + msq, + re, + ); + trchek_hold = Some(tr_full); + transition_detected_this_iter = true; + tran = true; + flow_type = FlowType::Turbulent; + if itbl == 0 { + stations[i].ctau = 0.03; + cti = 0.03; + } + } } } diff --git a/crates/rustfoil-python/src/lib.rs b/crates/rustfoil-python/src/lib.rs index 6dad2e12..e94605ca 100644 --- a/crates/rustfoil-python/src/lib.rs +++ b/crates/rustfoil-python/src/lib.rs @@ -43,7 +43,7 @@ fn flat_from_points(pts: &[Point]) -> Vec<(f64, f64)> { /// Returns a dict with keys: cl, cd, cm, converged, iterations, residual, /// x_tr_upper, x_tr_lower, cd_friction, cd_pressure, alpha_deg, success, error. #[pyfunction] -#[pyo3(signature = (coords, alpha_deg, reynolds=1.0e6, mach=0.0, ncrit=9.0, max_iterations=100, re_type=1))] +#[pyo3(signature = (coords, alpha_deg, reynolds=1.0e6, mach=0.0, ncrit=9.0, max_iterations=100, re_type=1, xstrip_upper=1.0, xstrip_lower=1.0))] fn analyze_faithful( py: Python<'_>, coords: Vec, @@ -53,6 +53,8 @@ fn analyze_faithful( ncrit: f64, max_iterations: usize, re_type: u8, + xstrip_upper: f64, + xstrip_lower: f64, ) -> PyResult> { if coords.len() < 6 || coords.len() % 2 != 0 { let d = PyDict::new(py); @@ -78,6 +80,8 @@ fn analyze_faithful( ncrit, max_iterations, re_type: re_type_from_int(re_type), + xstrip_upper, + xstrip_lower, ..Default::default() }; @@ -291,7 +295,7 @@ fn faithful_result_to_pydict(py: Python<'_>, r: &FaithfulResult) -> PyResult, coords: Vec, @@ -301,6 +305,8 @@ fn analyze_faithful_batch( ncrit: f64, max_iterations: usize, re_type: u8, + xstrip_upper: f64, + xstrip_lower: f64, ) -> PyResult>> { let err_msg = if coords.len() < 6 || coords.len() % 2 != 0 { Some("Invalid coordinates".to_string()) @@ -343,6 +349,7 @@ fn analyze_faithful_batch( let options = XfoilOptions { reynolds, mach, ncrit, max_iterations, re_type: re_type_from_int(re_type), + xstrip_upper, xstrip_lower, ..Default::default() }; @@ -440,7 +447,7 @@ fn analyze_inviscid_batch( /// ue_upper, ue_lower, x_tr_upper, x_tr_lower, converged, iterations, /// residual, success, error. #[pyfunction] -#[pyo3(signature = (coords, alpha_deg, reynolds=1.0e6, mach=0.0, ncrit=9.0, max_iterations=100, re_type=1))] +#[pyo3(signature = (coords, alpha_deg, reynolds=1.0e6, mach=0.0, ncrit=9.0, max_iterations=100, re_type=1, xstrip_upper=1.0, xstrip_lower=1.0))] fn get_bl_distribution( py: Python<'_>, coords: Vec, @@ -450,6 +457,8 @@ fn get_bl_distribution( ncrit: f64, max_iterations: usize, re_type: u8, + xstrip_upper: f64, + xstrip_lower: f64, ) -> PyResult> { use rustfoil_xfoil::oper::{build_state_from_coords, solve_operating_point_from_state, coords_from_body, AlphaSpec}; use rustfoil_xfoil::XfoilOptions; @@ -476,6 +485,7 @@ fn get_bl_distribution( let options = XfoilOptions { reynolds, mach, ncrit, max_iterations, re_type: re_type_from_int(re_type), + xstrip_upper, xstrip_lower, ..Default::default() }; diff --git a/crates/rustfoil-solver/src/viscous/config.rs b/crates/rustfoil-solver/src/viscous/config.rs index 6545b786..543436d5 100644 --- a/crates/rustfoil-solver/src/viscous/config.rs +++ b/crates/rustfoil-solver/src/viscous/config.rs @@ -76,6 +76,14 @@ pub struct ViscousSolverConfig { /// Global operating-variable mode used by the Newton/update loop. pub operating_mode: OperatingMode, + + /// Forced transition location on upper surface as x/c fraction (0–1). + /// Default 1.0 means no forced transition (use e^N method only). + pub xstrip_upper: f64, + + /// Forced transition location on lower surface as x/c fraction (0–1). + /// Default 1.0 means no forced transition (use e^N method only). + pub xstrip_lower: f64, } impl Default for ViscousSolverConfig { @@ -91,6 +99,8 @@ impl Default for ViscousSolverConfig { hk_max_turbulent: 2.5, allow_separation: true, operating_mode: OperatingMode::PrescribedAlpha, + xstrip_upper: 1.0, + xstrip_lower: 1.0, } } } diff --git a/crates/rustfoil-solver/src/viscous/viscal.rs b/crates/rustfoil-solver/src/viscous/viscal.rs index 9abd7acc..87feafd3 100644 --- a/crates/rustfoil-solver/src/viscous/viscal.rs +++ b/crates/rustfoil-solver/src/viscous/viscal.rs @@ -22,7 +22,7 @@ use nalgebra::DMatrix; use rayon::prelude::*; use rustfoil_bl::equations::{blvar, FlowType}; use rustfoil_bl::state::BlStation; -use rustfoil_coupling::march::{march_fixed_ue, MarchConfig, MarchResult}; +use rustfoil_coupling::march::{march_fixed_ue, xstrip_to_xiforc, MarchConfig, MarchResult}; use rustfoil_coupling::global_newton::{ apply_global_updates_from_view, preview_global_update_ue_from_view, solve_global_system, GlobalNewtonSystem, }; @@ -891,7 +891,7 @@ pub fn solve_viscous_two_surfaces( )); } - let march_config = MarchConfig { + let march_config_base = MarchConfig { ncrit: config.ncrit, hlmax: config.hk_max_laminar, htmax: config.hk_max_turbulent, @@ -914,10 +914,24 @@ pub fn solve_viscous_two_surfaces( canonical_state.set_panel_inviscid_arrays(ue_inviscid_full, ue_inviscid_alpha_full); refresh_canonical_panel_arrays(&mut canonical_state); + // Compute per-surface forced transition arc-length from xstrip (XIFSET equivalent). + let xiforc_upper = if config.xstrip_upper < 1.0 - 1e-9 { + Some(xstrip_to_xiforc(upper_stations, config.xstrip_upper)) + } else { + None + }; + let xiforc_lower = if config.xstrip_lower < 1.0 - 1e-9 { + Some(xstrip_to_xiforc(lower_stations, config.xstrip_lower)) + } else { + None + }; + + let upper_march_config = MarchConfig { xiforc: xiforc_upper, ..march_config_base.clone() }; + let lower_march_config = MarchConfig { xiforc: xiforc_lower, ..march_config_base.clone() }; // March upper surface (side 1) through the canonical state. let upper_result = - canonical_state.march_surface(XfoilSurface::Upper, re, msq, &march_config); + canonical_state.march_surface(XfoilSurface::Upper, re, msq, &upper_march_config); // Debug: Print march results if rustfoil_bl::is_debug_active() { @@ -937,7 +951,7 @@ pub fn solve_viscous_two_surfaces( // March lower surface (side 2) through the canonical state. let lower_result = - canonical_state.march_surface(XfoilSurface::Lower, re, msq, &march_config); + canonical_state.march_surface(XfoilSurface::Lower, re, msq, &lower_march_config); // Debug: Print lower march results if rustfoil_bl::is_debug_active() { @@ -1030,6 +1044,8 @@ pub fn solve_viscous_two_surfaces( let iblte_lower = canonical_state.iblte_lower; let mut global_system = GlobalNewtonSystem::new(n_upper, n_lower, iblte_upper, iblte_lower); + global_system.xiforc_upper = xiforc_upper; + global_system.xiforc_lower = xiforc_lower; // ANTE: trailing edge gap thickness (XFOIL's WGAP(1) = DWTE). // For sharp TE airfoils this is zero; for blunt TE it's the gap. @@ -1113,8 +1129,8 @@ pub fn solve_viscous_two_surfaces( // bubble and N-factor history across Newton iterations, unlike the // previous march_surface call which re-initialised from Thwaites. - canonical_state.march_mixed_du(XfoilSurface::Upper, re, msq, &march_config); - canonical_state.march_mixed_du(XfoilSurface::Lower, re, msq, &march_config); + canonical_state.march_mixed_du(XfoilSurface::Upper, re, msq, &upper_march_config); + canonical_state.march_mixed_du(XfoilSurface::Lower, re, msq, &lower_march_config); let upper_debug = canonical_state.upper_station_view(); let lower_debug = canonical_state.lower_station_view(); diff --git a/crates/rustfoil-wasm/src/lib.rs b/crates/rustfoil-wasm/src/lib.rs index 0729f00a..270a04ca 100644 --- a/crates/rustfoil-wasm/src/lib.rs +++ b/crates/rustfoil-wasm/src/lib.rs @@ -941,6 +941,8 @@ fn faithful_snapshot( ncrit: f64, max_iterations: usize, re_type: u8, + xstrip_upper: f64, + xstrip_lower: f64, ) -> Result { if coords.len() < 6 || coords.len() % 2 != 0 { return Err("Invalid coordinates: need at least 3 points (6 values)".to_string()); @@ -1006,6 +1008,8 @@ fn faithful_snapshot( ncrit, max_iterations, re_type: re_type_from_int(re_type), + xstrip_upper, + xstrip_lower, ..Default::default() }; let oper_result = solve_operating_point_from_state(&mut state, &factorized, &options) @@ -1152,8 +1156,10 @@ pub fn analyze_airfoil_faithful( ncrit: f64, max_iterations: usize, re_type: u8, + xstrip_upper: f64, + xstrip_lower: f64, ) -> JsValue { - let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, re_type) { + let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, re_type, xstrip_upper, xstrip_lower) { Ok(snapshot) => snapshot.result, Err(message) => analysis_error(message), }; @@ -1169,7 +1175,7 @@ pub fn get_bl_distribution_faithful( ncrit: f64, max_iterations: usize, ) -> JsValue { - let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1) { + let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1, 1.0, 1.0) { Ok(snapshot) => rows_to_bl_distribution( &snapshot.upper_rows, &snapshot.lower_rows, @@ -1275,7 +1281,7 @@ pub fn get_bl_visualization_faithful( ncrit: f64, max_iterations: usize, ) -> JsValue { - let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1) { + let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1, 1.0, 1.0) { Ok(snapshot) => build_bl_visualization(&snapshot), Err(message) => bl_vis_error(message), }; @@ -1453,7 +1459,7 @@ pub fn compute_streamlines_faithful( seed_count: u32, bounds: &[f64], ) -> JsValue { - let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1) { + let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1, 1.0, 1.0) { Ok(snapshot) => { if bounds.len() != 4 { StreamlineResult { @@ -1596,7 +1602,7 @@ pub fn compute_dividing_streamline_faithful( max_iterations: usize, bounds: &[f64], ) -> JsValue { - let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1) { + let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1, 1.0, 1.0) { Ok(snapshot) => { if bounds.len() != 4 { DividingStreamlineResult { @@ -1820,7 +1826,7 @@ pub fn compute_psi_grid_faithful( bounds: &[f64], resolution: &[u32], ) -> JsValue { - let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1) { + let result = match faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1, 1.0, 1.0) { Ok(snapshot) => { if bounds.len() != 4 { PsiGridResult { @@ -1963,7 +1969,7 @@ impl WasmSmokeSystem { ncrit: f64, max_iterations: usize, ) { - if let Ok(snapshot) = faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1) { + if let Ok(snapshot) = faithful_snapshot(coords, alpha_deg, reynolds, mach, ncrit, max_iterations, 1, 1.0, 1.0) { let field = FaithfulFlowField::from_snapshot(&snapshot, alpha_deg); self.coords = field.nodes.clone(); self.gamma = field.gamma.clone(); @@ -3514,7 +3520,7 @@ mod tests { let paneled_flat = repanel_xfoil(&buffer_flat, 160); let alpha_deg = 15.0; - let snapshot = faithful_snapshot(&paneled_flat, alpha_deg, 1.0e6, 0.0, 9.0, 100, 1) + let snapshot = faithful_snapshot(&paneled_flat, alpha_deg, 1.0e6, 0.0, 9.0, 100, 1, 1.0, 1.0) .expect("faithful snapshot should solve"); let field = FaithfulFlowField::from_snapshot(&snapshot, alpha_deg); let options = StreamlineOptions { diff --git a/crates/rustfoil-xfoil/src/assembly.rs b/crates/rustfoil-xfoil/src/assembly.rs index ac4cd594..867ae7d2 100644 --- a/crates/rustfoil-xfoil/src/assembly.rs +++ b/crates/rustfoil-xfoil/src/assembly.rs @@ -20,9 +20,11 @@ pub fn setbl( ncrit: f64, mach: f64, iteration: usize, + xstrip_upper: f64, + xstrip_lower: f64, ) -> AssemblyState { if !state.lblini { - mrchue(state, reynolds, ncrit, iteration); + mrchue(state, reynolds, ncrit, iteration, xstrip_upper, xstrip_lower); if std::env::var("RUSTFOIL_SETBL_HANDOFF_DEBUG").is_ok() { let start = state.nbl_upper.saturating_sub(3); for ibl in start..state.nbl_upper { @@ -42,7 +44,7 @@ pub fn setbl( } } } - mrchdu(state, reynolds, ncrit, iteration); + mrchdu(state, reynolds, ncrit, iteration, xstrip_upper, xstrip_lower); // XFOIL re-enters SETBL with MASS coherent to the current BL state. // Keep the row cache aligned before building canonical stations or the @@ -215,6 +217,17 @@ pub fn setbl( state.iblte_lower, ); system.ante = state.ante; + // Set user-specified forced transition locations (converted from xstrip to arc-length) + system.xiforc_upper = if xstrip_upper < 1.0 - 1e-9 { + Some(crate::march::xstrip_to_xiforc_rows(upper_rows, state.iblte_upper, xstrip_upper)) + } else { + None + }; + system.xiforc_lower = if xstrip_lower < 1.0 - 1e-9 { + Some(crate::march::xstrip_to_xiforc_rows(lower_rows, state.iblte_lower, xstrip_lower)) + } else { + None + }; system.set_stagnation_derivs(state.sst_go, state.sst_gp); system.build_global_system( &upper_stations, diff --git a/crates/rustfoil-xfoil/src/config.rs b/crates/rustfoil-xfoil/src/config.rs index ca8d6ac8..64877782 100644 --- a/crates/rustfoil-xfoil/src/config.rs +++ b/crates/rustfoil-xfoil/src/config.rs @@ -36,6 +36,12 @@ pub struct XfoilOptions { pub wake_length_chords: f64, pub operating_mode: OperatingMode, pub re_type: ReType, + /// Forced transition location on upper surface as x/c fraction (0–1). + /// Default 1.0 means no forced transition (use e^N method only). + pub xstrip_upper: f64, + /// Forced transition location on lower surface as x/c fraction (0–1). + /// Default 1.0 means no forced transition (use e^N method only). + pub xstrip_lower: f64, } impl Default for XfoilOptions { @@ -49,6 +55,8 @@ impl Default for XfoilOptions { wake_length_chords: 1.0, operating_mode: OperatingMode::PrescribedAlpha, re_type: ReType::Type1, + xstrip_upper: 1.0, + xstrip_lower: 1.0, } } } @@ -91,6 +99,12 @@ impl XfoilOptions { if self.tolerance <= 0.0 { return Err("tolerance must be positive"); } + if !(0.0..=1.0).contains(&self.xstrip_upper) { + return Err("xstrip_upper must be in [0, 1]"); + } + if !(0.0..=1.0).contains(&self.xstrip_lower) { + return Err("xstrip_lower must be in [0, 1]"); + } Ok(()) } } @@ -154,4 +168,26 @@ mod tests { let o = XfoilOptions::default(); assert_eq!(o.re_type, ReType::Type1); } + + #[test] + fn default_xstrip_is_1() { + let o = XfoilOptions::default(); + assert_eq!(o.xstrip_upper, 1.0); + assert_eq!(o.xstrip_lower, 1.0); + } + + #[test] + fn xstrip_validation() { + let mut o = XfoilOptions::default(); + o.xstrip_upper = 0.3; + o.xstrip_lower = 0.5; + assert!(o.validate().is_ok()); + + o.xstrip_upper = -0.1; + assert!(o.validate().is_err()); + + o.xstrip_upper = 0.0; + o.xstrip_lower = 1.1; + assert!(o.validate().is_err()); + } } diff --git a/crates/rustfoil-xfoil/src/march.rs b/crates/rustfoil-xfoil/src/march.rs index ebcb2ec5..ac4ae8e3 100644 --- a/crates/rustfoil-xfoil/src/march.rs +++ b/crates/rustfoil-xfoil/src/march.rs @@ -228,12 +228,46 @@ pub fn blpini(state: &mut XfoilState, _reynolds: f64) { let _ = state; } +// --------------------------------------------------------------------------- +// XSTRIP → XIFORC conversion for XfoilBlRow arrays +// --------------------------------------------------------------------------- + +/// Convert an XSTRIP value (x/c chord fraction) to an arc-length coordinate +/// by interpolating the BL row arrays. Returns TE arc length when +/// `xstrip >= 1.0` (no forced transition). +pub fn xstrip_to_xiforc_rows(rows: &[XfoilBlRow], iblte: usize, xstrip: f64) -> f64 { + let te_x = rows.get(iblte).map(|r| r.x).unwrap_or(1.0); + if xstrip >= 1.0 - 1e-9 || iblte < 2 { + return te_x; + } + for i in 1..=iblte.min(rows.len().saturating_sub(1)) { + let xc1 = rows[i - 1].x_coord; + let xc2 = rows[i].x_coord; + if (xc1 <= xstrip && xstrip <= xc2) || (xc2 <= xstrip && xstrip <= xc1) { + let dx = xc2 - xc1; + if dx.abs() < 1e-15 { + return rows[i - 1].x; + } + let frac = (xstrip - xc1) / dx; + return rows[i - 1].x + frac * (rows[i].x - rows[i - 1].x); + } + } + te_x +} + // --------------------------------------------------------------------------- // MRCHUE // --------------------------------------------------------------------------- /// MRCHUE: Direct-mode BL march with transition checking (xbl.f:584-986). -pub fn mrchue(state: &mut XfoilState, reynolds: f64, ncrit: f64, iteration: usize) { +pub fn mrchue( + state: &mut XfoilState, + reynolds: f64, + ncrit: f64, + iteration: usize, + xstrip_upper: f64, + xstrip_lower: f64, +) { let msq = 0.0_f64; // Upper surface (IS=1): other side = lower (stale values) @@ -241,7 +275,7 @@ pub fn mrchue(state: &mut XfoilState, reynolds: f64, ncrit: f64, iteration: usiz let (itran_u, xssitr_u) = march_ue_side( 1, &mut state.upper_rows, state.nbl_upper, state.iblte_upper, - other_te, state.ante, reynolds, ncrit, msq, iteration, + other_te, state.ante, reynolds, ncrit, msq, iteration, xstrip_upper, ); state.itran_upper = itran_u; state.xssitr_upper = xssitr_u; @@ -251,7 +285,7 @@ pub fn mrchue(state: &mut XfoilState, reynolds: f64, ncrit: f64, iteration: usiz let (itran_l, xssitr_l) = march_ue_side( 2, &mut state.lower_rows, state.nbl_lower, state.iblte_lower, - other_te, state.ante, reynolds, ncrit, msq, iteration, + other_te, state.ante, reynolds, ncrit, msq, iteration, xstrip_lower, ); state.itran_lower = itran_l; state.xssitr_lower = xssitr_l; @@ -285,6 +319,7 @@ fn march_ue_side( ncrit: f64, msq: f64, iteration: usize, + xstrip: f64, ) -> (usize, f64) { if nbl < 2 { return (iblte, 0.0); @@ -293,7 +328,7 @@ fn march_ue_side( let mut itran = xfoil_ibl(iblte); let mut xssitr = 0.0; let mut turb = false; - let x_forced = rows.get(iblte).map(|row| row.x); + let x_forced = Some(xstrip_to_xiforc_rows(rows, iblte, xstrip)); // ---- Similarity station init (IBL=2 = index 1) ---- let xsi0 = rows[1].x.max(1e-12); @@ -732,7 +767,14 @@ fn march_ue_side( // --------------------------------------------------------------------------- /// MRCHDU: Mixed Ue-Hk march with transition checking (xbl.f:990-1312). -pub fn mrchdu(state: &mut XfoilState, reynolds: f64, ncrit: f64, iteration: usize) { +pub fn mrchdu( + state: &mut XfoilState, + reynolds: f64, + ncrit: f64, + iteration: usize, + xstrip_upper: f64, + xstrip_lower: f64, +) { let msq = 0.0_f64; // Upper surface @@ -741,7 +783,7 @@ pub fn mrchdu(state: &mut XfoilState, reynolds: f64, ncrit: f64, iteration: usiz let (itran_u, xssitr_u) = march_du_side( 1, &mut state.upper_rows, state.nbl_upper, state.iblte_upper, - other_te, state.ante, itrold_u, reynolds, ncrit, msq, iteration, + other_te, state.ante, itrold_u, reynolds, ncrit, msq, iteration, xstrip_upper, ); state.itran_upper = itran_u; state.xssitr_upper = xssitr_u; @@ -752,7 +794,7 @@ pub fn mrchdu(state: &mut XfoilState, reynolds: f64, ncrit: f64, iteration: usiz let (itran_l, xssitr_l) = march_du_side( 2, &mut state.lower_rows, state.nbl_lower, state.iblte_lower, - other_te, state.ante, itrold_l, reynolds, ncrit, msq, iteration, + other_te, state.ante, itrold_l, reynolds, ncrit, msq, iteration, xstrip_lower, ); state.itran_lower = itran_l; state.xssitr_lower = xssitr_l; @@ -772,6 +814,7 @@ fn march_du_side( ncrit: f64, msq: f64, iteration: usize, + xstrip: f64, ) -> (usize, f64) { if nbl < 2 { return (iblte, 0.0); @@ -785,7 +828,7 @@ fn march_du_side( let mut turb = false; let mut sens = 0.0_f64; let mut ami = rows[1].ctau; - let x_forced = rows.get(iblte).map(|row| row.x); + let x_forced = Some(xstrip_to_xiforc_rows(rows, iblte, xstrip)); // Build station_1 from stored similarity station (index 1) let ft0 = flow_type_for(1, itrold, false); diff --git a/crates/rustfoil-xfoil/src/oper.rs b/crates/rustfoil-xfoil/src/oper.rs index b7978714..44552315 100644 --- a/crates/rustfoil-xfoil/src/oper.rs +++ b/crates/rustfoil-xfoil/src/oper.rs @@ -174,7 +174,7 @@ pub fn solve_operating_point_from_state( for iter in 1..=options.max_iterations { let re_eff = options.effective_reynolds(state.cl); - let mut assembly = setbl(state, re_eff, options.ncrit, options.mach, iter); + let mut assembly = setbl(state, re_eff, options.ncrit, options.mach, iter, options.xstrip_upper, options.xstrip_lower); let solve = blsolv(state, &mut assembly, iter); update( state, diff --git a/crates/rustfoil-xfoil/tests/march_parity.rs b/crates/rustfoil-xfoil/tests/march_parity.rs index cb1ed677..6a75e4ab 100644 --- a/crates/rustfoil-xfoil/tests/march_parity.rs +++ b/crates/rustfoil-xfoil/tests/march_parity.rs @@ -120,8 +120,8 @@ fn perf_marching_path_smoke() { || { state.clone_from(&base); blpini(&mut state, 1.0e6); - mrchue(&mut state, 1.0e6, 9.0, 0); - mrchdu(&mut state, 1.0e6, 9.0, 0); + mrchue(&mut state, 1.0e6, 9.0, 0, 1.0, 1.0); + mrchdu(&mut state, 1.0e6, 9.0, 0, 1.0, 1.0); }, ); diff --git a/crates/rustfoil-xfoil/tests/support.rs b/crates/rustfoil-xfoil/tests/support.rs index 076743ad..10e1857f 100644 --- a/crates/rustfoil-xfoil/tests/support.rs +++ b/crates/rustfoil-xfoil/tests/support.rs @@ -385,7 +385,7 @@ pub fn rust_iteration_debug(alpha_deg: f64) -> serde_json::Value { let (mut state, _) = build_march_state(alpha_deg); blpini(&mut state, 1.0e6); init_debug(&debug_path); - let assembly = setbl(&mut state, 1.0e6, 9.0, 0.0, 1); + let assembly = setbl(&mut state, 1.0e6, 9.0, 0.0, 1, 1.0, 1.0); let mut assembly = assembly; let solve = blsolv(&mut state, &mut assembly, 1); update(&mut state, &assembly, &solve, 0.0, 1.0e6, 1); @@ -608,14 +608,14 @@ pub fn build_march_state(alpha_deg: f64) -> (XfoilState, rustfoil_inviscid::Fact pub fn rust_setbl_state(alpha_deg: f64) -> FortranMarchState { let (mut state, _) = build_march_state(alpha_deg); blpini(&mut state, 1.0e6); - let _ = setbl(&mut state, 1.0e6, 9.0, 0.0, 1); + let _ = setbl(&mut state, 1.0e6, 9.0, 0.0, 1, 1.0, 1.0); march_state_snapshot(&state) } pub fn rust_mrchue_state(alpha_deg: f64) -> FortranMarchState { let (mut state, _) = build_march_state(alpha_deg); blpini(&mut state, 1.0e6); - mrchue(&mut state, 1.0e6, 9.0, 0); + mrchue(&mut state, 1.0e6, 9.0, 0, 1.0, 1.0); march_state_snapshot(&state) } diff --git a/docs-site/docs/python-api.mdx b/docs-site/docs/python-api.mdx index 1c167ee6..85e75446 100644 --- a/docs-site/docs/python-api.mdx +++ b/docs-site/docs/python-api.mdx @@ -157,6 +157,25 @@ print(result.reynolds_eff) # effective Re after Mode 2 adjustment polar = foil.polar(alpha=(-5, 15, 0.5), Re=1e6, re_type=2) ``` +### Forced transition (trip strips) + +Force the boundary layer to transition at a specific chord location, +regardless of the e^N criterion. This is equivalent to XFOIL's `XSTRIP` +command (VPAR menu). Default is `1.0` (no forced transition — use e^N only). + +```python +# Force transition at 5% chord on both surfaces (simulating trip strips) +result = foil.solve(alpha=5.0, Re=1e6, xstrip_upper=0.05, xstrip_lower=0.05) +print(result.x_tr_upper) # ≈ 0.05 +print(result.x_tr_lower) # ≈ 0.05 + +# Force upper only — lower uses natural (e^N) transition +polar = foil.polar(alpha=(-5, 15, 0.5), Re=1e6, xstrip_upper=0.1) +``` + +If the natural (e^N) transition occurs before the forced location, natural +transition wins — matching XFOIL's behaviour. + ### Inviscid analysis ```python @@ -260,9 +279,9 @@ are byte-identical. | `.panel_coords` | Repaneled coordinates `list[(x, y)]` | | `.hash` | SHA-256 hash of panel coords (cache key) | | `.with_flap(hinge_x, deflection, hinge_y_frac, n_panels)` | Return new Airfoil with flap deflected | -| `.solve(alpha, Re, mach, ncrit, max_iter, viscous, store, re_type)` | Single-point analysis (`re_type`: 1, 2, or 3) | -| `.polar(alpha, Re, mach, ncrit, max_iter, viscous, store, parallel, re_type)` | Sweep over alpha range (parallel by default) | -| `.bl_distribution(alpha, Re, mach, ncrit, max_iter, re_type)` | Boundary-layer distribution at a single alpha | +| `.solve(alpha, Re, mach, ncrit, max_iter, viscous, store, re_type, xstrip_upper, xstrip_lower)` | Single-point analysis | +| `.polar(alpha, Re, mach, ncrit, max_iter, viscous, store, parallel, re_type, xstrip_upper, xstrip_lower)` | Sweep over alpha range (parallel by default) | +| `.bl_distribution(alpha, Re, mach, ncrit, max_iter, re_type, xstrip_upper, xstrip_lower)` | Boundary-layer distribution at a single alpha | ### SolveResult diff --git a/flexfoil-ui/src/components/panels/SolvePanel.tsx b/flexfoil-ui/src/components/panels/SolvePanel.tsx index c1d8d316..c9860e98 100644 --- a/flexfoil-ui/src/components/panels/SolvePanel.tsx +++ b/flexfoil-ui/src/components/panels/SolvePanel.tsx @@ -63,6 +63,10 @@ export function SolvePanel() { setSolverMode, reType, setReType, + xstripUpper, + xstripLower, + setXstripUpper, + setXstripLower, upsertPolar, clearAllPolars, } = useAirfoilStore(); @@ -180,10 +184,10 @@ export function SolvePanel() { const runSolver = useCallback((coords: { x: number; y: number }[], alpha: number) => { if (isViscous) { - return analyzeAirfoil(coords, alpha, reynolds, mach, ncrit, maxIterations, reType); + return analyzeAirfoil(coords, alpha, reynolds, mach, ncrit, maxIterations, reType, xstripUpper, xstripLower); } return analyzeAirfoilInviscid(coords, alpha); - }, [isViscous, reynolds, mach, ncrit, maxIterations, reType]); + }, [isViscous, reynolds, mach, ncrit, maxIterations, reType, xstripUpper, xstripLower]); /** Cache key uses Re=0/Mach=0/Ncrit=0/maxIter=0 for inviscid to separate caches. */ const cacheRe = isViscous ? reynolds : 0; @@ -967,6 +971,32 @@ export function SolvePanel() { max={14} /> +
+ + setXstripUpper(Math.max(0, Math.min(1, parseFloat(e.target.value) || 1)))} + step={0.05} + min={0} + max={1} + /> +
+
+ + setXstripLower(Math.max(0, Math.min(1, parseFloat(e.target.value) || 1)))} + step={0.05} + min={0} + max={1} + /> +
void; setSolverMode: (mode: SolverMode) => void; setReType: (reType: ReType) => void; + setXstripUpper: (xstripUpper: number) => void; + setXstripLower: (xstripLower: number) => void; // Point manipulation (legacy, kept for compatibility) updatePoint: (index: number, point: AirfoilPoint) => void; @@ -291,6 +293,8 @@ export const useAirfoilStore = create()( maxIterations: 100, solverMode: 'viscous', reType: 1 as ReType, + xstripUpper: 1.0, + xstripLower: 1.0, polarData: [], spacingPanelMode: 'simple', // Default to simple curvature-based sspInterpolation: 'linear', // Default to linear (Mark Drela's original) @@ -320,6 +324,8 @@ export const useAirfoilStore = create()( setMaxIterations: (maxIterations) => set({ maxIterations: Math.max(10, Math.min(500, maxIterations)) }), setSolverMode: (solverMode) => set({ solverMode }), setReType: (reType) => set({ reType }), + setXstripUpper: (xstripUpper) => set({ xstripUpper: Math.max(0, Math.min(1, xstripUpper)) }), + setXstripLower: (xstripLower) => set({ xstripLower: Math.max(0, Math.min(1, xstripLower)) }), updatePoint: (index, point) => set((state) => { const newCoords = [...state.coordinates]; diff --git a/packages/flexfoil-python/README.md b/packages/flexfoil-python/README.md index 7727bafb..d37162e7 100644 --- a/packages/flexfoil-python/README.md +++ b/packages/flexfoil-python/README.md @@ -131,6 +131,20 @@ for defl in [0, 5, 10, 15]: polar = f.polar(alpha=(-4, 14, 1.0), Re=1e6) ``` +### Forced transition (trip strips) + +Force the boundary layer to transition at a specific chord location, +equivalent to XFOIL's `XSTRIP` (VPAR menu). Default `1.0` = free transition. + +```python +# Trip strips at 5% chord on both surfaces +result = foil.solve(alpha=5.0, Re=1e6, xstrip_upper=0.05, xstrip_lower=0.05) +print(result.x_tr_upper) # ≈ 0.05 + +# Force upper only — lower uses natural (e^N) transition +polar = foil.polar(alpha=(-5, 15, 0.5), Re=1e6, xstrip_upper=0.1) +``` + ### Inviscid analysis ```python @@ -201,8 +215,8 @@ flexfoil info # show config and DB location | `.panel_coords` | Repaneled coordinates `list[(x, y)]` | | `.hash` | SHA-256 hash of panel coords (cache key) | | `.with_flap(hinge_x, deflection, hinge_y_frac, n_panels)` | Return new Airfoil with flap deflected | -| `.solve(alpha, Re, mach, ncrit, max_iter, viscous, store)` | Single-point analysis | -| `.polar(alpha, Re, mach, ncrit, max_iter, viscous, store, parallel)` | Sweep over alpha range (parallel by default) | +| `.solve(alpha, Re, mach, ncrit, max_iter, viscous, store, re_type, xstrip_upper, xstrip_lower)` | Single-point analysis | +| `.polar(alpha, Re, mach, ncrit, max_iter, viscous, store, parallel, re_type, xstrip_upper, xstrip_lower)` | Sweep over alpha range (parallel by default) | ### `SolveResult` diff --git a/packages/flexfoil-python/src/flexfoil/_rustfoil.pyi b/packages/flexfoil-python/src/flexfoil/_rustfoil.pyi index dba89f37..ccb3821f 100644 --- a/packages/flexfoil-python/src/flexfoil/_rustfoil.pyi +++ b/packages/flexfoil-python/src/flexfoil/_rustfoil.pyi @@ -8,6 +8,8 @@ def analyze_faithful( ncrit: float = 9.0, max_iterations: int = 100, re_type: int = 1, + xstrip_upper: float = 1.0, + xstrip_lower: float = 1.0, ) -> dict[str, object]: """Viscous (XFOIL-faithful) analysis at a single operating point.""" ... @@ -53,6 +55,8 @@ def analyze_faithful_batch( ncrit: float = 9.0, max_iterations: int = 100, re_type: int = 1, + xstrip_upper: float = 1.0, + xstrip_lower: float = 1.0, ) -> list[dict[str, object]]: """Batch viscous analysis: solve multiple alphas in parallel via rayon.""" ... @@ -76,6 +80,8 @@ def get_bl_distribution( ncrit: float = 9.0, max_iterations: int = 100, re_type: int = 1, + xstrip_upper: float = 1.0, + xstrip_lower: float = 1.0, ) -> dict[str, object]: """Compute boundary-layer distributions for a viscous operating point.""" ... diff --git a/packages/flexfoil-python/src/flexfoil/airfoil.py b/packages/flexfoil-python/src/flexfoil/airfoil.py index fc2ef407..10a5b1b2 100644 --- a/packages/flexfoil-python/src/flexfoil/airfoil.py +++ b/packages/flexfoil-python/src/flexfoil/airfoil.py @@ -208,6 +208,8 @@ def solve( viscous: bool = True, store: bool = True, re_type: int = 1, + xstrip_upper: float = 1.0, + xstrip_lower: float = 1.0, ) -> SolveResult: """Run aerodynamic analysis at a single operating point. @@ -220,11 +222,13 @@ def solve( max_iter : maximum Newton iterations viscous : if False, run inviscid panel method only store : if True, cache result in the local database + xstrip_upper : forced transition x/c on upper surface (1.0 = free) + xstrip_lower : forced transition x/c on lower surface (1.0 = free) """ coords = self._flat_panels() if viscous: - raw = analyze_faithful(coords, alpha, Re, mach, ncrit, max_iter, re_type) + raw = analyze_faithful(coords, alpha, Re, mach, ncrit, max_iter, re_type, xstrip_upper, xstrip_lower) # For viscous, also get Cp from inviscid pass raw_inv = analyze_inviscid(coords, alpha) result = SolveResult( @@ -284,6 +288,8 @@ def polar( store: bool = True, parallel: bool = True, re_type: int = 1, + xstrip_upper: float = 1.0, + xstrip_lower: float = 1.0, ) -> PolarResult | list[PolarResult]: """Run a polar sweep over angles of attack, optionally with matrix sweeps. @@ -297,6 +303,8 @@ def polar( viscous : if False, inviscid only store : if True, cache each point in the local database parallel : if True (default), solve all alphas in parallel via rayon + xstrip_upper : forced transition x/c on upper surface (1.0 = free) + xstrip_lower : forced transition x/c on lower surface (1.0 = free) Returns ------- @@ -329,6 +337,7 @@ def polar( Re=re_val, mach=mach_val, ncrit=ncrit_val, max_iter=max_iter, viscous=viscous, store=store, re_type=re_type, + xstrip_upper=xstrip_upper, xstrip_lower=xstrip_lower, ) else: results = [ @@ -336,6 +345,7 @@ def polar( a, Re=re_val, mach=mach_val, ncrit=ncrit_val, max_iter=max_iter, viscous=viscous, store=store, re_type=re_type, + xstrip_upper=xstrip_upper, xstrip_lower=xstrip_lower, ) for a in alphas ] @@ -361,6 +371,8 @@ def _polar_batch( viscous: bool, store: bool, re_type: int = 1, + xstrip_upper: float = 1.0, + xstrip_lower: float = 1.0, ) -> list[SolveResult]: """Solve all alphas in a single parallel Rust call.""" coords = self._flat_panels() @@ -368,7 +380,7 @@ def _polar_batch( if viscous: from flexfoil._rustfoil import analyze_faithful_batch - raw_list = analyze_faithful_batch(coords, alphas, Re, mach, ncrit, max_iter, re_type) + raw_list = analyze_faithful_batch(coords, alphas, Re, mach, ncrit, max_iter, re_type, xstrip_upper, xstrip_lower) results = [ SolveResult( cl=raw.get("cl", 0.0), @@ -429,6 +441,8 @@ def bl_distribution( ncrit: float = 9.0, max_iter: int = 100, re_type: int = 1, + xstrip_upper: float = 1.0, + xstrip_lower: float = 1.0, ) -> BLResult: """Compute boundary-layer distributions (viscous only). @@ -436,7 +450,7 @@ def bl_distribution( theta, H, and ue for upper and lower surfaces. """ coords = self._flat_panels() - raw = get_bl_distribution(coords, alpha, Re, mach, ncrit, max_iter, re_type) + raw = get_bl_distribution(coords, alpha, Re, mach, ncrit, max_iter, re_type, xstrip_upper, xstrip_lower) if not raw.get("success", False): return BLResult( diff --git a/packages/flexfoil-python/tests/test_api.py b/packages/flexfoil-python/tests/test_api.py index 78c3f7ef..bb9dca1a 100644 --- a/packages/flexfoil-python/tests/test_api.py +++ b/packages/flexfoil-python/tests/test_api.py @@ -406,6 +406,71 @@ def test_bl_distribution_type2(self): assert bl.converged +# --------------------------------------------------------------------------- +# Forced transition (XSTRIP) +# --------------------------------------------------------------------------- + +class TestForcedTransition: + def test_solve_forced_upper(self): + import flexfoil + foil = flexfoil.naca("2412") + r = foil.solve(5.0, Re=1e6, xstrip_upper=0.3, store=False) + assert r.success + assert r.x_tr_upper < 0.35, f"Expected x_tr_upper near 0.3, got {r.x_tr_upper}" + + def test_solve_forced_lower(self): + import flexfoil + foil = flexfoil.naca("2412") + r = foil.solve(5.0, Re=1e6, xstrip_lower=0.4, store=False) + assert r.success + assert r.x_tr_lower < 0.45, f"Expected x_tr_lower near 0.4, got {r.x_tr_lower}" + + def test_forced_moves_transition_earlier(self): + import flexfoil + foil = flexfoil.naca("2412") + r_free = foil.solve(5.0, Re=1e6, store=False) + r_forced = foil.solve(5.0, Re=1e6, xstrip_upper=0.1, store=False) + assert r_free.success and r_forced.success + assert r_forced.x_tr_upper < r_free.x_tr_upper, \ + f"Forced transition ({r_forced.x_tr_upper}) should be earlier than free ({r_free.x_tr_upper})" + + def test_default_xstrip_no_regression(self): + import flexfoil + foil = flexfoil.naca("2412") + r1 = foil.solve(5.0, Re=1e6, store=False) + r2 = foil.solve(5.0, Re=1e6, xstrip_upper=1.0, xstrip_lower=1.0, store=False) + assert abs(r1.cl - r2.cl) < 1e-10 + assert abs(r1.cd - r2.cd) < 1e-10 + + def test_polar_forced_transition(self): + import flexfoil + foil = flexfoil.naca("2412") + polar = foil.polar(alpha=[3.0, 5.0], Re=1e6, xstrip_upper=0.2, store=False) + for r in polar.converged: + assert r.x_tr_upper < 0.25, f"Expected forced xtr < 0.25, got {r.x_tr_upper}" + + def test_batch_forced_matches_sequential(self): + import flexfoil + foil = flexfoil.naca("2412") + par = foil.polar(alpha=[3.0, 5.0], Re=1e6, xstrip_upper=0.3, + parallel=True, store=False) + seq = foil.polar(alpha=[3.0, 5.0], Re=1e6, xstrip_upper=0.3, + parallel=False, store=False) + for p, s in zip(par.converged, seq.converged): + assert abs(p.cl - s.cl) < 1e-10 + assert abs(p.x_tr_upper - s.x_tr_upper) < 1e-10 + + def test_forced_changes_cd(self): + """Forcing transition earlier should increase drag.""" + import flexfoil + foil = flexfoil.naca("2412") + r_free = foil.solve(5.0, Re=1e6, store=False) + r_forced = foil.solve(5.0, Re=1e6, xstrip_upper=0.05, xstrip_lower=0.05, store=False) + assert r_free.success and r_forced.success + assert r_forced.cd > r_free.cd, \ + f"Forced early transition cd ({r_forced.cd}) should exceed free ({r_free.cd})" + + # --------------------------------------------------------------------------- # Matrix sweep (multi-Re) # --------------------------------------------------------------------------- diff --git a/packages/flexfoil-python/tests/test_solver.py b/packages/flexfoil-python/tests/test_solver.py index b5990977..12aa60d9 100644 --- a/packages/flexfoil-python/tests/test_solver.py +++ b/packages/flexfoil-python/tests/test_solver.py @@ -461,6 +461,54 @@ def test_bl_distribution_type2(self): assert result["converged"] +class TestForcedTransition: + """Tests for forced transition (XSTRIP) at the Rust binding level.""" + + def _make_coords(self): + raw = generate_naca4(2412) + flat = [v for x, y in raw for v in (x, y)] + paneled = repanel_xfoil(flat, 160) + return [v for x, y in paneled for v in (x, y)] + + def test_forced_upper_moves_transition(self): + coords = self._make_coords() + r_free = analyze_faithful(coords, 5.0, 1e6, 0.0, 9.0, 100) + r_forced = analyze_faithful(coords, 5.0, 1e6, 0.0, 9.0, 100, 1, 0.1, 1.0) + assert r_free["success"] and r_forced["success"] + assert r_forced["x_tr_upper"] < r_free["x_tr_upper"] + + def test_forced_lower_moves_transition(self): + coords = self._make_coords() + r_free = analyze_faithful(coords, 5.0, 1e6, 0.0, 9.0, 100) + r_forced = analyze_faithful(coords, 5.0, 1e6, 0.0, 9.0, 100, 1, 1.0, 0.2) + assert r_free["success"] and r_forced["success"] + assert r_forced["x_tr_lower"] < r_free["x_tr_lower"] + + def test_default_xstrip_unchanged(self): + coords = self._make_coords() + r1 = analyze_faithful(coords, 5.0, 1e6, 0.0, 9.0, 100) + r2 = analyze_faithful(coords, 5.0, 1e6, 0.0, 9.0, 100, 1, 1.0, 1.0) + assert abs(r1["cl"] - r2["cl"]) < 1e-10 + assert abs(r1["cd"] - r2["cd"]) < 1e-10 + + def test_forced_transition_increases_drag(self): + coords = self._make_coords() + r_free = analyze_faithful(coords, 5.0, 1e6, 0.0, 9.0, 100) + r_forced = analyze_faithful(coords, 5.0, 1e6, 0.0, 9.0, 100, 1, 0.05, 0.05) + assert r_free["success"] and r_forced["success"] + assert r_forced["cd"] > r_free["cd"] + + def test_batch_forced_matches_single(self): + from flexfoil._rustfoil import analyze_faithful_batch + coords = self._make_coords() + alphas = [3.0, 5.0] + batch = analyze_faithful_batch(coords, alphas, 1e6, 0.0, 9.0, 100, 1, 0.3, 1.0) + singles = [analyze_faithful(coords, a, 1e6, 0.0, 9.0, 100, 1, 0.3, 1.0) for a in alphas] + for b, s in zip(batch, singles): + assert abs(b["cl"] - s["cl"]) < 1e-10 + assert abs(b["x_tr_upper"] - s["x_tr_upper"]) < 1e-10 + + class TestInviscidBatch: def _make_coords(self): raw = generate_naca4(2412) From 9f893d2640db390b487cee925a8d327d63936306 Mon Sep 17 00:00:00 2001 From: aeronauty Date: Wed, 8 Apr 2026 21:50:21 +0100 Subject: [PATCH 2/2] chore: bump versions to 1.1.4/1.1.6 and add What's New for forced transition MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit UI: 1.1.3 → 1.1.4 Python: 1.1.5 → 1.1.6 Merge to main triggers pypi-publish.yml which builds wheels for all platforms and publishes to PyPI automatically. Co-Authored-By: Claude Opus 4.6 (1M context) --- flexfoil-ui/package.json | 2 +- flexfoil-ui/src/lib/version.ts | 23 +++++++++++++++++++++++ packages/flexfoil-python/pyproject.toml | 2 +- 3 files changed, 25 insertions(+), 2 deletions(-) diff --git a/flexfoil-ui/package.json b/flexfoil-ui/package.json index a9be3116..6fea9a6a 100644 --- a/flexfoil-ui/package.json +++ b/flexfoil-ui/package.json @@ -2,7 +2,7 @@ "name": "flexfoil-ui", "private": true, "license": "MIT", - "version": "1.1.3", + "version": "1.1.4", "type": "module", "scripts": { "dev": "vite", diff --git a/flexfoil-ui/src/lib/version.ts b/flexfoil-ui/src/lib/version.ts index a38abef4..a989e272 100644 --- a/flexfoil-ui/src/lib/version.ts +++ b/flexfoil-ui/src/lib/version.ts @@ -28,6 +28,29 @@ export interface ChangelogEntry { } export const CHANGELOG: ChangelogEntry[] = [ + { + version: '1.1.4', + date: '2026-04-08', + items: [ + { category: 'added', text: 'Forced transition (trip strips) — set Trip Upper and Trip Lower in the Solve panel to force laminar-to-turbulent transition at a specific chord location (XFOIL\'s XSTRIP equivalent)' }, + { category: 'added', text: 'Python API: xstrip_upper and xstrip_lower parameters on solve(), polar(), and bl_distribution() — default 1.0 (free transition)' }, + { category: 'changed', text: 'If natural e^N transition occurs before the forced location, natural transition wins — matching XFOIL behaviour' }, + ], + tourSlides: [ + { + title: 'Forced Transition', + description: 'Simulate trip strips by forcing boundary layer transition at a specific chord location', + icon: '📍', + items: [ + 'Trip Upper / Trip Lower: set x/c position (0–1) where transition is forced', + 'Default 1.0 = free transition (e^N method only)', + 'If natural transition occurs earlier, it takes precedence', + 'Equivalent to XFOIL\'s XSTRIP command in the VPAR menu', + ], + goTo: { panel: 'solve', label: 'Open Solve panel' }, + }, + ], + }, { version: '1.1.3', date: '2026-03-23', diff --git a/packages/flexfoil-python/pyproject.toml b/packages/flexfoil-python/pyproject.toml index f112a1e1..7b57495c 100644 --- a/packages/flexfoil-python/pyproject.toml +++ b/packages/flexfoil-python/pyproject.toml @@ -4,7 +4,7 @@ build-backend = "maturin" [project] name = "flexfoil" -version = "1.1.5" +version = "1.1.6" description = "Airfoil analysis in Python — XFOIL-faithful solver with a local web UI" readme = "README.md" license = { text = "MIT" }