1   /* Copyright 2002-2025 CS GROUP
2    * Licensed to CS GROUP (CS) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * CS licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *   http://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  package org.orekit.propagation.numerical;
18  
19  import org.hamcrest.CoreMatchers;
20  import org.hamcrest.MatcherAssert;
21  import org.hipparchus.CalculusFieldElement;
22  import org.hipparchus.Field;
23  import org.hipparchus.exception.LocalizedCoreFormats;
24  import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
25  import org.hipparchus.geometry.euclidean.threed.Rotation;
26  import org.hipparchus.geometry.euclidean.threed.Vector3D;
27  import org.hipparchus.linear.MatrixUtils;
28  import org.hipparchus.ode.ODEIntegrator;
29  import org.hipparchus.ode.events.Action;
30  import org.hipparchus.ode.nonstiff.AdaptiveStepsizeIntegrator;
31  import org.hipparchus.ode.nonstiff.ClassicalRungeKuttaIntegrator;
32  import org.hipparchus.ode.nonstiff.DormandPrince853Integrator;
33  import org.hipparchus.ode.nonstiff.LutherIntegrator;
34  import org.hipparchus.util.FastMath;
35  import org.junit.jupiter.api.AfterEach;
36  import org.junit.jupiter.api.Assertions;
37  import org.junit.jupiter.api.BeforeEach;
38  import org.junit.jupiter.api.Test;
39  import org.junit.jupiter.params.ParameterizedTest;
40  import org.junit.jupiter.params.provider.EnumSource;
41  import org.junit.jupiter.params.provider.ValueSource;
42  import org.mockito.Mockito;
43  import org.orekit.OrekitMatchers;
44  import org.orekit.Utils;
45  import org.orekit.attitudes.*;
46  import org.orekit.bodies.CelestialBodyFactory;
47  import org.orekit.bodies.OneAxisEllipsoid;
48  import org.orekit.errors.OrekitException;
49  import org.orekit.errors.OrekitMessages;
50  import org.orekit.forces.ForceModel;
51  import org.orekit.forces.drag.DragForce;
52  import org.orekit.forces.drag.IsotropicDrag;
53  import org.orekit.forces.gravity.HolmesFeatherstoneAttractionModel;
54  import org.orekit.forces.gravity.ThirdBodyAttraction;
55  import org.orekit.forces.gravity.potential.GRGSFormatReader;
56  import org.orekit.forces.gravity.potential.GravityFieldFactory;
57  import org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider;
58  import org.orekit.forces.gravity.potential.SHMFormatReader;
59  import org.orekit.forces.maneuvers.ConstantThrustManeuver;
60  import org.orekit.forces.maneuvers.ImpulseManeuver;
61  import org.orekit.forces.maneuvers.Maneuver;
62  import org.orekit.forces.maneuvers.propulsion.BasicConstantThrustPropulsionModel;
63  import org.orekit.forces.maneuvers.trigger.ManeuverTriggers;
64  import org.orekit.forces.maneuvers.trigger.TimeIntervalsManeuverTrigger;
65  import org.orekit.forces.radiation.IsotropicRadiationSingleCoefficient;
66  import org.orekit.forces.radiation.SolarRadiationPressure;
67  import org.orekit.frames.Frame;
68  import org.orekit.frames.FramesFactory;
69  import org.orekit.frames.LOFType;
70  import org.orekit.models.earth.atmosphere.DTM2000;
71  import org.orekit.models.earth.atmosphere.data.MarshallSolarActivityFutureEstimation;
72  import org.orekit.orbits.*;
73  import org.orekit.propagation.*;
74  import org.orekit.propagation.conversion.DormandPrince853IntegratorBuilder;
75  import org.orekit.propagation.conversion.NumericalPropagatorBuilder;
76  import org.orekit.propagation.events.*;
77  import org.orekit.propagation.events.handlers.*;
78  import org.orekit.propagation.integration.AbstractIntegratedPropagator;
79  import org.orekit.propagation.integration.AdditionalDerivativesProvider;
80  import org.orekit.propagation.integration.CombinedDerivatives;
81  import org.orekit.propagation.sampling.OrekitFixedStepHandler;
82  import org.orekit.propagation.sampling.OrekitStepHandler;
83  import org.orekit.propagation.sampling.OrekitStepInterpolator;
84  import org.orekit.time.*;
85  import org.orekit.utils.*;
86  
87  import java.io.FileNotFoundException;
88  import java.io.IOException;
89  import java.text.ParseException;
90  import java.util.ArrayList;
91  import java.util.Collections;
92  import java.util.List;
93  import java.util.concurrent.ExecutionException;
94  import java.util.concurrent.atomic.AtomicInteger;
95  import java.util.function.Consumer;
96  import java.util.stream.IntStream;
97  import java.util.stream.Stream;
98  
99  class NumericalPropagatorTest {
100 
101     private double               mu;
102     private AbsoluteDate         initDate;
103     private SpacecraftState      initialState;
104     private NumericalPropagator  propagator;
105 
106     // This test highlights the fix of: https://gitlab.orekit.org/orekit/orekit/-/issues/1808
107     @Test
108     void testPropagationEndCoincidingWithManeuverThrustEndForIssue1808() {
109         final UTCScale utc = TimeScalesFactory.getUTC();
110         // Initialize spacecraft
111         final AbsoluteDate orbitEpoch = new AbsoluteDate(2025, 6, 25, 14, 30, 32.965, utc);
112         final CartesianOrbit cartesian = new CartesianOrbit(new PVCoordinates(new Vector3D(1897711.783316963, 4938498.102677712, -4414426.470208112),
113                                                                               new Vector3D(-525.3256319053179, -4937.568767036118, -5754.938592063509)),
114                                                             FramesFactory.getTOD(false), orbitEpoch, Constants.WGS84_EARTH_MU);
115         // Maneuver
116         final AbsoluteDate maneuverStart = new AbsoluteDate(2025, 6, 25, 14, 54, 53.247, utc);
117         final AbsoluteDate maneuverEnd = maneuverStart.shiftedBy(720.0);
118 
119         // Initialize propagator
120         final NumericalPropagator numerical = new NumericalPropagator(new ClassicalRungeKuttaIntegrator(60.0));
121         numerical.addForceModel(new Maneuver(null,
122                                              TimeIntervalsManeuverTrigger.of(TimeInterval.of(maneuverStart, maneuverEnd)),
123                                              new BasicConstantThrustPropulsionModel(0.02, 1160, Vector3D.PLUS_I, "tangential")));
124         numerical.setInitialState(new SpacecraftState(cartesian));
125 
126         // Propagate
127         final AbsoluteDate propagationStart = new AbsoluteDate(2025, 6, 25, 10, 23, 17.009, utc);
128         Assertions.assertNotNull(numerical.propagate(propagationStart, maneuverEnd));
129     }
130 
131     @Test
132     void testDependsOnTimeOnlyWrong() {
133         // GIVEN
134         final AdditionalDataProvider<Boolean> dummyDataProvider = new AdditionalDataProvider<Boolean>() {
135             @Override
136             public String getName() {
137                 return "dummy";
138             }
139             @Override
140             public Boolean getAdditionalData(SpacecraftState state) {
141                 return Boolean.TRUE;
142             }
143         };
144         propagator.addAdditionalDataProvider(dummyDataProvider);
145         propagator.setInitialState(initialState.addAdditionalData(dummyDataProvider.getName(), Boolean.TRUE));
146         final EventDetector detectorCallingData = new DateDetector() {
147             @Override
148             public double g(SpacecraftState s) {
149                 s.getAdditionalData(dummyDataProvider.getName());  // this will cause an exception since detector should not depend on additional variables
150                 return super.g(s);
151             }
152         };
153         // WHEN & THEN
154         propagator.addEventDetector(detectorCallingData);
155         final AbsoluteDate targetDate = initialState.getDate().shiftedBy(1);
156         Assertions.assertThrows(OrekitException.class, () -> propagator.propagate(targetDate));
157     }
158 
159     @Test
160     void testIssue1032() {
161         Assertions.assertEquals(PropagationType.OSCULATING, propagator.getPropagationType());
162     }
163 
164     @Test
165     void testPropagateWithNonResettableManeuver() {
166         final ManeuverTriggers triggers = new ManeuverTriggers() {
167             @Override
168             public boolean isFiring(AbsoluteDate date, double[] parameters) {
169                 return false;
170             }
171 
172             @Override
173             public <T extends CalculusFieldElement<T>> boolean isFiring(FieldAbsoluteDate<T> date, T[] parameters) {
174                 return false;
175             }
176 
177             @Override
178             public Stream<EventDetector> getEventDetectors() {
179                 return Stream.empty();
180             }
181 
182             @Override
183             public <T extends CalculusFieldElement<T>> Stream<FieldEventDetector<T>> getFieldEventDetectors(Field<T> field) {
184                 return Stream.empty();
185             }
186 
187             @Override
188             public List<ParameterDriver> getParametersDrivers() {
189                 return Collections.emptyList();
190             }
191         };
192         propagator.addForceModel(new Maneuver(null, triggers, new BasicConstantThrustPropulsionModel(0., 1., Vector3D.PLUS_I, "")));
193         propagator.setupMatricesComputation("stm", MatrixUtils.createRealIdentityMatrix(6), new DoubleArrayDictionary());
194         propagator.setUpStmAndJacobianGenerators();
195     }
196 
197     @Test
198     void testIssue879() {
199         // GIVEN
200         final EphemerisGenerator generator = propagator.getEphemerisGenerator();
201         // WHEN
202         propagator.clearEphemerisGenerators();
203         propagator.propagate(initDate.shiftedBy(1));
204         // THEN
205         Assertions.assertThrows(NullPointerException.class, generator::getGeneratedEphemeris);
206     }
207 
208     @Test
209     void testIssue1826() {
210         // GIVEN
211         propagator.setupMatricesComputation("stm", null, null);
212         propagator.setResetAtEnd(false);
213         final AbsoluteDate targetDate = initDate.shiftedBy(1);
214         propagator.propagate(targetDate);
215         // WHEN
216         propagator.clearMatricesComputation();
217         // THEN
218         Assertions.assertDoesNotThrow(() -> propagator.propagate(targetDate));
219     }
220 
221     @Test
222     void testEventsWithTimeRangePropagation() {
223         final AtomicInteger counter = new AtomicInteger(0);
224         final double dt = 60.0;
225         final EventDetector singleDetector = new DateDetector(initDate.shiftedBy(dt / 2)).
226                                              withHandler((state, detector, increasing) -> {
227                                                  counter.incrementAndGet();
228                                                  return Action.CONTINUE;
229                                              });
230         ForceModel force = new ForceModelAdapter() {
231             @Override
232             public Stream<EventDetector> getEventDetectors() {
233                 return Stream.of(singleDetector);
234             }
235         };
236 
237         // action
238         propagator.setOrbitType(OrbitType.CARTESIAN);
239         propagator.addForceModel(force);
240         AbsoluteDate target = initDate.shiftedBy(dt);
241         propagator.propagate(initDate.shiftedBy(1.0), target);
242 
243         Assertions.assertEquals(1, counter.intValue());
244 
245     }
246 
247     @Test
248     void testForceModelInitialized() {
249         // setup
250         // mutable holders
251         SpacecraftState[] actualState = new SpacecraftState[1];
252         AbsoluteDate[] actualDate = new AbsoluteDate[1];
253         ForceModel force = new ForceModelAdapter() {
254             @Override
255             public void init(SpacecraftState initialState, AbsoluteDate target) {
256                 actualState[0] = initialState;
257                 actualDate[0] = target;
258             }
259         };
260 
261         // action
262         propagator.setOrbitType(OrbitType.CARTESIAN);
263         propagator.addForceModel(force);
264         AbsoluteDate target = initDate.shiftedBy(60);
265         propagator.propagate(target);
266 
267         // verify
268         MatcherAssert.assertThat(actualDate[0], CoreMatchers.is(target));
269         MatcherAssert.assertThat(actualState[0].getDate().durationFrom(initDate),
270                 CoreMatchers.is(0.0));
271         MatcherAssert.assertThat(actualState[0].getPVCoordinates(),
272                 OrekitMatchers.pvIs(initialState.getPVCoordinates()));
273     }
274 
275     @Test
276     @Deprecated
277     void testTolerancesOrbitdV() {
278         // GIVEN
279         final double dP = 1e-3;
280         final Vector3D position = new Vector3D(7.0e6, 1.0e6, 4.0e6);
281         final Vector3D velocity = new Vector3D(-500.0, 8000.0, 1000.0);
282         final Orbit orbit = new CartesianOrbit(new TimeStampedPVCoordinates(AbsoluteDate.ARBITRARY_EPOCH, position,
283                 velocity, Vector3D.ZERO), FramesFactory.getGCRF(), Constants.EGM96_EARTH_MU);
284         // WHEN
285         final double[][] tolerancesWithDv = NumericalPropagator.tolerances(dP, 1e-6, orbit, OrbitType.CARTESIAN);
286         // THEN
287         final double[][] tolerances = ToleranceProvider.getDefaultToleranceProvider(dP).getTolerances(orbit, OrbitType.CARTESIAN);
288         for (int i = 0; i < 3; i++) {
289             Assertions.assertEquals(tolerances[0][i], tolerancesWithDv[0][i]);
290         }
291     }
292 
293     @Deprecated
294     @ParameterizedTest
295     @EnumSource(OrbitType.class)
296     void testTolerancesOrbit(final OrbitType orbitType) {
297         // GIVEN
298         final double dP = 1e-3;
299         final Vector3D position = new Vector3D(7.0e6, 1.0e6, 4.0e6);
300         final Vector3D velocity = new Vector3D(-500.0, 8000.0, 1000.0);
301         final Orbit orbit = new CartesianOrbit(new TimeStampedPVCoordinates(AbsoluteDate.ARBITRARY_EPOCH, position,
302                 velocity, Vector3D.ZERO), FramesFactory.getGCRF(), Constants.EGM96_EARTH_MU);
303         // WHEN
304         final double[][] actualTolerances = NumericalPropagator.tolerances(dP, orbit, orbitType);
305         // THEN
306         final double[][] expectedTolerances = ToleranceProvider.getDefaultToleranceProvider(dP).getTolerances(orbit, orbitType);
307         Assertions.assertArrayEquals(expectedTolerances[0], actualTolerances[0]);
308         Assertions.assertArrayEquals(expectedTolerances[1], actualTolerances[1]);
309     }
310 
311     @Deprecated
312     @Test
313     void testTolerances() {
314         // GIVEN
315         final double dP = 1e-3;
316         final Vector3D position = new Vector3D(7.0e6, 1.0e6, 4.0e6);
317         final Vector3D velocity = new Vector3D(-500.0, 8000.0, 1000.0);
318         final Orbit orbit = new CartesianOrbit(new TimeStampedPVCoordinates(AbsoluteDate.ARBITRARY_EPOCH, position,
319                 velocity, Vector3D.ZERO), FramesFactory.getGCRF(), Constants.EGM96_EARTH_MU);
320         // WHEN
321         final double[][] orbitTolerances = NumericalPropagator.tolerances(dP, orbit, OrbitType.CARTESIAN);
322         // THEN
323         final double[][] pvTolerances = ToleranceProvider.getDefaultToleranceProvider(dP).getTolerances(new AbsolutePVCoordinates(orbit.getFrame(),
324                 new TimeStampedPVCoordinates(orbit.getDate(), position, velocity)));
325         for (int i = 0; i < 3; i++) {
326             Assertions.assertEquals(pvTolerances[0][i], orbitTolerances[0][i]);
327             Assertions.assertEquals(pvTolerances[1][i], orbitTolerances[1][i]);
328         }
329     }
330 
331     @Test
332     void testEphemerisModeWithHandler() {
333         // setup
334         AbsoluteDate end = initDate.shiftedBy(90 * 60);
335 
336         // action
337         final List<SpacecraftState> states = new ArrayList<>();
338         EphemerisGenerator generator = propagator.getEphemerisGenerator();
339         propagator.setStepHandler(interpolator -> states.add(interpolator.getCurrentState()));
340         propagator.propagate(end);
341         final BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
342 
343         //verify
344         Assertions.assertTrue(states.size() > 10); // got some data
345         for (SpacecraftState state : states) {
346             PVCoordinates actual =
347                     ephemeris.propagate(state.getDate()).getPVCoordinates();
348             MatcherAssert.assertThat(actual, OrekitMatchers.pvIs(state.getPVCoordinates()));
349         }
350     }
351 
352     /** test for issue #238 */
353     @Test
354     void testEventAtEndOfEphemeris() {
355         // setup
356         // choose duration that will round up when expressed as a double
357         AbsoluteDate end = initDate.shiftedBy(100)
358                 .shiftedBy(3 * FastMath.ulp(100.0) / 4);
359         final EphemerisGenerator generator = propagator.getEphemerisGenerator();
360         propagator.propagate(end);
361         BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
362         CountAndContinue handler = new CountAndContinue();
363         DateDetector detector = new DateDetector(end).
364                                 withMaxCheck(10).
365                                 withThreshold(1e-9).
366                                 withMinGap(10).
367                                 withHandler(handler);
368         // propagation works fine w/o event detector, but breaks with it
369         ephemeris.addEventDetector(detector);
370 
371         //action
372         // fails when this throws an "out of range date for ephemerides"
373         SpacecraftState actual = ephemeris.propagate(end);
374 
375         //verify
376         Assertions.assertEquals(actual.getDate().durationFrom(end), 0.0, 0.0);
377         Assertions.assertEquals(1, handler.getCount());
378     }
379 
380     /** test for issue #238 */
381     @Test
382     void testEventAtBeginningOfEphemeris() {
383         // setup
384         // choose duration that will round up when expressed as a double
385         AbsoluteDate end = initDate.shiftedBy(100)
386                 .shiftedBy(3 * FastMath.ulp(100.0) / 4);
387         final EphemerisGenerator generator = propagator.getEphemerisGenerator();
388         propagator.propagate(end);
389         BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
390         CountAndContinue handler = new CountAndContinue();
391         // events directly on propagation start date are not triggered,
392         // so move the event date slightly after
393         AbsoluteDate eventDate = initDate.shiftedBy(FastMath.ulp(100.0) / 10.0);
394         DateDetector detector = new DateDetector(eventDate).
395                                 withMaxCheck(10).
396                                 withThreshold(1e-9).
397                                 withMinGap(10).
398                                 withHandler(handler);
399         // propagation works fine w/o event detector, but breaks with it
400         ephemeris.addEventDetector(detector);
401 
402         // action + verify
403         // propagate forward
404         Assertions.assertEquals(ephemeris.propagate(end).getDate().durationFrom(end), 0.0, 0.0);
405         // propagate backward
406         Assertions.assertEquals(ephemeris.propagate(initDate).getDate().durationFrom(initDate), 0.0, 0.0);
407         Assertions.assertEquals(2, handler.getCount());
408     }
409 
410     /**
411      * check propagation succeeds when two events are within the tolerance of
412      * each other.
413      */
414     @Test
415     void testCloseEventDates() {
416         // setup
417         DateDetector d1 = new DateDetector(initDate.shiftedBy(15)).
418                           withMaxCheck(10).
419                           withThreshold(1).
420                           withHandler(new ContinueOnEvent());
421         DateDetector d2 = new DateDetector(initDate.shiftedBy(15.5)).
422                           withMaxCheck(10).
423                           withThreshold(1).
424                           withHandler(new ContinueOnEvent());
425         propagator.addEventDetector(d1);
426         propagator.addEventDetector(d2);
427 
428         //action
429         AbsoluteDate end = initDate.shiftedBy(30);
430         SpacecraftState actual = propagator.propagate(end);
431 
432         //verify
433         Assertions.assertEquals(actual.getDate().durationFrom(end), 0.0, 0.0);
434     }
435 
436     @Test
437     void testEphemerisDates() {
438         //setup
439         TimeScale tai = TimeScalesFactory.getTAI();
440         AbsoluteDate initialDate = new AbsoluteDate("2015-07-01", tai);
441         AbsoluteDate startDate = new AbsoluteDate("2015-07-03", tai).shiftedBy(-0.1);
442         AbsoluteDate endDate = new AbsoluteDate("2015-07-04", tai);
443         Frame eci = FramesFactory.getGCRF();
444         KeplerianOrbit orbit = new KeplerianOrbit(
445                 600e3 + Constants.WGS84_EARTH_EQUATORIAL_RADIUS, 0, 0, 0, 0, 0,
446                 PositionAngleType.TRUE, eci, initialDate, mu);
447         OrbitType type = OrbitType.CARTESIAN;
448         double[][] tol = ToleranceProvider.of(CartesianToleranceProvider.of(1e-3)).getTolerances(orbit, type, NumericalPropagator.DEFAULT_POSITION_ANGLE_TYPE);
449         NumericalPropagator prop = new NumericalPropagator(
450                 new DormandPrince853Integrator(0.1, 500, tol[0], tol[1]));
451         prop.setOrbitType(type);
452         prop.resetInitialState(new SpacecraftState(new CartesianOrbit(orbit)));
453 
454         //action
455         final EphemerisGenerator generator = prop.getEphemerisGenerator();
456         prop.propagate(startDate, endDate);
457         BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
458 
459         //verify
460         TimeStampedPVCoordinates actualPV = ephemeris.getPVCoordinates(startDate, eci);
461         TimeStampedPVCoordinates expectedPV = orbit.getPVCoordinates(startDate, eci);
462         MatcherAssert.assertThat(actualPV.getPosition(),
463                 OrekitMatchers.vectorCloseTo(expectedPV.getPosition(), 1.0));
464         MatcherAssert.assertThat(actualPV.getVelocity(),
465                 OrekitMatchers.vectorCloseTo(expectedPV.getVelocity(), 1.0));
466         MatcherAssert.assertThat(ephemeris.getMinDate().durationFrom(startDate),
467                 OrekitMatchers.closeTo(0, 0));
468         MatcherAssert.assertThat(ephemeris.getMaxDate().durationFrom(endDate),
469                 OrekitMatchers.closeTo(0, 0));
470         //test date
471         AbsoluteDate date = endDate.shiftedBy(-0.11);
472         Assertions.assertEquals(
473                 ephemeris.propagate(date).getDate().durationFrom(date), 0, 0);
474     }
475 
476     @Test
477     void testEphemerisDatesBackward() {
478         //setup
479         TimeScale tai = TimeScalesFactory.getTAI();
480         AbsoluteDate initialDate = new AbsoluteDate("2015-07-05", tai);
481         AbsoluteDate startDate = new AbsoluteDate("2015-07-03", tai).shiftedBy(-0.1);
482         AbsoluteDate endDate = new AbsoluteDate("2015-07-04", tai);
483         Frame eci = FramesFactory.getGCRF();
484         KeplerianOrbit orbit = new KeplerianOrbit(
485                 600e3 + Constants.WGS84_EARTH_EQUATORIAL_RADIUS, 0, 0, 0, 0, 0,
486                 PositionAngleType.TRUE, eci, initialDate, mu);
487         OrbitType type = OrbitType.CARTESIAN;
488         double[][] tol = ToleranceProvider.of(CartesianToleranceProvider.of(1e-3)).getTolerances(orbit, type);
489         NumericalPropagator prop = new NumericalPropagator(
490                 new DormandPrince853Integrator(0.1, 500, tol[0], tol[1]));
491         prop.setOrbitType(type);
492         prop.resetInitialState(new SpacecraftState(new CartesianOrbit(orbit)));
493 
494         //action
495         final EphemerisGenerator generator = prop.getEphemerisGenerator();
496         prop.propagate(endDate, startDate);
497         BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
498 
499         //verify
500         TimeStampedPVCoordinates actualPV = ephemeris.getPVCoordinates(startDate, eci);
501         TimeStampedPVCoordinates expectedPV = orbit.getPVCoordinates(startDate, eci);
502         MatcherAssert.assertThat(actualPV.getPosition(),
503                 OrekitMatchers.vectorCloseTo(expectedPV.getPosition(), 1.0));
504         MatcherAssert.assertThat(actualPV.getVelocity(),
505                 OrekitMatchers.vectorCloseTo(expectedPV.getVelocity(), 1.0));
506         MatcherAssert.assertThat(ephemeris.getMinDate().durationFrom(startDate),
507                 OrekitMatchers.closeTo(0, 0));
508         MatcherAssert.assertThat(ephemeris.getMaxDate().durationFrom(endDate),
509                 OrekitMatchers.closeTo(0, 0));
510         //test date
511         AbsoluteDate date = endDate.shiftedBy(-0.11);
512         Assertions.assertEquals(
513                 ephemeris.propagate(date).getDate().durationFrom(date), 0, 0);
514     }
515 
516     @Test
517     void testNoExtrapolation() {
518 
519         // Propagate of the initial at the initial date
520         final SpacecraftState finalState = propagator.propagate(initDate);
521 
522         // Initial orbit definition
523         final Vector3D initialPosition = initialState.getPosition();
524         final Vector3D initialVelocity = initialState.getVelocity();
525 
526         // Final orbit definition
527         final Vector3D finalPosition   = finalState.getPosition();
528         final Vector3D finalVelocity   = finalState.getVelocity();
529 
530         // Check results
531         Assertions.assertEquals(initialPosition.getX(), finalPosition.getX(), 1.0e-10);
532         Assertions.assertEquals(initialPosition.getY(), finalPosition.getY(), 1.0e-10);
533         Assertions.assertEquals(initialPosition.getZ(), finalPosition.getZ(), 1.0e-10);
534         Assertions.assertEquals(initialVelocity.getX(), finalVelocity.getX(), 1.0e-10);
535         Assertions.assertEquals(initialVelocity.getY(), finalVelocity.getY(), 1.0e-10);
536         Assertions.assertEquals(initialVelocity.getZ(), finalVelocity.getZ(), 1.0e-10);
537 
538     }
539 
540     @Test
541     void testNotInitialised1() {
542         Assertions.assertThrows(OrekitException.class, () -> {
543             final AbstractIntegratedPropagator notInitialised =
544                     new NumericalPropagator(new ClassicalRungeKuttaIntegrator(10.0));
545             notInitialised.propagate(AbsoluteDate.J2000_EPOCH);
546         });
547     }
548 
549     @Test
550     void testNotInitialised2() {
551         Assertions.assertThrows(OrekitException.class, () -> {
552             final AbstractIntegratedPropagator notInitialised =
553                     new NumericalPropagator(new ClassicalRungeKuttaIntegrator(10.0));
554             notInitialised.propagate(AbsoluteDate.J2000_EPOCH, AbsoluteDate.J2000_EPOCH.shiftedBy(3600));
555         });
556     }
557 
558     @Test
559     void testKepler() {
560 
561         // Propagation of the initial at t + dt
562         final double dt = 3200;
563         final Orbit finalOrbit =
564             propagator.propagate(initDate.shiftedBy(-60), initDate.shiftedBy(dt)).getOrbit();
565 
566         // Check results
567         final Orbit initialOrbit = initialState.getOrbit();
568         final double n = FastMath.sqrt(initialOrbit.getMu() / initialOrbit.getA()) / initialOrbit.getA();
569         Assertions.assertEquals(initialOrbit.getA(),    finalOrbit.getA(),    1.0e-10);
570         Assertions.assertEquals(initialOrbit.getEquinoctialEx(),    finalOrbit.getEquinoctialEx(),    1.0e-10);
571         Assertions.assertEquals(initialOrbit.getEquinoctialEy(),    finalOrbit.getEquinoctialEy(),    1.0e-10);
572         Assertions.assertEquals(initialOrbit.getHx(),    finalOrbit.getHx(),    1.0e-10);
573         Assertions.assertEquals(initialOrbit.getHy(),    finalOrbit.getHy(),    1.0e-10);
574         Assertions.assertEquals(initialOrbit.getLM() + n * dt, finalOrbit.getLM(), 2.0e-9);
575 
576     }
577 
578     @Test
579     void testCartesian() {
580 
581         // Propagation of the initial at t + dt
582         final double dt = 3200;
583         propagator.setOrbitType(OrbitType.CARTESIAN);
584         final PVCoordinates finalState =
585             propagator.propagate(initDate.shiftedBy(dt)).getPVCoordinates();
586         final Vector3D pFin = finalState.getPosition();
587         final Vector3D vFin = finalState.getVelocity();
588 
589         // Check results
590         final PVCoordinates reference = initialState.shiftedBy(dt).getPVCoordinates();
591         final Vector3D pRef = reference.getPosition();
592         final Vector3D vRef = reference.getVelocity();
593         Assertions.assertEquals(0, pRef.subtract(pFin).getNorm(), 2e-4);
594         Assertions.assertEquals(0, vRef.subtract(vFin).getNorm(), 7e-8);
595 
596     }
597 
598     @Test
599     void testPropagationTypesElliptical() throws ParseException, IOException {
600      // setup
601         AbsoluteDate         initDate  = new AbsoluteDate();
602         SpacecraftState     initialState;
603         final Vector3D position = new Vector3D(7.0e6, 1.0e6, 4.0e6);
604         final Vector3D velocity = new Vector3D(-500.0, 8000.0, 1000.0);
605         initDate = AbsoluteDate.J2000_EPOCH;
606 
607         final Orbit orbit = new EquinoctialOrbit(new PVCoordinates(position,  velocity),
608                                                  FramesFactory.getEME2000(), initDate, mu);
609         initialState = new SpacecraftState(orbit);
610         OrbitType type = OrbitType.EQUINOCTIAL;
611         double[][] tolerance = ToleranceProvider.of(CartesianToleranceProvider.of(0.001)).getTolerances(orbit, type);
612         AdaptiveStepsizeIntegrator integrator =
613                 new DormandPrince853Integrator(0.001, 200, tolerance[0], tolerance[1]);
614         integrator.setInitialStepSize(60);
615         propagator = new NumericalPropagator(integrator);
616         propagator.setOrbitType(type);
617         propagator.setInitialState(initialState);
618 
619         ForceModel gravityField =
620             new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true),
621                                                   GravityFieldFactory.getNormalizedProvider(5, 5));
622         propagator.addForceModel(gravityField);
623 
624         // Propagation of the initial at t + dt
625         final PVCoordinates pv = initialState.getPVCoordinates();
626         final double dP = 0.001;
627         final double dV = initialState.getOrbit().getMu() * dP /
628                           (pv.getPosition().getNormSq() * pv.getVelocity().getNorm());
629 
630         final PVCoordinates pvcM = propagateInType(initialState, dP, OrbitType.CARTESIAN,   PositionAngleType.MEAN);
631         final PVCoordinates pviM = propagateInType(initialState, dP, OrbitType.CIRCULAR,    PositionAngleType.MEAN);
632         final PVCoordinates pveM = propagateInType(initialState, dP, OrbitType.EQUINOCTIAL, PositionAngleType.MEAN);
633         final PVCoordinates pvkM = propagateInType(initialState, dP, OrbitType.KEPLERIAN,   PositionAngleType.MEAN);
634 
635         final PVCoordinates pvcE = propagateInType(initialState, dP, OrbitType.CARTESIAN,   PositionAngleType.ECCENTRIC);
636         final PVCoordinates pviE = propagateInType(initialState, dP, OrbitType.CIRCULAR,    PositionAngleType.ECCENTRIC);
637         final PVCoordinates pveE = propagateInType(initialState, dP, OrbitType.EQUINOCTIAL, PositionAngleType.ECCENTRIC);
638         final PVCoordinates pvkE = propagateInType(initialState, dP, OrbitType.KEPLERIAN,   PositionAngleType.ECCENTRIC);
639 
640         final PVCoordinates pvcT = propagateInType(initialState, dP, OrbitType.CARTESIAN,   PositionAngleType.TRUE);
641         final PVCoordinates pviT = propagateInType(initialState, dP, OrbitType.CIRCULAR,    PositionAngleType.TRUE);
642         final PVCoordinates pveT = propagateInType(initialState, dP, OrbitType.EQUINOCTIAL, PositionAngleType.TRUE);
643         final PVCoordinates pvkT = propagateInType(initialState, dP, OrbitType.KEPLERIAN,   PositionAngleType.TRUE);
644 
645         Assertions.assertEquals(0, pvcM.getPosition().subtract(pveT.getPosition()).getNorm() / dP, 3.1);
646         Assertions.assertEquals(0, pvcM.getVelocity().subtract(pveT.getVelocity()).getNorm() / dV, 2.1);
647         Assertions.assertEquals(0, pviM.getPosition().subtract(pveT.getPosition()).getNorm() / dP, 0.8);
648         Assertions.assertEquals(0, pviM.getVelocity().subtract(pveT.getVelocity()).getNorm() / dV, 0.5);
649         Assertions.assertEquals(0, pvkM.getPosition().subtract(pveT.getPosition()).getNorm() / dP, 0.5);
650         Assertions.assertEquals(0, pvkM.getVelocity().subtract(pveT.getVelocity()).getNorm() / dV, 0.4);
651         Assertions.assertEquals(0, pveM.getPosition().subtract(pveT.getPosition()).getNorm() / dP, 0.3);
652         Assertions.assertEquals(0, pveM.getVelocity().subtract(pveT.getVelocity()).getNorm() / dV, 0.3);
653 
654         Assertions.assertEquals(0, pvcE.getPosition().subtract(pveT.getPosition()).getNorm() / dP, 3.1);
655         Assertions.assertEquals(0, pvcE.getVelocity().subtract(pveT.getVelocity()).getNorm() / dV, 2.1);
656         Assertions.assertEquals(0, pviE.getPosition().subtract(pveT.getPosition()).getNorm() / dP, 0.1);
657         Assertions.assertEquals(0, pviE.getVelocity().subtract(pveT.getVelocity()).getNorm() / dV, 0.09);
658         Assertions.assertEquals(0, pvkE.getPosition().subtract(pveT.getPosition()).getNorm() / dP, 0.5);
659         Assertions.assertEquals(0, pvkE.getVelocity().subtract(pveT.getVelocity()).getNorm() / dV, 0.3);
660         Assertions.assertEquals(0, pveE.getPosition().subtract(pveT.getPosition()).getNorm() / dP, 0.2);
661         Assertions.assertEquals(0, pveE.getVelocity().subtract(pveT.getVelocity()).getNorm() / dV, 0.2);
662 
663         Assertions.assertEquals(0, pvcT.getPosition().subtract(pveT.getPosition()).getNorm() / dP, 3.1);
664         Assertions.assertEquals(0, pvcT.getVelocity().subtract(pveT.getVelocity()).getNorm() / dV, 2.1);
665         Assertions.assertEquals(0, pviT.getPosition().subtract(pveT.getPosition()).getNorm() / dP, 0.4);
666         Assertions.assertEquals(0, pviT.getVelocity().subtract(pveT.getVelocity()).getNorm() / dV, 0.2);
667         Assertions.assertEquals(0, pvkT.getPosition().subtract(pveT.getPosition()).getNorm() / dP, 0.5);
668         Assertions.assertEquals(0, pvkT.getVelocity().subtract(pveT.getVelocity()).getNorm() / dV, 0.2);
669 
670     }
671 
672     @Test
673     void testPropagationTypesHyperbolic() throws ParseException, IOException {
674 
675         SpacecraftState state =
676             new SpacecraftState(new KeplerianOrbit(-10000000.0, 2.5, 0.3, 0, 0, 0.0,
677                                                    PositionAngleType.TRUE,
678                                                    FramesFactory.getEME2000(), initDate,
679                                                    mu));
680 
681         ForceModel gravityField =
682             new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true),
683                                                   GravityFieldFactory.getNormalizedProvider(5, 5));
684         propagator.addForceModel(gravityField);
685 
686         // Propagation of the initial at t + dt
687         final PVCoordinates pv = state.getPVCoordinates();
688         final double dP = 0.001;
689         final double dV = state.getOrbit().getMu() * dP /
690                           (pv.getPosition().getNormSq() * pv.getVelocity().getNorm());
691 
692         final PVCoordinates pvcM = propagateInType(state, dP, OrbitType.CARTESIAN, PositionAngleType.MEAN);
693         final PVCoordinates pvkM = propagateInType(state, dP, OrbitType.KEPLERIAN, PositionAngleType.MEAN);
694 
695         final PVCoordinates pvcE = propagateInType(state, dP, OrbitType.CARTESIAN, PositionAngleType.ECCENTRIC);
696         final PVCoordinates pvkE = propagateInType(state, dP, OrbitType.KEPLERIAN, PositionAngleType.ECCENTRIC);
697 
698         final PVCoordinates pvcT = propagateInType(state, dP, OrbitType.CARTESIAN, PositionAngleType.TRUE);
699         final PVCoordinates pvkT = propagateInType(state, dP, OrbitType.KEPLERIAN, PositionAngleType.TRUE);
700 
701         Assertions.assertEquals(0, pvcM.getPosition().subtract(pvkT.getPosition()).getNorm() / dP, 0.4);
702         Assertions.assertEquals(0, pvcM.getVelocity().subtract(pvkT.getVelocity()).getNorm() / dV, 0.6);
703         Assertions.assertEquals(0, pvkM.getPosition().subtract(pvkT.getPosition()).getNorm() / dP, 0.4);
704         Assertions.assertEquals(0, pvkM.getVelocity().subtract(pvkT.getVelocity()).getNorm() / dV, 0.6);
705 
706         Assertions.assertEquals(0, pvcE.getPosition().subtract(pvkT.getPosition()).getNorm() / dP, 0.3);
707         Assertions.assertEquals(0, pvcE.getVelocity().subtract(pvkT.getVelocity()).getNorm() / dV, 0.6);
708         Assertions.assertEquals(0, pvkE.getPosition().subtract(pvkT.getPosition()).getNorm() / dP, 0.1);
709         Assertions.assertEquals(0, pvkE.getVelocity().subtract(pvkT.getVelocity()).getNorm() / dV, 0.04);
710 
711         Assertions.assertEquals(0, pvcT.getPosition().subtract(pvkT.getPosition()).getNorm() / dP, 0.3);
712         Assertions.assertEquals(0, pvcT.getVelocity().subtract(pvkT.getVelocity()).getNorm() / dV, 0.6);
713 
714     }
715 
716     private PVCoordinates propagateInType(SpacecraftState state, double dP,
717                                           OrbitType type, PositionAngleType angle)
718         {
719 
720         final double dt = 3200;
721         final double minStep = 0.001;
722         final double maxStep = 1000;
723 
724         double[][] tol = ToleranceProvider.of(CartesianToleranceProvider.of(dP)).getTolerances(state.getOrbit(), type,
725                 angle);
726         AdaptiveStepsizeIntegrator integrator =
727                 new DormandPrince853Integrator(minStep, maxStep, tol[0], tol[1]);
728         NumericalPropagator newPropagator = new NumericalPropagator(integrator);
729         newPropagator.setOrbitType(type);
730         newPropagator.setPositionAngleType(angle);
731         newPropagator.setInitialState(state);
732         for (ForceModel force: propagator.getAllForceModels()) {
733             newPropagator.addForceModel(force);
734         }
735         return newPropagator.propagate(state.getDate().shiftedBy(dt)).getPVCoordinates();
736 
737     }
738 
739     @Test
740     void testException() {
741         Assertions.assertThrows(OrekitException.class, () -> {
742             propagator.setStepHandler(new OrekitStepHandler() {
743                 private int countDown = 3;
744                 private AbsoluteDate previousCall = null;
745                 public void init(SpacecraftState s0, AbsoluteDate t) {
746                 }
747                 public void handleStep(OrekitStepInterpolator interpolator) {
748                     if (previousCall != null) {
749                         Assertions.assertTrue(interpolator.getCurrentState().getDate().compareTo(previousCall) < 0);
750                     }
751                     if (--countDown == 0) {
752                         throw new OrekitException(LocalizedCoreFormats.SIMPLE_MESSAGE, "dummy error");
753                     }
754                 }
755             });
756             propagator.propagate(initDate.shiftedBy(-3600));
757         });
758     }
759 
760     @Test
761     void testStopEvent() {
762         final AbsoluteDate stopDate = initDate.shiftedBy(1000);
763         CheckingHandler checking = new CheckingHandler(Action.STOP);
764         propagator.addEventDetector(new DateDetector(stopDate).withHandler(checking));
765         Assertions.assertEquals(1, propagator.getEventDetectors().size());
766         checking.assertEvent(false);
767         final SpacecraftState finalState = propagator.propagate(initDate.shiftedBy(3200));
768         checking.assertEvent(true);
769         Assertions.assertEquals(0, finalState.getDate().durationFrom(stopDate), 1.0e-10);
770         propagator.clearEventsDetectors();
771         Assertions.assertEquals(0, propagator.getEventDetectors().size());
772         Assertions.assertTrue(checking.isFinished);
773     }
774 
775     @Test
776     void testResetStateEvent() {
777         final AbsoluteDate resetDate = initDate.shiftedBy(1000);
778         CheckingHandler checking = new CheckingHandler(Action.RESET_STATE) {
779             public SpacecraftState resetState(EventDetector detector, SpacecraftState oldState) {
780                 return new SpacecraftState(oldState.getOrbit(), oldState.getAttitude()).withMass(oldState.getMass() - 200.0);
781             }
782         };
783         propagator.addEventDetector(new DateDetector(resetDate).withHandler(checking));
784         checking.assertEvent(false);
785         final SpacecraftState finalState = propagator.propagate(initDate.shiftedBy(3200));
786         checking.assertEvent(true);
787         Assertions.assertEquals(initialState.getMass() - 200, finalState.getMass(), 1.0e-10);
788     }
789 
790     @Test
791     void testResetDerivativesEvent() {
792         final AbsoluteDate resetDate = initDate.shiftedBy(1000);
793         CheckingHandler checking = new CheckingHandler(Action.RESET_DERIVATIVES);
794         propagator.addEventDetector(new DateDetector(resetDate).withHandler(checking));
795         final double dt = 3200;
796         checking.assertEvent(false);
797         Assertions.assertEquals(0.0, propagator.getInitialState().getDate().durationFrom(initDate), 1.0e-10);
798         propagator.setResetAtEnd(true);
799         final Orbit finalOrbit =
800             propagator.propagate(initDate.shiftedBy(dt)).getOrbit();
801         Assertions.assertEquals(dt, propagator.getInitialState().getDate().durationFrom(initDate), 1.0e-10);
802         checking.assertEvent(true);
803         final Orbit initialOrbit = initialState.getOrbit();
804         final double n = FastMath.sqrt(initialOrbit.getMu() / initialOrbit.getA()) / initialOrbit.getA();
805         Assertions.assertEquals(initialOrbit.getA(),    finalOrbit.getA(),    1.0e-10);
806         Assertions.assertEquals(initialOrbit.getEquinoctialEx(),    finalOrbit.getEquinoctialEx(),    1.0e-10);
807         Assertions.assertEquals(initialOrbit.getEquinoctialEy(),    finalOrbit.getEquinoctialEy(),    1.0e-10);
808         Assertions.assertEquals(initialOrbit.getHx(),    finalOrbit.getHx(),    1.0e-10);
809         Assertions.assertEquals(initialOrbit.getHy(),    finalOrbit.getHy(),    1.0e-10);
810         Assertions.assertEquals(initialOrbit.getLM() + n * dt, finalOrbit.getLM(), 6.0e-10);
811     }
812 
813     @Test
814     void testContinueEvent() {
815         final AbsoluteDate resetDate = initDate.shiftedBy(1000);
816         CheckingHandler checking = new CheckingHandler(Action.CONTINUE);
817         propagator.addEventDetector(new DateDetector(resetDate).withHandler(checking));
818         final double dt = 3200;
819         checking.assertEvent(false);
820         Assertions.assertEquals(0.0, propagator.getInitialState().getDate().durationFrom(initDate), 1.0e-10);
821         propagator.setResetAtEnd(false);
822         final Orbit finalOrbit =
823             propagator.propagate(initDate.shiftedBy(dt)).getOrbit();
824         Assertions.assertEquals(0.0, propagator.getInitialState().getDate().durationFrom(initDate), 1.0e-10);
825         checking.assertEvent(true);
826         final Orbit initialOrbit = initialState.getOrbit();
827         final double n = FastMath.sqrt(initialOrbit.getMu() / initialOrbit.getA()) / initialOrbit.getA();
828         Assertions.assertEquals(initialOrbit.getA(),    finalOrbit.getA(),    1.0e-10);
829         Assertions.assertEquals(initialOrbit.getEquinoctialEx(),    finalOrbit.getEquinoctialEx(),    1.0e-10);
830         Assertions.assertEquals(initialOrbit.getEquinoctialEy(),    finalOrbit.getEquinoctialEy(),    1.0e-10);
831         Assertions.assertEquals(initialOrbit.getHx(),    finalOrbit.getHx(),    1.0e-10);
832         Assertions.assertEquals(initialOrbit.getHy(),    finalOrbit.getHy(),    1.0e-10);
833         Assertions.assertEquals(initialOrbit.getLM() + n * dt, finalOrbit.getLM(), 6.0e-10);
834     }
835 
836     @Test
837     void testAdditionalStateEvent() {
838         propagator.addAdditionalDerivativesProvider(new AdditionalDerivativesProvider() {
839 
840             public String getName() {
841                 return "linear";
842             }
843 
844             public int getDimension() {
845                 return 1;
846             }
847 
848             public CombinedDerivatives combinedDerivatives(SpacecraftState s) {
849                 return new CombinedDerivatives(new double[] { 1.0 }, null);
850             }
851 
852         });
853         try {
854             propagator.addAdditionalDerivativesProvider(new AdditionalDerivativesProvider() {
855 
856                 public String getName() {
857                     return "linear";
858                 }
859 
860                 public int getDimension() {
861                     return 1;
862                 }
863 
864                 public CombinedDerivatives combinedDerivatives(SpacecraftState s) {
865                     return new CombinedDerivatives(new double[] { 1.0 }, null);
866                 }
867 
868             });
869             Assertions.fail("an exception should have been thrown");
870         } catch (OrekitException oe) {
871             Assertions.assertEquals(oe.getSpecifier(), OrekitMessages.ADDITIONAL_STATE_NAME_ALREADY_IN_USE);
872         }
873         try {
874             propagator.addAdditionalDerivativesProvider(new AdditionalDerivativesProvider() {
875                public String getName() {
876                     return "linear";
877                 }
878 
879                public int getDimension() {
880                    return 1;
881                }
882 
883                public CombinedDerivatives combinedDerivatives(SpacecraftState s) {
884                    return new CombinedDerivatives(new double[1], null);
885                }
886             });
887             Assertions.fail("an exception should have been thrown");
888         } catch (OrekitException oe) {
889             Assertions.assertEquals(oe.getSpecifier(), OrekitMessages.ADDITIONAL_STATE_NAME_ALREADY_IN_USE);
890         }
891         propagator.addAdditionalDataProvider(new AdditionalDataProvider<Double>() {
892             public String getName() {
893                 return "constant";
894             }
895 
896             public Double getAdditionalData(SpacecraftState state) {
897                 return 1.0;
898             }
899         });
900         Assertions.assertTrue(propagator.isAdditionalDataManaged("linear"));
901         Assertions.assertTrue(propagator.isAdditionalDataManaged("constant"));
902         Assertions.assertFalse(propagator.isAdditionalDataManaged("non-managed"));
903         Assertions.assertEquals(2, propagator.getManagedAdditionalData().length);
904         propagator.setInitialState(propagator.getInitialState().addAdditionalData("linear", 1.5));
905 
906         CheckingHandler checking = new CheckingHandler(Action.STOP);
907         propagator.addEventDetector(new AdditionalStateLinearDetector(new EventDetectionSettings(10.0, 1.0e-8, EventDetectionSettings.DEFAULT_MAX_ITER),
908                 checking));
909 
910         final double dt = 3200;
911         checking.assertEvent(false);
912         final SpacecraftState finalState =
913             propagator.propagate(initDate.shiftedBy(dt));
914         checking.assertEvent(true);
915         Assertions.assertEquals(3.0, finalState.getAdditionalState("linear")[0], 1.0e-8);
916         Assertions.assertEquals(1.5, finalState.getDate().durationFrom(initDate), 1.0e-8);
917 
918     }
919 
920     private static class AdditionalStateLinearDetector implements EventDetector {
921 
922         private final EventHandler eventHandler;
923         private final EventDetectionSettings detectionSettings;
924 
925         AdditionalStateLinearDetector(final EventDetectionSettings detectionSettings,
926                                       final EventHandler eventHandler) {
927             this.detectionSettings = detectionSettings;
928             this.eventHandler = eventHandler;
929         }
930 
931         public double g(SpacecraftState s) {
932             return s.getAdditionalState("linear")[0] - 3.0;
933         }
934 
935         @Override
936         public EventDetectionSettings getDetectionSettings() {
937             return detectionSettings;
938         }
939 
940         @Override
941         public EventHandler getHandler() {
942             return eventHandler;
943         }
944 
945     }
946 
947     @Test
948     void testResetAdditionalStateEvent() {
949         propagator.addAdditionalDerivativesProvider(new AdditionalDerivativesProvider() {
950 
951             public String getName() {
952                 return "linear";
953             }
954 
955             public int getDimension() {
956                 return 1;
957             }
958 
959             public CombinedDerivatives combinedDerivatives(SpacecraftState s) {
960                 return new CombinedDerivatives(new double[] { 1.0 }, null);
961             }
962         });
963         propagator.setInitialState(propagator.getInitialState().addAdditionalData("linear", 1.5));
964 
965         CheckingHandler checking = new CheckingHandler(Action.RESET_STATE) {
966             public SpacecraftState resetState(EventDetector detector, SpacecraftState oldState)
967                 {
968                 return oldState.addAdditionalData("linear", oldState.getAdditionalState("linear")[0] * 2);
969             }
970         };
971 
972         propagator.addEventDetector(new AdditionalStateLinearDetector(new EventDetectionSettings(10.0, 1.0e-8, EventDetectionSettings.DEFAULT_MAX_ITER),
973                 checking));
974 
975         final double dt = 3200;
976         checking.assertEvent(false);
977         final SpacecraftState finalState = propagator.propagate(initDate.shiftedBy(dt));
978         checking.assertEvent(true);
979         Assertions.assertEquals(dt + 4.5, finalState.getAdditionalState("linear")[0], 1.0e-8);
980         Assertions.assertEquals(dt, finalState.getDate().durationFrom(initDate), 1.0e-8);
981 
982     }
983 
984     @Test
985     void testEventDetectionBug() throws IOException, ParseException {
986 
987         TimeScale utc = TimeScalesFactory.getUTC();
988         AbsoluteDate initialDate = new AbsoluteDate(2005, 1, 1, 0, 0, 0.0, utc);
989         double duration = 100000.;
990         AbsoluteDate endDate = new AbsoluteDate(initialDate, duration);
991 
992         // Initialization of the frame EME2000
993         Frame EME2000 = FramesFactory.getEME2000();
994 
995 
996         // Initial orbit
997         double a = 35786000. + 6378137.0;
998         double e = 0.70;
999         double rApogee = a*(1+e);
1000         double vApogee = FastMath.sqrt(mu*(1-e)/(a*(1+e)));
1001         Orbit geo = new CartesianOrbit(new PVCoordinates(new Vector3D(rApogee, 0., 0.),
1002                                                          new Vector3D(0., vApogee, 0.)), EME2000,
1003                                                          initialDate, mu);
1004 
1005 
1006         duration = geo.getKeplerianPeriod();
1007         endDate = new AbsoluteDate(initialDate, duration);
1008 
1009         // Numerical Integration
1010         final double minStep  = 0.001;
1011         final double maxStep  = 1000;
1012         final double initStep = 60;
1013         final double[] absTolerance = {
1014             0.001, 1.0e-9, 1.0e-9, 1.0e-6, 1.0e-6, 1.0e-6, 0.001};
1015         final double[] relTolerance = {
1016             1.0e-7, 1.0e-4, 1.0e-4, 1.0e-7, 1.0e-7, 1.0e-7, 1.0e-7};
1017 
1018         AdaptiveStepsizeIntegrator integrator =
1019             new DormandPrince853Integrator(minStep, maxStep, absTolerance, relTolerance);
1020         integrator.setInitialStepSize(initStep);
1021 
1022         // Numerical propagator based on the integrator
1023         propagator = new NumericalPropagator(integrator);
1024         double mass = 1000.;
1025         SpacecraftState initialState = new SpacecraftState(geo).withMass(mass);
1026         propagator.setInitialState(initialState);
1027         propagator.setOrbitType(OrbitType.CARTESIAN);
1028 
1029 
1030         // Set the events Detectors
1031         ApsideDetector event1 = new ApsideDetector(geo);
1032         propagator.addEventDetector(event1);
1033 
1034         // Set the propagation mode
1035         propagator.clearStepHandlers();
1036 
1037         // Propagate
1038         SpacecraftState finalState = propagator.propagate(endDate);
1039 
1040         // we should stop long before endDate
1041         Assertions.assertTrue(endDate.durationFrom(finalState.getDate()) > 40000.0);
1042     }
1043 
1044     @Test
1045     void testEphemerisGenerationIssue14() throws IOException {
1046 
1047         // Propagation of the initial at t + dt
1048         final double dt = 3200;
1049         propagator.getInitialState();
1050 
1051         propagator.setOrbitType(OrbitType.CARTESIAN);
1052         final EphemerisGenerator generator = propagator.getEphemerisGenerator();
1053         propagator.propagate(initDate.shiftedBy(dt));
1054         final BoundedPropagator ephemeris1 = generator.getGeneratedEphemeris();
1055         Assertions.assertEquals(initDate, ephemeris1.getMinDate());
1056         Assertions.assertEquals(initDate.shiftedBy(dt), ephemeris1.getMaxDate());
1057 
1058         propagator.getPVCoordinates(initDate.shiftedBy( 2 * dt), FramesFactory.getEME2000());
1059         propagator.getPVCoordinates(initDate.shiftedBy(-2 * dt), FramesFactory.getEME2000());
1060 
1061         // the new propagations should not have changed ephemeris1
1062         Assertions.assertEquals(initDate, ephemeris1.getMinDate());
1063         Assertions.assertEquals(initDate.shiftedBy(dt), ephemeris1.getMaxDate());
1064 
1065         final BoundedPropagator ephemeris2 = generator.getGeneratedEphemeris();
1066         Assertions.assertEquals(initDate.shiftedBy(-2 * dt), ephemeris2.getMinDate());
1067         Assertions.assertEquals(initDate.shiftedBy( 2 * dt), ephemeris2.getMaxDate());
1068 
1069         // generating ephemeris2 should not have changed ephemeris1
1070         Assertions.assertEquals(initDate, ephemeris1.getMinDate());
1071         Assertions.assertEquals(initDate.shiftedBy(dt), ephemeris1.getMaxDate());
1072 
1073     }
1074 
1075     @Test
1076     void testEphemerisAdditionalState() throws IOException {
1077 
1078         // Propagation of the initial at t + dt
1079         final double dt = -3200;
1080         final double rate = 2.0;
1081 
1082         propagator.addAdditionalDataProvider(new AdditionalDataProvider<Double>() {
1083             public String getName() {
1084                 return "squaredA";
1085             }
1086             public Double getAdditionalData(SpacecraftState state) {
1087                 return state.getOrbit().getA() * state.getOrbit().getA();
1088             }
1089         });
1090         propagator.addAdditionalDerivativesProvider(new AdditionalDerivativesProvider() {
1091             public String getName() {
1092                 return "extra";
1093             }
1094             public int getDimension() {
1095                 return 1;
1096             }
1097             public CombinedDerivatives combinedDerivatives(SpacecraftState s) {
1098                 return new CombinedDerivatives(new double[] { rate }, null);
1099             }
1100         });
1101         propagator.setInitialState(propagator.getInitialState().addAdditionalData("extra", 1.5));
1102 
1103         propagator.setOrbitType(OrbitType.CARTESIAN);
1104         final EphemerisGenerator generator = propagator.getEphemerisGenerator();
1105         propagator.propagate(initDate.shiftedBy(dt));
1106         final BoundedPropagator ephemeris1 = generator.getGeneratedEphemeris();
1107         Assertions.assertEquals(initDate.shiftedBy(dt), ephemeris1.getMinDate());
1108         Assertions.assertEquals(initDate, ephemeris1.getMaxDate());
1109         try {
1110             ephemeris1.propagate(ephemeris1.getMinDate().shiftedBy(-10.0));
1111             Assertions.fail("an exception should have been thrown");
1112         } catch (OrekitException pe) {
1113             Assertions.assertEquals(OrekitMessages.OUT_OF_RANGE_EPHEMERIDES_DATE_BEFORE, pe.getSpecifier());
1114         }
1115         try {
1116             ephemeris1.propagate(ephemeris1.getMaxDate().shiftedBy(+10.0));
1117             Assertions.fail("an exception should have been thrown");
1118         } catch (OrekitException pe) {
1119             Assertions.assertEquals(OrekitMessages.OUT_OF_RANGE_EPHEMERIDES_DATE_AFTER, pe.getSpecifier());
1120         }
1121 
1122         double shift = -60;
1123         SpacecraftState s = ephemeris1.propagate(initDate.shiftedBy(shift));
1124         Assertions.assertEquals(2, s.getAdditionalDataValues().size());
1125         Assertions.assertTrue(s.hasAdditionalData("squaredA"));
1126         Assertions.assertTrue(s.hasAdditionalData("extra"));
1127         Assertions.assertEquals(s.getOrbit().getA() * s.getOrbit().getA(), s.getAdditionalState("squaredA")[0], 1.0e-10);
1128         Assertions.assertEquals(1.5 + shift * rate, s.getAdditionalState("extra")[0], 1.0e-10);
1129 
1130         try {
1131             ephemeris1.resetInitialState(s);
1132             Assertions.fail("an exception should have been thrown");
1133         } catch (OrekitException oe) {
1134             Assertions.assertEquals(OrekitMessages.NON_RESETABLE_STATE, oe.getSpecifier());
1135         }
1136 
1137     }
1138 
1139     @Test
1140     void testIssue157() {
1141         try {
1142             Orbit orbit = new KeplerianOrbit(13378000, 0.05, 0, 0, FastMath.PI, 0, PositionAngleType.MEAN,
1143                                              FramesFactory.getTOD(false),
1144                                              new AbsoluteDate(2003, 5, 6, TimeScalesFactory.getUTC()),
1145                                              Constants.EIGEN5C_EARTH_MU);
1146             ToleranceProvider.of(CartesianToleranceProvider.of(1)).getTolerances(orbit, OrbitType.KEPLERIAN);
1147             Assertions.fail("an exception should have been thrown");
1148         } catch (OrekitException pe) {
1149             Assertions.assertEquals(OrekitMessages.SINGULAR_JACOBIAN_FOR_ORBIT_TYPE, pe.getSpecifier());
1150         }
1151     }
1152 
1153     @Test
1154     void testInternalEventDetectorsFromAttitudeProvider() {
1155         // GIVEN
1156         propagator.resetInitialState(initialState);
1157         final AbsoluteDate epoch = propagator.getInitialState().getDate();
1158         final AbsoluteDate interruptingDate = epoch.shiftedBy(1);
1159         propagator.setAttitudeProvider(new InterruptingAttitudeProvider(interruptingDate));
1160         // WHEN
1161         final SpacecraftState state = propagator.propagate(interruptingDate.shiftedBy(10));
1162         // THEN
1163         Assertions.assertEquals(state.getDate(), interruptingDate);
1164     }
1165 
1166     private static class InterruptingAttitudeProvider extends FrameAlignedProvider {
1167 
1168         private final AbsoluteDate interruptingDate;
1169 
1170         public InterruptingAttitudeProvider(final AbsoluteDate interruptingDate) {
1171             super(Rotation.IDENTITY);
1172             this.interruptingDate = interruptingDate;
1173         }
1174 
1175         @Override
1176         public Stream<EventDetector> getEventDetectors() {
1177             final DateDetector detector = new DateDetector(interruptingDate).withHandler(new StopOnEvent());
1178             return Stream.of(detector);
1179         }
1180     }
1181 
1182     private static class CheckingHandler implements EventHandler {
1183 
1184         private final Action actionOnEvent;
1185         private boolean gotHere;
1186         private boolean isFinished = false;
1187 
1188         public CheckingHandler(final Action actionOnEvent) {
1189             this.actionOnEvent = actionOnEvent;
1190             this.gotHere       = false;
1191         }
1192 
1193         public void assertEvent(boolean expected) {
1194             Assertions.assertEquals(expected, gotHere);
1195         }
1196 
1197         public Action eventOccurred(SpacecraftState s, EventDetector detector, boolean increasing) {
1198             gotHere = true;
1199             return actionOnEvent;
1200         }
1201 
1202         @Override
1203         public void finish(SpacecraftState finalState, EventDetector detector) {
1204             isFinished = true;
1205         }
1206     }
1207 
1208     @Test
1209     void testParallelismIssue258() throws InterruptedException, ExecutionException, FileNotFoundException {
1210 
1211         Utils.setDataRoot("regular-data:atmosphere:potential/grgs-format");
1212         GravityFieldFactory.addPotentialCoefficientsReader(new GRGSFormatReader("grim4s4_gr", true));
1213         final double mu = GravityFieldFactory.getNormalizedProvider(2, 2).getMu();
1214 
1215         // Geostationary transfer orbit
1216         final double a = 24396159; // semi major axis in meters
1217         final double e = 0.72831215; // eccentricity
1218         final double i = FastMath.toRadians(7); // inclination
1219         final double omega = FastMath.toRadians(180); // perigee argument
1220         final double raan = FastMath.toRadians(261); // right ascension of ascending node
1221         final double lM = 0; // mean anomaly
1222         final Frame inertialFrame = FramesFactory.getEME2000();
1223         final TimeScale utc = TimeScalesFactory.getUTC();
1224         final AbsoluteDate initialDate = new AbsoluteDate(2003, 1, 1, 00, 00, 00.000, utc);
1225         final Orbit initialOrbit = new CartesianOrbit( new KeplerianOrbit(a, e, i, omega, raan, lM, PositionAngleType.MEAN,
1226                                                                           inertialFrame, initialDate, mu));
1227         final SpacecraftState initialState = new SpacecraftState(initialOrbit).withMass( 1000);
1228 
1229         // initialize the testing points
1230         final List<SpacecraftState> states = new ArrayList<SpacecraftState>();
1231         final NumericalPropagator propagator = createPropagator(initialState, OrbitType.CARTESIAN, PositionAngleType.TRUE);
1232         final double samplingStep = 10000.0;
1233         propagator.setStepHandler(samplingStep, state -> states.add(state));
1234         propagator.propagate(initialDate.shiftedBy(5 * samplingStep));
1235 
1236         // compute reference errors, using serial computation in a for loop
1237         final double[][] referenceErrors = new double[states.size() - 1][];
1238         for (int startIndex = 0; startIndex < states.size() - 1; ++startIndex) {
1239             referenceErrors[startIndex] = recomputeFollowing(startIndex, states);
1240         }
1241 
1242         final Consumer<SpacecraftState> checker = point -> {
1243             try {
1244                 final int startIndex = states.indexOf(point);
1245                 double[] errors = recomputeFollowing(startIndex, states);
1246                 for (int k = 0; k < errors.length; ++k) {
1247                     Assertions.assertEquals(
1248                                         referenceErrors[startIndex][k], errors[k],
1249                                         1.0e-9,startIndex + " to " + (startIndex + k + 1));
1250                 }
1251             } catch (OrekitException oe) {
1252                 Assertions.fail(oe.getLocalizedMessage());
1253             }
1254         };
1255 
1256         // serial propagation using Stream
1257         states.stream().forEach(checker);
1258 
1259         // parallel propagation using parallelStream
1260         states.parallelStream().forEach(checker);
1261 
1262     }
1263 
1264     @Test
1265     void testShiftKeplerianEllipticTrueWithoutDerivatives() {
1266         doTestShift(createEllipticOrbit(), OrbitType.KEPLERIAN, PositionAngleType.TRUE, false,
1267                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1268     }
1269 
1270     @Test
1271     void testShiftKeplerianEllipticTrueWithDerivatives() {
1272         doTestShift(createEllipticOrbit(), OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
1273                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1274     }
1275 
1276     @Test
1277     void testShiftKeplerianEllipticEccentricWithoutDerivatives() {
1278         doTestShift(createEllipticOrbit(), OrbitType.KEPLERIAN, PositionAngleType.ECCENTRIC, false,
1279                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1280     }
1281 
1282     @Test
1283     void testShiftKeplerianEllipticEcentricWithDerivatives() {
1284         doTestShift(createEllipticOrbit(), OrbitType.KEPLERIAN, PositionAngleType.ECCENTRIC, true,
1285                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1286     }
1287 
1288     @Test
1289     void testShiftKeplerianEllipticMeanWithoutDerivatives() {
1290         doTestShift(createEllipticOrbit(), OrbitType.KEPLERIAN, PositionAngleType.MEAN, false,
1291                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1292     }
1293 
1294     @Test
1295     void testShiftKeplerianEllipticMeanWithDerivatives() {
1296         doTestShift(createEllipticOrbit(), OrbitType.KEPLERIAN, PositionAngleType.MEAN, true,
1297                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1298     }
1299 
1300     @Test
1301     void testShiftKeplerianHyperbolicTrueWithoutDerivatives() {
1302         doTestShift(createHyperbolicOrbit(), OrbitType.KEPLERIAN, PositionAngleType.TRUE, false,
1303                     0.484, 1.94, 12.1, 48.3, 108.5);
1304     }
1305 
1306     @Test
1307     void testShiftKeplerianHyperbolicTrueWithDerivatives() {
1308         doTestShift(createHyperbolicOrbit(), OrbitType.KEPLERIAN, PositionAngleType.TRUE, true,
1309                     1.38e-4, 1.10e-3, 1.72e-2, 1.37e-1, 4.62e-1);
1310     }
1311 
1312     @Test
1313     void testShiftKeplerianHyperbolicEccentricWithoutDerivatives() {
1314         doTestShift(createHyperbolicOrbit(), OrbitType.KEPLERIAN, PositionAngleType.ECCENTRIC, false,
1315                     0.484, 1.94, 12.1, 48.3, 108.5);
1316     }
1317 
1318     @Test
1319     void testShiftKeplerianHyperbolicEcentricWithDerivatives() {
1320         doTestShift(createHyperbolicOrbit(), OrbitType.KEPLERIAN, PositionAngleType.ECCENTRIC, true,
1321                     1.38e-4, 1.10e-3, 1.72e-2, 1.37e-1, 4.62e-1);
1322     }
1323 
1324     @Test
1325     void testShiftKeplerianHyperbolicMeanWithoutDerivatives() {
1326         doTestShift(createHyperbolicOrbit(), OrbitType.KEPLERIAN, PositionAngleType.MEAN, false,
1327                     0.484, 1.94, 12.1, 48.3, 108.5);
1328     }
1329 
1330     @Test
1331     void testShiftKeplerianHyperbolicMeanWithDerivatives() {
1332         doTestShift(createHyperbolicOrbit(), OrbitType.KEPLERIAN, PositionAngleType.MEAN, true,
1333                     1.38e-4, 1.10e-3, 1.72e-2, 1.37e-1, 4.62e-1);
1334     }
1335 
1336     @Test
1337     void testShiftCartesianEllipticTrueWithoutDerivatives() {
1338         doTestShift(createEllipticOrbit(), OrbitType.CARTESIAN, PositionAngleType.TRUE, false,
1339                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1340     }
1341 
1342     @Test
1343     void testShiftCartesianEllipticTrueWithDerivatives() {
1344         doTestShift(createEllipticOrbit(), OrbitType.CARTESIAN, PositionAngleType.TRUE, true,
1345                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1346     }
1347 
1348     @Test
1349     void testShiftCartesianEllipticEccentricWithoutDerivatives() {
1350         doTestShift(createEllipticOrbit(), OrbitType.CARTESIAN, PositionAngleType.ECCENTRIC, false,
1351                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1352     }
1353 
1354     @Test
1355     void testShiftCartesianEllipticEcentricWithDerivatives() {
1356         doTestShift(createEllipticOrbit(), OrbitType.CARTESIAN, PositionAngleType.ECCENTRIC, true,
1357                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1358     }
1359 
1360     @Test
1361     void testShiftCartesianEllipticMeanWithoutDerivatives() {
1362         doTestShift(createEllipticOrbit(), OrbitType.CARTESIAN, PositionAngleType.MEAN, false,
1363                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1364     }
1365 
1366     @Test
1367     void testShiftCartesianEllipticMeanWithDerivatives() {
1368         doTestShift(createEllipticOrbit(), OrbitType.CARTESIAN, PositionAngleType.MEAN, true,
1369                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1370     }
1371 
1372     @Test
1373     void testShiftCartesianHyperbolicTrueWithoutDerivatives() {
1374         doTestShift(createHyperbolicOrbit(), OrbitType.CARTESIAN, PositionAngleType.TRUE, false,
1375                     0.48, 1.93, 12.1, 48.3, 108.5);
1376     }
1377 
1378     @Test
1379     void testShiftCartesianHyperbolicTrueWithDerivatives() {
1380         doTestShift(createHyperbolicOrbit(), OrbitType.CARTESIAN, PositionAngleType.TRUE, true,
1381                     1.38e-4, 1.10e-3, 1.72e-2, 1.37e-1, 4.62e-1);
1382     }
1383 
1384     @Test
1385     void testShiftCartesianHyperbolicEccentricWithoutDerivatives() {
1386         doTestShift(createHyperbolicOrbit(), OrbitType.CARTESIAN, PositionAngleType.ECCENTRIC, false,
1387                     0.48, 1.93, 12.1, 48.3, 108.5);
1388     }
1389 
1390     @Test
1391     void testShiftCartesianHyperbolicEcentricWithDerivatives() {
1392         doTestShift(createHyperbolicOrbit(), OrbitType.CARTESIAN, PositionAngleType.ECCENTRIC, true,
1393                     1.38e-4, 1.10e-3, 1.72e-2, 1.37e-1, 4.62e-1);
1394     }
1395 
1396     @Test
1397     void testShiftCartesianHyperbolicMeanWithoutDerivatives() {
1398         doTestShift(createHyperbolicOrbit(), OrbitType.CARTESIAN, PositionAngleType.MEAN, false,
1399                     0.48, 1.93, 12.1, 48.3, 108.5);
1400     }
1401 
1402     @Test
1403     void testShiftCartesianHyperbolicMeanWithDerivatives() {
1404         doTestShift(createHyperbolicOrbit(), OrbitType.CARTESIAN, PositionAngleType.MEAN, true,
1405                     1.38e-4, 1.10e-3, 1.72e-2, 1.37e-1, 4.62e-1);
1406     }
1407 
1408     @Test
1409     void testShiftCircularTrueWithoutDerivatives() {
1410         doTestShift(createEllipticOrbit(), OrbitType.CIRCULAR, PositionAngleType.TRUE, false,
1411                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1412     }
1413 
1414     @Test
1415     void testShiftCircularTrueWithDerivatives() {
1416         doTestShift(createEllipticOrbit(), OrbitType.CIRCULAR, PositionAngleType.TRUE, true,
1417                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1418     }
1419 
1420     @Test
1421     void testShiftCircularEccentricWithoutDerivatives() {
1422         doTestShift(createEllipticOrbit(), OrbitType.CIRCULAR, PositionAngleType.ECCENTRIC, false,
1423                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1424     }
1425 
1426     @Test
1427     void testShiftCircularEcentricWithDerivatives() {
1428         doTestShift(createEllipticOrbit(), OrbitType.CIRCULAR, PositionAngleType.ECCENTRIC, true,
1429                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1430     }
1431 
1432     @Test
1433     void testShiftCircularMeanWithoutDerivatives() {
1434         doTestShift(createEllipticOrbit(), OrbitType.CIRCULAR, PositionAngleType.MEAN, false,
1435                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1436     }
1437 
1438     @Test
1439     void testShiftCircularMeanWithDerivatives() {
1440         doTestShift(createEllipticOrbit(), OrbitType.CIRCULAR, PositionAngleType.MEAN, true,
1441                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1442     }
1443 
1444     @Test
1445     void testShiftEquinoctialTrueWithoutDerivatives() {
1446         doTestShift(createEllipticOrbit(), OrbitType.EQUINOCTIAL, PositionAngleType.TRUE, false,
1447                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1448     }
1449 
1450     @Test
1451     void testShiftEquinoctialTrueWithDerivatives() {
1452         doTestShift(createEllipticOrbit(), OrbitType.EQUINOCTIAL, PositionAngleType.TRUE, true,
1453                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1454     }
1455 
1456     @Test
1457     void testShiftEquinoctialEccentricWithoutDerivatives() {
1458         doTestShift(createEllipticOrbit(), OrbitType.EQUINOCTIAL, PositionAngleType.ECCENTRIC, false,
1459                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1460     }
1461 
1462     @Test
1463     void testShiftEquinoctialEcentricWithDerivatives() {
1464         doTestShift(createEllipticOrbit(), OrbitType.EQUINOCTIAL, PositionAngleType.ECCENTRIC, true,
1465                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1466     }
1467 
1468     @Test
1469     void testShiftEquinoctialMeanWithoutDerivatives() {
1470         doTestShift(createEllipticOrbit(), OrbitType.EQUINOCTIAL, PositionAngleType.MEAN, false,
1471                     18.1, 72.0, 437.3, 1601.1, 3141.8);
1472     }
1473 
1474     @Test
1475     void testShiftEquinoctialMeanWithDerivatives() {
1476         doTestShift(createEllipticOrbit(), OrbitType.EQUINOCTIAL, PositionAngleType.MEAN, true,
1477                     1.14, 9.1, 140.3, 1066.7, 3306.9);
1478     }
1479 
1480     private static void doTestShift(final CartesianOrbit orbit, final OrbitType orbitType,
1481                                     final PositionAngleType angleType, final boolean withDerivatives,
1482                                     final double error60s, final double error120s,
1483                                     final double error300s, final double error600s,
1484                                     final double error900s)
1485         {
1486 
1487         Utils.setDataRoot("regular-data:atmosphere:potential/grgs-format");
1488         GravityFieldFactory.addPotentialCoefficientsReader(new GRGSFormatReader("grim4s4_gr", true));
1489         final NumericalPropagator np = createPropagator(new SpacecraftState(orbit), orbitType, angleType);
1490 
1491         // the reference date for shifts is set at 60s, so the propagator can provide derivatives if needed
1492         // (derivatives are not available in the initial orbit)
1493         final AbsoluteDate reference = orbit.getDate().shiftedBy(60.0);
1494         final ShiftChecker checker   = new ShiftChecker(withDerivatives, orbitType, angleType,
1495                                                         error60s,
1496                                                         error120s, error300s,
1497                                                         error600s, error900s);
1498         np.addEventDetector(new DateDetector(reference,
1499                                              reference.shiftedBy( 60.0),
1500                                              reference.shiftedBy(120.0),
1501                                              reference.shiftedBy(300.0),
1502                                              reference.shiftedBy(600.0),
1503                                              reference.shiftedBy(900.0)).
1504                             withMaxCheck(30.0).
1505                             withThreshold(1.0e-9).
1506                             withHandler(checker));
1507         np.propagate(reference.shiftedBy(1000.0));
1508     }
1509 
1510     private static class ShiftChecker implements EventHandler {
1511 
1512         private final boolean       withDerivatives;
1513         private final OrbitType     orbitType;
1514         private final PositionAngleType angleType;
1515         private final double        error60s;
1516         private final double        error120s;
1517         private final double        error300s;
1518         private final double        error600s;
1519         private final double        error900s;
1520         private SpacecraftState     referenceState;
1521 
1522         ShiftChecker(final boolean withDerivatives, final OrbitType orbitType,
1523                      final PositionAngleType angleType, final double error60s,
1524                      final double error120s, final double error300s,
1525                      final double error600s, final double error900s) {
1526             this.withDerivatives = withDerivatives;
1527             this.orbitType       = orbitType;
1528             this.angleType       = angleType;
1529             this.error60s        = error60s;
1530             this.error120s       = error120s;
1531             this.error300s       = error300s;
1532             this.error600s       = error600s;
1533             this.error900s       = error900s;
1534             this.referenceState  = null;
1535         }
1536 
1537         @Override
1538         public Action eventOccurred(final SpacecraftState s, final EventDetector detector,
1539                                     final boolean increasing) {
1540             if (referenceState == null) {
1541                 // first event, we retrieve the reference state for later use
1542                 if (withDerivatives) {
1543                     referenceState = s;
1544                 } else {
1545                     // remove derivatives, to check accuracy of the shiftedBy method decreases without them
1546                     final double[] stateVector = new double[6];
1547                     final Orbit o = s.getOrbit();
1548                     orbitType.mapOrbitToArray(o, angleType, stateVector, null);
1549                     final Orbit fixedOrbit = orbitType.mapArrayToOrbit(stateVector, null, angleType,
1550                                                                        o.getDate(), o.getMu(), o.getFrame());
1551                     referenceState = new SpacecraftState(fixedOrbit, s.getAttitude()).withMass(s.getMass());
1552                 }
1553             } else {
1554                 // recurring event, we compare with the shifted reference state
1555                 final double dt = s.getDate().durationFrom(referenceState.getDate());
1556                 final SpacecraftState shifted = referenceState.shiftedBy(dt);
1557                 final double error = Vector3D.distance(shifted.getPosition(),
1558                                                        s.getPosition());
1559                 switch ((int) FastMath.rint(dt)) {
1560                     case 60 :
1561                         Assertions.assertEquals(error60s,  error, 0.01 * error60s);
1562                         break;
1563                     case 120 :
1564                         Assertions.assertEquals(error120s, error, 0.01 * error120s);
1565                         break;
1566                     case 300 :
1567                         Assertions.assertEquals(error300s, error, 0.01 * error300s);
1568                         break;
1569                     case 600 :
1570                         Assertions.assertEquals(error600s, error, 0.01 * error600s);
1571                         break;
1572                     case 900 :
1573                         Assertions.assertEquals(error900s, error, 0.01 * error900s);
1574                         break;
1575                     default :
1576                         // this should never happen
1577                         Assertions.fail("no error set for dt = " + dt);
1578                         break;
1579                 }
1580             }
1581             return Action.CONTINUE;
1582         }
1583 
1584     }
1585 
1586     /** Test de-activation of event detection and step handling.
1587      *  When propagating out of start and target date in propagate(startDate, targetDate)
1588      *  <p>See issue 449 in Orekit forge and
1589      *  {@link org.orekit.propagation.Propagator#propagate(AbsoluteDate, AbsoluteDate)}.
1590      *  </p>
1591      */
1592     @Test
1593     void testEventAndStepHandlerDeactivationIssue449() {
1594 
1595         // Setup
1596         RecordAndContinue recordAndContinue = new RecordAndContinue();
1597         DateDetector dateDetector = new DateDetector(initDate.shiftedBy(10.),
1598                                                      initDate.shiftedBy(15.),
1599                                                      initDate.shiftedBy(20.)).
1600                                     withMaxCheck(1).
1601                                     withThreshold(1e-1).
1602                                     withHandler(recordAndContinue);
1603 
1604         propagator.addEventDetector(dateDetector);
1605 
1606         final AbsoluteDate startDate = initDate.shiftedBy(30.);
1607         final AbsoluteDate finalDate = initDate.shiftedBy(40.);
1608 
1609         final DateRecorderHandler dateRecorderHandler = new DateRecorderHandler(startDate, finalDate);
1610         propagator.setStepHandler(1.0, dateRecorderHandler);
1611 
1612         // Action
1613         propagator.propagate(startDate, finalDate);
1614 
1615         // Verify
1616         // No event is detected
1617         Assertions.assertEquals(0, recordAndContinue.getEvents().size());
1618 
1619         // Handler is deactivated (no dates recorded between start and stop date)
1620         Assertions.assertEquals(0, dateRecorderHandler.handledDatesOutOfInterval.size());
1621     }
1622 
1623     @Test
1624     void testResetStateForward() {
1625         final Frame eme2000 = FramesFactory.getEME2000();
1626         final AbsoluteDate date = new AbsoluteDate(new DateComponents(2008, 6, 23),
1627                                                    new TimeComponents(14, 0, 0),
1628                                                    TimeScalesFactory.getUTC());
1629         final Orbit orbit = new KeplerianOrbit(8000000.0, 0.01, 0.87, 2.44, 0.21, -1.05, PositionAngleType.MEAN,
1630                                            eme2000,
1631                                            date, Constants.EIGEN5C_EARTH_MU);
1632         final NumericalPropagator propagator =
1633                         new NumericalPropagator(new LutherIntegrator(300.0),
1634                                                 new LofOffset(eme2000, LOFType.LVLH));
1635         propagator.resetInitialState(new SpacecraftState(orbit));
1636 
1637         // maneuver along Z in attitude aligned with LVLH will change orbital plane
1638         final AbsoluteDate maneuverDate = date.shiftedBy(1000.0);
1639         propagator.addEventDetector(new ImpulseManeuver(new DateDetector(maneuverDate),
1640                                                         new Vector3D(0.0, 0.0, -100.0),
1641                                                         350.0));
1642 
1643         final Vector3D initialNormal = orbit.getPVCoordinates().getMomentum();
1644         propagator.setStepHandler(60.0, state -> {
1645             final Vector3D currentNormal = state.getPVCoordinates().getMomentum();
1646             if (state.getDate().isBefore(maneuverDate)) {
1647                 Assertions.assertEquals(0.000, Vector3D.angle(initialNormal, currentNormal), 1.0e-3);
1648             } else {
1649                 Assertions.assertEquals(0.014, Vector3D.angle(initialNormal, currentNormal), 1.0e-3);
1650             }
1651         });
1652 
1653         propagator.propagate(orbit.getDate().shiftedBy(1500.0));
1654 
1655     }
1656 
1657     @Test
1658     void testResetStateBackward() {
1659         final Frame eme2000 = FramesFactory.getEME2000();
1660         final AbsoluteDate date = new AbsoluteDate(new DateComponents(2008, 6, 23),
1661                                                    new TimeComponents(14, 0, 0),
1662                                                    TimeScalesFactory.getUTC());
1663         final Orbit orbit = new KeplerianOrbit(8000000.0, 0.01, 0.87, 2.44, 0.21, -1.05, PositionAngleType.MEAN,
1664                                            eme2000,
1665                                            date, Constants.EIGEN5C_EARTH_MU);
1666         final NumericalPropagator propagator =
1667                         new NumericalPropagator(new LutherIntegrator(300.0),
1668                                                 new LofOffset(eme2000, LOFType.LVLH));
1669         propagator.resetInitialState(new SpacecraftState(orbit));
1670 
1671         // maneuver along Z in attitude aligned with LVLH will change orbital plane
1672         final AbsoluteDate maneuverDate = date.shiftedBy(-1000.0);
1673         propagator.addEventDetector(new ImpulseManeuver(new DateDetector(maneuverDate),
1674                                                         new Vector3D(0.0, 0.0, -100.0),
1675                                                         350.0));
1676 
1677         final Vector3D initialNormal = orbit.getPVCoordinates().getMomentum();
1678         propagator.setStepHandler(60.0, state -> {
1679             final Vector3D currentNormal = state.getPVCoordinates().getMomentum();
1680             if (state.getDate().isAfter(maneuverDate)) {
1681                 Assertions.assertEquals(0.000, Vector3D.angle(initialNormal, currentNormal), 1.0e-3);
1682             } else {
1683                 Assertions.assertEquals(0.014, Vector3D.angle(initialNormal, currentNormal), 1.0e-3);
1684             }
1685         });
1686 
1687         propagator.propagate(orbit.getDate().shiftedBy(-1500.0));
1688 
1689     }
1690 
1691     @Test
1692     void testAdditionalDerivatives() {
1693 
1694         final Frame eme2000 = FramesFactory.getEME2000();
1695         final AbsoluteDate date = new AbsoluteDate(new DateComponents(2008, 6, 23),
1696                                                    new TimeComponents(14, 0, 0),
1697                                                    TimeScalesFactory.getUTC());
1698         final Orbit initialOrbit = new KeplerianOrbit(8000000.0, 0.01, 0.87, 2.44, 0.21, -1.05, PositionAngleType.MEAN,
1699                                            eme2000,
1700                                            date, Constants.EIGEN5C_EARTH_MU);
1701         NumericalPropagatorBuilder builder = new NumericalPropagatorBuilder(initialOrbit,
1702                 new DormandPrince853IntegratorBuilder(0.02, 0.2, 1.), PositionAngleType.TRUE, 10);
1703         NumericalPropagator propagator = (NumericalPropagator) builder.buildPropagator();
1704 
1705         IntStream.
1706         range(0, 2).
1707         mapToObj(i -> new EmptyDerivativeProvider("test_provider_" + i, new double[] { 10 * i, 20 * i })).
1708         forEach(provider -> addDerivativeProvider(propagator, provider));
1709 
1710         EphemerisGenerator generator = propagator.getEphemerisGenerator();
1711         propagator.propagate(initialOrbit.getDate().shiftedBy(600));
1712         BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
1713         final SpacecraftState finalState = ephemeris.propagate(initialOrbit.getDate().shiftedBy(300));
1714         Assertions.assertEquals(2,    finalState.getAdditionalDataValues().size());
1715         Assertions.assertEquals(2,    finalState.getAdditionalState("test_provider_0").length);
1716         Assertions.assertEquals(0.0,  finalState.getAdditionalState("test_provider_0")[0], 1.0e-15);
1717         Assertions.assertEquals(0.0,  finalState.getAdditionalState("test_provider_0")[1], 1.0e-15);
1718         Assertions.assertEquals(2,    finalState.getAdditionalState("test_provider_1").length);
1719         Assertions.assertEquals(10.0, finalState.getAdditionalState("test_provider_1")[0], 1.0e-15);
1720         Assertions.assertEquals(20.0, finalState.getAdditionalState("test_provider_1")[1], 1.0e-15);
1721     }
1722 
1723     private void addDerivativeProvider(NumericalPropagator propagator, EmptyDerivativeProvider provider) {
1724         SpacecraftState initialState = propagator.getInitialState();
1725         propagator.setInitialState(initialState.addAdditionalData(provider.getName(), provider.getInitialState()));
1726         propagator.addAdditionalDerivativesProvider(provider);
1727     }
1728 
1729     private static class EmptyDerivativeProvider implements AdditionalDerivativesProvider {
1730 
1731         private final String name;
1732         private final double[] state;
1733 
1734         public EmptyDerivativeProvider(String name, double[] state) {
1735             this.name = name;
1736             this.state = state;
1737         }
1738 
1739         @Override
1740         public CombinedDerivatives combinedDerivatives(SpacecraftState s) {
1741             return new CombinedDerivatives(new double[getDimension()], null);
1742         }
1743 
1744         @Override
1745         public String getName() {
1746             return name;
1747         }
1748 
1749         @Override
1750         public int getDimension() {
1751             return state.length;
1752         }
1753 
1754         public double[] getInitialState() {
1755             return state;
1756         }
1757     }
1758 
1759     @Test
1760     void testInfinitePropagation() {
1761 
1762         Utils.setDataRoot("regular-data:atmosphere:potential/grgs-format");
1763         GravityFieldFactory.addPotentialCoefficientsReader(new GRGSFormatReader("grim4s4_gr", true));
1764 
1765         final NumericalPropagator propag = createPropagator(initialState, OrbitType.KEPLERIAN, PositionAngleType.TRUE);
1766 
1767         // Stop condition
1768         final double convergenceThreshold = 1e-9;
1769         propag.addEventDetector(new DateDetector(initialState.getDate().shiftedBy(60)).
1770                                 withMaxCheck(1e10).
1771                                 withThreshold(convergenceThreshold));
1772 
1773         // Propagate until the stop condition is reached
1774         final SpacecraftState finalState =  propag.propagate(AbsoluteDate.FUTURE_INFINITY);
1775 
1776         // Check that the expected final state was reached
1777         Assertions.assertEquals(60, finalState.getDate().durationFrom(initialState.getDate()), convergenceThreshold);
1778 
1779     }
1780 
1781     @Test
1782     void getIntegratorNameTest() {
1783         // GIVEN
1784         final String expectedName = "Name";
1785         final ODEIntegrator mockedIntegrator = Mockito.mock(ODEIntegrator.class);
1786         Mockito.when(mockedIntegrator.getName()).thenReturn(expectedName);
1787         // WHEN
1788         final NumericalPropagator numericalPropagator = new NumericalPropagator(mockedIntegrator);
1789         final String actualName = numericalPropagator.getIntegratorName();
1790         // THEN
1791         Assertions.assertEquals(expectedName, actualName);
1792     }
1793 
1794     @Test
1795     void testIssue1395() {
1796         // GIVEN
1797         final Orbit initialOrbit = createEllipticOrbit();
1798         final ClassicalRungeKuttaIntegrator rungeKuttaIntegrator = new ClassicalRungeKuttaIntegrator(10);
1799         final NumericalPropagator numericalPropagator = new NumericalPropagator(rungeKuttaIntegrator);
1800         final SpacecraftState state = new SpacecraftState(initialOrbit);
1801         final String name = "test";
1802         numericalPropagator.setInitialState(state.addAdditionalData(name, 0.));
1803         numericalPropagator.addAdditionalDerivativesProvider(mockDerivativeProvider(name));
1804         numericalPropagator.addForceModel(createForceModelBasedOnAdditionalState(name));
1805         // WHEN & THEN
1806         final AbsoluteDate epoch = initialOrbit.getDate();
1807         final SpacecraftState propagateState = Assertions.assertDoesNotThrow(() ->
1808                 numericalPropagator.propagate(epoch.shiftedBy(10.)));
1809         Assertions.assertNotEquals(epoch, propagateState.getDate());
1810     }
1811 
1812     private static AdditionalDerivativesProvider mockDerivativeProvider(final String name) {
1813         final AdditionalDerivativesProvider mockedProvider = Mockito.mock(AdditionalDerivativesProvider.class);
1814         final int dimension = 1;
1815         Mockito.when(mockedProvider.getDimension()).thenReturn(dimension);
1816         Mockito.when(mockedProvider.getName()).thenReturn(name);
1817         final double[] yDot = new double[dimension];
1818         final CombinedDerivatives combinedDerivatives = new CombinedDerivatives(yDot, null);
1819         Mockito.when(mockedProvider.combinedDerivatives(Mockito.any(SpacecraftState.class)))
1820                 .thenReturn(combinedDerivatives);
1821         return mockedProvider;
1822     }
1823 
1824     private static ForceModel createForceModelBasedOnAdditionalState(final String name) {
1825         return new ForceModel() {
1826 
1827             @Override
1828             public void init(SpacecraftState initialState, AbsoluteDate target) {
1829                 ForceModel.super.init(initialState, target);
1830                 initialState.getAdditionalState(name);
1831             }
1832 
1833             @Override
1834             public boolean dependsOnPositionOnly() {
1835                 return false;
1836             }
1837 
1838             @Override
1839             public Vector3D acceleration(SpacecraftState s, double[] parameters) {
1840                 return Vector3D.ZERO;
1841             }
1842 
1843             @Override
1844             public <T extends CalculusFieldElement<T>> FieldVector3D<T> acceleration(FieldSpacecraftState<T> s, T[] parameters) {
1845                 return null; // not used
1846             }
1847 
1848             @Override
1849             public List<ParameterDriver> getParametersDrivers() {
1850                 return new ArrayList<>();
1851             }
1852         };
1853     }
1854 
1855     /** Record the dates treated by the handler.
1856      *  If they are out of an interval defined by a start and final date.
1857      */
1858     private static class DateRecorderHandler implements OrekitFixedStepHandler {
1859 
1860         /** Start date of the propagation. */
1861         private final AbsoluteDate startDate;
1862 
1863         /** Final date of the propagation. */
1864         private final AbsoluteDate finalDate;
1865 
1866         /** List of handled date. Recorded only if they are out of the propagation interval. */
1867         public final List<AbsoluteDate> handledDatesOutOfInterval;
1868 
1869         DateRecorderHandler(final AbsoluteDate startDate, final AbsoluteDate finalDate) {
1870           this.startDate = startDate;
1871           this.finalDate = finalDate;
1872           this.handledDatesOutOfInterval = new ArrayList<>();
1873         }
1874 
1875         @Override
1876         public void handleStep(SpacecraftState currentState)
1877                 {
1878           final AbsoluteDate date = currentState.getDate();
1879           if (date.compareTo(startDate) < 0 || date.compareTo(finalDate) > 0) {
1880             handledDatesOutOfInterval.add(currentState.getDate());
1881           }
1882         }
1883       }
1884 
1885     /**
1886      * Assume we have 5 epochs, we will propagate from the input epoch to all the following epochs.
1887      *   If we have [0, 1, 2, 3, 4], and input is 2, then we will do 2->3, 2->4.
1888      * @param startIndex index of start state
1889      * @return position error for recomputed following points
1890      */
1891     private static double[] recomputeFollowing(final int startIndex, List<SpacecraftState> allPoints)
1892         {
1893         SpacecraftState startState = allPoints.get(startIndex);
1894         NumericalPropagator innerPropagator = createPropagator(startState, OrbitType.CARTESIAN, PositionAngleType.TRUE);
1895         double[] errors = new double[allPoints.size() - startIndex - 1];
1896         for (int endIndex = startIndex + 1; endIndex < allPoints.size(); ++endIndex) {
1897             final TimeStampedPVCoordinates reference  = allPoints.get(endIndex).getPVCoordinates();
1898             final TimeStampedPVCoordinates recomputed = innerPropagator.propagate(reference.getDate()).getPVCoordinates();
1899             errors[endIndex - startIndex - 1] = Vector3D.distance(recomputed.getPosition(), reference.getPosition());
1900         }
1901         return errors;
1902     }
1903 
1904     private synchronized static NumericalPropagator createPropagator(SpacecraftState spacecraftState,
1905                                                                      OrbitType orbitType, PositionAngleType angleType)
1906         {
1907 
1908         final double minStep                         = 0.001;
1909         final double maxStep                         = 120.0;
1910         final double positionTolerance               = 0.1;
1911         final int    degree                          = 20;
1912         final int    order                           = 20;
1913         final double spacecraftArea                  = 1.0;
1914         final double spacecraftDragCoefficient       = 2.0;
1915         final double spacecraftReflectionCoefficient = 2.0;
1916 
1917         // propagator main configuration
1918         final double[][] tol           = ToleranceProvider.of(CartesianToleranceProvider.of(positionTolerance))
1919                 .getTolerances(spacecraftState.getOrbit(), orbitType);
1920         final ODEIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, tol[0], tol[1]);
1921         final NumericalPropagator np   = new NumericalPropagator(integrator);
1922         np.setOrbitType(orbitType);
1923         np.setPositionAngleType(angleType);
1924         np.setInitialState(spacecraftState);
1925 
1926         // Earth gravity field
1927         final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
1928                                                             Constants.WGS84_EARTH_FLATTENING,
1929                                                             FramesFactory.getITRF(IERSConventions.IERS_2010, true));
1930         final NormalizedSphericalHarmonicsProvider harmonicsGravityProvider = GravityFieldFactory.getNormalizedProvider(degree, order);
1931         np.addForceModel(new HolmesFeatherstoneAttractionModel(earth.getBodyFrame(), harmonicsGravityProvider));
1932 
1933         // Sun and Moon attraction
1934         np.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getSun()));
1935         np.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getMoon()));
1936 
1937         // atmospheric drag
1938         MarshallSolarActivityFutureEstimation msafe =
1939                         new MarshallSolarActivityFutureEstimation("Jan2000F10-edited-data\\.txt",
1940                                                                   MarshallSolarActivityFutureEstimation.StrengthLevel.AVERAGE);
1941         DTM2000 atmosphere = new DTM2000(msafe, CelestialBodyFactory.getSun(), earth);
1942         np.addForceModel(new DragForce(atmosphere, new IsotropicDrag(spacecraftArea, spacecraftDragCoefficient)));
1943 
1944         // solar radiation pressure
1945         np.addForceModel(new SolarRadiationPressure(CelestialBodyFactory.getSun(), earth,
1946                                                     new IsotropicRadiationSingleCoefficient(spacecraftArea, spacecraftReflectionCoefficient)));
1947 
1948         return np;
1949 
1950     }
1951 
1952     private CartesianOrbit createEllipticOrbit() {
1953         final AbsoluteDate date         = new AbsoluteDate("2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
1954         final Vector3D     position     = new Vector3D(6896874.444705,  1956581.072644,  -147476.245054);
1955         final Vector3D     velocity     = new Vector3D(169.816407662, -1126.783301861, -7332.745712770);
1956         final TimeStampedPVCoordinates pv = new TimeStampedPVCoordinates(date, position, velocity);
1957         final Frame frame = FramesFactory.getEME2000();
1958         final double mu   = Constants.EIGEN5C_EARTH_MU;
1959         return new CartesianOrbit(pv, frame, mu);
1960     }
1961 
1962     private CartesianOrbit createHyperbolicOrbit() {
1963         final AbsoluteDate date         = new AbsoluteDate("2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
1964         final Vector3D     position     = new Vector3D(224267911.905821, 290251613.109399, 45534292.777492);
1965         final Vector3D     velocity     = new Vector3D(-1494.068165293, 1124.771027677, 526.915286134);
1966         final TimeStampedPVCoordinates pv = new TimeStampedPVCoordinates(date, position, velocity);
1967         final Frame frame = FramesFactory.getEME2000();
1968         final double mu   = Constants.EIGEN5C_EARTH_MU;
1969         return new CartesianOrbit(pv, frame, mu);
1970     }
1971 
1972     @ParameterizedTest
1973     @ValueSource(ints = {6, 7})
1974     void testClearMatricesComputation(final int stateDimension) {
1975         // GIVEN
1976         propagator.setupMatricesComputation("stm", MatrixUtils.createRealIdentityMatrix(stateDimension), null);
1977         propagator.setInitialState(initialState);
1978         propagator.setResetAtEnd(false);
1979         final Maneuver maneuver = new ConstantThrustManeuver(initDate, 10, 1., 100, Vector3D.MINUS_I);
1980         propagator.addForceModel(maneuver);
1981         maneuver.getParameterDriver(ParameterDrivenDateIntervalDetector.MEDIAN_SUFFIX).setSelected(true);
1982         maneuver.getParameterDriver(ParameterDrivenDateIntervalDetector.DURATION_SUFFIX).setSelected(true);
1983         final AbsoluteDate targetDate = initialState.getDate().shiftedBy(1);
1984         propagator.propagate(targetDate);
1985         // WHEN
1986         propagator.clearMatricesComputation();
1987         // THEN
1988         Assertions.assertTrue(propagator.getAdditionalDataProviders().isEmpty());
1989         Assertions.assertTrue(propagator.getAdditionalDerivativesProviders().isEmpty());
1990     }
1991 
1992     @Test
1993     void testFinalAttitudeWithForcesNeedingRates() {
1994         final ForceModel dummyForceDependingOnRates = new ForceModel() {
1995             @Override
1996             public boolean dependsOnAttitudeRate() {
1997                 return true;
1998             }
1999 
2000             @Override
2001             public boolean dependsOnPositionOnly() {
2002                 return false;
2003             }
2004 
2005             @Override
2006             public Vector3D acceleration(SpacecraftState s, double[] parameters) {
2007                 return Vector3D.ZERO;
2008             }
2009 
2010             @Override
2011             public <T extends CalculusFieldElement<T>> FieldVector3D<T> acceleration(FieldSpacecraftState<T> s, T[] parameters) {
2012                 return null;
2013             }
2014 
2015             @Override
2016             public List<ParameterDriver> getParametersDrivers() {
2017                 return Collections.emptyList();
2018             }
2019         };
2020         final List<ForceModel> forceModels = new ArrayList<>();
2021         forceModels.add(dummyForceDependingOnRates);
2022         testTemplateFinalAttitudeWithoutForcesNeedingRates(forceModels);
2023     }
2024 
2025     @Test
2026     void testFinalAttitudeWithoutForcesNeedingRates() {
2027         testTemplateFinalAttitudeWithoutForcesNeedingRates(new ArrayList<>());
2028     }
2029 
2030     private void testTemplateFinalAttitudeWithoutForcesNeedingRates(final List<ForceModel> forceModels) {
2031         // GIVEN
2032         final AttitudeProvider attitudeProvider = createAttitudeProviderWithNonZeroRates();
2033         propagator.setAttitudeProvider(attitudeProvider);
2034         for (final ForceModel forceModel : forceModels) {
2035             propagator.addForceModel(forceModel);
2036         }
2037         // WHEN
2038         final SpacecraftState state = propagator.propagate(propagator.getInitialState().getDate().shiftedBy(10.));
2039         // THEN
2040         final Attitude attitude = state.getAttitude();
2041         Assertions.assertNotEquals(Vector3D.ZERO, attitude.getSpin());
2042         Assertions.assertNotEquals(Vector3D.ZERO, attitude.getRotationAcceleration());
2043     }
2044 
2045     private AttitudeProvider createAttitudeProviderWithNonZeroRates() {
2046         return new AttitudeProvider() {
2047             @Override
2048             public Attitude getAttitude(PVCoordinatesProvider pvProv, AbsoluteDate date, Frame frame) {
2049                 return createAttitudeWithNonZeroRates(date, frame);
2050             }
2051 
2052             @Override
2053             public <T extends CalculusFieldElement<T>> FieldAttitude<T> getAttitude(FieldPVCoordinatesProvider<T> pvProv, FieldAbsoluteDate<T> date, Frame frame) {
2054                 return null;
2055             }
2056         };
2057     }
2058 
2059     private Attitude createAttitudeWithNonZeroRates(final AbsoluteDate date, final Frame frame) {
2060         final AngularCoordinates angularCoordinates = new AngularCoordinates(Rotation.IDENTITY,
2061                 Vector3D.PLUS_K, Vector3D.MINUS_I);
2062         return new Attitude(date, frame, angularCoordinates);
2063     }
2064 
2065     @BeforeEach
2066     public void setUp() {
2067         Utils.setDataRoot("regular-data:potential/shm-format");
2068         GravityFieldFactory.addPotentialCoefficientsReader(new SHMFormatReader("^eigen_cg03c_coef$", false));
2069         mu  = GravityFieldFactory.getUnnormalizedProvider(0, 0).getMu();
2070         final Vector3D position = new Vector3D(7.0e6, 1.0e6, 4.0e6);
2071         final Vector3D velocity = new Vector3D(-500.0, 8000.0, 1000.0);
2072         initDate = AbsoluteDate.J2000_EPOCH;
2073         final Orbit orbit = new EquinoctialOrbit(new PVCoordinates(position,  velocity),
2074                                                  FramesFactory.getEME2000(), initDate, mu);
2075         initialState = new SpacecraftState(orbit);
2076         double[][] tolerance = ToleranceProvider.of(CartesianToleranceProvider.of(0.001)).getTolerances(orbit,
2077                 OrbitType.EQUINOCTIAL, PositionAngleType.TRUE);
2078         AdaptiveStepsizeIntegrator integrator =
2079                 new DormandPrince853Integrator(0.001, 200, tolerance[0], tolerance[1]);
2080         integrator.setInitialStepSize(60);
2081         propagator = new NumericalPropagator(integrator);
2082         propagator.setInitialState(initialState);
2083     }
2084 
2085     @AfterEach
2086     public void tearDown() {
2087         initDate = null;
2088         initialState = null;
2089         propagator = null;
2090     }
2091 
2092     /**
2093      * Adapter class for {@link ForceModel} so that sub classes only have to implement the
2094      * methods they want.
2095      */
2096     private static class ForceModelAdapter implements ForceModel {
2097 
2098         @Override
2099         public boolean dependsOnPositionOnly() {
2100             return false;
2101         }
2102 
2103         @Override
2104         public boolean isSupported(String name) {
2105             return false;
2106         }
2107 
2108         @Override
2109         public void addContribution(SpacecraftState s, TimeDerivativesEquations adder) {
2110         }
2111 
2112         @Override
2113         public <T extends CalculusFieldElement<T>> void
2114         addContribution(FieldSpacecraftState<T> s,
2115                         FieldTimeDerivativesEquations<T> adder) {
2116         }
2117 
2118         /** {@inheritDoc} */
2119         @Override
2120         public Vector3D acceleration(final SpacecraftState s, final double[] parameters)
2121             {
2122             return Vector3D.ZERO;
2123         }
2124 
2125         /** {@inheritDoc} */
2126         @Override
2127         public <T extends CalculusFieldElement<T>> FieldVector3D<T> acceleration(final FieldSpacecraftState<T> s,
2128                                                                              final T[] parameters)
2129             {
2130             return FieldVector3D.getZero(s.getDate().getField());
2131         }
2132 
2133         @Override
2134         public Stream<EventDetector> getEventDetectors() {
2135             return Stream.empty();
2136         }
2137 
2138         @Override
2139         public <T extends CalculusFieldElement<T>> Stream<FieldEventDetector<T>> getFieldEventDetectors(final Field<T> field) {
2140             return Stream.empty();
2141         }
2142 
2143         @Override
2144         public List<ParameterDriver> getParametersDrivers() {
2145             return Collections.emptyList();
2146         }
2147     }
2148 
2149 }
2150