Understanding Apside Detector with maneuvers - running into math domain error

Orekit usage
1

Hi,
I have a working low thrust maneuver generation script. However i wanted to experiment with the Apside detector as I might possibly use it in my script.

I simply initialized the apside detector and added it to my propagator hoping to be able to plot the event times, to start with. However mid-way during the propagation i get ValueError: math domain error. To be noted that the only extra lines i added to my code were:

event_detector = ApsideDetector(initialOrbit).withMaxCheck(60.0).withHandler(ContinueOnEvent())
event_logger = EventsLogger()
event_detector = event_logger.monitorDetector(event_detector)
propagator_num.addEventDetector(event_detector)

I am NOT using a near-circular orbit in this case.

initialOrbit = KeplerianOrbit(7000000., 0.01, 1, 0., 0., 0., PositionAngleType.TRUE, FramesFactory.getEME2000(), initialEpoch, mu)

I also constructed an alternate apside detector extending AbstractDetector whose switching function is based on zero crossing of the spacecraft’s radial velocity as follows. However that also runs into the same error.

To be noted that, NodeDetector used the same way works as expected.

Can someone help me understand what might be going on? Does the MaxCheck interval for the event detector have to do anything with the propagation step size? Can applying low-thrust maneuvers to the propagation cause this?

For reference, i set up my propagator as follows:

    initialOrbit = KeplerianOrbit(7000000., 0.01, 1, 0., 0., 0., PositionAngleType.TRUE, 
    FramesFactory.getEME2000(), initialEpoch, mu)
    initialOrbit = CartesianOrbit(initialOrbit)
    minStep = 0.001
    maxstep = 60.0
    initStep = 60.0
    positionTolerance = .01
    tolerances = NumericalPropagator.tolerances(positionTolerance,
                                                initialOrbit,
                                                initialOrbit.getType())
    integrator = DormandPrince853Integrator(minStep, maxstep,
                                            JArray(JDouble)(tolerances[0]),
                                            JArray(JDouble)(tolerances[1]))
    integrator.setInitialStepSize(initStep)
    initialState = SpacecraftState(initialOrbit, satellite_mass)
    propagator_num = NumericalPropagator(integrator)
    propagator_num.setOrbitType(OrbitType.CARTESIAN)
    propagator_num.setInitialState(initialState)

Any help will be much appreciated.

2

Hi there,

you do have a problem. You compute your tolerances using Keplerian elements bu then your propagation orbit type is Cartesian. You need to bee consistent. I would actually suggest using Equinoctial for the integration variables.

I would also suggest using 0.001 for the position tolerance, or even less since you’re simulating maneuvers.

Finaly your max step seems a bit too large, but here it is pretty much superseded by your max check of 60s.

Cheers,
Romain.

3

@Serrof Thank you for your response.

After you pointed out the inconsistency, i realized that i missed to include the wrapping of the initial orbit to Cartesian orbit in the code i shared earlier. I have corrected the error now.

I changed my position tolerance to 0.001 as you suggested and changed the max-step for my integrator to 60s (same as apside detector max step) and also to 10s in another simulation run (lower than 60s apside detector tolerance). In both cases, i continue getting the math domain error during propagation.

P.S.: I understand your advice to use Equinoctial variable for integration to avoid singularity issues and will consider that moving forward. However, as I am NOT working with a circular orbit in this case, and since i am propagating in cartesian space, i don’t think orbit element singularity is the reason for the error. Please correct me if i am wrong.

4

I’m just saying that Cartesian variables are fast so more prompt to integration noise.

Since you’re still getting the exception, a reproductible example is what we need to investigate further I think.

Cheers,
Romain

5

Could you tell us where the error is triggered?
If you get an exception, you can use e.printStackTrace() to get the methods calls stack, this would help finding the culprit.

6

@Serrof and @luc
Turned out the issue stemmed from the event detector causing the propagator to use a variable but very small step size near the apsides. This led to numerical errors in another routine that I am using that relies on acos() .

As @Serrof had suggested, i modified the integration variables to equinoctial elements as opposed to using Cartesian. And including a check to ignore numerical errors of 1E-12 order or less in the other function helped solve the issue.