Add SIM_orbinit set03 (semi-latus rectum + true anomaly) Tier 3 coverage

crates/astrodyn_dynamics/src/body_init.rs

@@ -78,6 +78,75 @@ pub fn init_from_orbital_elements(
     TranslationalState { position, velocity }
 }
 
+/// Initialize translational state from Keplerian orbital elements using the
+/// semi-latus rectum (rather than semi-major axis) plus true anomaly.
+///
+/// Port of the `SlrEccIncAscnodeArgperTanom` branch of JEOD
+/// `DynBodyInitOrbit::apply()` from
+/// `models/dynamics/body_action/src/dyn_body_init_orbit.cc:196-200, 285-321`.
+///
+/// JEOD selects `shape = ShapeSemiLatusRectum`, which **skips** the
+/// `semi_latus_rectum = semi_major_axis * (1 - e²)` derivation (that block
+/// runs only for `ShapeSemiMajorAxis`). The deck-supplied semi-latus rectum
+/// is therefore used verbatim as `elem.semiparam`. To match JEOD bit-for-bit
+/// we set `semiparam = semi_latus_rectum` directly here — routing through
+/// `init_from_orbital_elements` (which recomputes `semiparam = a·(1-e²)` from
+/// `a = p/(1-e²)`) would introduce a round-trip that JEOD never performs.
+///
+/// # Arguments
+/// * `semi_latus_rectum` - Semi-latus rectum p (m)
+/// * `eccentricity` - Orbital eccentricity
+/// * `inclination` - Inclination (rad)
+/// * `raan` - Right ascension of ascending node (rad)
+/// * `arg_periapsis` - Argument of periapsis (rad)
+/// * `true_anomaly` - True anomaly (rad)
+/// * `mu` - Gravitational parameter of central body (m^3/s^2)
+pub fn init_from_semi_latus_rectum_true_anomaly(
+    semi_latus_rectum: f64,
+    eccentricity: f64,
+    inclination: f64,
+    raan: f64,
+    arg_periapsis: f64,
+    true_anomaly: f64,
+    mu: f64,
+) -> TranslationalState {
+    // JEOD_INV: BA.05 — orbit initializer requires a valid gravity source (mu > 0)
+    // JEOD dyn_body_init_orbit.cc:98-111: validate mu before use.
+    assert!(
+        mu > 0.0,
+        "init_from_semi_latus_rectum_true_anomaly: mu must be positive, got {mu}"
+    );
+    assert!(
+        semi_latus_rectum > 0.0 && semi_latus_rectum.is_finite(),
+        "init_from_semi_latus_rectum_true_anomaly: semi_latus_rectum must be positive and finite, \
+         got {semi_latus_rectum}"
+    );
+    assert!(
+        (0.0..1.0).contains(&eccentricity),
+        "init_from_semi_latus_rectum_true_anomaly: eccentricity must be in [0, 1), \
+         got {eccentricity}"
+    );
+
+    // JEOD dyn_body_init_orbit.cc: ShapeSemiLatusRectum leaves semiparam as
+    // the deck value, then sets the angles, true_anom, and calls
+    // nu_to_anomalies() followed by to_cartesian().
+    use astrodyn_quantities::frame::SelfPlanet;
+    let mut oe = OrbitalElements::<SelfPlanet>::default();
+    oe.semiparam = semi_latus_rectum;
+    oe.e_mag = eccentricity;
+    oe.inclination = inclination;
+    oe.long_asc_node = raan;
+    oe.arg_periapsis = arg_periapsis;
+    oe.true_anom = true_anomaly;
+    oe.nu_to_anomalies();
+
+    let (position, velocity) = oe
+        .to_cartesian(mu)
+        .expect("init_from_semi_latus_rectum_true_anomaly: to_cartesian failed");
+
+    TranslationalState { position, velocity }
+}
+
 /// Typed sibling of [`init_from_orbital_elements`].
 ///
 /// Returns a [`TranslationalStateTyped<RootInertial>`] — Phase 3 callers
