1   /* Copyright 2022-2025 Romain Serra
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.control.indirect.adjoint.cost;
18  
19  import org.hipparchus.Field;
20  import org.hipparchus.complex.Complex;
21  import org.hipparchus.complex.ComplexField;
22  import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
23  import org.hipparchus.geometry.euclidean.threed.Vector3D;
24  import org.hipparchus.ode.events.Action;
25  import org.hipparchus.util.Binary64;
26  import org.hipparchus.util.Binary64Field;
27  import org.hipparchus.util.MathArrays;
28  import org.junit.jupiter.api.Assertions;
29  import org.junit.jupiter.api.Test;
30  import org.junit.jupiter.params.ParameterizedTest;
31  import org.junit.jupiter.params.provider.ValueSource;
32  import org.orekit.propagation.events.FieldEventDetector;
33  
34  import java.util.List;
35  import java.util.stream.Collectors;
36  import java.util.stream.Stream;
37  
38  class FieldBoundedCartesianEnergyTest {
39  
40      @ParameterizedTest
41      @ValueSource(booleans = {false, true})
42      void testUpdateFieldAdjointDerivatives(final boolean withMass) {
43          // GIVEN
44          final Binary64 massFlowRateFactor = withMass ? Binary64.ONE : Binary64.ZERO;
45          final FieldBoundedCartesianEnergy<Binary64> cost = new FieldBoundedCartesianEnergy<>("adjoint", massFlowRateFactor, Binary64.PI);
46          final Binary64[] adjoint = MathArrays.buildArray(Binary64Field.getInstance(), withMass ? 7 : 6);
47          final Binary64[] derivatives = adjoint.clone();
48          adjoint[3] = Binary64.ONE;
49          // WHEN
50          cost.updateFieldAdjointDerivatives(adjoint, Binary64.ONE, derivatives);
51          // THEN
52          final Binary64 zero = Binary64.ZERO;
53          for (int i = 0; i < 6; ++i) {
54              Assertions.assertEquals(zero, derivatives[i]);
55          }
56          if (withMass) {
57              Assertions.assertNotEquals(zero, derivatives[derivatives.length - 1]);
58          } else {
59              Assertions.assertEquals(zero, derivatives[derivatives.length - 1]);
60          }
61      }
62  
63      @Test
64      void getFieldEventDetectorsSizeAndActionTest() {
65          // GIVEN
66          final Complex maximumThrustMagnitude = Complex.ONE;
67          final Complex massFlowRateFactor = new Complex(2);
68          final FieldBoundedCartesianEnergy<Complex> boundedCartesianEnergy = new FieldBoundedCartesianEnergy<>("", massFlowRateFactor,
69                  maximumThrustMagnitude);
70          final Field<Complex> field = ComplexField.getInstance();
71          // WHEN
72          final Stream<FieldEventDetector<Complex>> eventDetectorStream = boundedCartesianEnergy.getFieldEventDetectors(field);
73          // THEN
74          final List<FieldEventDetector<Complex>> eventDetectors = eventDetectorStream.collect(Collectors.toList());
75          Assertions.assertEquals(2, eventDetectors.size());
76          for (final FieldEventDetector<Complex> eventDetector : eventDetectors) {
77              Assertions.assertInstanceOf(FieldCartesianEnergyConsideringMass.FieldSingularityDetector.class, eventDetector);
78              final FieldCartesianEnergyConsideringMass<Complex>.FieldSingularityDetector singularityDetector =
79                      (FieldCartesianEnergyConsideringMass<Complex>.FieldSingularityDetector) eventDetector;
80              Assertions.assertEquals(Action.RESET_DERIVATIVES, singularityDetector.getHandler().eventOccurred(null, null, true));
81          }
82      }
83  
84      @ParameterizedTest
85      @ValueSource(doubles = {0, 1})
86      void testGetFieldThrustForceNorm(final double massFlowRateFactor) {
87          // GIVEN
88          final FieldBoundedCartesianEnergy<Complex> boundedCartesianEnergy = new FieldBoundedCartesianEnergy<>("",
89                  new Complex(massFlowRateFactor), new Complex(2));
90          final ComplexField field = ComplexField.getInstance();
91          final Complex[] fieldAdjoint = MathArrays.buildArray(field, 7);
92          fieldAdjoint[3] = new Complex(1.0, 0.0);
93          fieldAdjoint[4] = new Complex(2.0, 0.0);
94          fieldAdjoint[5] = new Complex(3.0, 0.0);
95          fieldAdjoint[6] = new Complex(4.0, 0.0);
96          final Complex mass = new Complex(1, 0.);
97          // WHEN
98          final Complex fieldThrustNorm = boundedCartesianEnergy.getFieldThrustForceNorm(fieldAdjoint, mass);
99          // THEN
100         final double[] adjoint = new double[7];
101         for (int i = 0; i < adjoint.length; i++) {
102             adjoint[i] = fieldAdjoint[i].getReal();
103         }
104         final double thrustNorm = boundedCartesianEnergy.toCartesianCost().getThrustForceNorm(adjoint, mass.getReal());
105         Assertions.assertEquals(thrustNorm, fieldThrustNorm.getReal());
106     }
107 
108     @ParameterizedTest
109     @ValueSource(doubles = {1e-3, 1e-1, 1e2})
110     void testGetFieldThrustAccelerationVectorFieldFactor(final double massReal) {
111         // GIVEN
112         final Complex massRateFactor = Complex.ONE;
113         final FieldBoundedCartesianEnergy<Complex> boundedCartesianEnergy = new FieldBoundedCartesianEnergy<>("", massRateFactor,
114                 new Complex(2));
115         final ComplexField field = ComplexField.getInstance();
116         final Complex[] fieldAdjoint = MathArrays.buildArray(field, 7);
117         fieldAdjoint[3] = new Complex(1.0, 0.0);
118         fieldAdjoint[4] = new Complex(2.0, 0.0);
119         fieldAdjoint[5] = new Complex(3.0, 0.0);
120         fieldAdjoint[6] = new Complex(4.0, 0.0);
121         final Complex mass = new Complex(massReal, 0.);
122         // WHEN
123         final FieldVector3D<Complex> fieldThrustVector = boundedCartesianEnergy.getFieldThrustAccelerationVector(fieldAdjoint, mass);
124         // THEN
125         final double[] adjoint = new double[7];
126         for (int i = 0; i < adjoint.length; i++) {
127             adjoint[i] = fieldAdjoint[i].getReal();
128         }
129         final Vector3D thrustVector = boundedCartesianEnergy.toCartesianCost().getThrustAccelerationVector(adjoint, mass.getReal());
130         Assertions.assertEquals(thrustVector, fieldThrustVector.toVector3D());
131     }
132 
133     @ParameterizedTest
134     @ValueSource(doubles = {1e-4, 1e2})
135     void testFieldUpdateDerivatives(final double mass) {
136         // GIVEN
137         final Complex massRateFactor = Complex.ONE;
138         final FieldBoundedCartesianEnergy<Complex> boundedCartesianEnergy = new FieldBoundedCartesianEnergy<>("", massRateFactor,
139                 new Complex(2));
140         final ComplexField field = ComplexField.getInstance();
141         final Complex[] fieldAdjoint = MathArrays.buildArray(field, 7);
142         fieldAdjoint[3] = new Complex(1.0e-3, 0.0);
143         fieldAdjoint[4] = new Complex(2.0e-3, 0.0);
144         fieldAdjoint[6] = Complex.ONE;
145         final Complex fieldMass = new Complex(mass, 0.);
146         final Complex[] fieldDerivatives = MathArrays.buildArray(field, 7);
147         // WHEN
148         boundedCartesianEnergy.updateFieldAdjointDerivatives(fieldAdjoint, fieldMass, fieldDerivatives);
149         // THEN
150         final double[] adjoint = new double[7];
151         for (int i = 0; i < fieldAdjoint.length; i++) {
152             adjoint[i] = fieldAdjoint[i].getReal();
153         }
154         final double[] derivatives = new double[7];
155         boundedCartesianEnergy.toCartesianCost().updateAdjointDerivatives(adjoint, fieldMass.getReal(), derivatives);
156         Assertions.assertEquals(derivatives[6], fieldDerivatives[6].getReal());
157     }
158 
159     @Test
160     void testGetFieldThrustAccelerationVectorFieldVersusUnbounded() {
161         // GIVEN
162         final Complex massRateFactor = Complex.ONE;
163         final Complex maximumThrustMagnitude = new Complex(1);
164         final FieldBoundedCartesianEnergy<Complex> boundedCartesianEnergy = new FieldBoundedCartesianEnergy<>("", massRateFactor,
165                 maximumThrustMagnitude);
166         final ComplexField field = ComplexField.getInstance();
167         final Complex[] fieldAdjoint = MathArrays.buildArray(field, 7);
168         fieldAdjoint[3] = new Complex(1.0e-3, 0.0);
169         fieldAdjoint[4] = new Complex(2.0e-3, 0.0);
170         final Complex mass = new Complex(3.e-3, 0.);
171         // WHEN
172         final FieldVector3D<Complex> fieldThrustVector = boundedCartesianEnergy.getFieldThrustAccelerationVector(fieldAdjoint, mass);
173         // THEN
174         final FieldVector3D<Complex> expectedThrustVector = new FieldUnboundedCartesianEnergy<>("", massRateFactor)
175                 .getFieldThrustAccelerationVector(fieldAdjoint, mass).scalarMultiply(maximumThrustMagnitude);
176         Assertions.assertEquals(expectedThrustVector, fieldThrustVector);
177     }
178 
179     @Test
180     void testToCartesianCost() {
181         // GIVEN
182         final Complex massRateFactor = Complex.ONE;
183         final FieldBoundedCartesianEnergy<Complex> fieldCartesianEnergy = new FieldBoundedCartesianEnergy<>("",
184                 massRateFactor, new Complex(2));
185         // WHEN
186         final BoundedCartesianEnergy cartesianEnergy = fieldCartesianEnergy.toCartesianCost();
187         // THEN
188         Assertions.assertEquals(cartesianEnergy.getAdjointName(), fieldCartesianEnergy.getAdjointName());
189         Assertions.assertEquals(cartesianEnergy.getMaximumThrustMagnitude(),
190                 fieldCartesianEnergy.getMaximumThrustMagnitude().getReal());
191     }
192 }