I’m having a numerical propagation issue in LEO when using the drag force model with the NRLMSISE00 atmospheric model. I use Orekit v11.3.3 but applying the latest fix addressing the issue #1365 that was reported in this thread.
To give you the context, only the drag perturbation is considered, with the CSSI space weather data feeding the NRLMSISE00 model. Here are the Keplerian elements of the initial orbit, which is propagated for 10 days starting from 01/01/2010 with a DP54 integrator:
a = 6687136.3 m
e = 0.
i = 30°
pa = raan = v = 0°
By comparing the behavior before/after the fix, it seems that in my case, the above fix didn’t solve the corresponding bug. More precisely, I initially had the same NAN_APPEARING_DURING_INTEGRATION exception as in the thread above. Inside the NRLMSISE00.Output.densu() method, the variable expl (and therefore densu) evaluated to NaN. Applying the fix still leads to the INFINITE_NRLMSISE00_DENSITY exception, since expl is NaN, so are many other variables in the method. The exception doesn’t occur for a higher initial orbit. So the problem is certainly located upstream and is related to the higher atmospheric density at lower altitudes. But I’m yet to pin down the exact cause of my problem. Any help would be greatly appreciated.
You’ve got an exception because your satellite has crashed on the ground
To avoid this, you can set up an altitude detector:
final AltitudeDetector altitudeDetector = new AltitudeDetector(0., earth).withHandler(new StopOnEvent());
propagator.addEventDetector(altitudeDetector);
You’ll see that the propagation stops before the 10 days limit.
With a 10kg mass, 1m² drag area and 2 for drag coeff, “landing” occurs around 2010-01-07T18:18:00.
That being said, I don’t know what is the proper behaviour to have in the library about this:
An exception when altitude is lower than 0 is not a good idea since it will break the 0m altitude detector,
We could force the density to 1000 kg/m^3 when the altitude is negative (as if crashing in the ocean). I’ve tried but this also breaks the 0m altitude detector with a “negative eccentricity” exception.
Or we could throw an exception in NRLMSIS00 when the altitude is lower than say -1000m
Thanks Maxime for your quick reply. I’ve also figured out that a “RUD” had happened since then
An altitude detector seems indeed a proper way to handle this kind of situation. I also think that negative altitudes should be handled upstream, that is before running into the INFINITE_NRLMSISE00_DENSITY exception. I’m not sure your fix was specifically meant to handle this case anyway.
