Feature/330 gravity effector refactor

[juan-g-bonilla]

• Tickets addressed: Closes Refactor GravityEffector's to make gravity models a polymorphic class #330
• Review: By commit
• Merge strategy: Merge (no squash)

Description

Computation of gravity has been moved away from GravityEffector and into an abstract class GravityModel. Three new classes have been created that inherit from GravityModel: PolyhedralGravityModel, SphericalHarmonicsGravityModel, and PointMassGravityModel.

Utility functions have been added to switch gravity models and load coefficients in a single call:

planet.useSphericalHarmParams = True
simIncludeGravBody.loadGravFromFile('GGM2BData.txt', planet.spherHarm, 100)

can now be written as:

planet.useSphericalHarmonicsGravityModel('GGM2BData.txt', 100)

Apart from being a slightly shorter syntax, we avoid invalid states by condensing these calls into one:

• Calling simIncludeGravBody.loadGravFromFile without calling useSphericalHarmParams = True would lead to loading the coefficients but not using them. Confusing.
• Calling useSphericalHarmParams = True without calling simIncludeGravBody.loadGravFromFile (or manually setting the coefficients) would lead to an error during propagation.
• After this PR: calling simIncludeGravBody.loadGravFromFile before useSphericalHarmParams = True would raise an error, as planet.spherHarm raises an error if the body is not already using Spherical Harmonics.

A similar syntax was added for usePolyhedralGravityModel and usePointMassGravityModel. All scenarios were updated to use this syntax, but not tests.

In order to set any gravity model (including user-defined ones) without using the utility functions, one can do:

planet.gravityModel = MyCustomGravityModel()
planet.gravityModel.myCustomAttribute = ...

Additionally, there has been some general clean-up in GravityEffector.

Verification

Existing tests all pass except for src\simulation\dynamics\RadiationPressure\_UnitTest\test_radiation_pressure_integrated.py. Does anybody know how trust-worthy the truth data for this test is? I think that we are testing to an accuracy that is smaller than the integrator error we are making. Using the develop branch, I computed the final position of propagation using the currently set time step of 10 s and using 1 s and the integrator position error is in the order of 1 km. This is much higher than the accuracy of 1 m currently set in the test.

Documentation

The documentation file for GravityEffector was updated. Outdated information was removed. The code flow diagram was removed, as it is invalidated by the changes but also simplified slightly. Could be updated and added back, but to me this seems like an implementation detail that is covered in the text anyway.

To deprecate or not to deprecate

It's possible to deprecate the following syntax:

planet.useSphericalHarmParams = True
simIncludeGravBody.loadGravFromFile('GGM2BData.txt', planet.spherHarm, 100)

in favor of:

planet.useSphericalHarmonicsGravityModel('GGM2BData.txt', 100)

If users want to set coefficients themselves, they could do this with:

planet.gravityModel = sphericalHarmonics.SphericalHarmonicsGravityModel()
planet.gravityModel.cBar = ...

Should we deprecate the older syntax? This would lead to deprecation warnings in many user scripts and in 20 tests (which would need to be updated).

Future work

Implement new GravityModel subclasses. Issue #332.

[juan-g-bonilla]

@schaubh src\simulation\dynamics\RadiationPressure\_UnitTest\test_radiation_pressure_integrated.py is the only non-passing test. As discussed some time ago, this test uses some data that I suspect was itself generated with Basilisk because the integrator error is greater than the accuracy that is used in the test. Using the develop branch, I computed the final position of propagation using the currently set time step of 10 s and using 1 s and the integrator position error is in the order of 1 km. This is much higher than the accuracy of 1 m currently set in the test. Should we rewrite the data of the test? Change the test is some other way?

[schaubh]

@juan-g-bonilla , thanks for the update. Are you saying the test is failing on develop, or just your new branch?

[juan-g-bonilla]

@juan-g-bonilla , thanks for the update. Are you saying the test is failing on develop, or just your new branch?

The test is passing in develop. However, after the gravity effector refactor it appears that the numerical error builds up different and this test that integrates fails. If you try to reduce the integrator step size in develop, you find that the current step up of 10 s of time step produces integrator numerical errors of around 1 km, which is much greater than the accuracy of the test of 1 m. To me, this suggests that the true data in the test was not obtained from some 3rd party true source, but rather was generated by Basilisk itself.

[schaubh]

@juan-g-bonilla , I'm ok to update the truth values as you did. I looked at the new results and the change is comparatively small. As you say, we are doing different algebra to implement the same math, so there will be small round off difference that accumulate over time. The difference we are seeing here is due to this. You can keep the surface area at 1 m^2 as this causes a position difference on the order of 10's of km, much larger than our test tolerance.

