Skip to main content

hpr_sim/
issues.rs

1//! Warnings for known errors in hpr's drag, each by its issue number: where a flight's speed and
2//! its design's shape meet the condition under which a number is known to read wrong
3//! ([ADR-163 §4][adr-163], [ADR-180][adr-180]).
4//!
5//! Drag that reads high makes three numbers flatter than they are: the apogee reads low against
6//! a waiver's ceiling; the top speed reads low, so the flutter margin and the largest dynamic
7//! pressure look better; the drift reads small. Each warning names them. Every flight still
8//! flies; a warning only says which numbers to trust less, and which way they lean. A flight
9//! flown on a drag of its own (a table or a model in place of hpr's) meets none of these.
10//!
11//! | issue | the error | the condition |
12//! |---|---|---|
13//! | [#67][i67] | the closed-form cone's pressure drag reads high, transonic | a nose or a shoulder whose shape takes it ([`takes_cone_formula`]), past Mach 0.8 |
14//! | [#68][i68] | base drag reads high from Mach 0.8 to 1.2, low from 1.5 | every rocket, past Mach 0.8 |
15//! | [#70][i70] | a sharp fin section's drag reads high, supersonic | an `Airfoil` fin set, past Mach 1 |
16//! | [#72][i72] | a steep boattail's drag reads high, supersonic | a boattail steeper than 10°, past Mach 1 |
17//! | [#222][i222] | supersonic pressure drag about twice OpenRocket's | an `Ogive` nose or an `Airfoil` fin set, past Mach 1 |
18//! | [#73][i73] | a boattail's drag by Niskanen's subsonic rule, up to +41.7% on a measured one | a boattail steeper than 9.46° ([`SUBSONIC_BOATTAIL_RULE_RAD`]), any flight |
19//! | [#18][i18] | skin friction taken as fully turbulent reads high on a surface smooth enough for a laminar run | every flight ([`friction_issue_warnings`]) |
20//!
21//! Others warn where the margin reads high ([`stability_issue_warnings`],
22//! [`layout_stability_issue_warnings`]):
23//! a forward-swept fin past Mach 1 ([#64][i64]), fin sets sharing a station ([#325][i325]), a
24//! freeform fin set at any speed ([#326][i326], [ADR-190][adr-190]), the
25//! supersonic switches (#87, #120, #121) and #172 on every flight ([ADR-189][adr-189]).
26//!
27//! Two more warn where a part flies on its own after a separation, with only its devices' drag
28//! ([`separated_part_issue_warnings`]): a booster ([#179][i179]), and a part that starts with
29//! nothing open ([#354][i354]). They read the flight's bodies, not its drag, so a flight on a drag
30//! of its own meets them too.
31//!
32//! Each Mach edge is exclusive, as the envelope's are ([`crate::envelope`]): a flight whose top
33//! Mach number is exactly at an edge raises nothing, and a NaN top Mach number raises every
34//! warning whose shape the design has, so a broken number can't hide one.
35//!
36//! [adr-163]: https://github.com/nrdptel/hpr-sim/blob/main/docs/decisions/0163-the-second-pass-products-and-priorities.md
37//! [adr-180]: https://github.com/nrdptel/hpr-sim/blob/main/docs/decisions/0180-the-drag-issue-warnings.md
38//! [i67]: https://github.com/nrdptel/hpr-sim/issues/67
39//! [i68]: https://github.com/nrdptel/hpr-sim/issues/68
40//! [i70]: https://github.com/nrdptel/hpr-sim/issues/70
41//! [i72]: https://github.com/nrdptel/hpr-sim/issues/72
42//! [i222]: https://github.com/nrdptel/hpr-sim/issues/222
43//! [i64]: https://github.com/nrdptel/hpr-sim/issues/64
44//! [i73]: https://github.com/nrdptel/hpr-sim/issues/73
45//! [i325]: https://github.com/nrdptel/hpr-sim/issues/325
46//! [i18]: https://github.com/nrdptel/hpr-sim/issues/18
47//! [i326]: https://github.com/nrdptel/hpr-sim/issues/326
48//! [i179]: https://github.com/nrdptel/hpr-sim/issues/179
49//! [i354]: https://github.com/nrdptel/hpr-sim/issues/354
50//! [adr-190]: https://github.com/nrdptel/hpr-sim/blob/main/docs/decisions/0190-the-friction-and-freeform-fin-warnings.md
51//! [adr-189]: https://github.com/nrdptel/hpr-sim/blob/main/docs/decisions/0189-the-safety-and-shape-warnings.md
52
53use hpr_aero::nose_drag::takes_cone_formula;
54use hpr_design::{FinCrossSection, FinPlanform, Layout, NoseShape, Part, PlacedComponent};
55use serde::{Deserialize, Serialize};
56
57use crate::FlightResult;
58use crate::metrics::Peak;
59
60/// Where the closed-form cone's transonic pressure drag starts to read high: Mach 0.8
61/// ([issue #67](https://github.com/nrdptel/hpr-sim/issues/67): +87% at Mach 0.8 and +105% at
62/// 0.85 against Stoney's measured 3:1 cone, NASA TR R-100; 2 to 9 times MIL-HDBK-762's 3:1
63/// tangent ogive from Mach 0.9 to 1.2; still +15% at Mach 1.5).
64pub const TRANSONIC_NOSE_DRAG_MACH: f64 = 0.8;
65
66/// Where base drag starts to read high: Mach 0.8 ([issue #68](https://github.com/nrdptel/hpr-sim/issues/68):
67/// Fleeman's formula gives 0.225 at Mach 0.9 where MIL-HDBK-762's measured data give 0.156).
68pub const BASE_DRAG_HIGH_MACH: f64 = 0.8;
69
70/// The last Mach number where base drag is measured to read high: 1.2 (0.208 against
71/// MIL-HDBK-762's 0.194, [issue #68](https://github.com/nrdptel/hpr-sim/issues/68)). Past it a
72/// flight's warning adds the low side ([`BASE_DRAG_LOW_MACH`]); between the two the sign is
73/// unmeasured.
74pub const BASE_DRAG_HIGH_TOP_MACH: f64 = 1.2;
75
76/// Where base drag is measured to read low: from Mach 1.5 (5% under Love's correlation, NACA
77/// TN 3819, at Mach 1.5, 17% at 2.0, 14% at 3.0; [issue #68](https://github.com/nrdptel/hpr-sim/issues/68)).
78/// There the apogee reads high.
79pub const BASE_DRAG_LOW_MACH: f64 = 1.5;
80
81/// Where the supersonic drag warnings start: Mach 1 (issues
82/// [#70](https://github.com/nrdptel/hpr-sim/issues/70),
83/// [#72](https://github.com/nrdptel/hpr-sim/issues/72) and
84/// [#222](https://github.com/nrdptel/hpr-sim/issues/222)). Their measured gaps start at Mach 1.0
85/// to 1.5; the warning starts at the lowest.
86pub const SUPERSONIC_DRAG_MACH: f64 = 1.0;
87
88/// A boattail steeper than this, in radians, warns past [`SUPERSONIC_DRAG_MACH`]: 10°, the
89/// steepest that [issue #72](https://github.com/nrdptel/hpr-sim/issues/72) measured inside its
90/// spread (−6.3% to +17.4% at 10° and below). It measured 15° at +13.5% to +50.8% and 16° at
91/// +26.4% to +54.1%; between 10° and 15° is unmeasured, so the warning takes it in rather than
92/// leave a high drag unwarned.
93pub const STEEP_BOATTAIL_RAD: f64 = 10.0_f64.to_radians();
94
95/// Where the supersonic body switches start to matter: Mach 1.2, the earliest the
96/// shock-expansion method joins slender-body theory ([`hpr_aero::SUPERSONIC_JOIN_START_MACH`]).
97pub const SUPERSONIC_SWITCH_MACH: f64 = hpr_aero::SUPERSONIC_JOIN_START_MACH;
98
99/// How far the static margin read high against OpenRocket 24.12, calibres: the largest gap on
100/// the private designs, 0.1108 (four designs read 0.0350 to 0.1108 calibres more stable in hpr,
101/// none with a measured cause; *How far to trust the margin* on
102/// [the stability page](https://nrdptel.github.io/hpr-sim/stability-for-certification.html#how-far-to-trust-the-margin),
103/// [issue #172](https://github.com/nrdptel/hpr-sim/issues/172)).
104pub const MARGIN_HIGH_BOUND_CAL: f64 = 0.1108;
105
106/// Where a forward-swept fin's warning starts: Mach 1 ([ADR-188][adr-188]). Its measured error is
107/// a normal-force slope that rises by up to 7.7% past linear theory's start, Mach 1.2
108/// ([`hpr_aero::fins::SUPERSONIC_START_MACH`]), peaking up to 0.37 past it
109/// ([issue #64](https://github.com/nrdptel/hpr-sim/issues/64)); the warning starts at the lowest
110/// supersonic Mach number, below that rise. The transonic join into linear theory starts at Mach
111/// 0.8 ([`hpr_aero::fins::TRANSONIC_START_MACH`]); below Mach 1 the error isn't measured.
112///
113/// [adr-188]: https://github.com/nrdptel/hpr-sim/blob/main/docs/decisions/0188-release-notes-and-the-flattering-side-survey.md
114pub const FORWARD_SWEEP_MACH: f64 = 1.0;
115
116/// The steepest boattail that Niskanen's subsonic rule gives no drag, rad: where its length
117/// ratio `γ = l/(d₁ − d₂)` is 3, so `atan(1/6)`, 9.46° (Niskanen 2009 eq. 3.88,
118/// [`hpr_aero::drag::boattail_factor`]). A steeper one takes a share of base drag below Mach 0.8,
119/// which read from −40.8% to +41.7% on Cubbage's measured 16° boattails and high on the Arcas
120/// Robin's 15° one ([issue #73](https://github.com/nrdptel/hpr-sim/issues/73)); between 9.46° and
121/// 15° is unmeasured, so the warning takes it in.
122pub const SUBSONIC_BOATTAIL_RULE_RAD: f64 = 0.165_148_677_414_626_83;
123
124/// The most fins that one station can hold before the fin–fin interference factor drops below
125/// 1: four (Niskanen 2009 table 3.3, [`hpr_aero::fins::fin_count_factor`]). Fin sets at one station
126/// with more than this between them are counted apart by hpr and together by OpenRocket
127/// ([issue #325](https://github.com/nrdptel/hpr-sim/issues/325)).
128pub const FINS_WITHOUT_INTERFERENCE: u32 = 4;
129
130/// How much of the margin gap on OpenRocket's *Pods--airframes and winglets* example counting
131/// fin sets apart explains, calibres: about 0.029 of 0.071 to 0.076
132/// ([issue #325](https://github.com/nrdptel/hpr-sim/issues/325)).
133pub const FIN_SETS_AT_ONE_STATION_CAL: f64 = 0.029;
134
135/// How much of the friction a laminar run takes off on RocketPy's Calisto at Mach 0.3, where
136/// hpr takes it fully turbulent ([issue #18](https://github.com/nrdptel/hpr-sim/issues/18)):
137/// Barrowman 1967's transitional term `1700/R` (eq. 4-6) at `R = 1.8e7`, against Niskanen 2009's
138/// smooth turbulent coefficient `1/(1.50 ln R − 5.6)²` (eq. 3.78), in percent. Barrowman gives no
139/// rule for when a surface is too rough to stay laminar, so every flight on hpr's drag warns.
140pub const LAMINAR_FRICTION_PERCENT: f64 = 3.6;
141
142/// The share of the static margin, calibres, that a kinked freeform fin's center of pressure
143/// takes on OpenRocket's *Pods--airframes and winglets* example, where hpr's sits 1.6 mm aft of
144/// OpenRocket 24.12's ([issue #326](https://github.com/nrdptel/hpr-sim/issues/326)). Unprobed on
145/// other outlines, so every freeform fin set warns.
146pub const FREEFORM_FIN_MARGIN_CAL: f64 = 0.047;
147
148/// A known error in hpr's drag or stability that a flight can meet.
149#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
150#[serde(rename_all = "snake_case")]
151#[non_exhaustive]
152pub enum KnownIssue {
153    /// [#67](https://github.com/nrdptel/hpr-sim/issues/67): the closed-form cone's pressure
154    /// drag, which cones, ogives and some power and parabolic series take, reads high from
155    /// Mach 0.8.
156    TransonicNoseDrag,
157    /// [#68](https://github.com/nrdptel/hpr-sim/issues/68): base drag reads high from Mach 0.8
158    /// to 1.2 and low from 1.5.
159    BaseDrag,
160    /// [#70](https://github.com/nrdptel/hpr-sim/issues/70): a sharp fin section's drag reads
161    /// high faster than sound. hpr has no sharp section, so the warning goes with the airfoil.
162    SharpFinDrag,
163    /// [#72](https://github.com/nrdptel/hpr-sim/issues/72): a steep boattail's drag reads high
164    /// faster than sound.
165    SteepBoattailDrag,
166    /// [#222](https://github.com/nrdptel/hpr-sim/issues/222): supersonic pressure drag about twice
167    /// OpenRocket's on an ogive nose and airfoil fins.
168    SupersonicPressureDrag,
169    /// [#87](https://github.com/nrdptel/hpr-sim/issues/87): a step in radius takes the whole
170    /// body off the shock-expansion method faster than Mach 1.2.
171    RadiusStepFallback,
172    /// [#120](https://github.com/nrdptel/hpr-sim/issues/120): a lip longer than its boattail's
173    /// drop in diameter does the same.
174    LongLipFallback,
175    /// [#121](https://github.com/nrdptel/hpr-sim/issues/121): a pointed tip steeper than the cone
176    /// tables' 30° does the same.
177    SteepTipFallback,
178    /// [#172](https://github.com/nrdptel/hpr-sim/issues/172): the static margin read up to
179    /// 0.1108 calibres higher than OpenRocket's on four private designs, cause unknown.
180    MarginReadsHigh,
181    /// [#64](https://github.com/nrdptel/hpr-sim/issues/64): a fin whose leading edge sweeps
182    /// forward takes a supersonic normal-force slope that rises past linear theory's start,
183    /// moving the center of pressure aft.
184    ForwardSweptFinSlope,
185    /// [#73](https://github.com/nrdptel/hpr-sim/issues/73): below Mach 0.8, Niskanen's rule gives
186    /// a boattail steeper than 9.46° a share of base drag that read up to +41.7% on a measured
187    /// one.
188    SubsonicBoattailDrag,
189    /// [#325](https://github.com/nrdptel/hpr-sim/issues/325): fin sets at one station are counted
190    /// apart for fin–fin interference, where OpenRocket counts them together.
191    FinSetsAtOneStation,
192    /// [#18](https://github.com/nrdptel/hpr-sim/issues/18): skin friction is taken as fully
193    /// turbulent, so on a surface smooth enough to keep its boundary layer laminar for a way it
194    /// reads high.
195    TurbulentFriction,
196    /// [#326](https://github.com/nrdptel/hpr-sim/issues/326): a kinked freeform fin's center of
197    /// pressure sat 1.6 mm aft of OpenRocket's; other freeform outlines are unprobed.
198    FreeformFinCenterOfPressure,
199    /// [#179](https://github.com/nrdptel/hpr-sim/issues/179): a booster dropped at a separation
200    /// flies on as a point with only its devices' drag, tumbling side-on from the split where a
201    /// real one first coasts nose-first on its airframe's drag.
202    BoosterAirframeDrag,
203    /// [#354](https://github.com/nrdptel/hpr-sim/issues/354): a separated part has no drag at all
204    /// until its first device opens, where a real one has its airframe's.
205    DragFreeSeparatedPart,
206}
207
208impl KnownIssue {
209    /// The issue's number on GitHub.
210    #[must_use]
211    pub fn number(self) -> u32 {
212        match self {
213            Self::TransonicNoseDrag => 67,
214            Self::BaseDrag => 68,
215            Self::SharpFinDrag => 70,
216            Self::SteepBoattailDrag => 72,
217            Self::SupersonicPressureDrag => 222,
218            Self::RadiusStepFallback => 87,
219            Self::LongLipFallback => 120,
220            Self::SteepTipFallback => 121,
221            Self::MarginReadsHigh => 172,
222            Self::ForwardSweptFinSlope => 64,
223            Self::SubsonicBoattailDrag => 73,
224            Self::FinSetsAtOneStation => 325,
225            Self::TurbulentFriction => 18,
226            Self::FreeformFinCenterOfPressure => 326,
227            Self::BoosterAirframeDrag => 179,
228            Self::DragFreeSeparatedPart => 354,
229        }
230    }
231
232    /// Whether the issue is about the margin and center of pressure (#64, #87, #120, #121, #172,
233    /// #325, #326) rather than drag.
234    #[must_use]
235    pub fn is_stability(self) -> bool {
236        matches!(
237            self,
238            Self::RadiusStepFallback
239                | Self::LongLipFallback
240                | Self::SteepTipFallback
241                | Self::MarginReadsHigh
242                | Self::ForwardSweptFinSlope
243                | Self::FinSetsAtOneStation
244                | Self::FreeformFinCenterOfPressure
245        )
246    }
247
248    /// Whether the issue has a Mach condition, so its message quotes the flight's top speed:
249    /// every one but #73 (any flight is subsonic for a time), #18, #172, #325 and #326, which
250    /// hold at any speed, and #179 and #354, which are a separated part's.
251    #[must_use]
252    pub fn has_mach_condition(self) -> bool {
253        !matches!(
254            self,
255            Self::MarginReadsHigh
256                | Self::SubsonicBoattailDrag
257                | Self::FinSetsAtOneStation
258                | Self::TurbulentFriction
259                | Self::FreeformFinCenterOfPressure
260                | Self::BoosterAirframeDrag
261                | Self::DragFreeSeparatedPart
262        )
263    }
264
265    /// The issue's address.
266    #[must_use]
267    pub fn url(self) -> String {
268        format!(
269            "https://github.com/nrdptel/hpr-sim/issues/{}",
270            self.number()
271        )
272    }
273}
274
275/// A known issue a flight meets, with the flight's top Mach number and the parts whose shape
276/// meets its condition.
277#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
278#[non_exhaustive]
279pub struct IssueWarning {
280    /// Which issue.
281    pub issue: KnownIssue,
282    /// The flight's top Mach number.
283    pub max_mach: Peak,
284    /// The ids of the parts whose shape meets the issue's condition, in the layout's order;
285    /// empty for base drag ([`KnownIssue::BaseDrag`]), which every rocket has.
286    pub parts: Vec<String>,
287}
288
289impl IssueWarning {
290    /// The issue's number on GitHub.
291    #[must_use]
292    pub fn number(&self) -> u32 {
293        self.issue.number()
294    }
295
296    /// What the issue means for this flight, in a sentence or two, with which way its numbers
297    /// lean.
298    #[must_use]
299    pub fn message(&self) -> String {
300        let mach = self.max_mach.value;
301        let low = "so the apogee, the top speed and the drift read low, and the flutter margin \
302                   and the largest dynamic pressure look better than they are";
303        let what = match self.issue {
304            KnownIssue::TransonicNoseDrag => format!(
305                "a cone-like nose's or shoulder's pressure drag reads high from Mach \
306                 {TRANSONIC_NOSE_DRAG_MACH} (about twice a measured cone's at Mach 0.85, still \
307                 +15% at 1.5), {low}"
308            ),
309            KnownIssue::BaseDrag => {
310                let mut what = format!(
311                    "base drag reads high from Mach {BASE_DRAG_HIGH_MACH} to \
312                     {BASE_DRAG_HIGH_TOP_MACH} (0.225 against a measured 0.156 at Mach 0.9), \
313                     {low}"
314                );
315                if past(mach, BASE_DRAG_HIGH_TOP_MACH) {
316                    what.push_str(&format!(
317                        "; from Mach {BASE_DRAG_LOW_MACH} it reads low (17% at Mach 2), so the \
318                         apogee and the top speed there read high"
319                    ));
320                }
321                what
322            }
323            KnownIssue::SharpFinDrag => format!(
324                "hpr has no sharp-edged (double-wedge) fin section, and drawn as an airfoil such \
325                 a fin's drag reads high faster than sound: if these fins are sharp, {low}"
326            ),
327            KnownIssue::SteepBoattailDrag => format!(
328                "a boattail steeper than {:.0}° drags high faster than sound (+26% to +54% at \
329                 16°), {low}",
330                STEEP_BOATTAIL_RAD.to_degrees()
331            ),
332            KnownIssue::SupersonicPressureDrag => format!(
333                "supersonic pressure drag on an ogive nose or airfoil fins is about twice \
334                 OpenRocket's, and which is right is unresolved: if hpr's is high, {low}"
335            ),
336            KnownIssue::RadiusStepFallback
337            | KnownIssue::LongLipFallback
338            | KnownIssue::SteepTipFallback => {
339                let switch = match self.issue {
340                    KnownIssue::RadiusStepFallback => "a step in radius",
341                    KnownIssue::LongLipFallback => {
342                        "a lip longer than its boattail's drop in diameter"
343                    }
344                    _ => "a pointed tip steeper than the cone tables' 30°",
345                };
346                format!(
347                    "{switch} takes the whole body off the shock-expansion method from Mach \
348                     {SUPERSONIC_SWITCH_MACH}, so slender-body theory puts the center of pressure \
349                     aft there and the margin reads high"
350                )
351            }
352            KnownIssue::MarginReadsHigh => format!(
353                "the static margin read up to {MARGIN_HIGH_BOUND_CAL} calibres higher than \
354                 OpenRocket's on four private designs, for a reason not yet found, so this \
355                 flight's margin may read high by as much"
356            ),
357            KnownIssue::ForwardSweptFinSlope => format!(
358                "a fin whose leading edge sweeps forward takes a normal-force slope that rises by \
359                 up to 7.7% past Mach {}, where it should fall, so faster than Mach \
360                 {FORWARD_SWEEP_MACH} the center of pressure may sit aft and the margin read high",
361                hpr_aero::fins::SUPERSONIC_START_MACH
362            ),
363            KnownIssue::SubsonicBoattailDrag => format!(
364                "below Mach {} a boattail steeper than {:.1}° takes a share of base drag by a rule \
365                 that read from −40.8% to +41.7% on measured 16° boattails; where its drag reads \
366                 high, {low}",
367                hpr_aero::drag::SUBSONIC_MACH_LIMIT,
368                SUBSONIC_BOATTAIL_RULE_RAD.to_degrees()
369            ),
370            KnownIssue::FinSetsAtOneStation => format!(
371                "fin sets that share a station, more than {FINS_WITHOUT_INTERFERENCE} fins between \
372                 them, are each counted alone for fin–fin interference, where OpenRocket counts \
373                 them together (about {FIN_SETS_AT_ONE_STATION_CAL} calibres of margin on its pods \
374                 example), so the margin may read high"
375            ),
376            KnownIssue::TurbulentFriction => format!(
377                "skin friction is taken as fully turbulent, but on a smooth surface the flow stays \
378                 laminar near the nose, where friction is lower: hpr's reads high by \
379                 {LAMINAR_FRICTION_PERCENT}% on RocketPy's Calisto at Mach 0.3; if this surface \
380                 is that smooth, the drag reads high, {low}"
381            ),
382            KnownIssue::FreeformFinCenterOfPressure => format!(
383                "a kinked freeform fin's center of pressure sat 1.6 mm aft of OpenRocket's \
384                 (about {FREEFORM_FIN_MARGIN_CAL} calibres of margin on its pods example) and \
385                 other freeform outlines are unprobed, so the margin may read high"
386            ),
387            KnownIssue::BoosterAirframeDrag => "a booster dropped at a separation flies on as a \
388                 point with only its devices' drag, where a real one first coasts nose-first on \
389                 its airframe's drag: tumbling side-on from the split, as hpr tumbles a .ork \
390                 file's booster, the booster's peak and drift likely read short (with nothing \
391                 open at the split, see #354); its flight is not validated"
392                .to_owned(),
393            KnownIssue::DragFreeSeparatedPart => "a separated part flies on as a point with only \
394                 its open devices' drag, so it has no drag at all until its first device opens, \
395                 where a real one has its airframe's: falling, it likely lands early and the wind \
396                 carries it less far than it would, so its drift likely reads short; climbing, its \
397                 peak reads high"
398                .to_owned(),
399        };
400        // An issue with no Mach condition doesn't quote the flight's top speed.
401        if !self.issue.has_mach_condition() {
402            return format!("issue #{}: {what} ({})", self.number(), self.issue.url());
403        }
404        format!(
405            "issue #{}: {what}; this flight reaches Mach {mach:.2} at {:.1} s ({})",
406            self.number(),
407            self.max_mach.time_s,
408            self.issue.url()
409        )
410    }
411}
412
413/// Whether `value` is past `edge`, a NaN counting as past it.
414fn past(value: f64, edge: f64) -> bool {
415    value.is_nan() || value > edge
416}
417
418/// The half-angle of a transition that narrows aft, rad: `atan((r_fore − r_aft)/l)`, its mean
419/// angle (a conical boattail's own). `None` for one that doesn't narrow, and for one of no
420/// length, which the drag buildup takes as a step, not a boattail.
421fn boattail_angle_rad(fore_radius_m: f64, aft_radius_m: f64, length_m: f64) -> Option<f64> {
422    let drop_m = fore_radius_m - aft_radius_m;
423    // Written out so a NaN radius or length counts as narrowing: a broken number can't hide one.
424    if drop_m <= 0.0 || length_m == 0.0 {
425        return None;
426    }
427    Some(drop_m.atan2(length_m))
428}
429
430/// The issues a flight meets: one whose top Mach number is `max_mach`, of a rocket whose parts
431/// are `parts`, each an id and its part, flown on hpr's own drag. In [`KnownIssue`]'s order, each
432/// listing the parts that meet its condition. A flight with no top Mach number meets none.
433#[must_use]
434pub fn issue_warnings<'a>(
435    parts: impl IntoIterator<Item = (&'a str, &'a Part)>,
436    max_mach: Option<Peak>,
437) -> Vec<IssueWarning> {
438    let Some(max_mach) = max_mach else {
439        return Vec::new();
440    };
441    let mut cone_like = Vec::new();
442    let mut ogive_or_airfoil = Vec::new();
443    let mut airfoil = Vec::new();
444    let mut steep = Vec::new();
445    let mut ruled = Vec::new();
446    for (id, part) in parts {
447        match part {
448            Part::NoseCone(cone) => {
449                if takes_cone_formula(cone.shape) {
450                    cone_like.push(id.to_owned());
451                }
452                if matches!(cone.shape, NoseShape::Ogive { .. }) {
453                    ogive_or_airfoil.push(id.to_owned());
454                }
455            }
456            Part::FinSet(fins) if fins.cross_section == FinCrossSection::Airfoil => {
457                airfoil.push(id.to_owned());
458                ogive_or_airfoil.push(id.to_owned());
459            }
460            Part::Transition(transition) => {
461                let (fore, aft) = (transition.fore_radius_m, transition.aft_radius_m);
462                // A shoulder (widening aft) drags as a nose does; a NaN radius counts.
463                let widening = aft.is_nan() || fore.is_nan() || aft > fore;
464                if widening && takes_cone_formula(transition.shape) {
465                    cone_like.push(id.to_owned());
466                }
467                let angle = boattail_angle_rad(fore, aft, transition.length_m);
468                if angle.is_some_and(|angle| past(angle, STEEP_BOATTAIL_RAD)) {
469                    steep.push(id.to_owned());
470                }
471                if angle.is_some_and(|angle| past(angle, SUBSONIC_BOATTAIL_RULE_RAD)) {
472                    ruled.push(id.to_owned());
473                }
474            }
475            _ => {}
476        }
477    }
478    let mach = max_mach.value;
479    let supersonic = past(mach, SUPERSONIC_DRAG_MACH);
480    let mut warnings = Vec::new();
481    let mut warn = |issue, parts: Vec<String>| {
482        warnings.push(IssueWarning {
483            issue,
484            max_mach,
485            parts,
486        });
487    };
488    if past(mach, TRANSONIC_NOSE_DRAG_MACH) && !cone_like.is_empty() {
489        warn(KnownIssue::TransonicNoseDrag, cone_like);
490    }
491    if past(mach, BASE_DRAG_HIGH_MACH) {
492        warn(KnownIssue::BaseDrag, Vec::new());
493    }
494    if supersonic && !airfoil.is_empty() {
495        warn(KnownIssue::SharpFinDrag, airfoil);
496    }
497    if supersonic && !steep.is_empty() {
498        warn(KnownIssue::SteepBoattailDrag, steep);
499    }
500    if supersonic && !ogive_or_airfoil.is_empty() {
501        warn(KnownIssue::SupersonicPressureDrag, ogive_or_airfoil);
502    }
503    // Every flight is subsonic for a time, so the rule's boattails warn at any speed.
504    if !ruled.is_empty() {
505        warn(KnownIssue::SubsonicBoattailDrag, ruled);
506    }
507    warnings
508}
509
510/// The friction issue a flight on hpr's own drag meets: [#18][i18] whenever it has a top Mach
511/// number, at any speed and on any finish, since Barrowman 1967 gives no rule for when a surface
512/// is too rough to stay laminar ([`LAMINAR_FRICTION_PERCENT`]). Apart from [`issue_warnings`],
513/// which reads the design's shape.
514///
515/// [i18]: https://github.com/nrdptel/hpr-sim/issues/18
516#[must_use]
517pub fn friction_issue_warnings(max_mach: Option<Peak>) -> Vec<IssueWarning> {
518    max_mach
519        .map(|max_mach| IssueWarning {
520            issue: KnownIssue::TurbulentFriction,
521            max_mach,
522            parts: Vec::new(),
523        })
524        .into_iter()
525        .collect()
526}
527
528/// The separated parts' issues a flight meets ([ADR-200][adr-200]), read from the bodies that
529/// flew on their own ([`FlightResult::bodies`]), whatever the rocket's drag: each part's descent
530/// is a point under its open devices ([`crate::recovery`]), never the rocket's drag.
531///
532/// - [#179][i179] when a body aft of a stage boundary flew on its own: a booster dropped at a
533///   separation, or a parallel stage. Its first stage is past the nose's (stage 0), whatever
534///   devices it carries; the error is unsized, so no speed at the split is excluded.
535/// - [#354][i354] when a body started its own flight with no drag area open (zero, or not a
536///   number, so a broken one can't hide it), and no device opened at that same instant (an event
537///   at its start time, as a device fired by the split with no lag is): it flies with no drag
538///   until its first device opens. A canopy that fills from nothing (a fill time) opened on the
539///   stack at the split also starts at zero, so it warns too: an over-warning, on the safe side.
540///   The flight refuses such a part while it climbs through a fired device's lag, and refuses one
541///   that lands with nothing open ([`crate::recovery`]), so what is left is mostly a part already
542///   falling.
543///
544/// Each warns once, whatever the number of bodies that meet it, and only for a flight with a top
545/// Mach number (`max_mach`), as the others.
546///
547/// [adr-200]: https://github.com/nrdptel/hpr-sim/blob/main/docs/decisions/0200-the-separated-parts-warnings.md
548/// [i179]: https://github.com/nrdptel/hpr-sim/issues/179
549/// [i354]: https://github.com/nrdptel/hpr-sim/issues/354
550#[must_use]
551pub fn separated_part_issue_warnings(
552    result: &FlightResult,
553    max_mach: Option<Peak>,
554) -> Vec<IssueWarning> {
555    let Some(max_mach) = max_mach else {
556        return Vec::new();
557    };
558    let booster = result.bodies.iter().any(|body| body.stages.0 > 0);
559    let drag_free = result.bodies.iter().any(|body| {
560        // A device that opens at the start's own instant is an event at its time, after the
561        // start sample. A NaN area or time is never open.
562        let start = &body.start_sample;
563        let open = |area_m2: f64| area_m2 > 0.0;
564        !open(start.recovery_drag_area_m2)
565            && !body
566                .events
567                .iter()
568                .take_while(|event| event.sample.time_s <= start.time_s)
569                .any(|event| open(event.sample.recovery_drag_area_m2))
570    });
571    [
572        (booster, KnownIssue::BoosterAirframeDrag),
573        (drag_free, KnownIssue::DragFreeSeparatedPart),
574    ]
575    .into_iter()
576    .filter(|(meets, _)| *meets)
577    .map(|(_, issue)| IssueWarning {
578        issue,
579        max_mach,
580        parts: Vec::new(),
581    })
582    .collect()
583}
584
585/// Whether a fin's leading edge sweeps forward anywhere: a trapezoid's tip ahead of its root's
586/// leading edge; a freeform outline with a point ahead of it (the root's leading edge is the
587/// outline's origin), or with an edge that runs outward and forward, as a kinked leading edge
588/// does past its kink. An ellipse's edge never does. A NaN sweep or point counts, so a broken
589/// number can't hide one.
590fn sweeps_forward(planform: &FinPlanform) -> bool {
591    match planform {
592        FinPlanform::Trapezoidal { sweep_m, .. } => sweep_m.is_nan() || *sweep_m < 0.0,
593        FinPlanform::Elliptical { .. } => false,
594        FinPlanform::Freeform { points_m, .. } => {
595            points_m.iter().any(|[x, _]| x.is_nan() || *x < 0.0)
596                || points_m.windows(2).any(|pair| {
597                    let ([x0, h0], [x1, h1]) = (pair[0], pair[1]);
598                    h1 > h0 && x1 < x0
599                })
600        }
601        // A planform added later warns until it is read: the cautious side.
602        _ => true,
603    }
604}
605
606/// The fin sets of `layout` that share a station with others, more than
607/// [`FINS_WITHOUT_INTERFERENCE`] fins between them, in the layout's order
608/// ([issue #325](https://github.com/nrdptel/hpr-sim/issues/325)). Two sets share a station when
609/// their roots overlap or touch along the axis on the same stage and the same airframe or pod
610/// set. OpenRocket's rule is unsized between overlapping roots and a common aft station; the
611/// overlap is the wider of the two, so it can't miss a case the other would warn.
612fn fin_sets_at_one_station(layout: &Layout) -> Vec<String> {
613    let fin_sets: Vec<(usize, &PlacedComponent, u32)> = layout
614        .components
615        .iter()
616        .enumerate()
617        .filter_map(|(index, component)| match &component.part {
618            Part::FinSet(fins) => Some((index, component, fins.count)),
619            _ => None,
620        })
621        .collect();
622    let shares = |(i, a, _): &(usize, &PlacedComponent, u32),
623                  (j, b, _): &(usize, &PlacedComponent, u32)| {
624        let (a_fore, a_aft) = (a.fore_station_m, a.fore_station_m + a.length_m);
625        let (b_fore, b_aft) = (b.fore_station_m, b.fore_station_m + b.length_m);
626        a.stage == b.stage
627            && layout.pod_set_of(*i) == layout.pod_set_of(*j)
628            // Written out so a NaN station counts as shared.
629            && !(a_aft < b_fore || b_aft < a_fore)
630    };
631    fin_sets
632        .iter()
633        .filter(|set| {
634            let mut sharing = false;
635            // Summed wide, so a layout built by hand with huge counts can't overflow.
636            let fins: u64 = fin_sets
637                .iter()
638                .filter(|other| shares(set, other))
639                .map(|other| {
640                    sharing |= other.0 != set.0;
641                    u64::from(other.2)
642                })
643                .sum();
644            sharing && fins > u64::from(FINS_WITHOUT_INTERFERENCE)
645        })
646        .map(|(_, component, _)| component.id.clone())
647        .collect()
648}
649
650/// [`issue_warnings`] over every part of `layout`, pods' included.
651#[must_use]
652pub fn layout_issue_warnings(layout: &Layout, max_mach: Option<Peak>) -> Vec<IssueWarning> {
653    issue_warnings(
654        layout
655            .components
656            .iter()
657            .map(|component| (component.id.as_str(), &component.part)),
658        max_mach,
659    )
660}
661
662/// The stability issues a design's shape meets, on a flight whose top Mach number is
663/// `max_mach` ([ADR-189][adr-189]): [#64](https://github.com/nrdptel/hpr-sim/issues/64) for its
664/// forward-swept fin sets past [`FORWARD_SWEEP_MACH`], then
665/// [#325](https://github.com/nrdptel/hpr-sim/issues/325) for its fin sets that share a station,
666/// then [#326](https://github.com/nrdptel/hpr-sim/issues/326) for its freeform fin sets at any
667/// speed ([`FREEFORM_FIN_MARGIN_CAL`]). They are errors in hpr's normal force, not its drag, so
668/// a flight on a drag of its own meets them and one on a normal-force table of its own doesn't,
669/// as with [`stability_issue_warnings`]. A flight with no top Mach number meets none.
670///
671/// [adr-189]: https://github.com/nrdptel/hpr-sim/blob/main/docs/decisions/0189-the-safety-and-shape-warnings.md
672#[must_use]
673pub fn layout_stability_issue_warnings(
674    layout: &Layout,
675    max_mach: Option<Peak>,
676) -> Vec<IssueWarning> {
677    let Some(max_mach) = max_mach else {
678        return Vec::new();
679    };
680    let mut warnings = Vec::new();
681    let forward_swept: Vec<String> = layout
682        .components
683        .iter()
684        .filter(|component| {
685            matches!(&component.part, Part::FinSet(fins) if sweeps_forward(&fins.planform))
686        })
687        .map(|component| component.id.clone())
688        .collect();
689    if past(max_mach.value, FORWARD_SWEEP_MACH) && !forward_swept.is_empty() {
690        warnings.push(IssueWarning {
691            issue: KnownIssue::ForwardSweptFinSlope,
692            max_mach,
693            parts: forward_swept,
694        });
695    }
696    let parts = fin_sets_at_one_station(layout);
697    if !parts.is_empty() {
698        warnings.push(IssueWarning {
699            issue: KnownIssue::FinSetsAtOneStation,
700            max_mach,
701            parts,
702        });
703    }
704    let freeform: Vec<String> = layout
705        .components
706        .iter()
707        .filter(|component| {
708            matches!(&component.part, Part::FinSet(fins)
709                if matches!(fins.planform, FinPlanform::Freeform { .. }))
710        })
711        .map(|component| component.id.clone())
712        .collect();
713    if !freeform.is_empty() {
714        warnings.push(IssueWarning {
715            issue: KnownIssue::FreeformFinCenterOfPressure,
716            max_mach,
717            parts: freeform,
718        });
719    }
720    warnings
721}
722
723/// The stability issues a flight meets ([ADR-181][adr-181]): one whose top Mach number is
724/// `max_mach`, whose aerodynamics stopped the supersonic body run for `fallback`
725/// ([`hpr_aero::AeroModel::supersonic_fallback`]), and whose smallest static margin is
726/// `min_static_margin`. Issues #87, #120 and #121 when their switch stopped the run and the flight
727/// is past [`SUPERSONIC_SWITCH_MACH`], naming the component; #172 whenever the flight has a
728/// smallest static margin. A flight with no top Mach number meets none.
729///
730/// [adr-181]: https://github.com/nrdptel/hpr-sim/blob/main/docs/decisions/0181-the-stability-issue-warnings.md
731#[must_use]
732pub fn stability_issue_warnings(
733    fallback: Option<&hpr_aero::SupersonicFallback>,
734    min_static_margin: Option<Peak>,
735    max_mach: Option<Peak>,
736) -> Vec<IssueWarning> {
737    use hpr_aero::SupersonicFallback as Fallback;
738    let Some(max_mach) = max_mach else {
739        return Vec::new();
740    };
741    let mut warnings = Vec::new();
742    if past(max_mach.value, SUPERSONIC_SWITCH_MACH) {
743        let switch = match fallback {
744            Some(Fallback::RadiusStep { component }) => {
745                Some((KnownIssue::RadiusStepFallback, vec![component.clone()]))
746            }
747            Some(Fallback::LongLip { component }) => {
748                Some((KnownIssue::LongLipFallback, vec![component.clone()]))
749            }
750            Some(Fallback::SteepTip) => Some((KnownIssue::SteepTipFallback, Vec::new())),
751            _ => None,
752        };
753        if let Some((issue, parts)) = switch {
754            warnings.push(IssueWarning {
755                issue,
756                max_mach,
757                parts,
758            });
759        }
760    }
761    if min_static_margin.is_some() {
762        warnings.push(IssueWarning {
763            issue: KnownIssue::MarginReadsHigh,
764            max_mach,
765            parts: Vec::new(),
766        });
767    }
768    warnings
769}
770
771#[cfg(test)]
772mod tests {
773    use hpr_design::{FinPlanform, FinSet, Material, NoseCone, Transition, Wall};
774
775    use super::*;
776
777    fn wood() -> Material {
778        Material::bulk("wood", 600.0)
779    }
780
781    fn nose(shape: NoseShape) -> Part {
782        Part::NoseCone(NoseCone {
783            shape,
784            length_m: 0.3,
785            base_radius_m: 0.05,
786            wall: Wall::Filled {},
787            shoulder: None,
788            material: wood(),
789        })
790    }
791
792    fn fins(cross_section: FinCrossSection) -> Part {
793        Part::FinSet(FinSet {
794            count: 3,
795            planform: FinPlanform::Trapezoidal {
796                root_chord_m: 0.2,
797                tip_chord_m: 0.1,
798                span_m: 0.1,
799                sweep_m: 0.05,
800            },
801            thickness_m: 0.004,
802            cross_section,
803            tab: None,
804            fillet: None,
805            cant_rad: 0.0,
806            base_angle_rad: 0.0,
807            material: wood(),
808        })
809    }
810
811    /// A conical transition from radius 0.05 m narrowing at `angle_rad` over 0.1 m.
812    fn boattail(angle_rad: f64) -> Part {
813        let length_m = 0.1;
814        Part::Transition(Transition {
815            shape: NoseShape::Conical {},
816            clipped: false,
817            length_m,
818            fore_radius_m: 0.05,
819            aft_radius_m: 0.05 - length_m * angle_rad.tan(),
820            wall: Wall::Filled {},
821            fore_shoulder: None,
822            aft_shoulder: None,
823            material: wood(),
824        })
825    }
826
827    fn at(mach: f64) -> Option<Peak> {
828        Some(Peak {
829            value: mach,
830            time_s: 2.5,
831            height_above_ground_m: 300.0,
832        })
833    }
834
835    fn numbers(parts: &[Part], mach: f64) -> Vec<u32> {
836        let ids = ["a", "b", "c", "d", "e"];
837        issue_warnings(ids.into_iter().zip(parts), at(mach))
838            .iter()
839            .map(IssueWarning::number)
840            .collect()
841    }
842
843    /// The next `f64` above `x`.
844    fn above(x: f64) -> f64 {
845        f64::from_bits(x.to_bits() + 1)
846    }
847
848    #[test]
849    fn a_cone_or_an_ogive_nose_warns_just_past_mach_0_8() {
850        for shape in [
851            NoseShape::Conical {},
852            NoseShape::TANGENT_OGIVE,
853            NoseShape::Ogive { radius_ratio: 2.0 },
854        ] {
855            let parts = [nose(shape)];
856            assert_eq!(numbers(&parts, TRANSONIC_NOSE_DRAG_MACH), Vec::<u32>::new());
857            assert_eq!(numbers(&parts, above(TRANSONIC_NOSE_DRAG_MACH)), [67, 68]);
858        }
859    }
860
861    #[test]
862    fn other_noses_never_meet_67() {
863        for shape in [
864            NoseShape::VON_KARMAN,
865            NoseShape::Haack {
866                parameter: 1.0 / 3.0,
867            },
868            NoseShape::Elliptical {},
869            NoseShape::PowerSeries { exponent: 0.5 },
870            NoseShape::ParabolicSeries { parameter: 1.0 },
871        ] {
872            assert_eq!(numbers(&[nose(shape)], 0.95), [68], "{shape:?}");
873        }
874    }
875
876    #[test]
877    fn base_drag_warns_past_mach_0_8_and_adds_its_low_side_only_past_1_2() {
878        assert_eq!(numbers(&[], BASE_DRAG_HIGH_MACH), Vec::<u32>::new());
879        assert_eq!(numbers(&[], above(BASE_DRAG_HIGH_MACH)), [68]);
880        let message = |mach: f64| issue_warnings(std::iter::empty(), at(mach))[0].message();
881        let low_side = "from Mach 1.5 it reads low";
882        assert!(message(BASE_DRAG_HIGH_TOP_MACH).contains("read low"));
883        assert!(!message(BASE_DRAG_HIGH_TOP_MACH).contains(low_side));
884        let past_top = message(above(BASE_DRAG_HIGH_TOP_MACH));
885        assert!(past_top.contains(low_side), "{past_top}");
886        assert!(past_top.contains("there read high"), "{past_top}");
887    }
888
889    /// The edges are the issues' numbers: a test against each literal, so moving a constant
890    /// fails here as well as in the comparisons above.
891    #[test]
892    fn the_edges_are_the_issues_numbers() {
893        assert_eq!(TRANSONIC_NOSE_DRAG_MACH, 0.8);
894        assert_eq!(BASE_DRAG_HIGH_MACH, 0.8);
895        assert_eq!(BASE_DRAG_HIGH_TOP_MACH, 1.2);
896        assert_eq!(BASE_DRAG_LOW_MACH, 1.5);
897        assert_eq!(SUPERSONIC_DRAG_MACH, 1.0);
898        assert_eq!(STEEP_BOATTAIL_RAD, 10f64.to_radians());
899        // And the messages quote them.
900        let message = |issue: KnownIssue| {
901            IssueWarning {
902                issue,
903                max_mach: at(2.0).unwrap(),
904                parts: Vec::new(),
905            }
906            .message()
907        };
908        assert!(message(KnownIssue::TransonicNoseDrag).contains("from Mach 0.8 "));
909        assert!(message(KnownIssue::BaseDrag).contains("from Mach 0.8 to 1.2 "));
910        assert!(message(KnownIssue::BaseDrag).contains("from Mach 1.5 it reads low"));
911        assert!(message(KnownIssue::SteepBoattailDrag).contains("steeper than 10°"));
912    }
913
914    #[test]
915    fn a_cone_like_shoulder_meets_67_a_von_karman_one_or_a_boattail_doesnt() {
916        let transition = |shape: NoseShape, fore_radius_m: f64, aft_radius_m: f64| {
917            Part::Transition(Transition {
918                shape,
919                clipped: false,
920                length_m: 0.2,
921                fore_radius_m,
922                aft_radius_m,
923                wall: Wall::Filled {},
924                fore_shoulder: None,
925                aft_shoulder: None,
926                material: wood(),
927            })
928        };
929        let von_karman = nose(NoseShape::VON_KARMAN);
930        for shape in [NoseShape::Conical {}, NoseShape::TANGENT_OGIVE] {
931            let shoulder = transition(shape, 0.027, 0.049);
932            assert_eq!(numbers(&[von_karman.clone(), shoulder], 0.95), [67, 68]);
933            // Narrowing at 6°, the same shape is a gentle boattail: neither #67 nor #72.
934            let boattail = transition(shape, 0.049, 0.028);
935            assert_eq!(numbers(&[von_karman.clone(), boattail], 0.95), [68]);
936        }
937        let shoulder = transition(NoseShape::VON_KARMAN, 0.027, 0.049);
938        assert_eq!(numbers(&[von_karman, shoulder], 0.95), [68]);
939    }
940
941    #[test]
942    fn a_series_nose_that_blends_in_the_cone_meets_67() {
943        let below = |x: f64| f64::from_bits(x.to_bits() - 1);
944        for (shape, meets) in [
945            (NoseShape::PowerSeries { exponent: 1.0 }, true),
946            (
947                NoseShape::PowerSeries {
948                    exponent: above(0.75),
949                },
950                true,
951            ),
952            (NoseShape::PowerSeries { exponent: 0.75 }, false),
953            (NoseShape::ParabolicSeries { parameter: 0.0 }, true),
954            (
955                NoseShape::ParabolicSeries {
956                    parameter: below(0.5),
957                },
958                true,
959            ),
960            (NoseShape::ParabolicSeries { parameter: 0.5 }, false),
961        ] {
962            let want: &[u32] = if meets { &[67, 68] } else { &[68] };
963            assert_eq!(numbers(&[nose(shape)], 0.95), want, "{shape:?}");
964        }
965    }
966
967    #[test]
968    fn airfoil_fins_warn_just_past_mach_1_other_sections_never() {
969        let parts = [fins(FinCrossSection::Airfoil)];
970        assert_eq!(numbers(&parts, SUPERSONIC_DRAG_MACH), [68]);
971        assert_eq!(numbers(&parts, above(SUPERSONIC_DRAG_MACH)), [68, 70, 222]);
972        for section in [FinCrossSection::Square, FinCrossSection::Rounded] {
973            assert_eq!(numbers(&[fins(section)], 3.0), [68], "{section:?}");
974        }
975    }
976
977    #[test]
978    fn a_boattail_warns_just_steeper_than_10_degrees_past_mach_1() {
979        // The angle's round trip through tan and atan2 lands within a few ulps of where it
980        // started, so the edge is probed a hair either side.
981        let edge = STEEP_BOATTAIL_RAD;
982        // Steeper than 9.46°, each also meets #73 (below).
983        assert_eq!(numbers(&[boattail(edge * (1.0 - 1e-12))], 2.0), [68, 73]);
984        assert_eq!(
985            numbers(&[boattail(edge * (1.0 + 1e-12))], 2.0),
986            [68, 72, 73]
987        );
988        assert_eq!(numbers(&[boattail(0.5)], SUPERSONIC_DRAG_MACH), [68, 73]);
989        assert_eq!(
990            numbers(&[boattail(0.5)], above(SUPERSONIC_DRAG_MACH)),
991            [68, 72, 73]
992        );
993    }
994
995    #[test]
996    fn a_boattails_angle_is_its_mean_and_a_flare_or_a_level_transition_is_none() {
997        assert_eq!(
998            boattail_angle_rad(0.5, 0.25, 0.25),
999            Some(45f64.to_radians())
1000        );
1001        assert_eq!(boattail_angle_rad(0.04, 0.05, 0.01), None);
1002        assert_eq!(boattail_angle_rad(0.05, 0.05, 0.01), None);
1003        // A narrowing of no length is a step, which the drag buildup doesn't take as a boattail.
1004        assert_eq!(boattail_angle_rad(0.05, 0.04, 0.0), None);
1005    }
1006
1007    #[test]
1008    fn an_ogive_nose_meets_222_supersonic_a_cone_doesnt() {
1009        assert_eq!(
1010            numbers(
1011                &[nose(NoseShape::TANGENT_OGIVE)],
1012                above(SUPERSONIC_DRAG_MACH)
1013            ),
1014            [67, 68, 222]
1015        );
1016        assert_eq!(numbers(&[nose(NoseShape::Conical {})], 2.0), [67, 68]);
1017    }
1018
1019    #[test]
1020    fn each_warning_lists_its_parts_by_id() {
1021        let parts = [
1022            nose(NoseShape::TANGENT_OGIVE),
1023            boattail(0.4),
1024            fins(FinCrossSection::Airfoil),
1025            fins(FinCrossSection::Square),
1026        ];
1027        let ids = ["nose", "tail", "fins", "square"];
1028        let warnings = issue_warnings(ids.into_iter().zip(&parts), at(1.5));
1029        let listed: Vec<(u32, Vec<&str>)> = warnings
1030            .iter()
1031            .map(|warning| {
1032                let parts = warning.parts.iter().map(String::as_str).collect();
1033                (warning.number(), parts)
1034            })
1035            .collect();
1036        assert_eq!(
1037            listed,
1038            [
1039                (67, vec!["nose"]),
1040                (68, vec![]),
1041                (70, vec!["fins"]),
1042                (72, vec!["tail"]),
1043                (222, vec!["nose", "fins"]),
1044                (73, vec!["tail"]),
1045            ]
1046        );
1047        for warning in &warnings {
1048            let message = warning.message();
1049            assert!(
1050                message.contains(&format!("issue #{}", warning.number())),
1051                "{message}"
1052            );
1053            assert!(message.contains(&warning.issue.url()), "{message}");
1054        }
1055    }
1056
1057    #[test]
1058    fn a_nan_top_mach_raises_every_warning_its_shape_allows_and_none_raises_none() {
1059        let parts = [nose(NoseShape::TANGENT_OGIVE), boattail(0.4)];
1060        assert_eq!(numbers(&parts, f64::NAN), [67, 68, 72, 222, 73]);
1061        let ids = ["a", "b"];
1062        assert!(issue_warnings(ids.into_iter().zip(&parts), None).is_empty());
1063    }
1064
1065    #[test]
1066    fn a_nan_radius_counts_as_a_steep_boattail() {
1067        let angle = boattail_angle_rad(f64::NAN, 0.04, 0.1);
1068        assert!(angle.is_some_and(|angle| past(angle, STEEP_BOATTAIL_RAD)));
1069    }
1070
1071    #[test]
1072    fn each_switch_warns_just_past_mach_1_2_naming_its_part() {
1073        use hpr_aero::SupersonicFallback as Fallback;
1074        let step = Fallback::RadiusStep {
1075            component: "tube".to_owned(),
1076        };
1077        let lip = Fallback::LongLip {
1078            component: "lip".to_owned(),
1079        };
1080        let numbers = |fallback: &Fallback, mach: f64| -> Vec<(u32, Vec<String>)> {
1081            stability_issue_warnings(Some(fallback), None, at(mach))
1082                .into_iter()
1083                .map(|w| (w.number(), w.parts))
1084                .collect()
1085        };
1086        for (fallback, number, parts) in [
1087            (&step, 87, vec!["tube".to_owned()]),
1088            (&lip, 120, vec!["lip".to_owned()]),
1089            (&Fallback::SteepTip, 121, Vec::new()),
1090        ] {
1091            assert!(numbers(fallback, SUPERSONIC_SWITCH_MACH).is_empty());
1092            assert_eq!(
1093                numbers(fallback, above(SUPERSONIC_SWITCH_MACH)),
1094                vec![(number, parts.clone())]
1095            );
1096            assert_eq!(numbers(fallback, f64::NAN), vec![(number, parts)]);
1097        }
1098        assert!(numbers(&Fallback::Other, 3.0).is_empty());
1099        assert!(stability_issue_warnings(None, None, at(3.0)).is_empty());
1100        assert!(stability_issue_warnings(Some(&step), None, None).is_empty());
1101    }
1102
1103    #[test]
1104    fn a_margin_meets_172_on_any_flight_that_has_one() {
1105        let margin = at(1.5);
1106        let numbers: Vec<u32> = stability_issue_warnings(None, margin, at(0.3))
1107            .iter()
1108            .map(IssueWarning::number)
1109            .collect();
1110        assert_eq!(numbers, vec![172]);
1111        assert!(stability_issue_warnings(None, None, at(0.3)).is_empty());
1112        let message = stability_issue_warnings(None, margin, at(0.3))[0].message();
1113        assert!(message.contains("0.1108 calibres"), "{message}");
1114        assert!(message.contains("four private designs"), "{message}");
1115    }
1116
1117    fn trapezoid(sweep_m: f64) -> FinPlanform {
1118        FinPlanform::Trapezoidal {
1119            root_chord_m: 0.2,
1120            tip_chord_m: 0.1,
1121            span_m: 0.1,
1122            sweep_m,
1123        }
1124    }
1125
1126    /// The 54 mm test design with its fin set's outline `planform`.
1127    fn swept(planform: FinPlanform) -> Layout {
1128        let mut rocket = crate::testing::design("synthetic-54mm-three-fin");
1129        for child in &mut rocket.stages[0].components[1].children {
1130            if let Part::FinSet(fins) = &mut child.part {
1131                fins.planform = planform.clone();
1132            }
1133        }
1134        rocket.layout().unwrap_or_else(|error| panic!("{error}"))
1135    }
1136
1137    fn stability_numbers(layout: &Layout, mach: f64) -> Vec<u32> {
1138        layout_stability_issue_warnings(layout, at(mach))
1139            .iter()
1140            .map(IssueWarning::number)
1141            .collect()
1142    }
1143
1144    #[test]
1145    fn a_forward_swept_fin_meets_64_just_past_mach_1() {
1146        // The tip's leading edge a hair ahead of the root's: swept forward.
1147        let forward = swept(trapezoid(-1e-9));
1148        assert_eq!(stability_numbers(&forward, above(FORWARD_SWEEP_MACH)), [64]);
1149        assert!(stability_numbers(&forward, FORWARD_SWEEP_MACH).is_empty());
1150        // It is a margin error, not a drag one: the drag warnings never list it.
1151        assert!(
1152            !layout_issue_warnings(&forward, at(3.0))
1153                .iter()
1154                .any(|warning| warning.number() == 64)
1155        );
1156        // A straight leading edge, or one swept aft, never.
1157        for sweep_m in [0.0, 0.05] {
1158            assert!(stability_numbers(&swept(trapezoid(sweep_m)), 3.0).is_empty());
1159        }
1160        let ellipse = FinPlanform::Elliptical {
1161            root_chord_m: 0.2,
1162            span_m: 0.1,
1163        };
1164        assert!(stability_numbers(&swept(ellipse), 3.0).is_empty());
1165        // A freeform outline with a point ahead of the root's leading edge, and one without.
1166        let freeform = |tip_x_m: f64| FinPlanform::Freeform {
1167            points_m: vec![[0.0, 0.0], [tip_x_m, 0.1], [0.15, 0.1], [0.2, 0.0]],
1168            root_m: Vec::new(),
1169        };
1170        // Every freeform set also meets #326.
1171        assert_eq!(stability_numbers(&swept(freeform(-1e-9)), 1.5), [64, 326]);
1172        assert_eq!(stability_numbers(&swept(freeform(0.0)), 1.5), [326]);
1173        // A kinked leading edge that turns forward past its kink, its tip still aft of the root's
1174        // leading edge; and the same kink turning straight out.
1175        let kinked = |tip_x_m: f64| FinPlanform::Freeform {
1176            points_m: vec![
1177                [0.0, 0.0],
1178                [0.1, 0.05],
1179                [tip_x_m, 0.1],
1180                [0.15, 0.1],
1181                [0.2, 0.0],
1182            ],
1183            root_m: Vec::new(),
1184        };
1185        assert_eq!(
1186            stability_numbers(&swept(kinked(0.1 - 1e-9)), 1.5),
1187            [64, 326]
1188        );
1189        assert_eq!(stability_numbers(&swept(kinked(0.1)), 1.5), [326]);
1190        // A NaN sweep or point counts (a layout refuses one, so the rule is asked directly).
1191        assert!(sweeps_forward(&trapezoid(f64::NAN)));
1192        assert!(sweeps_forward(&freeform(f64::NAN)));
1193        assert!(KnownIssue::ForwardSweptFinSlope.is_stability());
1194        assert_eq!(FORWARD_SWEEP_MACH, 1.0);
1195        assert_eq!(hpr_aero::fins::SUPERSONIC_START_MACH, 1.2);
1196    }
1197
1198    #[test]
1199    fn a_boattail_meets_73_steeper_than_where_the_rule_gives_drag_at_any_speed() {
1200        // The edge is where Niskanen's rule starts to give a boattail drag: γ = 3.
1201        assert!((SUBSONIC_BOATTAIL_RULE_RAD - (1.0_f64 / 6.0).atan()).abs() < 1e-16);
1202        let gamma = |angle_rad: f64| 1.0 / (2.0 * angle_rad.tan());
1203        let factor = |angle_rad: f64| {
1204            hpr_aero::drag::boattail_factor(gamma(angle_rad), 1.0, 0.0).unwrap_or(f64::NAN)
1205        };
1206        let edge = SUBSONIC_BOATTAIL_RULE_RAD;
1207        assert!(factor(edge * (1.0 + 1e-9)) > 0.0);
1208        assert_eq!(factor(edge * (1.0 - 1e-9)), 0.0);
1209        // A hair either side of it, at a walking pace and past Mach 1.
1210        assert_eq!(
1211            numbers(&[boattail(edge * (1.0 - 1e-12))], 0.3),
1212            Vec::<u32>::new()
1213        );
1214        assert_eq!(numbers(&[boattail(edge * (1.0 + 1e-12))], 0.3), [73]);
1215        assert_eq!(numbers(&[boattail(edge * (1.0 - 1e-12))], 0.9), [68]);
1216        assert_eq!(numbers(&[boattail(edge * (1.0 + 1e-12))], 0.9), [68, 73]);
1217        // A level transition or a flare never.
1218        assert_eq!(numbers(&[boattail(0.0)], 0.3), Vec::<u32>::new());
1219        assert_eq!(numbers(&[boattail(-0.3)], 0.3), Vec::<u32>::new());
1220        assert!(!KnownIssue::SubsonicBoattailDrag.is_stability());
1221        let warning = &issue_warnings([("tail", &boattail(0.3))], at(0.3))[0];
1222        let message = warning.message();
1223        assert!(message.contains("steeper than 9.5°"), "{message}");
1224        assert!(message.contains("below Mach 0.8 "), "{message}");
1225        assert!(!message.contains("this flight reaches"), "{message}");
1226    }
1227
1228    /// The 54 mm test design with its fin set replaced by `sets`, each a count of fins and how
1229    /// far its root's trailing edge sits forward of the airframe's aft end, m. Root chords are
1230    /// 0.15 m.
1231    fn finned(sets: &[(u32, f64)]) -> Layout {
1232        use hpr_design::Position;
1233        let mut rocket = crate::testing::design("synthetic-54mm-three-fin");
1234        let airframe = &mut rocket.stages[0].components[1];
1235        let index = airframe
1236            .children
1237            .iter()
1238            .position(|child| child.id == "sustainer-fins")
1239            .unwrap_or_else(|| panic!("the test design has a fin set"));
1240        let original = airframe.children.remove(index);
1241        for (number, &(count, from_aft_m)) in sets.iter().enumerate() {
1242            let mut set = original.clone();
1243            set.id = format!("fins{number}");
1244            set.position = Some(Position::Bottom {
1245                aft_offset_m: -from_aft_m,
1246            });
1247            if let Part::FinSet(fins) = &mut set.part {
1248                fins.count = count;
1249            }
1250            airframe.children.push(set);
1251        }
1252        rocket.layout().unwrap_or_else(|error| panic!("{error}"))
1253    }
1254
1255    fn station_numbers(layout: &Layout) -> Vec<(u32, Vec<String>)> {
1256        layout_stability_issue_warnings(layout, at(0.3))
1257            .into_iter()
1258            .map(|warning| (warning.number(), warning.parts))
1259            .collect()
1260    }
1261
1262    #[test]
1263    fn fin_sets_sharing_a_station_meet_325_past_four_fins() {
1264        let both = vec![(325, vec!["fins0".to_owned(), "fins1".to_owned()])];
1265        // Three and three at one station, overlapping, and touching end to end.
1266        assert_eq!(station_numbers(&finned(&[(3, 0.0), (3, 0.0)])), both);
1267        assert_eq!(station_numbers(&finned(&[(3, 0.0), (3, 0.1)])), both);
1268        assert_eq!(station_numbers(&finned(&[(3, 0.0), (3, 0.15)])), both);
1269        // A millimeter apart, never.
1270        assert!(station_numbers(&finned(&[(3, 0.0), (3, 0.151)])).is_empty());
1271        // Four fins between them have no interference to miss: two and two, three and one.
1272        assert!(station_numbers(&finned(&[(2, 0.0), (2, 0.0)])).is_empty());
1273        assert!(station_numbers(&finned(&[(3, 0.0), (1, 0.0)])).is_empty());
1274        // Five fins, three and two: both listed.
1275        assert_eq!(station_numbers(&finned(&[(3, 0.0), (2, 0.0)])), both);
1276        // One set of six alone is counted as six already.
1277        assert!(station_numbers(&finned(&[(6, 0.0)])).is_empty());
1278        // A third set apart is not listed.
1279        assert_eq!(
1280            station_numbers(&finned(&[(3, 0.0), (3, 0.0), (3, 0.5)])),
1281            both
1282        );
1283        assert!(KnownIssue::FinSetsAtOneStation.is_stability());
1284        let warning = &layout_stability_issue_warnings(&finned(&[(3, 0.0), (3, 0.0)]), at(0.3))[0];
1285        let message = warning.message();
1286        assert!(message.contains("issue #325"), "{message}");
1287        assert!(message.contains("more than 4 fins"), "{message}");
1288        assert!(!message.contains("this flight reaches"), "{message}");
1289        // No top Mach number, no flight: nothing.
1290        assert!(layout_stability_issue_warnings(&finned(&[(3, 0.0), (3, 0.0)]), None).is_empty());
1291        // Nor among the drag warnings.
1292        assert!(layout_issue_warnings(&finned(&[(3, 0.0), (3, 0.0)]), at(0.3)).is_empty());
1293    }
1294
1295    #[test]
1296    fn every_flight_on_hprs_drag_meets_18_and_none_without_a_top_mach() {
1297        // At any speed, and apart from the shape's warnings.
1298        for mach in [0.05, 0.3, 3.0] {
1299            let warnings = friction_issue_warnings(at(mach));
1300            let numbers: Vec<u32> = warnings.iter().map(IssueWarning::number).collect();
1301            assert_eq!(numbers, [18]);
1302            assert!(!self::numbers(&[nose(NoseShape::VON_KARMAN)], mach).contains(&18));
1303        }
1304        assert!(friction_issue_warnings(None).is_empty());
1305        assert!(!KnownIssue::TurbulentFriction.is_stability());
1306        let warning = &friction_issue_warnings(at(0.3))[0];
1307        let message = warning.message();
1308        assert!(message.contains("issue #18"), "{message}");
1309        assert!(
1310            message.contains("by 3.6% on RocketPy's Calisto"),
1311            "{message}"
1312        );
1313        assert!(
1314            message.contains("the apogee, the top speed and the drift read low"),
1315            "{message}"
1316        );
1317        assert!(!message.contains("this flight reaches"), "{message}");
1318        assert!(warning.parts.is_empty());
1319    }
1320
1321    #[test]
1322    fn the_laminar_share_is_barrowmans_term_over_niskanens_turbulent_friction() {
1323        let reynolds: f64 = 1.8e7;
1324        let turbulent = 1.0 / (1.50 * reynolds.ln() - 5.6).powi(2);
1325        let percent = 100.0 * (1700.0 / reynolds) / turbulent;
1326        assert!(
1327            (percent - LAMINAR_FRICTION_PERCENT).abs() < 0.05,
1328            "{percent}"
1329        );
1330    }
1331
1332    #[test]
1333    fn a_freeform_fin_set_meets_326_at_any_speed_a_trapezoid_or_an_ellipse_never() {
1334        let outline = FinPlanform::Freeform {
1335            points_m: vec![[0.0, 0.0], [0.05, 0.05], [0.1, 0.1], [0.2, 0.0]],
1336            root_m: Vec::new(),
1337        };
1338        let freeform = swept(outline);
1339        for mach in [0.05, 0.3, 3.0] {
1340            assert!(stability_numbers(&freeform, mach).contains(&326), "{mach}");
1341        }
1342        let warning = layout_stability_issue_warnings(&freeform, at(0.3))
1343            .into_iter()
1344            .find(|warning| warning.number() == 326)
1345            .unwrap_or_else(|| panic!("no #326"));
1346        assert_eq!(warning.parts.len(), 1);
1347        let message = warning.message();
1348        assert!(message.contains("issue #326"), "{message}");
1349        assert!(message.contains("0.047 calibres"), "{message}");
1350        assert!(!message.contains("this flight reaches"), "{message}");
1351        assert!(KnownIssue::FreeformFinCenterOfPressure.is_stability());
1352        // The same design's trapezoid and an ellipse never.
1353        assert!(!stability_numbers(&swept(trapezoid(0.05)), 3.0).contains(&326));
1354        let ellipse = FinPlanform::Elliptical {
1355            root_chord_m: 0.2,
1356            span_m: 0.1,
1357        };
1358        assert!(!stability_numbers(&swept(ellipse), 3.0).contains(&326));
1359        // Not among the drag warnings, and no top Mach number, no flight: nothing.
1360        assert!(
1361            !layout_issue_warnings(&freeform, at(0.3))
1362                .iter()
1363                .any(|warning| warning.number() == 326)
1364        );
1365        assert!(layout_stability_issue_warnings(&freeform, None).is_empty());
1366    }
1367}