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;
18  
19  import org.hipparchus.ode.ODEIntegrator;
20  import org.hipparchus.ode.nonstiff.DormandPrince853Integrator;
21  import org.hipparchus.util.FastMath;
22  import org.junit.jupiter.api.AfterEach;
23  import org.junit.jupiter.api.Assertions;
24  import org.junit.jupiter.api.BeforeEach;
25  import org.junit.jupiter.api.Test;
26  import org.orekit.Utils;
27  import org.orekit.attitudes.AttitudeProvider;
28  import org.orekit.attitudes.BodyCenterPointing;
29  import org.orekit.bodies.OneAxisEllipsoid;
30  import org.orekit.errors.OrekitException;
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.events.DateDetector;
44  import org.orekit.propagation.numerical.NumericalPropagator;
45  import org.orekit.propagation.semianalytical.dsst.DSSTPropagator;
46  import org.orekit.propagation.semianalytical.dsst.forces.DSSTForceModel;
47  import org.orekit.propagation.semianalytical.dsst.forces.DSSTZonal;
48  import org.orekit.time.AbsoluteDate;
49  import org.orekit.time.DateComponents;
50  import org.orekit.time.TimeComponents;
51  import org.orekit.time.TimeScalesFactory;
52  import org.orekit.utils.Constants;
53  import org.orekit.utils.IERSConventions;
54  
55  import java.util.Arrays;
56  import java.util.List;
57  import java.util.stream.Collectors;
58  
59  
60  public class PropagatorsParallelizerEphemerisTest {
61  
62      private double mass;
63      private Orbit orbit;
64      private AttitudeProvider attitudeLaw;
65      private UnnormalizedSphericalHarmonicsProvider unnormalizedGravityField;
66      private NormalizedSphericalHarmonicsProvider normalizedGravityField;
67  
68      @BeforeEach
69      public void setUp() {
70          try {
71              Utils.setDataRoot("regular-data:potential/icgem-format");
72              unnormalizedGravityField = GravityFieldFactory.getUnnormalizedProvider(6, 0);
73              normalizedGravityField   = GravityFieldFactory.getNormalizedProvider(6, 0);
74  
75              mass = 2500;
76              double a = 7187990.1979844316;
77              double e = 0.5e-4;
78              double i = 1.7105407051081795;
79              double omega = 1.9674147913622104;
80              double OMEGA = FastMath.toRadians(261);
81              double lv = 0;
82  
83              AbsoluteDate date = new AbsoluteDate(new DateComponents(2004, 01, 01),
84                      TimeComponents.H00,
85                      TimeScalesFactory.getUTC());
86              orbit = new KeplerianOrbit(a, e, i, omega, OMEGA, lv, PositionAngleType.TRUE,
87                      FramesFactory.getEME2000(), date, normalizedGravityField.getMu());
88  
89              OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
90                      Constants.WGS84_EARTH_FLATTENING,
91                      FramesFactory.getITRF(IERSConventions.IERS_2010, true));
92              attitudeLaw = new BodyCenterPointing(orbit.getFrame(), earth);
93  
94          } catch (OrekitException oe) {
95              Assertions.fail(oe.getLocalizedMessage());
96          }
97      }
98  
99      @AfterEach
100     public void tearDown() {
101         mass                     = Double.NaN;
102         orbit                    = null;
103         attitudeLaw              = null;
104         unnormalizedGravityField = null;
105         normalizedGravityField   = null;
106     }
107 
108 
109 
110     /*
111      * Test set to check that ephemeris are produced at the end of the propagation parallelizer
112      * in the case of NuericalPropagators.
113      * Check existence of all ephemeris.
114      * Check end time of all ephemeris for varying step times.
115      *
116      * Should validate the modification, except for the retrieveNextParameters added check.
117      *
118      * Tests are based on Anne-Laure LUGAN's proposed test, and on PropagatorsParallelizerTest.
119      */
120 
121     @Test
122     public void testSeveralEphemeris() {
123         /*
124          * The closing behaviour is checked by verifying the presence of generated ephemeris as a result of the
125          * following process, using PropagatorsParallelizer.
126          */
127         final AbsoluteDate startDate = orbit.getDate();
128         final AbsoluteDate endDate = startDate.shiftedBy(3600);
129 
130         List<Propagator> propagators = Arrays.asList(buildNumerical(), buildNumerical(), buildDSST(), buildEcksteinHechler());
131         List<EphemerisGenerator> generators = propagators.stream().map(Propagator::getEphemerisGenerator).collect(Collectors.toList());
132         PropagatorsParallelizer parallelizer = new PropagatorsParallelizer(propagators, interpolators -> {
133             // Do nothing
134         });
135 
136         parallelizer.propagate(startDate, endDate);
137         for ( EphemerisGenerator generator : generators ) {
138             Assertions.assertNotNull(generator.getGeneratedEphemeris());
139             Assertions.assertEquals(endDate, generator.getGeneratedEphemeris().getMaxDate());
140             Assertions.assertEquals(startDate, generator.getGeneratedEphemeris().getMinDate());
141         }
142     }
143 
144     @Test
145     public void testSeveralEphemerisDateDetector() {
146         /*
147          * The closing behaviour is checked for a stop event occuring during the propagation.
148          * The test isn't applied to analytical propagators as their behaviour differs.
149          * (Analytical propagator's ephemerides are the analytical propagators.)
150          */
151         final AbsoluteDate startDate = orbit.getDate();
152         final AbsoluteDate endDate = startDate.shiftedBy(3600.015);
153         List<Propagator> propagators = Arrays.asList(buildNumerical(1,300), buildNumerical(0.001,300), buildDSST());
154 
155         // Add new instance of event with same date. DateDetector behaviour at event is stop.
156         AbsoluteDate detectorDate = startDate.shiftedBy(1800);
157         propagators.stream().forEach(propagator -> propagator.addEventDetector(new DateDetector(detectorDate)));
158 
159         List<EphemerisGenerator> generators = propagators.stream().map(Propagator::getEphemerisGenerator).collect(Collectors.toList());
160         PropagatorsParallelizer parallelizer = new PropagatorsParallelizer(propagators, interpolators -> {
161             // Do nothing
162         });
163 
164         parallelizer.propagate(startDate, endDate);
165 
166         // Check for all generators
167         for ( EphemerisGenerator generator : generators ) {
168             Assertions.assertNotNull(generator.getGeneratedEphemeris());
169             Assertions.assertEquals(startDate, generator.getGeneratedEphemeris().getMinDate());
170             Assertions.assertEquals(detectorDate, generator.getGeneratedEphemeris().getMaxDate());
171         }
172     }
173 
174 
175 
176     private EcksteinHechlerPropagator buildEcksteinHechler() {
177         return new EcksteinHechlerPropagator(orbit, attitudeLaw, mass, unnormalizedGravityField);
178     }
179 
180     private DSSTPropagator buildDSST(final double minStep, final double maxStep) {
181         // Gravity
182         Utils.setDataRoot("regular-data:potential/icgem-format");
183         GravityFieldFactory.addPotentialCoefficientsReader(new ICGEMFormatReader("^eigen-6s-truncated$", false));
184         UnnormalizedSphericalHarmonicsProvider gravity = GravityFieldFactory.getUnnormalizedProvider(8, 8);
185 
186         // Propagator
187         final double[][] tol = ToleranceProvider.getDefaultToleranceProvider(0.01).getTolerances(orbit, OrbitType.EQUINOCTIAL);
188         final DSSTPropagator propagator = new DSSTPropagator(new DormandPrince853Integrator(minStep, maxStep, tol[0], tol[1]), PropagationType.MEAN);
189 
190         // Force models
191         final DSSTForceModel zonal = new DSSTZonal(gravity, 4, 3, 9);
192         propagator.addForceModel(zonal);
193 
194         propagator.setInitialState(new SpacecraftState(orbit));
195 
196         return propagator;
197     }
198 
199     private DSSTPropagator buildDSST() {
200         return buildDSST(0.01,300);
201     }
202 
203     private NumericalPropagator buildNumerical() {
204         return buildNumerical(0.001,300);
205     }
206 
207     private NumericalPropagator buildNumerical(double minStep, double maxStep) {
208         NumericalPropagator numericalPropagator = buildNotInitializedNumerical(minStep, maxStep);
209         numericalPropagator.setInitialState(new SpacecraftState(orbit,
210                 attitudeLaw.getAttitude(orbit,
211                         orbit.getDate(),
212                         orbit.getFrame()),
213                 mass));
214         return numericalPropagator;
215     }
216 
217 
218     private NumericalPropagator buildNotInitializedNumerical(double minStep, double maxStep) {
219         OrbitType type = OrbitType.CARTESIAN;
220         double[][] tolerances = ToleranceProvider.getDefaultToleranceProvider(10.).getTolerances(orbit, type);
221         ODEIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, tolerances[0], tolerances[1]);
222         NumericalPropagator numericalPropagator = new NumericalPropagator(integrator);
223         ForceModel gravity = new HolmesFeatherstoneAttractionModel(FramesFactory.getITRF(IERSConventions.IERS_2010, true),
224                 normalizedGravityField);
225         numericalPropagator.addForceModel(gravity);
226         return numericalPropagator;
227     }
228 
229 
230 }