@@ -569,7 +638,8 @@ mod tests {
             .time_periapsis
             .expect("ISS set01 should have time_periapsis");
         let state = init_from_time_periapsis(
-            init.semi_major_axis,
+            init.semi_major_axis
+                .expect("ISS set01 should have semi_major_axis"),
             init.eccentricity,
             init.inclination,
             init.ascending_node,
@@ -1323,6 +1393,63 @@ mod tests {
     // assert and leave the others intact.
     // =======================================================================
 
+    #[test]
+    fn slr_true_anomaly_matches_orbital_elements_within_roundoff() {
+        // The slr+true-anomaly converter sets semiparam = p directly (JEOD's
+        // SlrEccIncAscnodeArgperTanom path), whereas init_from_orbital_elements
+        // takes sma and recomputes semiparam = a·(1-e²). Feeding the
+        // algebraically-equivalent sma = p/(1-e²) into the sma path must agree
+        // to within the round-trip roundoff (a few ULP × radius).
+        let p = 6_732_889.984_55;
+        let e = 0.00129073350;
+        let inc = 51.670450765_f64.to_radians();
+        let raan = 49.708417385_f64.to_radians();
+        let argp = 100.582445989_f64.to_radians();
+        let nu = 299.884499026_f64.to_radians();
+
+        let slr = init_from_semi_latus_rectum_true_anomaly(p, e, inc, raan, argp, nu, EARTH_MU);
+        let a = p / (1.0 - e * e);
+        let sma = init_from_orbital_elements(a, e, inc, raan, argp, nu, EARTH_MU);
+
+        let pos_err = (slr.position - sma.position).length();
+        let vel_err = (slr.velocity - sma.velocity).length();
+        // ~7e6 m radius × ~1e-15 relative ULP ≈ 1e-8 m; allow generous margin.
+        assert!(
+            pos_err < 1e-6,
+            "slr vs sma position roundoff too large: {pos_err} m"
+        );
+        assert!(
+            vel_err < 1e-9,
+            "slr vs sma velocity roundoff too large: {vel_err} m/s"
+        );
+    }
+
+    #[test]
+    fn slr_true_anomaly_position_magnitude_matches_conic() {
+        // r = p / (1 + e·cos ν) — verify the converter reproduces the conic
+        // radius for a non-trivial true anomaly.
+        let p = 6_700_000.0;
+        let e = 0.01;
+        let nu = 1.3_f64; // rad
+        let state = init_from_semi_latus_rectum_true_anomaly(p, e, 0.5, 0.3, 0.7, nu, EARTH_MU);
+        let r_expected = p / (1.0 + e * nu.cos());
+        let r_actual = state.position.length();
+        assert!(
+            (r_actual - r_expected).abs() < 1e-6,
+            "conic radius: expected {r_expected}, got {r_actual}"
+        );
+    }
+
+    #[test]
+    #[should_panic(expected = "mu must be positive")]
+    fn ba_05_panics_on_zero_mu_in_slr_true_anomaly_init() {
+        // JEOD_INV: BA.05 — `init_from_semi_latus_rectum_true_anomaly` shares
+        // the mu>0 guard so the set03 path can't slip a misconfigured gravity
+        // source past the others.
+        let _ =
+            init_from_semi_latus_rectum_true_anomaly(6_700_000.0, 0.01, 0.0, 0.0, 0.0, 0.0, 0.0);
+    }
+
     #[test]
     #[should_panic(expected = "mu must be positive")]
     fn ba_05_panics_on_zero_mu_in_orbital_elements_init() {

crates/astrodyn_dynamics/src/lib.rs

@@ -83,7 +83,7 @@ pub use abm4::{abm4_translational_step, Abm4State};
 pub use attach::{combine_states_at_attach, AttachCombineInputs, AttachCombineOutputs};
 pub use body_init::{
     compute_ned_rotation, init_from_lvlh, init_from_mean_anomaly, init_from_ned,
-    init_from_orbital_elements, init_from_time_periapsis,
+    init_from_orbital_elements, init_from_semi_latus_rectum_true_anomaly, init_from_time_periapsis,
 };
 pub use constraints::{apply_constraint, BaumgarteSolver, HolonomicConstraint, PendulumConstraint};
 pub use forces::{

crates/astrodyn_dynamics/tests/tier2_body_init.rs

@@ -95,7 +95,8 @@ fn iss_set01_time_periapsis() {
         .expect("ISS set01 must have time_periapsis");
 
     let computed = init_from_time_periapsis(
-        init.semi_major_axis,
+        init.semi_major_axis
+            .expect("ISS set01 must have semi_major_axis"),
         init.eccentricity,
         init.inclination,
         init.ascending_node,
@@ -144,7 +145,8 @@ fn iss_set02_mean_anomaly() {
     let mean_anomaly = init.mean_anomaly.expect("ISS set02 must have mean_anomaly");
 
     let computed = init_from_mean_anomaly(
-        init.semi_major_axis,
+        init.semi_major_axis
+            .expect("ISS set02 must have semi_major_axis"),
         init.eccentricity,
         init.inclination,
         init.ascending_node,
@@ -193,7 +195,8 @@ fn iss_set10_true_anomaly() {
     let true_anomaly = init.true_anomaly.expect("ISS set10 must have true_anomaly");
 
     let computed = init_from_orbital_elements(
-        init.semi_major_axis,
+        init.semi_major_axis
+            .expect("ISS set10 must have semi_major_axis"),
         init.eccentricity,
         init.inclination,
         init.ascending_node,
@@ -240,7 +243,9 @@ fn iss_element_sets_cross_consistent() {
         "trans_Orbit_inertial_body_set01",
     );
     let state01 = init_from_time_periapsis(
-        init01.semi_major_axis,
+        init01
+            .semi_major_axis
+            .expect("ISS set01 must have semi_major_axis"),
         init01.eccentricity,
         init01.inclination,
         init01.ascending_node,
@@ -255,7 +260,9 @@ fn iss_element_sets_cross_consistent() {
         "trans_Orbit_inertial_body_set02",
     );
     let state02 = init_from_mean_anomaly(
-        init02.semi_major_axis,
+        init02
+            .semi_major_axis
+            .expect("ISS set02 must have semi_major_axis"),
         init02.eccentricity,
         init02.inclination,
         init02.ascending_node,
@@ -270,7 +277,9 @@ fn iss_element_sets_cross_consistent() {
         "trans_Orbit_inertial_body_set10",
     );
     let state10 = init_from_orbital_elements(
-        init10.semi_major_axis,
+        init10
+            .semi_major_axis
+            .expect("ISS set10 must have semi_major_axis"),
         init10.eccentricity,
         init10.inclination,
         init10.ascending_node,

crates/astrodyn_math/tests/jeod_validation.rs

@@ -19,12 +19,14 @@ fn validate_iss_orbital_elements_to_cartesian() {
     // not `$JEOD_HOME` at runtime.
     let init = orbital_init::load_orbital_init("ISS", "trans_Orbit_inertial_body_set01");
     let expected = reference_state::load_reference_state("ISS", "inertial");
+    let sma = init
+        .semi_major_axis
+        .expect("ISS set01 should have semi_major_axis");
 
     // Verify parsed values are sensible
     assert!(
-        init.semi_major_axis > 6_000_000.0 && init.semi_major_axis < 7_000_000.0,
-        "ISS semi-major axis should be ~6732 km, got {} m",
-        init.semi_major_axis
+        sma > 6_000_000.0 && sma < 7_000_000.0,
+        "ISS semi-major axis should be ~6732 km, got {sma} m",
     );
     assert!(
         init.eccentricity < 0.01,
@@ -45,15 +47,15 @@ fn validate_iss_orbital_elements_to_cartesian() {
     //   mean_anomaly = time_periapsis * sqrt(mu / a) / a
     // which is algebraically identical to M = n·t_peri with n = sqrt(mu/a^3)
     // but matches JEOD's arithmetic order for bit-parity with the port.
-    let a = init.semi_major_axis;
+    let a = sma;
     let t_peri = init
         .time_periapsis
         .expect("ISS set01 should have time_periapsis");
     let mean_anomaly = t_peri * (MU_EARTH / a).sqrt() / a;
     let n = (MU_EARTH / (a * a * a)).sqrt();
 
     let mut oe = OrbitalElements::<astrodyn_quantities::frame::SelfPlanet>::default();
-    oe.semi_major_axis = init.semi_major_axis;
+    oe.semi_major_axis = a;
     oe.e_mag = init.eccentricity;
     oe.inclination = init.inclination;
     oe.long_asc_node = init.ascending_node;
@@ -368,10 +370,12 @@ fn validate_orbital_init_parser() {
     // Cross-check parsed values against known file contents
     let deg2rad = std::f64::consts::PI / 180.0;
 
+    let sma = init
+        .semi_major_axis
+        .expect("ISS set01 should have semi_major_axis");
     assert!(
-        (init.semi_major_axis - 6_732_901.201_52).abs() < 0.01,
-        "semi_major_axis: expected 6732901.20152 m, got {}",
-        init.semi_major_axis
+        (sma - 6_732_901.201_52).abs() < 0.01,
+        "semi_major_axis: expected 6732901.20152 m, got {sma}",
     );
     assert!(
         (init.eccentricity - 0.00129073350).abs() < 1e-12,

crates/astrodyn_verif_jeod/src/bin/extract_body_init.rs

@@ -27,6 +27,7 @@
 //!   `reference_inertial_trans_state.py`        — ECI reference state
 //!   `trans_Orbit_inertial_body_set01.py`       — orbit (sma/ecc/inc/raan/argp/t_peri)
 //!   `trans_Orbit_inertial_body_set02.py`       — orbit (mean anomaly)
+//!   `trans_Orbit_inertial_body_set03.py`       — orbit (semi-latus rectum + true anomaly)
 //!   `trans_Orbit_inertial_body_set10.py`       — orbit (true anomaly)
 //!   `trans_Orbit_pfix_body_set01.py`           — pfix orbit (set01 form)
 //!   `trans_TransState_inertial_body.py`        — direct Cartesian (STS_114 only)
@@ -63,6 +64,7 @@ const SCENARIOS: &[Scenario] = &[
         orbit_inits: &[
             "trans_Orbit_inertial_body_set01",
             "trans_Orbit_inertial_body_set02",
+            "trans_Orbit_inertial_body_set03",
             "trans_Orbit_inertial_body_set10",
             "trans_Orbit_pfix_body_set01",
         ],
@@ -74,6 +76,7 @@ const SCENARIOS: &[Scenario] = &[
         orbit_inits: &[
             "trans_Orbit_inertial_body_set01",
             "trans_Orbit_inertial_body_set02",
+            "trans_Orbit_inertial_body_set03",
             "trans_Orbit_pfix_body_set01",
         ],
         trans_states: &["trans_TransState_inertial_body"],
@@ -313,7 +316,13 @@ fn write_bundle(
         writeln!(
             out,
             "      \"semi_major_axis\": {},",
-            fmt(init.semi_major_axis)
+            fmt_opt(init.semi_major_axis)
+        )
+        .unwrap();
+        writeln!(
+            out,
+            "      \"semi_latus_rectum\": {},",
+            fmt_opt(init.semi_latus_rectum)
         )
         .unwrap();
         writeln!(out, "      \"eccentricity\": {},", fmt(init.eccentricity)).unwrap();