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