Skip to content

Commit

Permalink
Solar eclipsing in ANISE is significantly better than it was in Nyx
Browse files Browse the repository at this point in the history
  • Loading branch information
ChristopherRabotin committed Sep 7, 2024
1 parent d2f9497 commit ff4130f
Show file tree
Hide file tree
Showing 2 changed files with 16 additions and 10 deletions.
14 changes: 10 additions & 4 deletions src/dynamics/solarpressure.rs
Original file line number Diff line number Diff line change
Expand Up @@ -128,15 +128,18 @@ impl ForceModel for SolarPressure {

let r_sun_unit = r_sun / r_sun.norm();

// Compute the shaddowing factor.
let k: f64 = self
// ANISE returns the occultation percentage (or factor), which is the opposite as the illumination factor.
let occult = self
.e_loc
.compute(osc, almanac)
.context(DynamicsAlmanacSnafu {
action: "solar radiation pressure computation",
})?
.factor();

// Compute the illumination factor.
let k: f64 = (occult - 1.0).abs();

let r_sun_au = r_sun.norm() / AU;
// in N/(m^2)
let flux_pressure = (k * self.phi / SPEED_OF_LIGHT_M_S) * (1.0 / r_sun_au).powi(2);
Expand All @@ -163,15 +166,18 @@ impl ForceModel for SolarPressure {
let r_sun_d: Vector3<OHyperdual<f64, Const<9>>> = hyperspace_from_vector(&r_sun);
let r_sun_unit = r_sun_d / norm(&r_sun_d);

// Compute the shadowing factor.
let k: f64 = self
// ANISE returns the occultation percentage (or factor), which is the opposite as the illumination factor.
let occult = self
.e_loc
.compute(osc, almanac.clone())
.context(DynamicsAlmanacSnafu {
action: "solar radiation pressure computation",
})?
.factor();

// Compute the illumination factor.
let k: f64 = (occult - 1.0).abs();

let r_sun_au = norm(&r_sun_d) / AU;
let inv_r_sun_au = OHyperdual::<f64, Const<9>>::from_real(1.0) / (r_sun_au);
let inv_r_sun_au_p2 = inv_r_sun_au.powi(2);
Expand Down
12 changes: 6 additions & 6 deletions tests/mission_design/force_models.rs
Original file line number Diff line number Diff line change
Expand Up @@ -67,8 +67,8 @@ fn srp_earth_full_vis(almanac: Arc<Almanac>) {
err_r * 1e3,
err_v * 1e3
);
assert!(err_r < 5e-3, "position error too large for SRP");
assert!(err_v < 1e-7, "velocity error too large for SRP");
assert!(err_r < 5e-4, "position error too large for SRP");
assert!(err_v < 9e-8, "velocity error too large for SRP");
}

#[rstest]
Expand Down Expand Up @@ -119,8 +119,8 @@ fn srp_earth_leo(almanac: Arc<Almanac>) {
err_r * 1e3,
err_v * 1e3
);
assert!(err_r < 3e-1, "position error too large for SRP");
assert!(err_v < 5e-4, "velocity error too large for SRP");
assert!(err_r < 6e-3, "position error too large for SRP");
assert!(err_v < 7e-6, "velocity error too large for SRP");
}

#[rstest]
Expand Down Expand Up @@ -173,8 +173,8 @@ fn srp_earth_meo_ecc_inc(almanac: Arc<Almanac>) {
err_r * 1e3,
err_v * 1e3
);
assert!(err_r < 5e-2, "position error too large for SRP");
assert!(err_v < 2e-5, "velocity error too large for SRP");
assert!(err_r < 2e-3, "position error too large for SRP");
assert!(err_v < 1e-6, "velocity error too large for SRP");

// Compare the case with the hyperdual EOMs (computation uses another part of the code)
let mut prop = setup.with(sc.with_stm(), almanac);
Expand Down

0 comments on commit ff4130f

Please sign in to comment.