1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.propagation;
18
19 import org.hipparchus.exception.LocalizedCoreFormats;
20 import org.hipparchus.geometry.euclidean.threed.Vector3D;
21 import org.hipparchus.ode.ODEIntegrator;
22 import org.hipparchus.ode.nonstiff.DormandPrince853Integrator;
23 import org.hipparchus.util.FastMath;
24 import org.junit.jupiter.api.*;
25 import org.orekit.Utils;
26 import org.orekit.attitudes.AttitudeProvider;
27 import org.orekit.attitudes.BodyCenterPointing;
28 import org.orekit.bodies.OneAxisEllipsoid;
29 import org.orekit.errors.OrekitException;
30 import org.orekit.errors.OrekitMessages;
31 import org.orekit.forces.ForceModel;
32 import org.orekit.forces.gravity.HolmesFeatherstoneAttractionModel;
33 import org.orekit.forces.gravity.potential.GravityFieldFactory;
34 import org.orekit.forces.gravity.potential.ICGEMFormatReader;
35 import org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider;
36 import org.orekit.forces.gravity.potential.UnnormalizedSphericalHarmonicsProvider;
37 import org.orekit.frames.FramesFactory;
38 import org.orekit.orbits.KeplerianOrbit;
39 import org.orekit.orbits.Orbit;
40 import org.orekit.orbits.OrbitType;
41 import org.orekit.orbits.PositionAngleType;
42 import org.orekit.propagation.analytical.EcksteinHechlerPropagator;
43 import org.orekit.propagation.analytical.KeplerianPropagator;
44 import org.orekit.propagation.events.DateDetector;
45 import org.orekit.propagation.events.handlers.StopOnEvent;
46 import org.orekit.propagation.integration.AdditionalDerivativesProvider;
47 import org.orekit.propagation.integration.CombinedDerivatives;
48 import org.orekit.propagation.numerical.NumericalPropagator;
49 import org.orekit.propagation.semianalytical.dsst.DSSTPropagator;
50 import org.orekit.propagation.semianalytical.dsst.forces.DSSTForceModel;
51 import org.orekit.propagation.semianalytical.dsst.forces.DSSTZonal;
52 import org.orekit.time.AbsoluteDate;
53 import org.orekit.time.DateComponents;
54 import org.orekit.time.TimeComponents;
55 import org.orekit.time.TimeScalesFactory;
56 import org.orekit.utils.Constants;
57 import org.orekit.utils.IERSConventions;
58
59 import java.util.ArrayList;
60 import java.util.Arrays;
61 import java.util.List;
62 import java.util.concurrent.TimeUnit;
63 import java.util.concurrent.atomic.AtomicInteger;
64
65
66 public class PropagatorsParallelizerTest {
67
68
69
70
71
72
73 @Test
74 public void testIssue717() {
75
76
77 Utils.setDataRoot("regular-data:potential/icgem-format");
78 GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("^eigen-6s-truncated$", false));
79 UnnormalizedSphericalHarmonicsProvider gravity = GravityFieldFactory.getUnnormalizedProvider(8, 8);
80
81
82 Orbit orbit = new KeplerianOrbit(15000000.0, 0.125, 1.25,
83 0.250, 1.375, 0.0625, PositionAngleType.MEAN,
84 FramesFactory.getEME2000(),
85 new AbsoluteDate(2000, 2, 24, 11, 35, 47.0, TimeScalesFactory.getUTC()),
86 gravity.getMu());
87
88
89
90 final double[][] tol = ToleranceProvider.getDefaultToleranceProvider(0.01).getTolerances(orbit, OrbitType.EQUINOCTIAL);
91 final DSSTPropagator propagator = new DSSTPropagator(new DormandPrince853Integrator(0.01, 600.0, tol[0], tol[1]), PropagationType.OSCULATING);
92
93
94 final DSSTForceModel zonal = new DSSTZonal(gravity, 4, 3, 9);
95 propagator.addForceModel(zonal);
96
97 propagator.setInitialState(new SpacecraftState(orbit));
98
99
100 final double deltaT = 30.0;
101
102 final PropagatorsParallelizer parallelizer =
103 new PropagatorsParallelizer(Arrays.asList(propagator), interpolators -> {interpolators.get(0).getCurrentState().getDate();});
104
105 final SpacecraftState state = parallelizer.propagate(orbit.getDate().shiftedBy(deltaT).shiftedBy(+1.0), orbit.getDate().shiftedBy(-2.0 * deltaT).shiftedBy(-1.0)).get(0);
106
107
108 Assertions.assertNotNull(state);
109
110 }
111
112
113
114
115
116
117
118 @Test
119 @Timeout(value = 5, unit = TimeUnit.SECONDS)
120 public void testIssue1105() {
121
122
123
124
125
126 final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
127 final Orbit orbit = new KeplerianOrbit(Constants.IERS2003_EARTH_EQUATORIAL_RADIUS + 1.e6,
128 0.125, 1.25, 0.250, 1.375, 0.0625, PositionAngleType.MEAN,
129 FramesFactory.getGCRF(), date, Constants.IERS2010_EARTH_MU);
130
131
132
133 final List<Propagator> propagators = new ArrayList<>();
134 propagators.add(new KeplerianPropagator(orbit));
135 propagators.add(new KeplerianPropagator(orbit));
136
137
138 final PropagatorsParallelizer parallelizer = new PropagatorsParallelizer(propagators, interpolators -> {});
139
140
141
142
143
144 SpacecraftState state = parallelizer.propagate(date, date.shiftedBy(10.)).get(0);
145
146
147
148 state = parallelizer.propagate(date.shiftedBy(10.), date.shiftedBy(20.)).get(0);
149
150
151
152
153
154 Assertions.assertNotNull(state);
155
156 }
157
158 @Test
159 public void testNumericalNotInitialized() {
160
161 final AbsoluteDate startDate = orbit.getDate();
162 final AbsoluteDate endDate = startDate.shiftedBy(3600.0);
163 List<Propagator> propagators = Arrays.asList(buildEcksteinHechler(),
164 buildNotInitializedNumerical());
165
166 PropagatorsParallelizer parallelizer =
167 new PropagatorsParallelizer(propagators, interpolators -> Assertions.fail("should not be called"));
168 try {
169 parallelizer.propagate(startDate, endDate);
170 Assertions.fail("an exception should have been thrown");
171 } catch (OrekitException oe) {
172 Assertions.assertEquals(OrekitMessages.INITIAL_STATE_NOT_SPECIFIED_FOR_ORBIT_PROPAGATION,
173 oe.getSpecifier());
174 }
175
176 }
177
178 @Test
179 public void testAnalyticalAndNumericalSameOrbit() {
180
181 final AbsoluteDate startDate = orbit.getDate();
182 final AbsoluteDate endDate = startDate.shiftedBy(3600.0);
183 List<Propagator> propagators = Arrays.asList(buildEcksteinHechler(),
184 buildNumerical());
185
186 PropagatorsParallelizer parallelizer =
187 new PropagatorsParallelizer(propagators,
188 interpolators -> {
189 Assertions.assertEquals(2, interpolators.size());
190 AbsoluteDate aPrev = interpolators.get(0).getPreviousState().getDate();
191 AbsoluteDate aCurr = interpolators.get(0).getCurrentState().getDate();
192 AbsoluteDate nPrev = interpolators.get(1).getPreviousState().getDate();
193 AbsoluteDate nCurr = interpolators.get(1).getCurrentState().getDate();
194 Assertions.assertEquals(0.0, aPrev.durationFrom(nPrev), 3.0e-13);
195 Assertions.assertEquals(0.0, aCurr.durationFrom(nCurr), 3.0e-13);
196 Vector3D aPos = interpolators.get(0).getCurrentState().getPosition();
197 Vector3D nPos = interpolators.get(1).getCurrentState().getPosition();
198 Assertions.assertTrue(Vector3D.distance(aPos, nPos) < 111.0);
199 });
200 List<SpacecraftState> results = parallelizer.propagate(startDate, endDate);
201
202 Assertions.assertEquals(2, results.size());
203 for (final SpacecraftState state : results) {
204 Assertions.assertEquals(0.0, state.getDate().durationFrom(endDate), 1.0e-15);
205 }
206
207 }
208
209 @Test
210 public void testVsAnalyticalMonoSat() {
211
212 final AbsoluteDate startDate = orbit.getDate();
213 final AbsoluteDate endDate = startDate.shiftedBy(3600.0);
214 Propagator mono = buildEcksteinHechler();
215 final EphemerisGenerator generator = mono.getEphemerisGenerator();
216 mono.propagate(startDate, endDate);
217 final BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
218
219 List<Propagator> propagators = Arrays.asList(buildEcksteinHechler(),
220 buildNumerical());
221
222 PropagatorsParallelizer parallelizer =
223 new PropagatorsParallelizer(propagators,
224 interpolators -> {
225 AbsoluteDate aCurr = interpolators.get(0).getCurrentState().getDate();
226 Vector3D aPos = interpolators.get(0).getCurrentState().getPosition();
227 Vector3D ePos = ephemeris.getPosition(aCurr, orbit.getFrame());
228 Assertions.assertEquals(0, Vector3D.distance(ePos, aPos), 1.0e-15);
229 });
230 List<SpacecraftState> results = parallelizer.propagate(startDate, endDate);
231
232 Assertions.assertEquals(2, parallelizer.getPropagators().size());
233 Assertions.assertEquals(2, results.size());
234 for (final SpacecraftState state : results) {
235 Assertions.assertEquals(0.0, state.getDate().durationFrom(endDate), 1.0e-15);
236 }
237
238 }
239
240 @Test
241 public void testVsNumericalMonoSat() {
242
243 final AbsoluteDate startDate = orbit.getDate();
244 final AbsoluteDate endDate = startDate.shiftedBy(3600.0);
245 Propagator mono = buildNumerical();
246 final EphemerisGenerator generator = mono.getEphemerisGenerator();
247 mono.propagate(startDate, endDate);
248 final BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
249
250 List<Propagator> propagators = Arrays.asList(buildEcksteinHechler(),
251 buildNumerical());
252
253 PropagatorsParallelizer parallelizer =
254 new PropagatorsParallelizer(propagators,
255 interpolators -> {
256 AbsoluteDate nCurr = interpolators.get(1).getCurrentState().getDate();
257 Vector3D nPos = interpolators.get(1).getCurrentState().getPosition();
258 Vector3D ePos = ephemeris.getPosition(nCurr, orbit.getFrame());
259 Assertions.assertEquals(0, Vector3D.distance(ePos, nPos), 1.0e-15);
260 });
261 List<SpacecraftState> results = parallelizer.propagate(startDate, endDate);
262
263 Assertions.assertEquals(2, results.size());
264 for (final SpacecraftState state : results) {
265 Assertions.assertEquals(0.0, state.getDate().durationFrom(endDate), 1.0e-15);
266 }
267
268 }
269
270 @Test
271 public void testOrekitException() {
272 final AbsoluteDate startDate = orbit.getDate();
273 final AbsoluteDate endDate = startDate.shiftedBy(3600.0);
274 List<Propagator> propagators = Arrays.asList(buildEcksteinHechler(),
275 buildNumerical());
276 propagators.get(0).addEventDetector(new DateDetector(startDate.shiftedBy(900.0)).
277 withHandler((state, detector, increasing) -> {
278 throw new OrekitException(LocalizedCoreFormats.SIMPLE_MESSAGE,
279 "inTest");
280 }));
281 try {
282 new PropagatorsParallelizer(propagators, interpolators -> {}).propagate(startDate, endDate);
283 Assertions.fail("an exception should have been thrown");
284 } catch (OrekitException oe) {
285 Assertions.assertNull(oe.getCause());
286 Assertions.assertEquals(LocalizedCoreFormats.SIMPLE_MESSAGE, oe.getSpecifier());
287 Assertions.assertEquals("inTest", (String) oe.getParts()[0]);
288 }
289 }
290
291 @Test
292 public void testEarlyStop() {
293 final AbsoluteDate startDate = orbit.getDate();
294 final AbsoluteDate endDate = startDate.shiftedBy(3600.0);
295 List<Propagator> propagators = Arrays.asList(buildEcksteinHechler(),
296 buildNumerical());
297 propagators.get(0).addEventDetector(new DateDetector(startDate.shiftedBy(900.0)).
298 withHandler((state, detector, increasing) -> {
299 throw new RuntimeException("boo!");
300 }));
301 try {
302 new PropagatorsParallelizer(propagators, interpolators -> {}).propagate(startDate, endDate);
303 Assertions.fail("an exception should have been thrown");
304 } catch (OrekitException oe) {
305 Assertions.assertNotNull(oe.getCause());
306 Assertions.assertTrue(oe.getCause() instanceof RuntimeException);
307 Assertions.assertEquals(LocalizedCoreFormats.SIMPLE_MESSAGE, oe.getSpecifier());
308 Assertions.assertTrue(((String) oe.getParts()[0]).endsWith("boo!"));
309 }
310 }
311
312 @Test
313 public void testStopOnEarlyEvent() {
314 final AbsoluteDate startDate = orbit.getDate();
315 final AbsoluteDate endDate = startDate.shiftedBy(3600.0);
316 final AbsoluteDate stopDate = startDate.shiftedBy(0.01);
317 List<Propagator> propagators = Arrays.asList(buildEcksteinHechler(),
318 buildNumerical());
319 propagators.get(0).addEventDetector(new DateDetector(stopDate).withHandler(new StopOnEvent()));
320 List<SpacecraftState> results = new PropagatorsParallelizer(propagators, interpolators -> {}).
321 propagate(startDate, endDate);
322 Assertions.assertEquals(2, results.size());
323 Assertions.assertEquals(0.0, results.get(0).getDate().durationFrom(stopDate), 1.0e-15);
324 Assertions.assertEquals(0.0, results.get(1).getDate().durationFrom(stopDate), 1.0e-15);
325 }
326
327 @Test
328 public void testStopOnLateEvent() {
329 final AbsoluteDate startDate = orbit.getDate();
330 final AbsoluteDate endDate = startDate.shiftedBy(3600.0);
331 final AbsoluteDate stopDate = startDate.shiftedBy(900.0);
332 List<Propagator> propagators = Arrays.asList(buildEcksteinHechler(),
333 buildNumerical());
334 propagators.get(0).addEventDetector(new DateDetector(stopDate).withHandler(new StopOnEvent()));
335 List<SpacecraftState> results = new PropagatorsParallelizer(propagators, interpolators -> {}).
336 propagate(startDate, endDate);
337 Assertions.assertEquals(2, results.size());
338 Assertions.assertEquals(0.0, results.get(0).getDate().durationFrom(stopDate), 1.0e-15);
339 Assertions.assertEquals(0.0, results.get(1).getDate().durationFrom(stopDate), 1.0e-15);
340 }
341
342 @Test
343 public void testInternalStepHandler() {
344 final AbsoluteDate startDate = orbit.getDate();
345 final AbsoluteDate endDate = startDate.shiftedBy(3600.0);
346 List<Propagator> propagators = Arrays.asList(buildEcksteinHechler(),
347 buildNumerical());
348 final AtomicInteger called0 = new AtomicInteger();
349 propagators.get(0).getMultiplexer().add(interpolator -> called0.set(11));
350 final AtomicInteger called1 = new AtomicInteger();
351 propagators.get(1).getMultiplexer().add(interpolator -> called1.set(22));
352 List<SpacecraftState> results = new PropagatorsParallelizer(propagators, interpolators -> {}).
353 propagate(startDate, endDate);
354 Assertions.assertEquals(2, results.size());
355 Assertions.assertEquals(0.0, results.get(0).getDate().durationFrom(endDate), 1.0e-15);
356 Assertions.assertEquals(0.0, results.get(1).getDate().durationFrom(endDate), 1.0e-15);
357 Assertions.assertEquals(11, called0.get());
358 Assertions.assertEquals(22, called1.get());
359 }
360
361 @Test
362 public void testFixedStepHandler() {
363 final AbsoluteDate startDate = orbit.getDate();
364 final AbsoluteDate endDate = startDate.shiftedBy(3600.0);
365 final double h = 60.0;
366 List<Propagator> propagators = Arrays.asList(buildEcksteinHechler(),
367 buildNumerical());
368 final AtomicInteger counter = new AtomicInteger();
369 new PropagatorsParallelizer(propagators,
370 h,
371 states -> {
372 for (final SpacecraftState state : states) {
373 Assertions.assertEquals(h * counter.get(),
374 state.getDate().durationFrom(startDate),
375 1.0e-10);
376 }
377 counter.addAndGet(1);
378 }).propagate(startDate, endDate);
379 Assertions.assertEquals(1 + (int) FastMath.rint(endDate.durationFrom(startDate) / h), counter.get());
380 }
381
382 @Test
383 public void testIntegrableGenerator() {
384 final AbsoluteDate startDate = orbit.getDate();
385 final AbsoluteDate endDate = startDate.shiftedBy(3600.0);
386 final NumericalPropagator p0 = buildNumerical();
387 final NumericalPropagator p1 = buildNumerical();
388 final String name = "generator";
389 final double base0 = 2.0e-3;
390 final double base1 = 2.5e-3;
391 p0.addAdditionalDerivativesProvider(new Exponential(name, base0));
392 p0.setInitialState(p0.getInitialState().addAdditionalData(name, 1.0));
393 p1.addAdditionalDerivativesProvider(new Exponential(name, base1));
394 p1.setInitialState(p1.getInitialState().addAdditionalData(name, 1.0));
395 List<SpacecraftState> results = new PropagatorsParallelizer(Arrays.asList(p0, p1), interpolators -> {}).
396 propagate(startDate, endDate);
397 double expected0 = FastMath.exp(base0 * endDate.durationFrom(startDate));
398 double expected1 = FastMath.exp(base1 * endDate.durationFrom(startDate));
399 Assertions.assertEquals(expected0, results.get(0).getAdditionalState(name)[0], 6.0e-9 * expected0);
400 Assertions.assertEquals(expected1, results.get(1).getAdditionalState(name)[0], 5.0e-8 * expected1);
401 }
402
403 private static class Exponential implements AdditionalDerivativesProvider {
404 final String name;
405 final double base;
406 Exponential(final String name, final double base) {
407 this.name = name;
408 this.base = base;
409 }
410 public String getName() {
411 return name;
412 }
413 public boolean yields(final SpacecraftState state) {
414 return !state.hasAdditionalData(name);
415 }
416 public CombinedDerivatives combinedDerivatives(SpacecraftState state) {
417 return new CombinedDerivatives(new double[] { base * state.getAdditionalState(name)[0] },
418 null);
419 }
420 public int getDimension() {
421 return 1;
422 }
423 }
424
425 private EcksteinHechlerPropagator buildEcksteinHechler() {
426 return new EcksteinHechlerPropagator(orbit, attitudeLaw, mass, unnormalizedGravityField);
427 }
428
429 private NumericalPropagator buildNumerical() {
430 NumericalPropagator numericalPropagator = buildNotInitializedNumerical();
431 numericalPropagator.setInitialState(new SpacecraftState(orbit,
432 attitudeLaw.getAttitude(orbit,
433 orbit.getDate(),
434 orbit.getFrame()),
435 mass));
436 return numericalPropagator;
437 }
438
439 private NumericalPropagator buildNotInitializedNumerical() {
440 OrbitType type = OrbitType.CARTESIAN;
441 double minStep = 0.001;
442 double maxStep = 300;
443 double[][] tolerances = ToleranceProvider.getDefaultToleranceProvider(10.).getTolerances(orbit, type);
444 ODEIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, tolerances[0], tolerances[1]);
445 NumericalPropagator numericalPropagator = new NumericalPropagator(integrator);
446 ForceModel gravity = new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true),
447 normalizedGravityField);
448 numericalPropagator.addForceModel(gravity);
449 return numericalPropagator;
450 }
451
452 @BeforeEach
453 public void setUp() {
454 try {
455 Utils.setDataRoot("regular-data:potential/icgem-format");
456 unnormalizedGravityField = GravityFieldFactory.getUnnormalizedProvider(6, 0);
457 normalizedGravityField = GravityFieldFactory.getNormalizedProvider(6, 0);
458
459 mass = 2500;
460 double a = 7187990.1979844316;
461 double e = 0.5e-4;
462 double i = 1.7105407051081795;
463 double omega = 1.9674147913622104;
464 double OMEGA = FastMath.toRadians(261);
465 double lv = 0;
466
467 AbsoluteDate date = new AbsoluteDate(new DateComponents(2004, 01, 01),
468 TimeComponents.H00,
469 TimeScalesFactory.getUTC());
470 orbit = new KeplerianOrbit(a, e, i, omega, OMEGA, lv, PositionAngleType.TRUE,
471 FramesFactory.getEME2000(), date, normalizedGravityField.getMu());
472 OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
473 Constants.WGS84_EARTH_FLATTENING,
474 FramesFactory.getITRF(IERSConventions.IERS_2010, true));
475 attitudeLaw = new BodyCenterPointing(orbit.getFrame(), earth);
476
477 } catch (OrekitException oe) {
478 Assertions.fail(oe.getLocalizedMessage());
479 }
480 }
481
482 @AfterEach
483 public void tearDown() {
484 mass = Double.NaN;
485 orbit = null;
486 attitudeLaw = null;
487 unnormalizedGravityField = null;
488 normalizedGravityField = null;
489 }
490
491 private double mass;
492 private Orbit orbit;
493 private AttitudeProvider attitudeLaw;
494 private UnnormalizedSphericalHarmonicsProvider unnormalizedGravityField;
495 private NormalizedSphericalHarmonicsProvider normalizedGravityField;
496
497 }