marss 0.0.3

Mars celestial simulation crate for the MilkyWay SolarSystem workspace
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71

use sciforge::hub::domain::common::constants::DEGREE_TO_RAD;

pub const SIDEREAL_DAY_S: f64 = crate::SIDEREAL_DAY_S;
pub const SOL_S: f64 = crate::SOL_S;
pub const AXIAL_TILT_DEG: f64 = crate::AXIAL_TILT_DEG;
pub const AXIAL_TILT_RAD: f64 = AXIAL_TILT_DEG * DEGREE_TO_RAD;
pub const PRECESSION_PERIOD_YEARS: f64 = 170_000.0;

pub struct MarsRotation {
    pub angular_velocity_rad_s: f64,
    pub axial_tilt_rad: f64,
}

impl Default for MarsRotation {
    fn default() -> Self {
        Self::new()
    }
}

impl MarsRotation {
    pub fn new() -> Self {
        Self {
            angular_velocity_rad_s: 2.0 * std::f64::consts::PI / SIDEREAL_DAY_S,
            axial_tilt_rad: AXIAL_TILT_RAD,
        }
    }

    pub fn surface_velocity_at_latitude(&self, latitude_deg: f64) -> f64 {
        let lat_rad = latitude_deg.to_radians();
        self.angular_velocity_rad_s * crate::MARS_EQUATORIAL_RADIUS * lat_rad.cos()
    }

    pub fn coriolis_parameter(&self, latitude_deg: f64) -> f64 {
        let lat_rad = latitude_deg.to_radians();
        2.0 * self.angular_velocity_rad_s * lat_rad.sin()
    }

    pub fn moment_of_inertia(&self) -> f64 {
        crate::MOI_FACTOR * crate::MARS_MASS * crate::MARS_RADIUS * crate::MARS_RADIUS
    }

    pub fn rotational_kinetic_energy(&self) -> f64 {
        0.5 * self.moment_of_inertia() * self.angular_velocity_rad_s.powi(2)
    }

    pub fn angular_momentum(&self) -> f64 {
        self.moment_of_inertia() * self.angular_velocity_rad_s
    }

    pub fn precession_rate_rad_per_year(&self) -> f64 {
        2.0 * std::f64::consts::PI / PRECESSION_PERIOD_YEARS
    }

    pub fn solar_declination(&self, areocentric_longitude_deg: f64) -> f64 {
        let ls_rad = areocentric_longitude_deg.to_radians();
        (self.axial_tilt_rad.sin() * ls_rad.sin()).asin()
    }

    pub fn day_length_hours(&self, latitude_deg: f64, areocentric_longitude_deg: f64) -> f64 {
        let lat_rad = latitude_deg.to_radians();
        let decl = self.solar_declination(areocentric_longitude_deg);
        let cos_ha = -(lat_rad.tan() * decl.tan());
        if cos_ha > 1.0 {
            return 0.0;
        }
        if cos_ha < -1.0 {
            return 24.0 * SOL_S / 86400.0;
        }
        2.0 * cos_ha.acos() / (2.0 * std::f64::consts::PI) * SOL_S / 3600.0
    }
}