The integrated test remains a consistency test and can keep the 1m tolerance for now. The compute values should not change unless we change how the math is implemented.

[juan-g-bonilla]

@schaubh Finally polished up this PR for merging. Some of the improvements here should also help fix some of the leaks Mark is experiencing. Lmk of any requested changes.

[schaubh]

Thanks. I’ll review ASAP. 😎 Blue Skies, Hanspeter Schaub [http://hanspeterschaub.info/pics/AVSlogo5-Email.png] Hanspeter Schaub, Ph.D. Schaden Leadership Chair Professor and Chair of Smead Aerospace Engineering Sciences Department LASP Senior Associate for Engineering Research and Education Director of the Autonomous Vehicle Systems (AVS) Laboratory Fellow of AIAA and AAS Address: Aerospace Engineering Sciences (ASEN) Department Colorado Center for Astrodynamics Research (CCAR) 3775 Discovery Drive, 429 UCB University of Colorado Boulder, CO 80303-0429 Office: AERO 224N Phone: +1 (303) 492-2767 Email: ***@***.******@***.***> Web: http://hanspeterschaub.info<http://hanspeterschaub.info/> On Oct 25, 2023, at 19:00, Juan Garcia Bonilla ***@***.***> wrote:  @schaubh<https://github.com/schaubh> Finally polished up this PR for merging. Some of the improvements here should also help fix some of the leaks Mark is experiencing. Lmk of any requested changes. — Reply to this email directly, view it on GitHub<#331 (comment)>, or unsubscribe<https://github.com/notifications/unsubscribe-auth/AVGSPDJGT4DVNYQ4APH25ATYBGYZBAVCNFSM6AAAAAAYFYSVDKVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMYTOOBQGI2TMNZSHE>. You are receiving this because you were mentioned.Message ID: ***@***.***>

[schaubh]

Nice refactor.

• examples/OpNavScenarios/scenariosOpNav/scenario_faultDetOpNav.py is failing now? Isn't this gravity refactor supposed to be backwards compatible? I updated the script to use the new way to set spherical harmonics.
• if this is not backwards compatible, then we should have notes in the release notes warning users to update these spherical harmonics setting part of their code?

[juan-g-bonilla]

Nice refactor.

* `examples/OpNavScenarios/scenariosOpNav/scenario_faultDetOpNav.py` is failing now?  Isn't this gravity refactor supposed to be backwards compatible?  I updated the script to use the new way to set spherical harmonics.

* if this is not backwards compatible, then we should have notes in the release notes warning users to update these spherical harmonics setting part of their code?

The issue was that when we called:

simIncludeGravBody.loadGravFromFile(bskPath + '/supportData/LocalGravData/GGM2BData.txt'
                                            , gravBodies['mars barycenter'].spherHarm
                                            , 2
                                            )

we were retrieving spherHarm from gravBodies['mars barycenter'], but at no point we had specified that gravBodies['mars barycenter'] should use spherical harmonics. After this update, this raises an error. What BSK_OpNavDynamics was doing here was loading spherical harmonics data, but never using it (it was using a point mass). This will now raise an error, since this is likely a user error (who provided data but forgot to indicate that it should be used).

With the change you published, we are saying that we should use spherical harmonics, and the loading them. So after this update, scripts that use BSK_OpNavDynamics will generate different results since now spherical harmonics are being used.

It is possible that this wasnt a design error, but that BSK_OpNavDynamics loaded the data always and then users could enable spherical harmonics and the data would be pre-loaded. If this is the case, then this is no longer supported. If we want to retain the behaviour of defaulting to a point mass, we need to not call useSphericalHarmonicsGravityModel and let users do so manually (and specying the SH data they want to use, which I believe is best anyway). Alternatively, a utility function could be added to BSK_OpNavDynamics to do this for the user.

Let me know of what you think is the best way forward. Personally, I believe this was a mistake and BSK_OpNavDynamics was always meant to use spherical harmonics, so this update will just fix that mistake. Thus, the changes while being backwards incompatible, will only be so for cases where user error is the most likely explanation.

[schaubh]

Let me know of what you think is the best way forward. Personally, I believe this was a mistake and BSK_OpNavDynamics was always meant to use spherical harmonics, so this update will just fix that mistake. Thus, the changes while being backwards incompatible, will only be so for cases where user error is the most likely explanation.

I agree that this was a mistake in the original script. Should be good to go then.

[schaubh]

Ok, good to go with this.