Commit b8c08b67 authored by Luc Maisonobe's avatar Luc Maisonobe
Browse files

Fixed ephemeris generation with several derivatives providers.

Fixes #919
parent d8f94d68
......@@ -21,6 +21,9 @@
</properties>
<body>
<release version="11.2" date="TBD" description="TBD">
<action dev="luc" type="fix" issue="919">
Fixed ephemeris generation with several derivatives providers.
</action>
<action dev="luc" type="add" issue="918">
Added support for ITRF-2020.
</action>
......
......@@ -1115,15 +1115,18 @@ public abstract class AbstractIntegratedPropagator extends AbstractPropagator {
}
// get the names of additional states managed by differential equations
final String[] names = new String[additionalDerivativesProviders.size()];
final String[] names = new String[additionalDerivativesProviders.size()];
final int[] dimensions = new int[additionalDerivativesProviders.size()];
for (int i = 0; i < names.length; ++i) {
names[i] = additionalDerivativesProviders.get(i).getName();
dimensions[i] = additionalDerivativesProviders.get(i).getDimension();
}
// create the ephemeris
ephemeris = new IntegratedEphemeris(startDate, minDate, maxDate,
stateMapper, propagationType, model,
unmanaged, getAdditionalStateProviders(), names);
unmanaged, getAdditionalStateProviders(),
names, dimensions);
}
......
......@@ -16,6 +16,7 @@
*/
package org.orekit.propagation.integration;
import java.util.Arrays;
import java.util.List;
import java.util.Map;
......@@ -23,6 +24,7 @@ import org.hipparchus.ode.DenseOutputModel;
import org.hipparchus.ode.ODEStateAndDerivative;
import org.orekit.attitudes.AttitudeProvider;
import org.orekit.errors.OrekitException;
import org.orekit.errors.OrekitInternalError;
import org.orekit.errors.OrekitMessages;
import org.orekit.frames.Frame;
import org.orekit.orbits.Orbit;
......@@ -107,9 +109,9 @@ public class IntegratedEphemeris
* @param unmanaged unmanaged additional states that must be simply copied
* @param providers providers for pre-integrated states
* @param equations names of additional equations
* @deprecated as of 11.1, replaced by {@link #IntegratedEphemeris(AbsoluteDate,
* @deprecated as of 11.1.2, replaced by {@link #IntegratedEphemeris(AbsoluteDate,
* AbsoluteDate, AbsoluteDate, StateMapper, PropagationType, DenseOutputModel,
* DoubleArrayDictionary, List, String[])}
* DoubleArrayDictionary, List, String[], int[])}
*/
@Deprecated
public IntegratedEphemeris(final AbsoluteDate startDate,
......@@ -134,7 +136,11 @@ public class IntegratedEphemeris
* @param providers providers for pre-integrated states
* @param equations names of additional equations
* @since 11.1
* @deprecated as of 11.1.2, replaced by {@link #IntegratedEphemeris(AbsoluteDate,
* AbsoluteDate, AbsoluteDate, StateMapper, PropagationType, DenseOutputModel,
* DoubleArrayDictionary, List, String[], int[])}
*/
@Deprecated
public IntegratedEphemeris(final AbsoluteDate startDate,
final AbsoluteDate minDate, final AbsoluteDate maxDate,
final StateMapper mapper, final PropagationType type,
......@@ -142,6 +148,31 @@ public class IntegratedEphemeris
final DoubleArrayDictionary unmanaged,
final List<AdditionalStateProvider> providers,
final String[] equations) {
this(startDate, minDate, maxDate, mapper, type, model,
unmanaged, providers, equations,
remainingDimensions(model, unmanaged, providers, equations));
}
/** Creates a new instance of IntegratedEphemeris.
* @param startDate Start date of the integration (can be minDate or maxDate)
* @param minDate first date of the range
* @param maxDate last date of the range
* @param mapper mapper between raw double components and spacecraft state
* @param type type of orbit to output (mean or osculating)
* @param model underlying raw mathematical model
* @param unmanaged unmanaged additional states that must be simply copied
* @param providers providers for pre-integrated states
* @param equations names of additional equations
* @param dimensions dimensions of additional equations
* @since 11.1.2
*/
public IntegratedEphemeris(final AbsoluteDate startDate,
final AbsoluteDate minDate, final AbsoluteDate maxDate,
final StateMapper mapper, final PropagationType type,
final DenseOutputModel model,
final DoubleArrayDictionary unmanaged,
final List<AdditionalStateProvider> providers,
final String[] equations, final int[] dimensions) {
super(mapper.getAttitudeProvider());
......@@ -160,11 +191,37 @@ public class IntegratedEphemeris
// set up providers to map the final elements of the model array to additional states
for (int i = 0; i < equations.length; ++i) {
addAdditionalStateProvider(new LocalGenerator(equations[i], i));
addAdditionalStateProvider(new LocalGenerator(equations[i], i, dimensions[i]));
}
}
/** Compute remaining dimensions for additional equations.
* @param model underlying raw mathematical model
* @param unmanaged unmanaged additional states that must be simply copied
* @param providers providers for pre-integrated states
* @param equations names of additional equations
* @return dimensions of additional equations
* @deprecated as of 11.1.2 this method is temporary and should be removed
* when the calling constructors are removed
* @since 11.1.2
*/
@Deprecated
private static int[] remainingDimensions(final DenseOutputModel model,
final DoubleArrayDictionary unmanaged,
final List<AdditionalStateProvider> providers,
final String[] equations) {
final ODEStateAndDerivative osd = model.getInterpolatedState(model.getInitialTime());
if (equations.length != osd.getNumberOfSecondaryStates()) {
throw new OrekitInternalError(null);
}
final int[] dimensions = new int[equations.length];
for (int i = 0; i < dimensions.length; ++i) {
dimensions[i] = osd.getSecondaryStateDimension(i + 1);
}
return dimensions;
}
/** Interpolate the model at some date.
* @param date desired interpolation date
* @return state interpolated at date
......@@ -270,13 +327,18 @@ public class IntegratedEphemeris
/** Index of the additional state. */
private final int index;
/** Dimension of the additional state. */
private final int dimension;
/** Simple constructor.
* @param name name of the additional state
* @param index index of the additional state
* @param dimension dimension of the additional state
*/
LocalGenerator(final String name, final int index) {
this.name = name;
this.index = index;
LocalGenerator(final String name, final int index, final int dimension) {
this.name = name;
this.index = index;
this.dimension = dimension;
}
/** {@inheritDoc} */
......@@ -288,7 +350,8 @@ public class IntegratedEphemeris
public double[] getAdditionalState(final SpacecraftState state) {
// extract the part of the interpolated array corresponding to the additional state
return getInterpolatedState(state.getDate()).getSecondaryState(index + 1);
final double[] combined = getInterpolatedState(state.getDate()).getSecondaryState(1);
return Arrays.copyOfRange(combined, index, index + dimension);
}
......
......@@ -16,11 +16,19 @@
*/
package org.orekit.propagation.integration;
import java.util.Collections;
import org.hipparchus.geometry.euclidean.threed.Vector3D;
import org.hipparchus.linear.MatrixUtils;
import org.hipparchus.linear.RealMatrix;
import org.hipparchus.ode.DenseOutputModel;
import org.hipparchus.ode.ExpandableODE;
import org.hipparchus.ode.ODEState;
import org.hipparchus.ode.OrdinaryDifferentialEquation;
import org.hipparchus.ode.SecondaryODE;
import org.hipparchus.ode.nonstiff.AdaptiveStepsizeIntegrator;
import org.hipparchus.ode.nonstiff.DormandPrince853Integrator;
import org.hipparchus.ode.nonstiff.EulerIntegrator;
import org.junit.Assert;
import org.junit.Before;
import org.junit.Test;
......@@ -28,15 +36,18 @@ import org.orekit.Utils;
import org.orekit.attitudes.CelestialBodyPointed;
import org.orekit.attitudes.InertialProvider;
import org.orekit.bodies.CelestialBodyFactory;
import org.orekit.errors.OrekitInternalError;
import org.orekit.forces.gravity.potential.GravityFieldFactory;
import org.orekit.forces.gravity.potential.ICGEMFormatReader;
import org.orekit.frames.FramesFactory;
import org.orekit.orbits.EquinoctialOrbit;
import org.orekit.orbits.Orbit;
import org.orekit.orbits.OrbitType;
import org.orekit.orbits.PositionAngle;
import org.orekit.propagation.BoundedPropagator;
import org.orekit.propagation.EphemerisGenerator;
import org.orekit.propagation.MatricesHarvester;
import org.orekit.propagation.PropagationType;
import org.orekit.propagation.SpacecraftState;
import org.orekit.propagation.analytical.KeplerianPropagator;
import org.orekit.propagation.numerical.NumericalPropagator;
......@@ -155,6 +166,54 @@ public class IntegratedEphemerisTest {
}
@Deprecated
@Test
public void testDeprecated() {
EulerIntegrator integ = new EulerIntegrator(1.0);
DenseOutputModel dom = new DenseOutputModel();
integ.addStepHandler(dom);
ExpandableODE eode = new ExpandableODE(new OrdinaryDifferentialEquation() {
public int getDimension() { return 1; }
public double[] computeDerivatives(double t, double[] y) { return y; }
});
eode.addSecondaryEquations(new SecondaryODE() {
public int getDimension() { return 1; }
public double[] computeDerivatives(double t, double[] primary,
double[] primaryDot, double[] secondary) { return secondary; }
});
integ.integrate(eode, new ODEState(0.0, new double[1], new double[1][1]), 1.0);
StateMapper mapper = new StateMapper(AbsoluteDate.ARBITRARY_EPOCH, Constants.EIGEN5C_EARTH_MU,
OrbitType.CARTESIAN, PositionAngle.TRUE,
new InertialProvider(FramesFactory.getEME2000()),
FramesFactory.getEME2000()) {
public void mapStateToArray(SpacecraftState state, double[] y, double[] yDot) {}
public SpacecraftState mapArrayToState(AbsoluteDate date, double[] y, double[] yDot, PropagationType type) {
return null;
}
};
try {
new IntegratedEphemeris(AbsoluteDate.ARBITRARY_EPOCH,
AbsoluteDate.PAST_INFINITY, AbsoluteDate.FUTURE_INFINITY,
mapper, PropagationType.OSCULATING,
dom, Collections.emptyMap(), Collections.emptyList(),
new String[] { "equation-1", "equation-2" });
Assert.fail("an exception should have been thrown");
} catch (OrekitInternalError oie) {
// expected as only one equation could be handled properly by this deprecated constructor
}
IntegratedEphemeris ie = new IntegratedEphemeris(AbsoluteDate.ARBITRARY_EPOCH,
AbsoluteDate.PAST_INFINITY, AbsoluteDate.FUTURE_INFINITY,
mapper, PropagationType.OSCULATING,
dom, Collections.emptyMap(), Collections.emptyList(),
new String[] { "equation-1" });
Assert.assertNotNull(ie);
}
@Before
public void setUp() {
......
......@@ -25,6 +25,7 @@ import java.util.List;
import java.util.concurrent.ExecutionException;
import java.util.concurrent.atomic.AtomicInteger;
import java.util.function.Consumer;
import java.util.stream.IntStream;
import java.util.stream.Stream;
import org.hamcrest.CoreMatchers;
......@@ -81,6 +82,8 @@ import org.orekit.propagation.BoundedPropagator;
import org.orekit.propagation.EphemerisGenerator;
import org.orekit.propagation.FieldSpacecraftState;
import org.orekit.propagation.SpacecraftState;
import org.orekit.propagation.conversion.DormandPrince853IntegratorBuilder;
import org.orekit.propagation.conversion.NumericalPropagatorBuilder;
import org.orekit.propagation.events.AbstractDetector;
import org.orekit.propagation.events.ApsideDetector;
import org.orekit.propagation.events.DateDetector;
......@@ -1523,6 +1526,71 @@ public class NumericalPropagatorTest {
}
@Test
public void testAdditionalDerivatives() {
final Frame eme2000 = FramesFactory.getEME2000();
final AbsoluteDate date = new AbsoluteDate(new DateComponents(2008, 6, 23),
new TimeComponents(14, 0, 0),
TimeScalesFactory.getUTC());
final Orbit initialOrbit = new KeplerianOrbit(8000000.0, 0.01, 0.87, 2.44, 0.21, -1.05, PositionAngle.MEAN,
eme2000,
date, Constants.EIGEN5C_EARTH_MU);
NumericalPropagatorBuilder builder = new NumericalPropagatorBuilder(initialOrbit, new DormandPrince853IntegratorBuilder(0.02, 0.2, 1.), PositionAngle.TRUE, 10);
NumericalPropagator propagator = builder.buildPropagator(builder.getSelectedNormalizedParameters());
IntStream.
range(0, 2).
mapToObj(i -> new EmptyDerivativeProvider("test_provider_" + i, new double[] { 10 * i })).
forEach(provider -> addDerivativeProvider(propagator, provider));
EphemerisGenerator generator = propagator.getEphemerisGenerator();
propagator.propagate(initialOrbit.getDate().shiftedBy(600));
BoundedPropagator ephemeris = generator.getGeneratedEphemeris();
final SpacecraftState finalState = ephemeris.propagate(initialOrbit.getDate().shiftedBy(300));
Assert.assertEquals(2, finalState.getAdditionalStatesValues().size());
Assert.assertEquals(1, finalState.getAdditionalState("test_provider_0").length);
Assert.assertEquals(0.0, finalState.getAdditionalState("test_provider_0")[0], 1.0e-15);
Assert.assertEquals(1, finalState.getAdditionalState("test_provider_1").length);
Assert.assertEquals(10.0, finalState.getAdditionalState("test_provider_1")[0], 1.0e-15);
}
private void addDerivativeProvider(NumericalPropagator propagator, EmptyDerivativeProvider provider) {
SpacecraftState initialState = propagator.getInitialState();
propagator.setInitialState(initialState.addAdditionalState(provider.getName(), provider.getInitialState()));
propagator.addAdditionalDerivativesProvider(provider);
}
private static class EmptyDerivativeProvider implements AdditionalDerivativesProvider {
private final String name;
private final double[] state;
public EmptyDerivativeProvider(String name, double[] state) {
this.name = name;
this.state = state;
}
@Override
public double[] derivatives(SpacecraftState s) {
return new double[getDimension()];
}
@Override
public String getName() {
return name;
}
@Override
public int getDimension() {
return state.length;
}
public double[] getInitialState() {
return state;
}
}
/** Record the dates treated by the handler.
* If they are out of an interval defined by a start and final date.
*/
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment