1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.control.indirect.shooting;
18
19 import org.hipparchus.analysis.differentiation.Gradient;
20 import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
21 import org.hipparchus.geometry.euclidean.threed.Vector3D;
22 import org.hipparchus.ode.FieldOrdinaryDifferentialEquation;
23 import org.hipparchus.util.Binary64;
24 import org.hipparchus.util.Binary64Field;
25 import org.hipparchus.util.MathArrays;
26 import org.junit.jupiter.api.Assertions;
27 import org.junit.jupiter.api.Test;
28 import org.orekit.TestUtils;
29 import org.orekit.control.indirect.adjoint.CartesianAdjointDerivativesProvider;
30 import org.orekit.control.indirect.adjoint.CartesianAdjointKeplerianTerm;
31 import org.orekit.control.indirect.shooting.propagation.*;
32 import org.orekit.forces.ForceModel;
33 import org.orekit.forces.gravity.NewtonianAttraction;
34 import org.orekit.propagation.FieldSpacecraftState;
35 import org.orekit.propagation.SpacecraftState;
36 import org.orekit.propagation.integration.CombinedDerivatives;
37 import org.orekit.propagation.numerical.NumericalPropagator;
38 import org.orekit.time.AbsoluteDate;
39 import org.orekit.time.FieldAbsoluteDate;
40 import org.orekit.utils.Constants;
41
42 import java.util.ArrayList;
43 import java.util.List;
44
45 class AbstractFixedInitialCartesianSingleShootingTest {
46
47 @Test
48 void createFieldInitialStateWithMassAndAdjointTest() {
49
50 final ShootingIntegrationSettings integrationSettings = ShootingIntegrationSettingsFactory.getClassicalRungeKuttaIntegratorSettings(1);
51 final ShootingPropagationSettings propagationSettings = new ShootingPropagationSettings(new ArrayList<>(),
52 CartesianAdjointDynamicsProviderFactory.buildUnboundedEnergyProviderNeglectingMass("adjoint"), integrationSettings);
53 final SpacecraftState state = new SpacecraftState(TestUtils.getDefaultOrbit(AbsoluteDate.ARBITRARY_EPOCH));
54 final TestSingleShooting testSingleShooting = new TestSingleShooting(propagationSettings, state);
55 final double expectedMass = 1.;
56 final double[] expectedAdjoint = new double[1];
57
58 final FieldSpacecraftState<Gradient> fieldState = testSingleShooting.createFieldInitialStateWithMassAndAdjoint(expectedMass,
59 expectedAdjoint);
60
61 Assertions.assertEquals(state.getDate(), fieldState.getDate().toAbsoluteDate());
62 Assertions.assertEquals(expectedMass, fieldState.getMass().getReal());
63 final FieldVector3D<Gradient> fieldPosition = fieldState.getPosition();
64 final FieldVector3D<Gradient> fieldVelocity = fieldState.getPVCoordinates().getVelocity();
65 Assertions.assertEquals(state.getPosition(), fieldPosition.toVector3D());
66 Assertions.assertEquals(state.getPVCoordinates().getVelocity(), fieldVelocity.toVector3D());
67 Assertions.assertEquals(1., fieldState.getAdditionalState("adjoint")[0].getGradient()[0]);
68 }
69
70 @Test
71 void buildPropagatorTest() {
72
73 final ShootingIntegrationSettings integrationSettings = ShootingIntegrationSettingsFactory.getClassicalRungeKuttaIntegratorSettings(1);
74 final ShootingPropagationSettings propagationSettings = new ShootingPropagationSettings(new ArrayList<>(),
75 CartesianAdjointDynamicsProviderFactory.buildUnboundedEnergyProviderNeglectingMass("adjoint"), integrationSettings);
76 final SpacecraftState state = new SpacecraftState(TestUtils.getDefaultOrbit(AbsoluteDate.ARBITRARY_EPOCH));
77 final TestSingleShooting testSingleShooting = new TestSingleShooting(propagationSettings, state);
78
79 final NumericalPropagator propagator = testSingleShooting.buildPropagator(state);
80
81 Assertions.assertEquals(1, propagator.getMultiplexer().getHandlers().size());
82 }
83
84 @Test
85 void buildFieldODETest() {
86
87 final double mu = Constants.EGM96_EARTH_MU;
88 final NewtonianAttraction attraction = new NewtonianAttraction(mu);
89 final List<ForceModel> forces = new ArrayList<>();
90 forces.add(attraction);
91 final CartesianAdjointKeplerianTerm adjointKeplerianTerm = new CartesianAdjointKeplerianTerm(mu);
92 final ShootingIntegrationSettings integrationSettings = ShootingIntegrationSettingsFactory.getClassicalRungeKuttaIntegratorSettings(1);
93 final String adjointName = "adjoint";
94 final CartesianAdjointDynamicsProvider adjointDynamicsProvider = CartesianAdjointDynamicsProviderFactory.buildUnboundedEnergyProviderNeglectingMass(adjointName,
95 adjointKeplerianTerm);
96 final ShootingPropagationSettings propagationSettings = new ShootingPropagationSettings(forces,
97 adjointDynamicsProvider, integrationSettings);
98 final double[] adjoint = new double[] {1, 2, 3, 4, 5, 6};
99 final SpacecraftState state = new SpacecraftState(TestUtils.getDefaultOrbit(AbsoluteDate.ARBITRARY_EPOCH))
100 .addAdditionalData(adjointName, adjoint);
101 final TestSingleShooting testSingleShooting = new TestSingleShooting(propagationSettings, state);
102 final FieldAbsoluteDate<Binary64> fieldDate = new FieldAbsoluteDate<>(Binary64Field.getInstance(), state.getDate());
103
104 final FieldOrdinaryDifferentialEquation<Binary64> fieldODE = testSingleShooting.buildFieldODE(fieldDate);
105
106 Assertions.assertEquals(7 + propagationSettings.getAdjointDynamicsProvider().getDimension(),
107 fieldODE.getDimension());
108 checkDerivatives(fieldODE, state, adjoint, attraction, adjointDynamicsProvider);
109 }
110
111 private static void checkDerivatives(final FieldOrdinaryDifferentialEquation<Binary64> fieldODE,
112 final SpacecraftState state, final double[] adjoint,
113 final NewtonianAttraction attraction,
114 final CartesianAdjointDynamicsProvider adjointDynamicsProvider) {
115 final Binary64[] fullState = MathArrays.buildArray(Binary64Field.getInstance(), fieldODE.getDimension());
116 fullState[0] = new Binary64(state.getPosition().getX());
117 fullState[1] = new Binary64(state.getPosition().getY());
118 fullState[2] = new Binary64(state.getPosition().getZ());
119 fullState[3] = new Binary64(state.getPVCoordinates().getVelocity().getX());
120 fullState[4] = new Binary64(state.getPVCoordinates().getVelocity().getY());
121 fullState[5] = new Binary64(state.getPVCoordinates().getVelocity().getZ());
122 fullState[6] = new Binary64(state.getMass());
123 for (int i = 0; i < adjoint.length; i++) {
124 fullState[7 + i] = new Binary64(adjoint[i]);
125 }
126 final Binary64[] derivatives = fieldODE.computeDerivatives(Binary64.ZERO, fullState);
127 Assertions.assertEquals(derivatives[0], fullState[3]);
128 Assertions.assertEquals(derivatives[1], fullState[4]);
129 Assertions.assertEquals(derivatives[2], fullState[5]);
130 final double mu = attraction.getMu(AbsoluteDate.ARBITRARY_EPOCH);
131 final Vector3D newtonianAcceleration = attraction.acceleration(state, new double[] {mu});
132 final CartesianAdjointDerivativesProvider derivativesProvider = adjointDynamicsProvider.buildAdditionalDerivativesProvider();
133 final CombinedDerivatives combinedDerivatives = derivativesProvider.combinedDerivatives(state);
134 final double[] mainDerivativesIncrement = combinedDerivatives.getMainStateDerivativesIncrements();
135 Assertions.assertEquals(derivatives[3].getReal(), newtonianAcceleration.getX() + mainDerivativesIncrement[3]);
136 Assertions.assertEquals(derivatives[4].getReal(), newtonianAcceleration.getY() + mainDerivativesIncrement[4]);
137 Assertions.assertEquals(derivatives[5].getReal(), newtonianAcceleration.getZ() + mainDerivativesIncrement[5]);
138 Assertions.assertEquals(derivatives[6].getReal(), mainDerivativesIncrement[6], 1e-20);
139 final double[] adjointDerivatives = combinedDerivatives.getAdditionalDerivatives();
140 for (int i = 0; i < adjoint.length; i++) {
141 Assertions.assertEquals(adjointDerivatives[i], derivatives[i + 7].getReal());
142 }
143 }
144
145 private static class TestSingleShooting extends AbstractFixedInitialCartesianSingleShooting {
146
147 protected TestSingleShooting(ShootingPropagationSettings propagationSettings,
148 SpacecraftState initialSpacecraftStateTemplate) {
149 super(propagationSettings, initialSpacecraftStateTemplate);
150 }
151
152 @Override
153 public int getMaximumIterationCount() {
154 return 0;
155 }
156
157 @Override
158 public ShootingBoundaryOutput computeCandidateSolution(SpacecraftState initialState, int iterationCount) {
159 return null;
160 }
161
162 @Override
163 protected double[] updateShootingVariables(double[] originalShootingVariables, FieldSpacecraftState<Gradient> fieldTerminalState) {
164 return new double[0];
165 }
166
167 @Override
168 protected double[] getScales() {
169 return new double[] {1., 1., 1., 1., 1., 1.};
170 }
171 }
172 }