brahe 1.4.0

Brahe is a modern satellite dynamics library for research and engineering applications designed to be easy-to-learn, high-performance, and quick-to-deploy. The north-star of the development is enabling users to solve meaningful problems and answer questions quickly, easily, and correctly.
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
72
73

/*!
 Acceleration due to non-Newtionian gravitational effects. Specifically, special and general
relativity.
 */

use nalgebra::{Vector3, Vector6};

use crate::constants::{C_LIGHT, GM_EARTH};

/// Calculate the acceleration due to special and general relativity for an Earth orbiting object.
///
/// # Arguments
///
/// - `x_object`: State vector of the object in the ECI frame.
///
/// # Returns
///
/// - `a_relativity` : Acceleration due to special and general relativity.
///
/// # Examples
///
/// ```
/// use brahe::coordinates::state_koe_to_eci;
/// use brahe::orbit_dynamics::accel_relativity;
/// use nalgebra::Vector6;
/// use brahe::R_EARTH;
///
/// let x_object = Vector6::new(R_EARTH + 500.0e3, 0.0, 0.0, 0.0, 0.0, 0.0);
/// let a_relativity = accel_relativity(x_object);
/// ```
pub fn accel_relativity(x_object: Vector6<f64>) -> Vector3<f64> {
    // Extract state variables
    let r = x_object.fixed_rows::<3>(0);
    let v = x_object.fixed_rows::<3>(3);

    // Intermediate computations
    let norm_r = r.norm();
    let r2 = norm_r.powi(2);
    let norm_v = v.norm();
    let v2 = norm_v.powi(2);
    let c2 = C_LIGHT.powi(2);

    // Compute unit vectors
    let er = r / norm_r;
    let ev = v / norm_v;

    // Compute perturbation acceleration and return
    GM_EARTH / r2
        * ((4.0 * GM_EARTH / (c2 * norm_r) - v2 / c2) * er + 4.0 * v2 / c2 * er.dot(&ev) * ev)
}

#[cfg(test)]
#[cfg_attr(coverage_nightly, coverage(off))]
mod tests {
    use crate::DEGREES;
    use crate::constants::R_EARTH;
    use crate::coordinates::*;

    #[test]
    fn test_accel_relativity() {
        use super::*;

        let oe = Vector6::new(R_EARTH + 500e3, 0.01, 97.3, 15.0, 30.0, 45.0);

        let x_object = state_koe_to_eci(oe, DEGREES);

        let a = accel_relativity(x_object);

        // According to Motenbruck and Gill this should be on ghd order of ~1e-8 for a satellite
        // around 500 km altitude.
        assert!(a.norm() < 1.0e-7);
    }
}