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.orbits;
18  
19  import org.hamcrest.MatcherAssert;
20  import org.hipparchus.CalculusFieldElement;
21  import org.hipparchus.Field;
22  import org.hipparchus.analysis.UnivariateFunction;
23  import org.hipparchus.analysis.differentiation.DSFactory;
24  import org.hipparchus.analysis.differentiation.FiniteDifferencesDifferentiator;
25  import org.hipparchus.analysis.differentiation.UnivariateDifferentiableFunction;
26  import org.hipparchus.complex.Complex;
27  import org.hipparchus.complex.ComplexField;
28  import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
29  import org.hipparchus.geometry.euclidean.threed.Vector3D;
30  import org.hipparchus.linear.FieldMatrixPreservingVisitor;
31  import org.hipparchus.linear.MatrixUtils;
32  import org.hipparchus.util.Binary64;
33  import org.hipparchus.util.Binary64Field;
34  import org.hipparchus.util.FastMath;
35  import org.hipparchus.util.MathArrays;
36  import org.hipparchus.util.MathUtils;
37  import org.junit.jupiter.api.Assertions;
38  import org.junit.jupiter.api.BeforeEach;
39  import org.junit.jupiter.api.Test;
40  import org.junit.jupiter.params.ParameterizedTest;
41  import org.junit.jupiter.params.provider.EnumSource;
42  import org.orekit.Utils;
43  import org.orekit.errors.OrekitIllegalArgumentException;
44  import org.orekit.errors.OrekitMessages;
45  import org.orekit.frames.Frame;
46  import org.orekit.frames.FramesFactory;
47  import org.orekit.frames.Transform;
48  import org.orekit.time.AbsoluteDate;
49  import org.orekit.time.FieldAbsoluteDate;
50  import org.orekit.time.TimeScalesFactory;
51  import org.orekit.utils.Constants;
52  import org.orekit.utils.FieldPVCoordinates;
53  import org.orekit.utils.PVCoordinates;
54  import org.orekit.utils.TimeStampedFieldPVCoordinates;
55  
56  import java.util.function.Function;
57  
58  import static org.orekit.OrekitMatchers.relativelyCloseTo;
59  
60  
61  class FieldCircularOrbitTest {
62  
63      // Body mu
64      private double mu;
65  
66      @BeforeEach
67      public void setUp() {
68  
69          Utils.setDataRoot("regular-data");
70  
71          // Body mu
72          mu = 3.9860047e14;
73  
74      }
75  
76      @ParameterizedTest
77      @EnumSource(PositionAngleType.class)
78      void testWithCachedPositionAngleType(final PositionAngleType positionAngleType) {
79          // GIVEN
80          final Vector3D position = new Vector3D(-29536113.0, 30329259.0, -100125.0);
81          final Vector3D velocity = new Vector3D(-2194.0, -2141.0, -8.0);
82          final PVCoordinates pvCoordinates = new PVCoordinates(position, velocity);
83          final double muEarth = 3.9860047e14;
84          final CartesianOrbit cartesianOrbit = new CartesianOrbit(pvCoordinates, FramesFactory.getEME2000(), AbsoluteDate.ARBITRARY_EPOCH, muEarth);
85          final Binary64Field field = Binary64Field.getInstance();
86          final FieldCircularOrbit<Binary64> fieldOrbit = new FieldCircularOrbit<>(field, cartesianOrbit);
87          // WHEN
88          final FieldCircularOrbit<Binary64> orbit = fieldOrbit.withCachedPositionAngleType(positionAngleType);
89          // THEN
90          Assertions.assertEquals(fieldOrbit.getFrame(), orbit.getFrame());
91          Assertions.assertEquals(fieldOrbit.getDate(), orbit.getDate());
92          Assertions.assertEquals(fieldOrbit.getMu(), orbit.getMu());
93          final Vector3D relativePosition = fieldOrbit.getPosition(orbit.getFrame()).subtract(
94                  orbit.getPosition()).toVector3D();
95          Assertions.assertEquals(0., relativePosition.getNorm(), 1e-6);
96          Assertions.assertEquals(fieldOrbit.hasNonKeplerianAcceleration(),
97                  orbit.hasNonKeplerianAcceleration());
98      }
99  
100     @Test
101     void testInFrameNonKeplerian() {
102         testTemplateInFrame(Vector3D.MINUS_J);
103     }
104 
105     @Test
106     void testInFrameKeplerian() {
107         testTemplateInFrame(Vector3D.ZERO);
108     }
109 
110     private void testTemplateInFrame(final Vector3D acceleration) {
111         // GIVEN
112         final Vector3D position = new Vector3D(-29536113.0, 30329259.0, -100125.0);
113         final Vector3D velocity = new Vector3D(-2194.0, -2141.0, -8.0);
114         final PVCoordinates pvCoordinates = new PVCoordinates(position, velocity, acceleration);
115         final double muEarth = 3.9860047e14;
116         final CartesianOrbit cartesianOrbit = new CartesianOrbit(pvCoordinates, FramesFactory.getEME2000(),
117                 AbsoluteDate.ARBITRARY_EPOCH, muEarth);
118         final KeplerianOrbit keplerianOrbit = new KeplerianOrbit(cartesianOrbit);
119         final FieldCircularOrbit<Binary64> fieldOrbit = new FieldCircularOrbit<>(Binary64Field.getInstance(), keplerianOrbit);
120         // WHEN
121         final FieldCircularOrbit<Binary64> fieldOrbitWithOtherFrame = fieldOrbit.inFrame(FramesFactory.getGCRF());
122         // THEN
123         Assertions.assertNotEquals(fieldOrbit.getFrame(), fieldOrbitWithOtherFrame.getFrame());
124         Assertions.assertEquals(fieldOrbit.getDate(), fieldOrbitWithOtherFrame.getDate());
125         Assertions.assertEquals(fieldOrbit.getMu(), fieldOrbitWithOtherFrame.getMu());
126         final FieldVector3D<Binary64> relativePosition = fieldOrbit.getPosition(fieldOrbitWithOtherFrame.getFrame()).subtract(
127                 fieldOrbitWithOtherFrame.getPosition());
128         Assertions.assertEquals(0., relativePosition.getNorm().getReal(), 1e-6);
129         Assertions.assertEquals(fieldOrbit.hasNonKeplerianAcceleration(),
130                 fieldOrbitWithOtherFrame.hasNonKeplerianAcceleration());
131     }
132 
133     @Test
134     void testCircularToEquinoc() {
135         doTestCircularToEquinoctialEll(Binary64Field.getInstance());
136     }
137 
138     @Test
139     void testCircToEquinoc() {
140         doTestCircularToEquinoctialCirc(Binary64Field.getInstance());
141     }
142 
143     @Test
144     void testAnomalyCirc() {
145         doTestAnomalyCirc(Binary64Field.getInstance());
146     }
147 
148     @Test
149     void testAnomalyEll() {
150         doTestAnomalyEll(Binary64Field.getInstance());
151     }
152 
153     @Test
154     void testCircToCart() {
155         doTestCircularToCartesian(Binary64Field.getInstance());
156     }
157 
158     @Test
159     void testCircToKepl() {
160         doTestCircularToKeplerian(Binary64Field.getInstance());
161     }
162 
163     @Test
164     void testGeometryCirc() {
165         doTestGeometryCirc(Binary64Field.getInstance());
166     }
167 
168     @Test
169     void testGeometryEll() {
170         doTestGeometryEll(Binary64Field.getInstance());
171     }
172 
173     @Test
174     void testJacobianFinited() {
175         doTestJacobianFinitedifferences(Binary64Field.getInstance());
176     }
177 
178     @Test
179     void testJacoabianReference() {
180         doTestJacobianReference(Binary64Field.getInstance());
181     }
182 
183     @Test
184     void testNumericalIssue25() {
185         doTestNumericalIssue25(Binary64Field.getInstance());
186     }
187 
188     @Test
189     void testPerfectlyEquatorial() {
190         doTestPerfectlyEquatorial(Binary64Field.getInstance());
191     }
192 
193     @Test
194     void testPositionVelocityNormsCirc() {
195         doTestPositionVelocityNormsCirc(Binary64Field.getInstance());
196     }
197 
198     @Test
199     void testPositionVelocity() {
200         doTestPositionVelocityNormsEll(Binary64Field.getInstance());
201     }
202 
203     @Test
204     void testSymmetryCir() {
205         doTestSymmetryCir(Binary64Field.getInstance());
206     }
207 
208     @Test
209     void testSymmetryEll() {
210         doTestSymmetryEll(Binary64Field.getInstance());
211     }
212 
213     @Test
214     void testErrors()  {
215         Assertions.assertThrows(IllegalArgumentException.class,
216                                 () -> doTestNonInertialFrame(Binary64Field.getInstance()));
217     }
218 
219     @Test
220     void testHyperbolic1() {
221         doTestHyperbolic1(Binary64Field.getInstance());
222     }
223 
224     @Test
225     void testHyperbolic2() {
226         doTestHyperbolic2(Binary64Field.getInstance());
227     }
228 
229     @Test
230     void testToOrbitWithoutDerivatives() {
231         doTestToOrbitWithoutDerivatives(Binary64Field.getInstance());
232     }
233 
234     @Test
235     void testToOrbitWithDerivatives() {
236         doTestToOrbitWithDerivatives(Binary64Field.getInstance());
237     }
238 
239     @Test
240     void testDerivativesConversionSymmetry() {
241         doTestDerivativesConversionSymmetry(Binary64Field.getInstance());
242     }
243 
244     @Test
245     void testToString() {
246         doTestToString(Binary64Field.getInstance());
247     }
248 
249     @Test
250     void testNonKeplerianDerivatives() {
251         doTestNonKeplerianDerivatives(Binary64Field.getInstance());
252     }
253 
254     @Test
255     void testPositionAngleDerivatives() {
256         doTestPositionAngleDerivatives(Binary64Field.getInstance());
257     }
258 
259     @Test
260     void testEquatorialRetrograde() {
261         doTestEquatorialRetrograde(Binary64Field.getInstance());
262     }
263 
264     @Test
265     void testCopyNonKeplerianAcceleration() {
266         doTestCopyNonKeplerianAcceleration(Binary64Field.getInstance());
267     }
268 
269     @Test
270     void testNormalize() {
271         doTestNormalize(Binary64Field.getInstance());
272     }
273 
274     @Test
275     void testWithKeplerianRates() {
276         // GIVEN
277         final ComplexField field = ComplexField.getInstance();
278         final CircularOrbit expectedOrbit = createOrbitTestFromCircularOrbit(true);
279         final FieldCircularOrbit<Complex> fieldOrbit = new FieldCircularOrbit<>(field, expectedOrbit);
280         // WHEN
281         final FieldCircularOrbit<Complex> actualFieldOrbit = fieldOrbit.withKeplerianRates();
282         // THEN
283         Assertions.assertFalse(actualFieldOrbit.hasNonKeplerianRates());
284         Assertions.assertEquals(fieldOrbit.getMu(), actualFieldOrbit.getMu());
285         Assertions.assertEquals(fieldOrbit.getDate(), actualFieldOrbit.getDate());
286         Assertions.assertEquals(fieldOrbit.getFrame(), actualFieldOrbit.getFrame());
287         Assertions.assertEquals(fieldOrbit.getPosition(), actualFieldOrbit.getPosition());
288     }
289 
290     @Test
291     void testFromCircularOrbitWithDerivatives() {
292         // GIVEN
293         final ComplexField field = ComplexField.getInstance();
294         final CircularOrbit expectedOrbit = createOrbitTestFromCircularOrbit(true);
295         // WHEN
296         final FieldCircularOrbit<Complex> fieldOrbit = new FieldCircularOrbit<>(field, expectedOrbit);
297         // THEN
298         compareFieldOrbitToOrbit(expectedOrbit, fieldOrbit);
299     }
300 
301     @Test
302     void testFromCircularOrbitWithoutDerivatives() {
303         // GIVEN
304         final ComplexField field = ComplexField.getInstance();
305         final CircularOrbit expectedOrbit = createOrbitTestFromCircularOrbit(false);
306         // WHEN
307         final FieldCircularOrbit<Complex> fieldOrbit = new FieldCircularOrbit<>(field, expectedOrbit);
308         // THEN
309         compareFieldOrbitToOrbit(expectedOrbit, fieldOrbit);
310     }
311 
312     private CircularOrbit createOrbitTestFromCircularOrbit(final boolean withDerivatives) {
313         final PositionAngleType positionAngleType = PositionAngleType.TRUE;
314         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
315         final Frame frame = FramesFactory.getEME2000();
316         final double a = 10000.e3;
317         final double ex = 1e-4;
318         final double ey = 1e-3;
319         final double i = 1.;
320         final double raan = -2.;
321         final double alphaV = 0.5;
322         if (withDerivatives) {
323             final double derivative = 0.;
324             return new CircularOrbit(a, ex, ey, i, raan, alphaV, derivative, derivative, derivative, derivative,
325                     derivative, derivative, positionAngleType, frame, date, mu);
326         } else {
327             return new CircularOrbit(a, ex, ey, i, raan, alphaV, positionAngleType, frame, date, mu);
328         }
329     }
330 
331     private <T extends CalculusFieldElement<T>> void compareFieldOrbitToOrbit(final CircularOrbit orbit,
332                                                                               final FieldCircularOrbit<T> fieldOrbit) {
333         Assertions.assertEquals(orbit.getFrame(), fieldOrbit.getFrame());
334         Assertions.assertEquals(orbit.getMu(), fieldOrbit.getMu().getReal());
335         Assertions.assertEquals(orbit.getDate(), fieldOrbit.getDate().toAbsoluteDate());
336         Assertions.assertEquals(orbit.getA(), fieldOrbit.getA().getReal());
337         Assertions.assertEquals(orbit.getCircularEx(), fieldOrbit.getCircularEx().getReal());
338         Assertions.assertEquals(orbit.getCircularEy(), fieldOrbit.getCircularEy().getReal());
339         Assertions.assertEquals(orbit.getRightAscensionOfAscendingNode(),
340                 fieldOrbit.getRightAscensionOfAscendingNode().getReal());
341         Assertions.assertEquals(orbit.getI(), fieldOrbit.getI().getReal());
342         Assertions.assertEquals(orbit.getAlphaV(), fieldOrbit.getAlphaV().getReal());
343         Assertions.assertEquals(orbit.hasNonKeplerianAcceleration(), fieldOrbit.hasNonKeplerianAcceleration());
344         Assertions.assertEquals(orbit.getADot(), fieldOrbit.getADot().getReal());
345         Assertions.assertEquals(orbit.getCircularExDot(), fieldOrbit.getCircularExDot().getReal());
346         Assertions.assertEquals(orbit.getCircularEyDot(), fieldOrbit.getCircularEyDot().getReal());
347         Assertions.assertEquals(orbit.getRightAscensionOfAscendingNodeDot(),
348                 fieldOrbit.getRightAscensionOfAscendingNodeDot().getReal());
349         Assertions.assertEquals(orbit.getIDot(), fieldOrbit.getIDot().getReal());
350         Assertions.assertEquals(orbit.getAlphaVDot(), fieldOrbit.getAlphaVDot().getReal());
351     }
352 
353     @Test
354     void testCoverageCachedPositionAngleTypeWithRates() {
355         // GIVEN
356         final double semiMajorAxis = 1e4;
357         final double eccentricity = 0.;
358         final double expectedAlpha = 0.;
359         final double expectedAlphaDot = 0.;
360         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
361         final Binary64Field field = Binary64Field.getInstance();
362         final Binary64 zero = field.getZero();
363         // WHEN & THEN
364         for (final PositionAngleType inputPositionAngleType: PositionAngleType.values()) {
365             for (final PositionAngleType cachedPositionAngleType: PositionAngleType.values()) {
366                 final FieldCircularOrbit<Binary64> fieldOrbit = new FieldCircularOrbit<>(
367                         zero.newInstance(semiMajorAxis), zero.newInstance(eccentricity), zero, zero, zero,
368                         zero.newInstance(expectedAlpha), zero, zero, zero, zero, zero,
369                         zero.newInstance(expectedAlphaDot), inputPositionAngleType, cachedPositionAngleType,
370                         FramesFactory.getGCRF(), new FieldAbsoluteDate<>(field, date), zero.newInstance(mu));
371                 Assertions.assertEquals(cachedPositionAngleType, fieldOrbit.getCachedPositionAngleType());
372                 Assertions.assertEquals(expectedAlpha, fieldOrbit.getAlphaV().getReal());
373                 Assertions.assertEquals(expectedAlpha, fieldOrbit.getAlphaM().getReal());
374                 Assertions.assertEquals(expectedAlpha, fieldOrbit.getAlphaE().getReal());
375                 Assertions.assertEquals(expectedAlphaDot, fieldOrbit.getAlphaVDot().getReal());
376                 Assertions.assertEquals(expectedAlphaDot, fieldOrbit.getAlphaMDot().getReal());
377                 Assertions.assertEquals(expectedAlphaDot, fieldOrbit.getAlphaEDot().getReal());
378             }
379         }
380     }
381 
382     @Test
383     void testGetAlphaVersusDouble() {
384         // GIVEN
385         final double semiMajorAxis = 1e7;
386         final double eccentricity = 1e-2;
387         final double expectedAlpha = 4.;
388         final AbsoluteDate date = AbsoluteDate.ARBITRARY_EPOCH;
389         final Binary64Field field = Binary64Field.getInstance();
390         final Binary64 zero = field.getZero();
391         // WHEN & THEN
392         for (final PositionAngleType inputPositionAngleType: PositionAngleType.values()) {
393             for (final PositionAngleType cachedPositionAngleType: PositionAngleType.values()) {
394                 final FieldCircularOrbit<Binary64> fieldOrbit = new FieldCircularOrbit<>(
395                         zero.newInstance(semiMajorAxis), zero.newInstance(eccentricity), zero, zero, zero,
396                         zero.newInstance(expectedAlpha), inputPositionAngleType, cachedPositionAngleType,
397                         FramesFactory.getGCRF(), new FieldAbsoluteDate<>(field, date), zero.newInstance(mu));
398                 final CircularOrbit circularOrbit = fieldOrbit.toOrbit();
399                 Assertions.assertEquals(circularOrbit.getAlphaV(), fieldOrbit.getAlphaV().getReal());
400                 Assertions.assertEquals(circularOrbit.getAlphaM(), fieldOrbit.getAlphaM().getReal());
401                 Assertions.assertEquals(circularOrbit.getAlphaE(), fieldOrbit.getAlphaE().getReal());
402             }
403         }
404     }
405 
406     private <T extends CalculusFieldElement<T>> void doTestCircularToEquinoctialEll(Field<T> field) {
407 
408         T zero =  field.getZero();
409         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
410 
411         T ix = zero.add(1.200e-04);
412         T iy = zero.add(-1.16e-04);
413         T i  = ix.multiply(ix).add(iy.multiply(iy)).divide(4).sqrt().asin().multiply(2);
414 
415 
416         T raan = iy.atan2(ix);
417 
418         // elliptic orbit
419         FieldCircularOrbit<T> circ =
420             new FieldCircularOrbit<>(zero.add(42166712.0), zero.add(0.5), zero.add(-0.5), i, raan,
421                                      zero.add(5.300).subtract(raan), PositionAngleType.MEAN,
422                                      FramesFactory.getEME2000(), date, zero.add(mu));
423         FieldVector3D<T> pos = circ.getPosition();
424         FieldVector3D<T> vit = circ.getPVCoordinates().getVelocity();
425 
426         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(pos, vit);
427 
428 
429 
430         FieldEquinoctialOrbit<T> param = new FieldEquinoctialOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
431 
432         Assertions.assertEquals(param.getA().getReal(),  circ.getA().getReal(), Utils.epsilonTest * circ.getA().getReal());
433         Assertions.assertEquals(param.getEquinoctialEx().getReal(), circ.getEquinoctialEx().getReal(), Utils.epsilonE * FastMath.abs(circ.getE().getReal()));
434         Assertions.assertEquals(param.getEquinoctialEy().getReal(), circ.getEquinoctialEy().getReal(), Utils.epsilonE * FastMath.abs(circ.getE().getReal()));
435         Assertions.assertEquals(param.getHx().getReal(), circ.getHx().getReal(), Utils.epsilonAngle * FastMath.abs(circ.getI().getReal()));
436         Assertions.assertEquals(param.getHy().getReal(), circ.getHy().getReal(), Utils.epsilonAngle * FastMath.abs(circ.getI().getReal()));
437 
438 
439         Assertions.assertEquals(MathUtils.normalizeAngle(param.getLv().getReal(), circ.getLv().getReal()),
440                                 circ.getLv().getReal(),
441                                 Utils.epsilonAngle * FastMath.abs(circ.getLv().getReal()));
442 
443         Assertions.assertFalse(circ.hasNonKeplerianAcceleration());
444     }
445 
446     private <T extends CalculusFieldElement<T>> void doTestToOrbitWithoutDerivatives(Field<T> field) {
447         T zero =  field.getZero();
448         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
449 
450         FieldVector3D<T> position = new FieldVector3D<>(zero.add(7.0e6), zero.add(1.0e6), zero.add(4.0e6));
451         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-500.0), zero.add(8000.0), zero.add(1000.0));
452         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity);
453         FieldCircularOrbit<T>  fieldOrbit = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
454         CircularOrbit orbit = fieldOrbit.toOrbit();
455         Assertions.assertFalse(orbit.hasNonKeplerianAcceleration());
456         MatcherAssert.assertThat(orbit.getA(),                             relativelyCloseTo(fieldOrbit.getA().getReal(),                             0));
457         MatcherAssert.assertThat(orbit.getCircularEx(),                    relativelyCloseTo(fieldOrbit.getCircularEx().getReal(),                    0));
458         MatcherAssert.assertThat(orbit.getCircularEy(),                    relativelyCloseTo(fieldOrbit.getCircularEy().getReal(),                    0));
459         MatcherAssert.assertThat(orbit.getI(),                             relativelyCloseTo(fieldOrbit.getI().getReal(),                             0));
460         MatcherAssert.assertThat(orbit.getRightAscensionOfAscendingNode(), relativelyCloseTo(fieldOrbit.getRightAscensionOfAscendingNode().getReal(), 0));
461         MatcherAssert.assertThat(orbit.getAlphaV(),                        relativelyCloseTo(fieldOrbit.getAlphaV().getReal(),                        0));
462     }
463 
464     private <T extends CalculusFieldElement<T>> void doTestToOrbitWithDerivatives(Field<T> field) {
465         T zero =  field.getZero();
466         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
467 
468         FieldVector3D<T> position = new FieldVector3D<>(zero.add(7.0e6), zero.add(1.0e6), zero.add(4.0e6));
469         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-500.0), zero.add(8000.0), zero.add(1000.0));
470         T r2 = position.getNormSq();
471         T r = r2.sqrt();
472         final FieldVector3D<T> keplerianAcceleration = new FieldVector3D<>(r.multiply(r2).reciprocal().multiply(zero.add(mu).negate()),
473                 position);
474         final FieldVector3D<T> nonKeplerianAcceleration = keplerianAcceleration.scalarMultiply(1.1);
475         final FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity, nonKeplerianAcceleration);
476         FieldCircularOrbit<T>  fieldOrbit = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
477         CircularOrbit orbit = fieldOrbit.toOrbit();
478         Assertions.assertTrue(orbit.hasNonKeplerianAcceleration());
479         MatcherAssert.assertThat(orbit.getA(),                                 relativelyCloseTo(fieldOrbit.getA().getReal(),                                 0));
480         MatcherAssert.assertThat(orbit.getCircularEx(),                        relativelyCloseTo(fieldOrbit.getCircularEx().getReal(),                        0));
481         MatcherAssert.assertThat(orbit.getCircularEy(),                        relativelyCloseTo(fieldOrbit.getCircularEy().getReal(),                        0));
482         MatcherAssert.assertThat(orbit.getI(),                                 relativelyCloseTo(fieldOrbit.getI().getReal(),                                 0));
483         MatcherAssert.assertThat(orbit.getRightAscensionOfAscendingNode(),     relativelyCloseTo(fieldOrbit.getRightAscensionOfAscendingNode().getReal(),     0));
484         MatcherAssert.assertThat(orbit.getAlphaV(),                            relativelyCloseTo(fieldOrbit.getAlphaV().getReal(),                            0));
485         MatcherAssert.assertThat(orbit.getADot(),                              relativelyCloseTo(fieldOrbit.getADot().getReal(),                              0));
486         MatcherAssert.assertThat(orbit.getEquinoctialExDot(),                  relativelyCloseTo(fieldOrbit.getEquinoctialExDot().getReal(),                  0));
487         MatcherAssert.assertThat(orbit.getEquinoctialEyDot(),                  relativelyCloseTo(fieldOrbit.getEquinoctialEyDot().getReal(),                  0));
488         MatcherAssert.assertThat(orbit.getHxDot(),                             relativelyCloseTo(fieldOrbit.getHxDot().getReal(),                             0));
489         MatcherAssert.assertThat(orbit.getHyDot(),                             relativelyCloseTo(fieldOrbit.getHyDot().getReal(),                             0));
490         MatcherAssert.assertThat(orbit.getLvDot(),                             relativelyCloseTo(fieldOrbit.getLvDot().getReal(),                             0));
491         MatcherAssert.assertThat(orbit.getLEDot(),                             relativelyCloseTo(fieldOrbit.getLEDot().getReal(),                             0));
492         MatcherAssert.assertThat(orbit.getLMDot(),                             relativelyCloseTo(fieldOrbit.getLMDot().getReal(),                             0));
493         MatcherAssert.assertThat(orbit.getEDot(),                              relativelyCloseTo(fieldOrbit.getEDot().getReal(),                              0));
494         MatcherAssert.assertThat(orbit.getIDot(),                              relativelyCloseTo(fieldOrbit.getIDot().getReal(),                              0));
495         MatcherAssert.assertThat(orbit.getCircularExDot(),                     relativelyCloseTo(fieldOrbit.getCircularExDot().getReal(),                     0));
496         MatcherAssert.assertThat(orbit.getCircularEyDot(),                     relativelyCloseTo(fieldOrbit.getCircularEyDot().getReal(),                     0));
497         MatcherAssert.assertThat(orbit.getRightAscensionOfAscendingNodeDot(),  relativelyCloseTo(fieldOrbit.getRightAscensionOfAscendingNodeDot().getReal(),  0));
498         MatcherAssert.assertThat(orbit.getAlphaVDot(),                         relativelyCloseTo(fieldOrbit.getAlphaVDot().getReal(),                         0));
499         MatcherAssert.assertThat(orbit.getAlphaEDot(),                         relativelyCloseTo(fieldOrbit.getAlphaEDot().getReal(),                         0));
500         MatcherAssert.assertThat(orbit.getAlphaMDot(),                         relativelyCloseTo(fieldOrbit.getAlphaMDot().getReal(),                         0));
501         MatcherAssert.assertThat(orbit.getAlphaDot(PositionAngleType.TRUE),        relativelyCloseTo(fieldOrbit.getAlphaDot(PositionAngleType.TRUE).getReal(),        0));
502         MatcherAssert.assertThat(orbit.getAlphaDot(PositionAngleType.ECCENTRIC),   relativelyCloseTo(fieldOrbit.getAlphaDot(PositionAngleType.ECCENTRIC).getReal(),   0));
503         MatcherAssert.assertThat(orbit.getAlphaDot(PositionAngleType.MEAN),        relativelyCloseTo(fieldOrbit.getAlphaDot(PositionAngleType.MEAN).getReal(),        0));
504     }
505 
506     private <T extends CalculusFieldElement<T>> void doTestCircularToEquinoctialCirc(Field<T> field) {
507 
508         T zero =  field.getZero();
509         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
510 
511         T ix = zero.add(1.200e-04);
512         T iy = zero.add(-1.16e-04);
513         T i  = ix.multiply(ix).add(iy.multiply(iy)).divide(4).sqrt().asin().multiply(2);
514         T raan = iy.atan2(ix);
515 
516         // circular orbit
517         FieldEquinoctialOrbit<T> circCir =
518             new FieldEquinoctialOrbit<>(zero.add(42166712.0), zero.add(0.1e-10), zero.add(-0.1e-10), i, raan,
519                                         raan.negate().add(5.300), PositionAngleType.MEAN,
520                                         FramesFactory.getEME2000(), date, zero.add(mu));
521         FieldVector3D<T> posCir = circCir.getPosition();
522         FieldVector3D<T> vitCir = circCir.getPVCoordinates().getVelocity();
523 
524         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>( posCir, vitCir);
525 
526         FieldEquinoctialOrbit<T> paramCir = new FieldEquinoctialOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
527         Assertions.assertEquals(paramCir.getA().getReal(), circCir.getA().getReal(), Utils.epsilonTest * circCir.getA().getReal());
528         Assertions.assertEquals(paramCir.getEquinoctialEx().getReal(), circCir.getEquinoctialEx().getReal(), Utils.epsilonEcir * FastMath.abs(circCir.getE().getReal()));
529         Assertions.assertEquals(paramCir.getEquinoctialEy().getReal(), circCir.getEquinoctialEy().getReal(), Utils.epsilonEcir * FastMath.abs(circCir.getE().getReal()));
530         Assertions.assertEquals(paramCir.getHx().getReal(), circCir.getHx().getReal(), Utils.epsilonAngle * FastMath.abs(circCir.getI().getReal()));
531         Assertions.assertEquals(paramCir.getHy().getReal(), circCir.getHy().getReal(), Utils.epsilonAngle * FastMath.abs(circCir.getI().getReal()));
532         Assertions.assertEquals(MathUtils.normalizeAngle(paramCir.getLv().getReal(), circCir.getLv().getReal()), circCir.getLv().getReal(), Utils.epsilonAngle * FastMath.abs(circCir.getLv().getReal()));
533 
534     }
535 
536     private <T extends CalculusFieldElement<T>> void doTestCircularToCartesian(Field<T> field) {
537         T zero =  field.getZero();
538         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
539 
540         T ix = zero.add(1.200e-04);
541         T iy = zero.add(-1.16e-04);
542         T i  = ix.multiply(ix).add(iy.multiply(iy)).divide(4).sqrt().asin().multiply(2);
543         T raan = iy.atan2(ix);
544         T cosRaan = raan.cos();
545         T sinRaan = raan.sin();
546         T exTilde = zero.add(-7.900e-6);
547         T eyTilde = zero.add(1.100e-4);
548         T ex = exTilde.multiply(cosRaan).add(eyTilde.multiply(sinRaan));
549         T ey = eyTilde.multiply(cosRaan).subtract(exTilde.multiply(sinRaan));
550 
551         FieldCircularOrbit<T> circ=
552             new FieldCircularOrbit<>(zero.add(42166712.0), ex, ey, i, raan,
553                                      raan.negate().add(5.300), PositionAngleType.MEAN,
554                                      FramesFactory.getEME2000(), date, zero.add(mu));
555         FieldVector3D<T> pos = circ.getPosition();
556         FieldVector3D<T> vel = circ.getPVCoordinates().getVelocity();
557 
558         // check 1/a = 2/r  - V2/mu
559         T r = pos.getNorm();
560         T v = vel.getNorm();
561         Assertions.assertEquals(2 / r.getReal() - v.getReal() * v.getReal() / mu, 1 / circ.getA().getReal(), 1.0e-7);
562 
563         Assertions.assertEquals( 0.233745668678733e+08, pos.getX().getReal(), Utils.epsilonTest * r.getReal());
564         Assertions.assertEquals(-0.350998914352669e+08, pos.getY().getReal(), Utils.epsilonTest * r.getReal());
565         Assertions.assertEquals(-0.150053723123334e+04, pos.getZ().getReal(), Utils.epsilonTest * r.getReal());
566 
567         Assertions.assertEquals(2558.7096558809967, vel.getX().getReal(), Utils.epsilonTest * v.getReal());
568         Assertions.assertEquals(1704.1586039092576, vel.getY().getReal(), Utils.epsilonTest * v.getReal());
569         Assertions.assertEquals(   0.5013093577879, vel.getZ().getReal(), Utils.epsilonTest * v.getReal());
570 
571     }
572 
573     private <T extends CalculusFieldElement<T>> void doTestCircularToKeplerian(Field<T> field) {
574         T zero =  field.getZero();
575         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
576 
577         T ix   = zero.add(1.20e-4);
578         T iy   = zero.add(-1.16e-4);
579         T i    = ix.multiply(ix).add(iy.multiply(iy)).divide(4).sqrt().asin().multiply(2);
580         T raan = iy.atan2(ix);
581         T cosRaan = raan.cos();
582         T sinRaan = raan.sin();
583         T exTilde = zero.add(-7.900e-6);
584         T eyTilde = zero.add(1.100e-4);
585         T ex = exTilde.multiply(cosRaan).add(eyTilde.multiply(sinRaan));
586         T ey = eyTilde.multiply(cosRaan).subtract(exTilde.multiply(sinRaan));
587 
588         FieldCircularOrbit<T> circ=
589             new FieldCircularOrbit<>(zero.add(42166712.0), ex, ey, i, raan,
590                                      raan.negate().add(5.300), PositionAngleType.MEAN,
591                                      FramesFactory.getEME2000(), date, zero.add(mu));
592         FieldKeplerianOrbit<T> kep = new FieldKeplerianOrbit<>(circ);
593 
594         Assertions.assertEquals(42166712.000, circ.getA().getReal(), Utils.epsilonTest * kep.getA().getReal());
595         Assertions.assertEquals(0.110283316961361e-03, kep.getE().getReal(), Utils.epsilonE * FastMath.abs(kep.getE().getReal()));
596         Assertions.assertEquals(0.166901168553917e-03, kep.getI().getReal(),
597                      Utils.epsilonAngle * FastMath.abs(kep.getI().getReal()));
598         Assertions.assertEquals(MathUtils.normalizeAngle(-3.87224326008837, kep.getPerigeeArgument().getReal()),
599                      kep.getPerigeeArgument().getReal(),
600                      Utils.epsilonTest * 6 * FastMath.abs(kep.getPerigeeArgument().getReal())); //numerical propagation we changed to 6 time the precision used
601         Assertions.assertEquals(MathUtils.normalizeAngle(5.51473467358854, kep.getRightAscensionOfAscendingNode().getReal()),
602                      kep.getRightAscensionOfAscendingNode().getReal(),
603                      Utils.epsilonTest * FastMath.abs(kep.getRightAscensionOfAscendingNode().getReal()));
604 
605         Assertions.assertEquals(MathUtils.normalizeAngle(zero.add(3.65750858649982), kep.getMeanAnomaly()).getReal(),
606                      kep.getMeanAnomaly().getReal(),
607                      Utils.epsilonTest * 5 * FastMath.abs(kep.getMeanAnomaly().getReal())); //numerical propagation we changed to 6 time the precision used
608 
609     }
610 
611     private <T extends CalculusFieldElement<T>> void doTestHyperbolic1(Field<T> field) {
612         T zero =  field.getZero();
613         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
614         try {
615             new FieldCircularOrbit<>(zero.add(-42166712.0), zero.add(1.9), zero.add(0.5), zero.add(0.01), zero.add(-0.02), zero.add( 5.300),
616                                      PositionAngleType.MEAN,  FramesFactory.getEME2000(), date, zero.add(mu));
617             Assertions.fail("an exception should have been thrown");
618         } catch (OrekitIllegalArgumentException oe) {
619             Assertions.assertEquals(OrekitMessages.HYPERBOLIC_ORBIT_NOT_HANDLED_AS, oe.getSpecifier());
620         }
621     }
622 
623     private <T extends CalculusFieldElement<T>>  void doTestHyperbolic2(Field<T> field) {
624         T zero =  field.getZero();
625        FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
626         FieldOrbit<T> orbit = new FieldKeplerianOrbit<>(zero.add(-42166712.0), zero.add(1.9), zero.add(0.5), zero.add(0.01), zero.add(-0.02), zero.add( 5.300),
627                                                         PositionAngleType.MEAN,  FramesFactory.getEME2000(), date, zero.add(mu));
628         try {
629             new FieldCircularOrbit<>(orbit.getPVCoordinates(), orbit.getFrame(), orbit.getMu());
630             Assertions.fail("an exception should have been thrown");
631         } catch (OrekitIllegalArgumentException oe) {
632             Assertions.assertEquals(OrekitMessages.HYPERBOLIC_ORBIT_NOT_HANDLED_AS, oe.getSpecifier());
633         }
634     }
635 
636     private <T extends CalculusFieldElement<T>> void doTestAnomalyEll(Field<T> field) {
637         T zero =  field.getZero();
638         T one =  field.getOne();
639         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
640         // elliptic orbit
641         FieldVector3D<T> position = new FieldVector3D<>(zero.add(7.0e6), zero.add(1.0e6), zero.add(4.0e6));
642         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-500.0), zero.add(8000.0), zero.add(1000.0));
643 
644         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>( position, velocity);
645 
646         FieldCircularOrbit<T>  p   = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
647         FieldKeplerianOrbit<T> kep = new FieldKeplerianOrbit<>(p);
648 
649         T e       = p.getE();
650         T eRatio  = one.subtract(e).divide(one.add(e)).sqrt();
651         T raan    = kep.getRightAscensionOfAscendingNode();
652         T paPraan = kep.getPerigeeArgument().add(raan);
653 
654         T lv = zero.add(1.1);
655         // formulations for elliptic case
656         T lE = lv.subtract(paPraan).divide(2).tan().multiply(eRatio).atan().multiply(2).add(paPraan);
657         T lM = lE.subtract(e.multiply(lE.subtract(paPraan).sin()));
658 
659         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
660                                      p.getRightAscensionOfAscendingNode(),
661                                      p.getAlphaV(), lv.subtract(raan), PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
662         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
663         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
664         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
665         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
666                                      p.getRightAscensionOfAscendingNode(),
667                                      p.getAlphaV(), zero, PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
668 
669 
670         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
671                                      p.getRightAscensionOfAscendingNode(),
672                                      p.getAlphaV(), lE.subtract(raan), PositionAngleType.ECCENTRIC, p.getFrame(), date, zero.add(mu));
673         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
674         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
675         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
676         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
677                                      p.getRightAscensionOfAscendingNode(),
678                                      p.getAlphaV(), zero, PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
679 
680         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
681                                      p.getRightAscensionOfAscendingNode(),
682                                      p.getAlphaV(), lM.subtract(raan), PositionAngleType.MEAN, p.getFrame(), date, zero.add(mu));
683         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
684         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
685         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
686 
687     }
688 
689     private <T extends CalculusFieldElement<T>> void doTestAnomalyCirc(Field<T> field) {
690         T zero =  field.getZero();
691         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
692 
693         FieldVector3D<T> position = new FieldVector3D<>(zero.add(7.0e6), zero.add(1.0e6), zero.add(4.0e6));
694         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-500.0), zero.add(8000.0), zero.add(1000.0));
695         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>( position, velocity);
696         FieldCircularOrbit<T>  p = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
697         T raan = p.getRightAscensionOfAscendingNode();
698 
699         // circular orbit
700         p = new FieldCircularOrbit<>(p.getA() , zero, zero, p.getRightAscensionOfAscendingNode(),
701                                      p.getAlphaV(), p.getAlphaV(), PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
702 
703         T lv = zero.add(1.1);
704         T lE = lv;
705         T lM = lE;
706 
707         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
708                                      p.getRightAscensionOfAscendingNode(),
709                                      p.getAlphaV(), lv.subtract(raan), PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
710         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
711         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
712         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
713         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
714                                      p.getRightAscensionOfAscendingNode(),
715                                      p.getAlphaV(), zero, PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
716 
717         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
718                                      p.getRightAscensionOfAscendingNode(),
719                                      p.getAlphaV(), lE.subtract(raan), PositionAngleType.ECCENTRIC, p.getFrame(), date, zero.add(mu));
720 
721         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
722         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
723         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
724         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
725                                      p.getRightAscensionOfAscendingNode(),
726                                      p.getAlphaV(), zero, PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
727 
728         p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(),
729                                      p.getRightAscensionOfAscendingNode(),
730                                      p.getAlphaV(), lM.subtract(raan), PositionAngleType.MEAN, p.getFrame(), date, zero.add(mu));
731         Assertions.assertEquals(p.getAlphaV().getReal() + raan.getReal(), lv.getReal(), Utils.epsilonAngle * FastMath.abs(lv.getReal()));
732         Assertions.assertEquals(p.getAlphaE().getReal() + raan.getReal(), lE.getReal(), Utils.epsilonAngle * FastMath.abs(lE.getReal()));
733         Assertions.assertEquals(p.getAlphaM().getReal() + raan.getReal(), lM.getReal(), Utils.epsilonAngle * FastMath.abs(lM.getReal()));
734 
735     }
736 
737     private <T extends CalculusFieldElement<T>> void doTestPositionVelocityNormsEll(Field<T> field) {
738         T zero =  field.getZero();
739         T one =  field.getOne();
740         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
741 
742         // elliptic and non equatorial (i retrograde) orbit
743         T hx =  zero.add(1.2);
744         T hy =  zero.add(2.1);
745         T i  = hx.multiply(hx).add(hy.multiply(hy)).sqrt().atan().multiply(2);
746         T raan = hy.atan2(hx);
747         FieldCircularOrbit<T> p =
748             new FieldCircularOrbit<>(zero.add(42166712.0), zero.add(0.5), zero.add(-0.5), i, raan,
749                                      raan.negate().add(0.67), PositionAngleType.TRUE,
750                                      FramesFactory.getEME2000(), date, zero.add(mu));
751 
752         T ex = p.getEquinoctialEx();
753         T ey = p.getEquinoctialEy();
754         T lv = p.getLv();
755         T ksi     = ex.multiply(lv.cos()).add(1).add(ey.multiply(lv.sin()));
756         T nu      = ex.multiply(lv.sin()).subtract(ey.multiply(lv.cos()));
757         T epsilon = one.subtract(ex.multiply(ex)).subtract(ey.multiply(ey)).sqrt();
758 
759         T a  = p.getA();
760         T na = a.reciprocal().multiply(mu).sqrt();
761 
762         Assertions.assertEquals(a.getReal() * epsilon.getReal() * epsilon.getReal() / ksi.getReal(),
763                      p.getPosition().getNorm().getReal(),
764                      Utils.epsilonTest * FastMath.abs(p.getPosition().getNorm().getReal()));
765         Assertions.assertEquals(na.getReal() * FastMath.sqrt(ksi.getReal() * ksi.getReal() + nu.getReal() * nu.getReal()) / epsilon.getReal(),
766                      p.getPVCoordinates().getVelocity().getNorm().getReal(),
767                      Utils.epsilonTest * FastMath.abs(p.getPVCoordinates().getVelocity().getNorm().getReal()));
768 
769     }
770 
771     private <T extends CalculusFieldElement<T>> void doTestNumericalIssue25(Field<T> field) {
772 
773         T zero =  field.getZero();
774         FieldVector3D<T> position = new FieldVector3D<>(zero.add(3782116.14107698), zero.add(416663.11924914), zero.add(5875541.62103057));
775         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-6349.7848910501), zero.add(288.4061811651), zero.add(4066.9366759691));
776         FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(new FieldPVCoordinates<>(position, velocity),
777                                                                FramesFactory.getEME2000(),
778                                                                new FieldAbsoluteDate<>(field, "2004-01-01T23:00:00.000",
779                                                                                        TimeScalesFactory.getUTC()),
780                                                                zero.add(3.986004415E14));
781         Assertions.assertEquals(0.0, orbit.getE().getReal(), 2.0e-14);
782     }
783 
784     private <T extends CalculusFieldElement<T>> void doTestPerfectlyEquatorial(Field<T> field) {
785         T zero =  field.getZero();
786         FieldVector3D<T> position = new FieldVector3D<>(zero.add(-7293947.695148368), zero.add( 5122184.668436634), zero.add(0.0));
787         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-3890.4029433398), zero.add( -5369.811285264604), zero.add(0.0));
788         FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(new FieldPVCoordinates<>(position, velocity),
789                                                                FramesFactory.getEME2000(),
790                                                                new FieldAbsoluteDate<>(field, "2004-01-01T23:00:00.000",
791                                                                                        TimeScalesFactory.getUTC()),
792                                                                zero.add(3.986004415E14));
793         Assertions.assertEquals(0.0, orbit.getI().getReal(), 2.0e-14);
794         Assertions.assertEquals(0.0, orbit.getRightAscensionOfAscendingNode().getReal(), 2.0e-14);
795     }
796 
797     private <T extends CalculusFieldElement<T>> void doTestPositionVelocityNormsCirc(Field<T> field) {
798         T zero =  field.getZero();
799         T one =  field.getOne();
800         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
801 
802         // elliptic and non equatorial (i retrograde) orbit
803         T hx =  zero.add(0.1e-8);
804         T hy =  zero.add(0.1e-8);
805         T i  = hx.multiply(hx).add(hy.multiply(hy)).sqrt().atan().multiply(2);
806         T raan = hy.atan2(hx);
807         FieldCircularOrbit<T> pCirEqua =
808             new FieldCircularOrbit<>(zero.add(42166712.0), zero.add(0.1e-8), zero.add(0.1e-8), i, raan,
809                                      raan.negate().add(0.67), PositionAngleType.TRUE,
810                                      FramesFactory.getEME2000(), date, zero.add(mu));
811 
812         T ex = pCirEqua.getEquinoctialEx();
813         T ey = pCirEqua.getEquinoctialEy();
814         T lv = pCirEqua.getLv();
815         T ksi     = ex.multiply(lv.cos()).add(1).add(ey.multiply(lv.sin()));
816         T nu      = ex.multiply(lv.sin()).subtract(ey.multiply(lv.cos()));
817         T epsilon = one.subtract(ex.multiply(ex)).subtract(ey.multiply(ey)).sqrt();
818 
819         T a  = pCirEqua.getA();
820         T na = a.reciprocal().multiply(mu).sqrt();
821 
822         Assertions.assertEquals(a.getReal() * epsilon.getReal() * epsilon.getReal() / ksi.getReal(),
823                             pCirEqua.getPosition().getNorm().getReal(),
824                             Utils.epsilonTest * FastMath.abs(pCirEqua.getPosition().getNorm().getReal()));
825         Assertions.assertEquals(na.getReal() * FastMath.sqrt(ksi.getReal() * ksi.getReal() + nu.getReal() * nu.getReal()) / epsilon.getReal(),
826                             pCirEqua.getPVCoordinates().getVelocity().getNorm().getReal(),
827                             Utils.epsilonTest * FastMath.abs(pCirEqua.getPVCoordinates().getVelocity().getNorm().getReal()));
828     }
829 
830     private <T extends CalculusFieldElement<T>> void doTestGeometryEll(Field<T> field) {
831 
832         T zero =  field.getZero();
833         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
834         // elliptic and non equatorial (i retrograde) orbit
835         T hx =  zero.add(1.2);
836         T hy =  zero.add(2.1);
837         T i  = hx.multiply(hx).add(hy.multiply(hy)).sqrt().atan().multiply(2);
838         T raan = hy.atan2(hx);
839         FieldCircularOrbit<T> p =
840             new FieldCircularOrbit<>(zero.add(42166712.0), zero.add(0.5), zero.add(-0.5), i, raan,
841                                      raan.negate().add(0.67), PositionAngleType.TRUE,
842                                      FramesFactory.getEME2000(), date, zero.add(mu));
843 
844         FieldVector3D<T> position = p.getPosition();
845         FieldVector3D<T> velocity = p.getPVCoordinates().getVelocity();
846         FieldVector3D<T> momentum = p.getPVCoordinates().getMomentum().normalize();
847 
848         T apogeeRadius  = p.getA().multiply( p.getE().add(1));
849         T perigeeRadius = p.getA().multiply(p.getE().negate().add(1));
850 
851         for (T alphaV = zero; alphaV.getReal() <= 2 * FastMath.PI; alphaV=alphaV.add(zero.add(2).multiply(FastMath.PI/100.))) {
852             p = new FieldCircularOrbit<>(p.getA() , p.getCircularEx(), p.getCircularEy(), p.getI(),
853                                          p.getRightAscensionOfAscendingNode(),
854                                          alphaV, PositionAngleType.TRUE, p.getFrame(), date, zero.add(mu));
855             position = p.getPosition();
856             // test if the norm of the position is in the range [perigee radius, apogee radius]
857             // Warning: these tests are without absolute value by choice
858             Assertions.assertTrue((position.getNorm().getReal() - apogeeRadius.getReal())  <= (  apogeeRadius.getReal() * Utils.epsilonTest));
859             Assertions.assertTrue((position.getNorm().getReal() - perigeeRadius.getReal()) >= (- perigeeRadius.getReal() * Utils.epsilonTest));
860 
861             position= position.normalize();
862             velocity = p.getPVCoordinates().getVelocity();
863             velocity= velocity.normalize();
864 
865             // at this stage of computation, all the vectors (position, velocity and momemtum) are normalized here
866 
867             // test of orthogonality between position and momentum
868             Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(position.toVector3D(), momentum.toVector3D())) < Utils.epsilonTest);
869             // test of orthogonality between velocity and momentum
870             Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(velocity.toVector3D(), momentum.toVector3D())) < Utils.epsilonTest);
871         }
872 
873     }
874 
875     private <T extends CalculusFieldElement<T>> void doTestGeometryCirc(Field<T> field) {
876 
877         T zero =  field.getZero();
878         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
879 
880         //  circular and equatorial orbit
881         T hx =  zero.add(0.1e-8);
882         T hy =  zero.add(0.1e-8);
883         T i  =  hx.multiply(hx).add(hy.multiply(hy)).sqrt().atan().multiply(2);
884         T raan = hy.atan2(hx);
885         FieldCircularOrbit<T> pCirEqua =
886             new FieldCircularOrbit<>(zero.add(42166712.0), zero.add(0.1e-8), zero.add(0.1e-8), i, raan,
887                                     raan.negate().add(0.67), PositionAngleType.TRUE,
888                                     FramesFactory.getEME2000(), date, zero.add(mu));
889 
890         FieldVector3D<T> position = pCirEqua.getPosition();
891         FieldVector3D<T> velocity = pCirEqua.getPVCoordinates().getVelocity();
892         FieldVector3D<T> momentum = pCirEqua.getPVCoordinates().getMomentum().normalize();
893 
894         T apogeeRadius  = pCirEqua.getA().multiply( pCirEqua.getE().add(1));
895         T perigeeRadius = pCirEqua.getA().multiply(pCirEqua.getE().negate().add(1));
896         // test if apogee equals perigee
897         Assertions.assertEquals(perigeeRadius.getReal(), apogeeRadius.getReal(), 1.e+4 * Utils.epsilonTest * apogeeRadius.getReal());
898 
899         for (T alphaV = zero; alphaV.getReal() <= 2 * FastMath.PI; alphaV = alphaV.add(zero.add(2 * FastMath.PI/100.))) {
900             pCirEqua = new FieldCircularOrbit<>(pCirEqua.getA() , pCirEqua.getCircularEx(), pCirEqua.getCircularEy(), pCirEqua.getI(),
901                                                pCirEqua.getRightAscensionOfAscendingNode(),
902                                                alphaV, PositionAngleType.TRUE, pCirEqua.getFrame(), date, zero.add(mu));
903             position = pCirEqua.getPosition();
904 
905             // test if the norm pf the position is in the range [perigee radius, apogee radius]
906             Assertions.assertTrue((position.getNorm().getReal() - apogeeRadius.getReal())  <= (  apogeeRadius.getReal() * Utils.epsilonTest));
907             Assertions.assertTrue((position.getNorm().getReal() - perigeeRadius.getReal()) >= (- perigeeRadius.getReal() * Utils.epsilonTest));
908 
909             position= position.normalize();
910             velocity = pCirEqua.getPVCoordinates().getVelocity();
911             velocity= velocity.normalize();
912 
913             // at this stage of computation, all the vectors (position, velocity and momemtum) are normalized here
914 
915             // test of orthogonality between position and momentum
916             Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(position.toVector3D(), momentum.toVector3D())) < Utils.epsilonTest);
917             // test of orthogonality between velocity and momentum
918             Assertions.assertTrue(FastMath.abs(Vector3D.dotProduct(velocity.toVector3D(), momentum.toVector3D())) < Utils.epsilonTest);
919         }
920     }
921 
922     private <T extends CalculusFieldElement<T>> void doTestSymmetryEll(Field<T> field) {
923 
924         T zero =  field.getZero();
925         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
926 
927         // elliptic and non equatorial orbit
928         FieldVector3D<T> position = new FieldVector3D<>(zero.add(4512.9), zero.add(18260.), zero.add(-5127.));
929         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(134664.6), zero.add(90066.8), zero.add(72047.6));
930         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity);
931 
932         FieldCircularOrbit<T> p = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
933 
934         FieldVector3D<T> positionOffset = p.getPosition();
935         FieldVector3D<T> velocityOffset = p.getPVCoordinates().getVelocity();
936 
937         positionOffset = positionOffset.subtract(position);
938         velocityOffset = velocityOffset.subtract(velocity);
939 
940         Assertions.assertEquals(0.0, positionOffset.getNorm().getReal(), position.getNorm().getReal() * Utils.epsilonTest);
941         Assertions.assertEquals(0.0, velocityOffset.getNorm().getReal(), velocity.getNorm().getReal() * Utils.epsilonTest);
942 
943     }
944 
945     private <T extends CalculusFieldElement<T>> void doTestSymmetryCir(Field<T> field) {
946         T zero =  field.getZero();
947         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
948 
949         // circular and equatorial orbit
950         FieldVector3D<T> position = new FieldVector3D<>(zero.add(33051.2), zero.add(26184.9), zero.add(-1.3E-5));
951         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-60376.2), zero.add(76208.), zero.add(2.7E-4));
952         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity);
953 
954         FieldCircularOrbit<T> p = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(), date, zero.add(mu));
955 
956         FieldVector3D<T> positionOffset = p.getPosition().subtract(position);
957         FieldVector3D<T> velocityOffset = p.getPVCoordinates().getVelocity().subtract(velocity);
958 
959         Assertions.assertEquals(0.0, positionOffset.getNorm().getReal(), position.getNorm().getReal() * Utils.epsilonTest);
960         Assertions.assertEquals(0.0, velocityOffset.getNorm().getReal(), velocity.getNorm().getReal() * Utils.epsilonTest);
961 
962     }
963 
964     private <T extends CalculusFieldElement<T>> void doTestNonInertialFrame(Field<T> field) throws IllegalArgumentException {
965         T zero =  field.getZero();
966         FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field);
967 
968         FieldVector3D<T> position = new FieldVector3D<>(zero.add(33051.2), zero.add(26184.9), zero.add(-1.3E-5));
969         FieldVector3D<T> velocity = new FieldVector3D<>(zero.add(-60376.2), zero.add(76208.), zero.add(2.7E-4));
970         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>( position, velocity);
971         new FieldCircularOrbit<>(pvCoordinates,
972                                  new Frame(FramesFactory.getEME2000(), Transform.IDENTITY, "non-inertial", false),
973                                  date, zero.add(mu));
974     }
975 
976     private <T extends CalculusFieldElement<T>> void doTestJacobianReference(Field<T> field) {
977         T zero =  field.getZero();
978         FieldAbsoluteDate<T> dateTca = new FieldAbsoluteDate<>(field, 2000, 4, 1, 0, 0, 0.000, TimeScalesFactory.getUTC());
979         double mu =  3.986004415e+14;
980         FieldCircularOrbit<T> orbCir = new FieldCircularOrbit<>(zero.add(7000000.0), zero.add(0.01), zero.add(-0.02), zero.add(1.2), zero.add(2.1),
981                                                                 zero.add(0.7), PositionAngleType.MEAN,
982                                                                 FramesFactory.getEME2000(), dateTca, zero.add(mu));
983 
984         // the following reference values have been computed using the free software
985         // version 6.2 of the MSLIB fortran library by the following program:
986         //        program cir_jacobian
987         //
988         //        use mslib
989         //        implicit none
990         //
991         //        integer, parameter :: nb = 11
992         //        integer :: i,j
993         //        type(tm_code_retour)      ::  code_retour
994         //
995         //        real(pm_reel), parameter :: mu= 3.986004415e+14_pm_reel
996         //        real(pm_reel),dimension(3)::vit_car,pos_car
997         //        type(tm_orb_cir)::cir
998         //        real(pm_reel), dimension(6,6)::jacob
999         //        real(pm_reel)::norme
1000         //
1001         //
1002         //        cir%a=7000000_pm_reel
1003         //        cir%ex=0.01_pm_reel
1004         //        cir%ey=-0.02_pm_reel
1005         //        cir%i=1.2_pm_reel
1006         //        cir%gom=2.1_pm_reel
1007         //        cir%pso_M=0.7_pm_reel
1008         //
1009         //        call mv_cir_car(mu,cir,pos_car,vit_car,code_retour)
1010         //        write(*,*)code_retour%valeur
1011         //        write(*,1000)pos_car,vit_car
1012         //
1013         //
1014         //        call mu_norme(pos_car,norme,code_retour)
1015         //        write(*,*)norme
1016         //
1017         //        call mv_car_cir (mu, pos_car, vit_car, cir, code_retour, jacob)
1018         //        write(*,*)code_retour%valeur
1019         //
1020         //        write(*,*)"circular = ", cir%a, cir%ex, cir%ey, cir%i, cir%gom, cir%pso_M
1021         //
1022         //        do i = 1,6
1023         //           write(*,*) " ",(jacob(i,j),j=1,6)
1024         //        end do
1025         //
1026         //        1000 format (6(f24.15,1x))
1027         //        end program cir_jacobian
1028         FieldVector3D<T> pRef = new FieldVector3D<>(zero.add(-4106905.105389204807580), zero.add( 3603162.539798960555345), zero.add(4439730.167038885876536));
1029         FieldVector3D<T> vRef = new FieldVector3D<>(zero.add(740.132407342422994), zero.add(-5308.773280141396754), zero.add( 5250.338353483879473));
1030         double[][] jRefR = {
1031             { -1.1535467596325562      ,        1.0120556393573172,        1.2470306024626943,        181.96913090864561,       -1305.2162699469984,        1290.8494448855752 },
1032             { -5.07367368325471104E-008, -1.27870567070456834E-008,  1.31544531338558113E-007, -3.09332106417043592E-005, -9.60781276304445404E-005,  1.91506964883791605E-004 },
1033             { -6.59428471712402018E-008,  1.24561703203882533E-007, -1.41907027322388158E-008,  7.63442601186485441E-005, -1.77446722746170009E-004,  5.99464401287846734E-005 },
1034             {  7.55079920652274275E-008,  4.41606835295069131E-008,  3.40079310688458225E-008,  7.89724635377817962E-005,  4.61868720707717372E-005,  3.55682891687782599E-005 },
1035             { -9.20788748896973282E-008, -5.38521280004949642E-008, -4.14712660805579618E-008,  7.78626692360739821E-005,  4.55378113077967091E-005,  3.50684505810897702E-005 },
1036             {  1.85082436324531617E-008,  1.20506219457886855E-007, -8.31277842285972640E-008,  1.27364008345789645E-004, -1.54770720974742483E-004, -1.78589436862677754E-004 }
1037         };
1038 
1039         T[][] jRef = MathArrays.buildArray(field, 6, 6);
1040 
1041         for (int ii = 0; ii<6 ; ii++){
1042             for (int jj = 0; jj<6 ; jj++){
1043                 jRef[ii][jj] = zero.add(jRefR[ii][jj]);
1044             }
1045         }
1046 
1047         FieldPVCoordinates<T> pv = orbCir.getPVCoordinates();
1048         Assertions.assertEquals(0, pv.getPosition().subtract(pRef).getNorm().getReal(), 3.0e-16 * pRef.getNorm().getReal());
1049         Assertions.assertEquals(0, pv.getVelocity().subtract(vRef).getNorm().getReal(), 2.0e-12 * vRef.getNorm().getReal());
1050 
1051         T[][] jacobian = MathArrays.buildArray(field, 6, 6);
1052 
1053         orbCir.getJacobianWrtCartesian(PositionAngleType.MEAN, jacobian);
1054 
1055         for (int i = 0; i < jacobian.length; i++) {
1056             T[] row    = jacobian[i];
1057             T[] rowRef = jRef[i];
1058             for (int j = 0; j < row.length; j++) {
1059                 Assertions.assertEquals(0, (row[j].getReal() - rowRef[j].getReal()) / rowRef[j].getReal(), 1e-14);
1060             }
1061         }
1062 
1063     }
1064 
1065     private <T extends CalculusFieldElement<T>> void doTestJacobianFinitedifferences(Field<T> field) {
1066         T zero =  field.getZero();
1067 
1068         FieldAbsoluteDate<T> dateTca = new FieldAbsoluteDate<>(field, 2000, 4, 1, 0, 0, 0.000, TimeScalesFactory.getUTC());
1069         double mu =  3.986004415e+14;
1070         FieldCircularOrbit<T> orbCir = new FieldCircularOrbit<>(zero.add(7000000.0), zero.add(0.01), zero.add(-0.02), zero.add(1.2), zero.add(2.1),
1071                                                                 zero.add(0.7), PositionAngleType.MEAN,
1072                                                                 FramesFactory.getEME2000(), dateTca, zero.add(mu));
1073 
1074         for (PositionAngleType type : PositionAngleType.values()) {
1075             T hP = zero.add(2.0);
1076             T[][] finiteDiffJacobian = finiteDifferencesJacobian(type, orbCir, hP);
1077             T[][] jacobian = MathArrays.buildArray(field, 6, 6);
1078             orbCir.getJacobianWrtCartesian(type, jacobian);
1079 
1080             for (int i = 0; i < jacobian.length; i++) {
1081                 T[] row    = jacobian[i];
1082                 T[] rowRef = finiteDiffJacobian[i];
1083 
1084                 for (int j = 0; j < row.length; j++) {
1085                     Assertions.assertEquals(0, (row[j].getReal() - rowRef[j].getReal()) / rowRef[j].getReal(), 8.0e-9);
1086                 }
1087 
1088            }
1089 
1090             T[][] invJacobian = MathArrays.buildArray(field, 6, 6);
1091             orbCir.getJacobianWrtParameters(type, invJacobian);
1092             MatrixUtils.createFieldMatrix(jacobian).
1093                             multiply(MatrixUtils.createFieldMatrix(invJacobian)).
1094             walkInRowOrder(new FieldMatrixPreservingVisitor<T>() {
1095                 public void start(int rows, int columns,
1096                                   int startRow, int endRow, int startColumn, int endColumn) {
1097                 }
1098 
1099                 public void visit(int row, int column, T value) {
1100                     Assertions.assertEquals(row == column ? 1.0 : 0.0, value.getReal(), 3.0e-9);
1101                 }
1102 
1103                 public T end() {
1104                     return null;
1105                 }
1106             });
1107 
1108         }
1109 
1110     }
1111 
1112     private <T extends CalculusFieldElement<T>> T[][] finiteDifferencesJacobian(PositionAngleType type, FieldCircularOrbit<T> orbit, T hP)
1113         {
1114         Field<T> field = hP.getField();
1115         T[][] jacobian = MathArrays.buildArray(field, 6, 6);
1116         for (int i = 0; i < 6; ++i) {
1117             fillColumn(type, i, orbit, hP, jacobian);
1118         }
1119         return jacobian;
1120     }
1121 
1122     private <T extends CalculusFieldElement<T>> void fillColumn(PositionAngleType type, int i, FieldCircularOrbit<T> orbit, T hP, T[][] jacobian) {
1123 
1124         T zero = hP.getField().getZero();
1125         // at constant energy (i.e. constant semi major axis), we have dV = -mu dP / (V * r^2)
1126         // we use this to compute a velocity step size from the position step size
1127         FieldVector3D<T> p = orbit.getPosition();
1128         FieldVector3D<T> v = orbit.getPVCoordinates().getVelocity();
1129         T hV =  hP.multiply(orbit.getMu()).divide(v.getNorm().multiply(p.getNormSq()));
1130 
1131         T h;
1132         FieldVector3D<T> dP = new FieldVector3D<>(zero, zero, zero);
1133         FieldVector3D<T> dV = new FieldVector3D<>(zero, zero, zero);
1134         switch (i) {
1135         case 0:
1136             h = hP;
1137             dP = new FieldVector3D<>(hP, zero, zero);
1138             break;
1139         case 1:
1140             h = hP;
1141             dP = new FieldVector3D<>(zero, hP, zero);
1142             break;
1143         case 2:
1144             h = hP;
1145             dP = new FieldVector3D<>(zero, zero, hP);
1146             break;
1147         case 3:
1148             h = hV;
1149             dV = new FieldVector3D<>(hV, zero, zero);
1150             break;
1151         case 4:
1152             h = hV;
1153             dV = new FieldVector3D<>(zero, hV, zero);
1154             break;
1155         default:
1156             h = hV;
1157             dV = new FieldVector3D<>(zero, zero, hV);
1158             break;
1159         }
1160 
1161         FieldCircularOrbit<T> oM4h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, -4, dP), new FieldVector3D<>(1, v, -4, dV)),
1162                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1163         FieldCircularOrbit<T> oM3h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, -3, dP), new FieldVector3D<>(1, v, -3, dV)),
1164                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1165         FieldCircularOrbit<T> oM2h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, -2, dP), new FieldVector3D<>(1, v, -2, dV)),
1166                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1167         FieldCircularOrbit<T> oM1h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, -1, dP), new FieldVector3D<>(1, v, -1, dV)),
1168                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1169         FieldCircularOrbit<T> oP1h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, +1, dP), new FieldVector3D<>(1, v, +1, dV)),
1170                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1171         FieldCircularOrbit<T> oP2h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, +2, dP), new FieldVector3D<>(1, v, +2, dV)),
1172                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1173         FieldCircularOrbit<T> oP3h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, +3, dP), new FieldVector3D<>(1, v, +3, dV)),
1174                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1175         FieldCircularOrbit<T> oP4h = new FieldCircularOrbit<>(new FieldPVCoordinates<>(new FieldVector3D<>(1, p, +4, dP), new FieldVector3D<>(1, v, +4, dV)),
1176                                                               orbit.getFrame(), orbit.getDate(), orbit.getMu());
1177 
1178 
1179         jacobian[0][i] =(zero.add( -3).multiply(oP4h.getA()                            .subtract( oM4h.getA()))  .add(
1180                          zero.add( 32).multiply(oP3h.getA()                            .subtract( oM3h.getA()))) .subtract(
1181                          zero.add(168).multiply(oP2h.getA()                            .subtract( oM2h.getA()))) .add(
1182                          zero.add(672).multiply(oP1h.getA()                            .subtract( oM1h.getA())))).divide(h.multiply(840));
1183         jacobian[1][i] =(zero.add( -3).multiply(oP4h.getCircularEx()                   .subtract( oM4h.getCircularEx()))  .add(
1184                          zero.add( 32).multiply(oP3h.getCircularEx()                   .subtract( oM3h.getCircularEx()))) .subtract(
1185                          zero.add(168).multiply(oP2h.getCircularEx()                   .subtract( oM2h.getCircularEx()))) .add(
1186                          zero.add(672).multiply(oP1h.getCircularEx()                   .subtract( oM1h.getCircularEx())))).divide(h.multiply(840));
1187         jacobian[2][i] =(zero.add( -3).multiply(oP4h.getCircularEy()                   .subtract( oM4h.getCircularEy()))  .add(
1188                          zero.add( 32).multiply(oP3h.getCircularEy()                   .subtract( oM3h.getCircularEy()))) .subtract(
1189                          zero.add(168).multiply(oP2h.getCircularEy()                   .subtract( oM2h.getCircularEy()))) .add(
1190                          zero.add(672).multiply(oP1h.getCircularEy()                   .subtract( oM1h.getCircularEy())))).divide(h.multiply(840));
1191         jacobian[3][i] =(zero.add( -3).multiply(oP4h.getI()                            .subtract( oM4h.getI()))   .add(
1192                          zero.add( 32).multiply(oP3h.getI()                            .subtract( oM3h.getI())) ) .subtract(
1193                          zero.add(168).multiply(oP2h.getI()                            .subtract( oM2h.getI())) ) .add(
1194                          zero.add(672).multiply(oP1h.getI()                            .subtract( oM1h.getI())))).divide(h.multiply(840));
1195         jacobian[4][i] =(zero.add( -3).multiply(oP4h.getRightAscensionOfAscendingNode().subtract( oM4h.getRightAscensionOfAscendingNode()))   .add(
1196                          zero.add( 32).multiply(oP3h.getRightAscensionOfAscendingNode().subtract( oM3h.getRightAscensionOfAscendingNode()))) .subtract(
1197                          zero.add(168).multiply(oP2h.getRightAscensionOfAscendingNode().subtract( oM2h.getRightAscensionOfAscendingNode()))) .add(
1198                          zero.add(672).multiply(oP1h.getRightAscensionOfAscendingNode().subtract( oM1h.getRightAscensionOfAscendingNode())))).divide(h.multiply(840));
1199         jacobian[5][i] =(zero.add( -3).multiply(oP4h.getAlpha(type)                    .subtract( oM4h.getAlpha(type)))   .add(
1200                          zero.add( 32).multiply(oP3h.getAlpha(type)                    .subtract( oM3h.getAlpha(type)))) .subtract(
1201                          zero.add(168).multiply(oP2h.getAlpha(type)                    .subtract( oM2h.getAlpha(type)))) .add(
1202                          zero.add(672).multiply(oP1h.getAlpha(type)                    .subtract( oM1h.getAlpha(type))))).divide(h.multiply(840));
1203 
1204     }
1205 
1206     private <T extends CalculusFieldElement<T>> void doTestNonKeplerianDerivatives(Field<T> field) {
1207         final T zero = field.getZero();
1208 
1209         final FieldAbsoluteDate<T> date         = new FieldAbsoluteDate<>(field, "2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
1210         final FieldVector3D<T>     position     = new FieldVector3D<>(field.getZero().add(6896874.444705),  field.getZero().add(1956581.072644),  field.getZero().add(-147476.245054));
1211         final FieldVector3D<T>     velocity     = new FieldVector3D<>(field.getZero().add(166.816407662), field.getZero().add(-1106.783301861), field.getZero().add(-7372.745712770));
1212         final FieldVector3D <T>    acceleration = new FieldVector3D<>(field.getZero().add(-7.466182457944), field.getZero().add(-2.118153357345),  field.getZero().add(0.160004048437));
1213         final TimeStampedFieldPVCoordinates<T> pv = new TimeStampedFieldPVCoordinates<>(date, position, velocity, acceleration);
1214         final Frame frame = FramesFactory.getEME2000();
1215         final T mu   = zero.add(Constants.EIGEN5C_EARTH_MU);
1216         final FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(pv, frame, mu);
1217 
1218         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getA),
1219                             orbit.getADot().getReal(),
1220                             4.3e-8);
1221         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getEquinoctialEx),
1222                             orbit.getEquinoctialExDot().getReal(),
1223                             2.1e-15);
1224         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getEquinoctialEy),
1225                             orbit.getEquinoctialEyDot().getReal(),
1226                             5.4e-16);
1227         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getHx),
1228                             orbit.getHxDot().getReal(),
1229                             1.6e-15);
1230         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getHy),
1231                             orbit.getHyDot().getReal(),
1232                             7.3e-17);
1233         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getLv),
1234                             orbit.getLvDot().getReal(),
1235                             3.4e-16);
1236         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getLE),
1237                             orbit.getLEDot().getReal(),
1238                             3.5e-15);
1239         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getLM),
1240                             orbit.getLMDot().getReal(),
1241                             5.3e-15);
1242         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getE),
1243                             orbit.getEDot().getReal(),
1244                             6.8e-16);
1245         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getI),
1246                             orbit.getIDot().getReal(),
1247                             5.7e-16);
1248         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getCircularEx),
1249                             orbit.getCircularExDot().getReal(),
1250                             2.2e-15);
1251         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getCircularEy),
1252                             orbit.getCircularEyDot().getReal(),
1253                             5.3e-17);
1254         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getAlphaV),
1255                             orbit.getAlphaVDot().getReal(),
1256                             4.3e-15);
1257         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getAlphaE),
1258                             orbit.getAlphaEDot().getReal(),
1259                             1.2e-15);
1260         Assertions.assertEquals(differentiate(pv, frame, mu, FieldCircularOrbit::getAlphaM),
1261                             orbit.getAlphaMDot().getReal(),
1262                             3.7e-15);
1263         Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAlpha(PositionAngleType.TRUE)),
1264                             orbit.getAlphaDot(PositionAngleType.TRUE).getReal(),
1265                             4.3e-15);
1266         Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAlpha(PositionAngleType.ECCENTRIC)),
1267                             orbit.getAlphaDot(PositionAngleType.ECCENTRIC).getReal(),
1268                             1.2e-15);
1269         Assertions.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getAlpha(PositionAngleType.MEAN)),
1270                             orbit.getAlphaDot(PositionAngleType.MEAN).getReal(),
1271                             3.7e-15);
1272 
1273     }
1274 
1275     private <T extends CalculusFieldElement<T>, S extends Function<FieldCircularOrbit<T>, T>>
1276     double differentiate(TimeStampedFieldPVCoordinates<T> pv, Frame frame, T mu, S picker) {
1277         final DSFactory factory = new DSFactory(1, 1);
1278         FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(8, 0.1);
1279         UnivariateDifferentiableFunction diff =
1280                 differentiator.differentiate((UnivariateFunction) dt -> picker.apply(new FieldCircularOrbit<>(pv.shiftedBy(dt),
1281                                                                                                               frame, mu)).getReal());
1282         return diff.value(factory.variable(0, 0.0)).getPartialDerivative(1);
1283      }
1284 
1285     private <T extends CalculusFieldElement<T>> void doTestPositionAngleDerivatives(final Field<T> field) {
1286         final T zero = field.getZero();
1287 
1288         final FieldAbsoluteDate<T> date         = new FieldAbsoluteDate<>(field, "2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
1289         final FieldVector3D<T>     position     = new FieldVector3D<>(field.getZero().add(6896874.444705),  field.getZero().add(1956581.072644),  field.getZero().add(-147476.245054));
1290         final FieldVector3D<T>     velocity     = new FieldVector3D<>(field.getZero().add(166.816407662), field.getZero().add(-1106.783301861), field.getZero().add(-7372.745712770));
1291         final FieldVector3D <T>    acceleration = new FieldVector3D<>(field.getZero().add(-7.466182457944), field.getZero().add(-2.118153357345),  field.getZero().add(0.160004048437));
1292         final TimeStampedFieldPVCoordinates<T> pv = new TimeStampedFieldPVCoordinates<>(date, position, velocity, acceleration);
1293         final Frame frame = FramesFactory.getEME2000();
1294         final double mu   = Constants.EIGEN5C_EARTH_MU;
1295         final FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(pv, frame, zero.add(mu));
1296 
1297         for (PositionAngleType type : PositionAngleType.values()) {
1298             final FieldCircularOrbit<T> rebuilt = new FieldCircularOrbit<>(orbit.getA(),
1299                                                                            orbit.getCircularEx(),
1300                                                                            orbit.getCircularEy(),
1301                                                                            orbit.getI(),
1302                                                                            orbit.getRightAscensionOfAscendingNode(),
1303                                                                            orbit.getAlpha(type),
1304                                                                            orbit.getADot(),
1305                                                                            orbit.getCircularExDot(),
1306                                                                            orbit.getCircularEyDot(),
1307                                                                            orbit.getIDot(),
1308                                                                            orbit.getRightAscensionOfAscendingNodeDot(),
1309                                                                            orbit.getAlphaDot(type),
1310                                                                            type, orbit.getFrame(), orbit.getDate(), orbit.getMu());
1311             MatcherAssert.assertThat(rebuilt.getA().getReal(),                                relativelyCloseTo(orbit.getA().getReal(),                                1));
1312             MatcherAssert.assertThat(rebuilt.getCircularEx().getReal(),                       relativelyCloseTo(orbit.getCircularEx().getReal(),                       1));
1313             MatcherAssert.assertThat(rebuilt.getCircularEy().getReal(),                       relativelyCloseTo(orbit.getCircularEy().getReal(),                       1));
1314             MatcherAssert.assertThat(rebuilt.getE().getReal(),                                relativelyCloseTo(orbit.getE().getReal(),                                1));
1315             MatcherAssert.assertThat(rebuilt.getI().getReal(),                                relativelyCloseTo(orbit.getI().getReal(),                                1));
1316             MatcherAssert.assertThat(rebuilt.getRightAscensionOfAscendingNode().getReal(),    relativelyCloseTo(orbit.getRightAscensionOfAscendingNode().getReal(),    1));
1317             MatcherAssert.assertThat(rebuilt.getADot().getReal(),                             relativelyCloseTo(orbit.getADot().getReal(),                             1));
1318             MatcherAssert.assertThat(rebuilt.getCircularExDot().getReal(),                    relativelyCloseTo(orbit.getCircularExDot().getReal(),                    1));
1319             MatcherAssert.assertThat(rebuilt.getCircularEyDot().getReal(),                    relativelyCloseTo(orbit.getCircularEyDot().getReal(),                    1));
1320             MatcherAssert.assertThat(rebuilt.getEDot().getReal(),                             relativelyCloseTo(orbit.getEDot().getReal(),                             1));
1321             MatcherAssert.assertThat(rebuilt.getIDot().getReal(),                             relativelyCloseTo(orbit.getIDot().getReal(),                             1));
1322             MatcherAssert.assertThat(rebuilt.getRightAscensionOfAscendingNodeDot().getReal(), relativelyCloseTo(orbit.getRightAscensionOfAscendingNodeDot().getReal(), 1));
1323             for (PositionAngleType type2 : PositionAngleType.values()) {
1324                 MatcherAssert.assertThat(rebuilt.getAlpha(type2).getReal(),    relativelyCloseTo(orbit.getAlpha(type2).getReal(),    1));
1325                 MatcherAssert.assertThat(rebuilt.getAlphaDot(type2).getReal(), relativelyCloseTo(orbit.getAlphaDot(type2).getReal(), 1));
1326             }
1327         }
1328 
1329     }
1330 
1331     private <T extends CalculusFieldElement<T>> void doTestEquatorialRetrograde(final Field<T> field) {
1332         T zero = field.getZero();
1333         FieldVector3D<T> position = new FieldVector3D<>(field.getZero().add(10000000.0), field.getZero(), field.getZero());
1334         FieldVector3D<T> velocity = new FieldVector3D<>(field.getZero(), field.getZero().add(-6500.0), field.getZero().add(1.0e-10));
1335         T r2 = position.getNormSq();
1336         T r  = r2.sqrt();
1337         FieldVector3D<T> acceleration = new FieldVector3D<>(r.multiply(r2.reciprocal().multiply(-mu)), position,
1338                                                             field.getOne(), new FieldVector3D<>(field.getZero().add(-0.1),
1339                                                                                                 field.getZero().add(0.2),
1340                                                                                                 field.getZero().add(0.3)));
1341         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity, acceleration);
1342         FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(),
1343                                                                FieldAbsoluteDate.getJ2000Epoch(field), zero.add(mu));
1344         Assertions.assertEquals(10637829.465, orbit.getA().getReal(), 1.0e-3);
1345         Assertions.assertEquals(-738.145, orbit.getADot().getReal(), 1.0e-3);
1346         Assertions.assertEquals(0.05995861, orbit.getE().getReal(), 1.0e-8);
1347         Assertions.assertEquals(-6.523e-5, orbit.getEDot().getReal(), 1.0e-8);
1348         Assertions.assertEquals(FastMath.PI, orbit.getI().getReal(), 2.0e-14);
1349         Assertions.assertEquals(-4.615e-5, orbit.getIDot().getReal(), 1.0e-8);
1350         Assertions.assertTrue(Double.isNaN(orbit.getHx().getReal()));
1351         Assertions.assertTrue(Double.isNaN(orbit.getHxDot().getReal()));
1352         Assertions.assertTrue(Double.isNaN(orbit.getHy().getReal()));
1353         Assertions.assertTrue(Double.isNaN(orbit.getHyDot().getReal()));
1354     }
1355 
1356     private <T extends CalculusFieldElement<T>> void doTestDerivativesConversionSymmetry(Field<T> field) {
1357         T zero = field.getZero();
1358         final FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field, "2003-05-01T00:01:20.000", TimeScalesFactory.getUTC());
1359         FieldVector3D<T> position     = new FieldVector3D<>(zero.add(6893443.400234382),
1360                                                             zero.add(1886406.1073757345),
1361                                                             zero.add(-589265.1150359757));
1362         FieldVector3D<T> velocity     = new FieldVector3D<>(zero.add(-281.1261461082365),
1363                                                             zero.add(-1231.6165642450928),
1364                                                             zero.add(-7348.756363469432));
1365         FieldVector3D<T> acceleration = new FieldVector3D<>(zero.add(-7.460341170581685),
1366                                                             zero.add(-2.0415957334584527),
1367                                                             zero.add(0.6393322823627762));
1368         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>( position, velocity, acceleration);
1369         FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(),
1370                                                                date, zero.add(Constants.EIGEN5C_EARTH_MU));
1371         Assertions.assertTrue(orbit.hasNonKeplerianAcceleration());
1372         T r2 = position.getNormSq();
1373         T r  = r2.sqrt();
1374         FieldVector3D<T> keplerianAcceleration = new FieldVector3D<>(r.multiply(r2).reciprocal().multiply(orbit.getMu().negate()),
1375                                                                      position);
1376         Assertions.assertEquals(0.0101, FieldVector3D.distance(keplerianAcceleration, acceleration).getReal(), 1.0e-4);
1377 
1378         for (OrbitType type : OrbitType.values()) {
1379             FieldOrbit<T> converted = type.convertType(orbit);
1380             Assertions.assertTrue(converted.hasNonKeplerianAcceleration());
1381             FieldCircularOrbit<T> rebuilt = (FieldCircularOrbit<T>) OrbitType.CIRCULAR.convertType(converted);
1382             Assertions.assertTrue(rebuilt.hasNonKeplerianAcceleration());
1383             Assertions.assertEquals(orbit.getADot().getReal(),                             rebuilt.getADot().getReal(),                             3.0e-13);
1384             Assertions.assertEquals(orbit.getCircularExDot().getReal(),                    rebuilt.getCircularExDot().getReal(),                    1.0e-15);
1385             Assertions.assertEquals(orbit.getCircularEyDot().getReal(),                    rebuilt.getCircularEyDot().getReal(),                    1.0e-15);
1386             Assertions.assertEquals(orbit.getIDot().getReal(),                             rebuilt.getIDot().getReal(),                             1.0e-15);
1387             Assertions.assertEquals(orbit.getRightAscensionOfAscendingNodeDot().getReal(), rebuilt.getRightAscensionOfAscendingNodeDot().getReal(), 1.0e-15);
1388             Assertions.assertEquals(orbit.getAlphaVDot().getReal(),                        rebuilt.getAlphaVDot().getReal(),                        1.0e-15);
1389         }
1390 
1391     }
1392 
1393     private <T extends CalculusFieldElement<T>> void doTestToString(Field<T> field) {
1394         T zero = field.getZero();
1395         FieldVector3D<T> position = new FieldVector3D<>(field.getZero().add(-29536113.0),
1396                                                         field.getZero().add(30329259.0),
1397                                                         field.getZero().add(-100125.0));
1398         FieldVector3D<T> velocity = new FieldVector3D<>(field.getZero().add(-2194.0),
1399                                                         field.getZero().add(-2141.0),
1400                                                         field.getZero().add(-8.0));
1401         FieldPVCoordinates<T> pvCoordinates = new FieldPVCoordinates<>(position, velocity);
1402         FieldCircularOrbit<T> orbit = new FieldCircularOrbit<>(pvCoordinates, FramesFactory.getEME2000(),
1403                                                                FieldAbsoluteDate.getJ2000Epoch(field), zero.add(mu));
1404         Assertions.assertEquals("circular parameters: {a: 4.225517000282565E7, ex: 0.002082917137146049, ey: 5.173980074371024E-4, i: 0.20189257051515358, raan: -87.91788415673473, alphaV: -137.84099636616548;}",
1405                             orbit.toString());
1406     }
1407 
1408     private <T extends CalculusFieldElement<T>> void doTestCopyNonKeplerianAcceleration(Field<T> field) {
1409 
1410         final T zero = field.getZero();
1411         final Frame eme2000     = FramesFactory.getEME2000();
1412 
1413         // Define GEO satellite position
1414         final FieldVector3D<T> position = new FieldVector3D<>(field.getZero().add(42164140),
1415                         field.getZero(),
1416                         field.getZero());
1417         // Build PVCoodrinates starting from its position and computing the corresponding circular velocity
1418         final FieldPVCoordinates<T> pv  =
1419                         new FieldPVCoordinates<>(position,
1420                                         new FieldVector3D<>(field.getZero(),
1421                                                         position.getNorm().reciprocal().multiply(mu).sqrt(),
1422                                                         field.getZero()));
1423         // Build a KeplerianOrbit in eme2000
1424         final FieldOrbit<T> orbit = new FieldCircularOrbit<>(pv, eme2000, FieldAbsoluteDate.getJ2000Epoch(field), zero.add(mu));
1425 
1426         // Build another KeplerianOrbit as a copy of the first one
1427         final FieldOrbit<T> orbitCopy = new FieldCircularOrbit<>(orbit);
1428 
1429         // Shift the orbit of a time-interval
1430         final FieldOrbit<T> shiftedOrbit     = orbit.shiftedBy(10); // This works good
1431         final FieldOrbit<T> shiftedOrbitCopy = orbitCopy.shiftedBy(10); // This does not work
1432 
1433         Assertions.assertEquals(0.0,
1434                             FieldVector3D.distance(shiftedOrbit.getPosition(),
1435                                                    shiftedOrbitCopy.getPosition()).getReal(),
1436                             1.0e-10);
1437         Assertions.assertEquals(0.0,
1438                             FieldVector3D.distance(shiftedOrbit.getPVCoordinates().getVelocity(),
1439                                                    shiftedOrbitCopy.getPVCoordinates().getVelocity()).getReal(),
1440                             1.0e-10);
1441 
1442     }
1443 
1444     private <T extends CalculusFieldElement<T>> void doTestNormalize(Field<T> field) {
1445         final T zero = field.getZero();
1446         FieldCircularOrbit<T> withoutDerivatives =
1447                         new FieldCircularOrbit<>(zero.newInstance(42166712.0), zero.newInstance(0.005),
1448                                                  zero.newInstance(-0.025), zero.newInstance(1.6),
1449                                                  zero.newInstance(1.25), zero.newInstance(0.4), PositionAngleType.MEAN,
1450                                                  FramesFactory.getEME2000(), FieldAbsoluteDate.getJ2000Epoch(field),
1451                                                  zero.newInstance(mu));
1452         FieldCircularOrbit<T> ref =
1453                         new FieldCircularOrbit<>(zero.newInstance(24000000.0), zero.newInstance(-0.012),
1454                                                  zero.newInstance(0.01), zero.newInstance(0.2),
1455                                                  zero.newInstance(-6.28), zero.newInstance(6.28), PositionAngleType.MEAN,
1456                                                  FramesFactory.getEME2000(), FieldAbsoluteDate.getJ2000Epoch(field),
1457                                                  zero.newInstance(mu));
1458 
1459         FieldCircularOrbit<T> normalized1 = (FieldCircularOrbit<T>) OrbitType.CIRCULAR.normalize(withoutDerivatives, ref);
1460         Assertions.assertFalse(normalized1.hasNonKeplerianAcceleration());
1461         Assertions.assertEquals(0.0, normalized1.getA().subtract(withoutDerivatives.getA()).getReal(),          1.0e-6);
1462         Assertions.assertEquals(0.0, normalized1.getCircularEx().subtract(withoutDerivatives.getCircularEx()).getReal(), 1.0e-10);
1463         Assertions.assertEquals(0.0, normalized1.getCircularEy().subtract(withoutDerivatives.getCircularEy()).getReal(), 1.0e-10);
1464         Assertions.assertEquals(0.0, normalized1.getI().subtract(withoutDerivatives.getI()).getReal(),          1.0e-10);
1465         Assertions.assertEquals(-MathUtils.TWO_PI, normalized1.getRightAscensionOfAscendingNode().subtract(withoutDerivatives.getRightAscensionOfAscendingNode()).getReal(), 1.0e-10);
1466         Assertions.assertEquals(+MathUtils.TWO_PI, normalized1.getAlphaV().subtract(withoutDerivatives.getAlphaV()).getReal(), 1.0e-10);
1467 
1468         T[] p    = MathArrays.buildArray(field, 6);
1469         T[] pDot = MathArrays.buildArray(field, 6);
1470         for (int i = 0; i < pDot.length; ++i) {
1471             pDot[i] = zero.newInstance(i);
1472         }
1473         OrbitType.CIRCULAR.mapOrbitToArray(withoutDerivatives, PositionAngleType.TRUE, p, null);
1474         FieldCircularOrbit<T> withDerivatives = (FieldCircularOrbit<T>) OrbitType.CIRCULAR.mapArrayToOrbit(p, pDot,
1475                                                                                                            PositionAngleType.TRUE,
1476                                                                                                            withoutDerivatives.getDate(),
1477                                                                                                            withoutDerivatives.getMu(),
1478                                                                                                            withoutDerivatives.getFrame());
1479         FieldCircularOrbit<T> normalized2 = (FieldCircularOrbit<T>) OrbitType.CIRCULAR.normalize(withDerivatives, ref);
1480         Assertions.assertTrue(normalized2.hasNonKeplerianAcceleration());
1481         Assertions.assertEquals(0.0, normalized2.getA().subtract(withDerivatives.getA()).getReal(),          1.0e-6);
1482         Assertions.assertEquals(0.0, normalized2.getCircularEx().subtract(withDerivatives.getCircularEx()).getReal(), 1.0e-10);
1483         Assertions.assertEquals(0.0, normalized2.getCircularEy().subtract(withDerivatives.getCircularEy()).getReal(), 1.0e-10);
1484         Assertions.assertEquals(0.0, normalized2.getI().subtract(withDerivatives.getI()).getReal(),          1.0e-10);
1485         Assertions.assertEquals(-MathUtils.TWO_PI, normalized2.getRightAscensionOfAscendingNode().subtract(withDerivatives.getRightAscensionOfAscendingNode()).getReal(), 1.0e-10);
1486         Assertions.assertEquals(+MathUtils.TWO_PI, normalized2.getAlphaV().subtract(withDerivatives.getAlphaV()).getReal(), 1.0e-10);
1487         Assertions.assertEquals(0.0, normalized2.getADot().subtract(withDerivatives.getADot()).getReal(),          1.0e-6);
1488         Assertions.assertEquals(0.0, normalized2.getCircularExDot().subtract(withDerivatives.getCircularExDot()).getReal(), 1.0e-10);
1489         Assertions.assertEquals(0.0, normalized2.getCircularEyDot().subtract(withDerivatives.getCircularEyDot()).getReal(), 1.0e-10);
1490         Assertions.assertEquals(0.0, normalized2.getIDot().subtract(withDerivatives.getIDot()).getReal(),          1.0e-10);
1491         Assertions.assertEquals(0.0, normalized2.getRightAscensionOfAscendingNodeDot().subtract(withDerivatives.getRightAscensionOfAscendingNodeDot()).getReal(), 1.0e-10);
1492         Assertions.assertEquals(0.0, normalized2.getAlphaVDot().subtract(withDerivatives.getAlphaVDot()).getReal(), 1.0e-10);
1493 
1494     }
1495 
1496 }