Add tolerance in DetectorBasedEventState check

Hipparchus development
1

Hi everyone,

I encounter a bug using AbstractIntegrator and DetectorBasedEventState in Hipparchus 4.0.1.

  • Inside AbstractIntegrator, the state is retrieved at detected date using an interpolator of ClassicalRungeKuttaStateInterpolator. As the first screenshot shows, the event date is 17016.237999999998 not rounded whereas after the interpolation, the eventState.time is truncated to 17016.238.
  • After that, the event is validated by the method doEvent of the DetectorBasedEventState, where a check of time’s equality is made (screenshot 3). This is where the snake bites its tail, as the interpolated state time is truncated and not the detected event time and the exception is raised.
  • When I look at how it is interpolated, it appears that the state is returned as it is when the time difference is lower than the least significant bit, LSB (screenshot 2).

The most straight-forward solution I think of, is to add the same tolerance of LSB, during the check, instead of the equality.
The deep correction would be to have the detection giving the right time but I suppose it is harder as the LSB is expected.

Can the first correction be foreseen ?

Thanks a lot !

Clement

Screenshot 1
Screenshot 11920×810 84 KB

2

screenshot 2
screenshot 21463×691 64.8 KB

3

pourOrekit2
pourOrekit21390×453 44.9 KB

(For some reason, since I’m a new user, I’m not allowed to upload more than one figure, so I just put in the comment section. Sorry for that)

4

Could you post a failing test case? Then we can add it to CloseEventsTest.

5

Hi @evan.ward and @clementyng

I reproduced the problem in an Orekit context. It shall be noted that in this case the propagation end epoch is equal to the maneuver thrust end. I thinks that’s an intersting starting point for investigation.

    @Test
    void testBugRk4() {
        GravityFieldFactory.addPotentialCoefficientsReader(new GRGSFormatReader("grim4s4_gr", true));
        UTCScale utc = TimeScalesFactory.getUTC();
        // Initialize spacecraft
        final CartesianOrbit cartesian = new CartesianOrbit(new PVCoordinates(new Vector3D(1897711.783316963, 4938498.102677712, -4414426.470208112),
                                                                              new Vector3D(-525.3256319053179, -4937.568767036118, -5754.938592063509)),
                                                            FramesFactory.getTOD(false),
                                                            new AbsoluteDate(2025, 6, 25, 14, 30, 32.965, utc),
                                                            Constants.WGS84_EARTH_MU);
        final double mass = 300.0;

        // Initialize propagator
        final NumericalPropagator numerical = new NumericalPropagator(new ClassicalRungeKuttaIntegrator(60.0));
        numerical.addForceModel(new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, false),
                                                                      GravityFieldFactory.getNormalizedProvider(10, 10)));
        numerical.addForceModel(new Maneuver(null,
                                             new DateBasedManeuverTriggers(new AbsoluteDate(2025, 6, 25, 14, 54, 53.247, utc),
                                                                           720.0),
                                             new BasicConstantThrustPropulsionModel(0.015, 1161, Vector3D.PLUS_I, "Cucchiethurst")));
        numerical.setInitialState(new SpacecraftState(cartesian).withMass(mass));

        // Propagate
        numerical.propagate(new AbsoluteDate(2025, 6, 25, 10, 23, 17.009, utc), new AbsoluteDate(2025, 6, 25, 15, 6, 53.247, utc));
    }

The stacktrace is:

Caused by: org.hipparchus.exception.MathRuntimeException: internal error, please fill a bug report at ``https://github.com/Hipparchus-Math/hipparchus/issues
at org.hipparchus.exception.MathRuntimeException.createInternalError(MathRuntimeException.java:69)
at org.hipparchus.ode.events.DetectorBasedEventState.check(DetectorBasedEventState.java:585)
at org.hipparchus.ode.events.DetectorBasedEventState.doEvent(DetectorBasedEventState.java:505)
at org.hipparchus.ode.AbstractIntegrator.acceptStep(AbstractIntegrator.java:369)
at org.hipparchus.ode.nonstiff.FixedStepRungeKuttaIntegrator.integrate(FixedStepRungeKuttaIntegrator.java:168)

An interesting observation is that using a Dormand-Prince integrator, there is no issue.

