Add universal gravitational constant to constants.rs

[philiplinden]

High level description

Universal gravitational constant, G, is surprisingly absent from anise::constants. I find it odd that speed of light is defined here but G is not.

Requirements

Add G in canonical units to anise/src/constants.rs as an f64 to the greatest reasonable precision.

Test plans

None. This does not impact anise as an application, just adds an endpoint for anise as a library.

Open Questions

• Should there be a Quantity numerical type that accounts for uncertainties in universal constants like c and G?
• What is "sufficient precision" for G?
• If G isn't used by anise, should it be in anise at all?

[ChristopherRabotin]

The gravitational constant is not used in ANISE anywhere because it doesn't attempt to solve the N-body problem. Instead, it relies on the N-body solutions published by NAIF via the developmental ephemerides (the "DE" files).

Do you have examples where it would make sense for ANISE to publish the gravitational constant?

[philiplinden]

I use ANISE to evaluate celestial body positions at a particular epoch. I need G to evaluate the total potential at a given coordinate.

ANISE provides an interface to import high fidelity frames used for gravity field and mass concentrations, but doesn't provide me any interface for calculating what the potential would actually be.

The speed of light is included in constants.rs and it looks to me like there are plenty of functions for calculations precisely the apparent position of bodies. Doesn't it make sense to observe the gravity with great precision too?

It seems reasonable to me to add a "canonical" definition of G in ANISE. The alternatives are:

1. I define G myself, using whatever precision I find in my reference book. If I want c I use anise::constants::SPEED_OF_LIGHT_KM_S and if I want G I use crate::my_extra_constants::GRAVITATIONAL_PARAMETER. I need to make sure my constant uses the correct units and data type to be compatible with ANISE.
2. I install a whole separate crate like uom just to import one constant. Now I have to write extra junk to make sure the data types from the other almanac of constants is compatible with ANISE definitions.

It's not a matter of whether ANISE needs the constant, it's a matter of consistency. ANISE is my almanac for everything related to celestial bodies. Oh except for G, that's over there.

I do not use SPICE, so I am likely not using ANISE in a way that is directly applicable there.

[ChristopherRabotin]

I understand your point. I'm OK adding that constants, and I would recommend using the value from NIST directly: https://physics.nist.gov/cgi-bin/cuu/Value?bg (that's the reference on Wikipedia).

As for the UOM crate, or another framework that allows defining quantities, I've long debated and considered developing my own crate for this (and there was a super minimalistic version of this in ANISE a few months ago). Three keys issues with uom to my knowledge are the following:

1. All quantities are in fact converted back to their base unit before any computation. Although this is probably fine in a lot of applications, converting to/from kilometers and meters repeatedly will lead to rounding errors. It isn't obvious in the case of a single query in ANISE, but in the case of a propagator (like in Nyx), the errors become visible: although the math does not change much it seriously affects the validation via other tools.
2. UOM isn't trivially compatible with Python and writing those bindings seems like a daunting task. At the moment, ANISE allows users to initialize classes with simple doubles that are expected to be in the correct units. This is definitely error-prone, but the units should be documented in the function signature and/or in the documentation of that function. I'm open to other ideas here.
3. When I tried UOM some years ago or so, it seemed like it was far from trivial to write function signatures whose input was "forced" to be in a given unit (and maybe that is not how one should be using the library). As such, it seemed to me that UOM was great for generalizing operations across different units and ensuring that the results were homogeneous, but in my experience, ensuring that the inputs and outputs of a given function are in the correct unit is more important. Again, I might have totally misunderstood how to use uom, and obviously, that crate is a used a whole lot.

When you mention uncertainties on the quantity, do you envision ANISE supporting the uncertainty in the calculation itself?

[philiplinden]

I don't see ANISE natively supporting uncertainties with quantities as a must have.

I agree that there are issues with uom, especially when it comes to orders of magnitude. Gravity calcs alone have G (e-11) planetary masses (e+24 or more) and astronomical distances (e+12 m right?)

[ChristopherRabotin]

Sounds like we're in agreement. Would you like to work on adding that constant to anise? I can work on it later this week if you wish.

[philiplinden]

I think I'll have time. If you don't see a PR by Friday assume I drowned in work and it fell off my desk.

[philiplinden]

This might not be necessary. The Almanac has a row Gravity param (km^3/s^2) for a given body, which satisfies the needs presented here. Even uses convenient units.

https://github.com/nyx-space/anise/blob/0ee9d7cd50c4346d98567ed749a30bc19e5ba84d/anise/src/almanac/planetary.rs#L86

It would be more convenient to use anise::constants::usual_planetary_constants but its not a whole lot better and I don't want to include duplicate information from more than one source.

I'll tinker some more before closing this ticket. I am going to see if it makes sense to add a public function to Almanac to make parameters easier to obtain besides parsing the output of .describe().

[philiplinden]

It's even easier than I thought. Mu can be pulled directly from a given Frame.

https://github.com/nyx-space/anise/blob/0ee9d7cd50c4346d98567ed749a30bc19e5ba84d/anise/src/frames/frame.rs#L35-L42

If I need the total gravitational force of a point in space, I can use mu from the current frame, distance to the body, and that's all you need to get "N / kg"

I have to double check my arithmetic in order to use an arbitrary frame, but it's easy enough to store mu of each body until it's needed and then combine them however it makes sense.

Before I close the ticket I'll try to write a test that demos this.

[philiplinden]

#[test]
#[should_panic]
fn get_mu_from_frame() {
    // mu is only defined for celestial frames, but this frame is not
    let earth_j2k = anise::constants::frames::EARTH_J2000;
    assert!(!earth_j2k.is_celestial());
    earth_j2k.mu_km3_s2().unwrap();
}

#[test]
fn get_mu_from_orbit() {
    use anise::constants::frames::EARTH_J2000;
    use anise::prelude::{Epoch, MetaAlmanac, Orbit};

    let eme2k = MetaAlmanac::latest().unwrap().frame_from_uid(EARTH_J2000).unwrap();
    let epoch = Epoch::from_unix_seconds(0.0);

    // Define the orbit from its Keplerian orbital elements.
    let test_orbit = Orbit::keplerian(
        8_191.93, 1e-6, 12.85, 306.614, 314.19, 99.887_7, epoch, eme2k,
    );
    assert!(test_orbit.frame.is_celestial());
    let test_gm = test_orbit.frame.mu_km3_s2().unwrap();
    assert_eq!(test_gm, 398_600.441_8)
}

assertion `left == right` failed
  left: 398600.435436096
 right: 398600.4418

Close enough. Closing as wontfix.

[ChristopherRabotin]

The value of the gravitational parameter is needed for a number of computations of the orbital elements, and it's stored in the Frame to avoid having to call the almanac every time such a computation is needed.

One limitation I've realized, as part of using ANISE in Nyx, is that it isn't easy to change the value. So the branch with these changes will include a function on each frame called with_mu_km3_s2 that can modify the value: https://github.com/nyx-space/anise/compare/master...191-nyx-usability-alterations#diff-aa290b0f071d2018651a0fc553a8726853fac7e89a910d92f9b9cf7bbe302e16R185 .

The GM values are copied from https://github.com/nyx-space/anise/blob/master/data/gm_de431.tpc#L84 into the pck11.pca binary file, which is a platform-independent serialization of the (PlanetaryDataSet)[https://docs.rs/anise/latest/anise/structure/type.PlanetaryDataSet.html]. However, it isn't currently easy to modify a data set, so I've created #200 to address this in the future.