Hi,
After discussing with @bcazabonne, we found a bug in DSSTPropagator. I was trying to propagate using this propagator with impulse maneuvers and I noticed that the change in semi major axis was way too large compared to the amplitude of the maneuvers.
After digging in orekit code, I think I understand why:
- when the DSST is in osculating mode, the propagation is in mean elements, but it outputs the state in osculating elements
- when a handler reset the state, the SpacecraftState used as an input is then in osculating elements
- the handler reset the sate, so it outputs a new state (so still in osculating elements)
- the state output by EventHandler.resetState is considered as the new starting point for the propagator. It should be in mean elements, but for now it is in osculating elements. So the main effect is actually to add a short period to the mean state in the integration, which explains the change I see in the semi major axis in my case. Actually each time resetState is called, new short periods are added to the osculating state.
This effect can be seen of the following graph (I add a NodeDetector in a free drift propagation, with a simple reset state handler)
I also created a short test that shows the issue. It can be seen that the effect of the handler increases the SMA of 1 km in 0.02 seconds of propagation.
import org.hipparchus.ode.events.Action;
import org.hipparchus.ode.nonstiff.ClassicalRungeKuttaIntegrator;
import org.hipparchus.util.FastMath;
import org.junit.jupiter.api.Test;
import org.orekit.data.LazyLoadedDataContext;
import org.orekit.forces.gravity.potential.UnnormalizedSphericalHarmonicsProvider;
import org.orekit.frames.FactoryManagedFrame;
import org.orekit.orbits.CircularOrbit;
import org.orekit.orbits.PositionAngleType;
import org.orekit.propagation.PropagationType;
import org.orekit.propagation.Propagator;
import org.orekit.propagation.SpacecraftState;
import org.orekit.propagation.events.DateDetector;
import org.orekit.propagation.semianalytical.dsst.DSSTPropagator;
import org.orekit.propagation.semianalytical.dsst.forces.DSSTZonal;
import org.orekit.time.AbsoluteDate;
class DsstPropagatorWithResetStateTest {
@Test
public void testDsstOsculating() {
LazyLoadedDataContext dataContext = new LazyLoadedDataContext();
UnnormalizedSphericalHarmonicsProvider provider = dataContext.getGravityFields().getUnnormalizedProvider(2, 2);
AbsoluteDate start = new AbsoluteDate(2020, 1, 1, 0, 0, 0.0, dataContext.getTimeScales().getUTC());
FactoryManagedFrame frame = dataContext.getFrames().getEME2000();
CircularOrbit initialOrbit = new CircularOrbit(7000.0e3, 0.0, 0.001, FastMath.toRadians(97.0), 0.0, 0.0, PositionAngleType.MEAN, frame, start, provider.getMu());
DSSTPropagator dsst = new DSSTPropagator(new ClassicalRungeKuttaIntegrator(60.0), PropagationType.OSCULATING, Propagator.getDefaultLaw(dataContext.getFrames()));
dsst.addForceModel(new DSSTZonal(dataContext.getFrames().getGTOD(true), provider));
dsst.setInitialState(new SpacecraftState(initialOrbit));
dsst.addEventDetector(new DateDetector(start.shiftedBy(10.0)).withHandler((s, detector, increasing) -> Action.RESET_STATE));
SpacecraftState beforeReset = dsst.propagate(start.shiftedBy(9.99));
SpacecraftState afterReset = dsst.propagate(start.shiftedBy(10.01));
System.out.println(beforeReset.getOrbit().getA());
System.out.println(afterReset.getOrbit().getA());
}
}