Thank you,
Bryan

6

Additionaly, I did this small test on epochs

      // Epochs
        final AbsoluteDate maneuverStart = new AbsoluteDate(2025, 6, 25, 14, 54, 53.247, utc);
        final AbsoluteDate maneuverEnd = maneuverStart.shiftedBy(720.0);
        final AbsoluteDate propagationStart = new AbsoluteDate(2025, 6, 25, 10, 23, 17.009, utc);
        final AbsoluteDate propagationEnd = new AbsoluteDate(2025, 6, 25, 15, 6, 53.247, utc);

        // Print some delta time
        System.out.println(propagationEnd.toString(utc));
        System.out.println(maneuverEnd.toString(utc));
        System.out.println(propagationEnd.durationFrom(maneuverEnd));
        System.out.println(propagationEnd.durationFrom(propagationStart));
        System.out.println(maneuverEnd.durationFrom(propagationStart));

The obtained results are

2025-06-25T15:06:53.246999999999999872
2025-06-25T15:06:53.246999999999999872
0.0
17016.238
17016.238
7

After investigation, the problem is in DetectorBasedEventState#doEvent. With an RK4, state.getTime() is equal to 17016.238 while pendingEventTime is 17016.237999999998.
By looking more closely at the origin of the parameter state and how it is calculated, we need to go back to AbstractODEStateInterpolator#getInterpolatedState(). In this method, two ‘state’ are used: (1) globalCurrentState whose getTime() is 17016.238 and (2) globalPreviousState whose getTime() is 17016.237999999998. The difference between the two is less than the least significant bit, so globalCurrentState is returned by the method, causing the problem mentioned above.

Given that doEvent undergoes this if done before. I think it should also apply it for safety, otherwise the current problem may occur. By looking at the Git history of the addition of the if condition between globalCurrentState and globalPreviousState, it turns out that it is a fix for Action.RESET_DERIVATIVES at the end of propagation. This is exactly what we are encountering in our case (i.e., for Action.CONTINUE there is no problem). Therefore, it could be a simple small oversight in the methods undergoing this choice between globalCurrentState and globalPreviousState (like doEvent).

Adding the safety of FastMath.ulp in doEvent fixes the problem and all Orekit + Hipparchus unit tests are green. The unit test above for propagation with the end of thrust coinciding with the end of propagation is also green.

What do you think?

I opened an issue: Can't propagate with maneuver thrust end coinciding with propagation end (#1808) · Issues · Orekit / Orekit · GitLab
We would be happy to fix it for 13.1.1 release

Best regards,
Bryan

1 Like
8

Hipparchus pull request: Fixed inconsistency between Abstract{Field}ODEStateInterpolator and {Field}DetectorBasedEventState for event time handling by BryanCazabonne · Pull Request #415 · Hipparchus-Math/hipparchus · GitHub

I currently have some difficulties to reproduce it in a 100% Hipparchus context… since it requires specific conditions. But it can be reproduced in an Orekit context with the test in the Gitlab issue description.

9

Hi Bryan,

Thanks for the thorough analysis and the solution.
If you want, in your pull request you can target the Hipparchus branch patch-4.0.2
I’ve already fixed an issue there so I think we could then release and use it as dependency for Orekit 13.1.1

Cheers,
Romain.

10

Hi Romain!

Thank you for your comment. That’s a very good idea, I’ll do it tomorrow morning.

Best regards,

Bryan

11

Hipparchus pull request to patch-4.0.2 is opened: Fixed inconsistency between Abstract{Field}ODEStateInterpolator and {Field}DetectorBasedEventState for event time handling. by BryanCazabonne · Pull Request #417 · Hipparchus-Math/hipparchus · GitHub

Orekit merge requrest to patch-13.1.1 is also opened (need Hipparchus 4.0.2 released before): Draft: Fixed inconsistency between Abstract{Field}ODEStateInterpolator and {Field}DetectorBasedEventState for event time handling. (!927) · Merge requests · Orekit / Orekit · GitLab

Do we expect new fixes to include in Hipparchus 4.0.2?

Best regards

Bryan

1 Like
12

No I don’t think so

Cheers,
Romain.

13

Great! I’ll do a post this week to propose a patch release.

3 Likes