NumericalPropagator.java

  1. /* Copyright 2002-2022 CS GROUP
  2.  * Licensed to CS GROUP (CS) under one or more
  3.  * contributor license agreements.  See the NOTICE file distributed with
  4.  * this work for additional information regarding copyright ownership.
  5.  * CS licenses this file to You under the Apache License, Version 2.0
  6.  * (the "License"); you may not use this file except in compliance with
  7.  * the License.  You may obtain a copy of the License at
  8.  *
  9.  *   http://www.apache.org/licenses/LICENSE-2.0
  10.  *
  11.  * Unless required by applicable law or agreed to in writing, software
  12.  * distributed under the License is distributed on an "AS IS" BASIS,
  13.  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  14.  * See the License for the specific language governing permissions and
  15.  * limitations under the License.
  16.  */
  17. package org.orekit.propagation.numerical;

  19. import java.util.ArrayList;
  20. import java.util.Arrays;
  21. import java.util.Collections;
  22. import java.util.List;

  24. import org.hipparchus.exception.LocalizedCoreFormats;
  25. import org.hipparchus.geometry.euclidean.threed.Vector3D;
  26. import org.hipparchus.linear.MatrixUtils;
  27. import org.hipparchus.linear.QRDecomposition;
  28. import org.hipparchus.linear.RealMatrix;
  29. import org.hipparchus.ode.ODEIntegrator;
  30. import org.hipparchus.util.FastMath;
  31. import org.hipparchus.util.Precision;
  32. import org.orekit.annotation.DefaultDataContext;
  33. import org.orekit.attitudes.Attitude;
  34. import org.orekit.attitudes.AttitudeProvider;
  35. import org.orekit.data.DataContext;
  36. import org.orekit.errors.OrekitException;
  37. import org.orekit.errors.OrekitIllegalArgumentException;
  38. import org.orekit.errors.OrekitInternalError;
  39. import org.orekit.errors.OrekitMessages;
  40. import org.orekit.forces.ForceModel;
  41. import org.orekit.forces.gravity.NewtonianAttraction;
  42. import org.orekit.forces.inertia.InertialForces;
  43. import org.orekit.forces.maneuvers.Maneuver;
  44. import org.orekit.forces.maneuvers.jacobians.Duration;
  45. import org.orekit.forces.maneuvers.jacobians.MedianDate;
  46. import org.orekit.forces.maneuvers.jacobians.TriggerDate;
  47. import org.orekit.forces.maneuvers.trigger.AbstractManeuverTriggers;
  48. import org.orekit.forces.maneuvers.trigger.ManeuverTriggers;
  49. import org.orekit.frames.Frame;
  50. import org.orekit.orbits.Orbit;
  51. import org.orekit.orbits.OrbitType;
  52. import org.orekit.orbits.PositionAngle;
  53. import org.orekit.propagation.AbstractMatricesHarvester;
  54. import org.orekit.propagation.AdditionalStateProvider;
  55. import org.orekit.propagation.MatricesHarvester;
  56. import org.orekit.propagation.PropagationType;
  57. import org.orekit.propagation.Propagator;
  58. import org.orekit.propagation.SpacecraftState;
  59. import org.orekit.propagation.events.EventDetector;
  60. import org.orekit.propagation.events.ParameterDrivenDateIntervalDetector;
  61. import org.orekit.propagation.integration.AbstractIntegratedPropagator;
  62. import org.orekit.propagation.integration.AdditionalDerivativesProvider;
  63. import org.orekit.propagation.integration.StateMapper;
  64. import org.orekit.time.AbsoluteDate;
  65. import org.orekit.utils.AbsolutePVCoordinates;
  66. import org.orekit.utils.DoubleArrayDictionary;
  67. import org.orekit.utils.PVCoordinates;
  68. import org.orekit.utils.ParameterDriver;
  69. import org.orekit.utils.ParameterDriversList;
  70. import org.orekit.utils.ParameterDriversList.DelegatingDriver;
  71. import org.orekit.utils.ParameterObserver;
  72. import org.orekit.utils.TimeStampedPVCoordinates;

  74. /** This class propagates {@link org.orekit.orbits.Orbit orbits} using
  75.  * numerical integration.
  76.  * <p>Numerical propagation is much more accurate than analytical propagation
  77.  * like for example {@link org.orekit.propagation.analytical.KeplerianPropagator
  78.  * Keplerian} or {@link org.orekit.propagation.analytical.EcksteinHechlerPropagator
  79.  * Eckstein-Hechler}, but requires a few more steps to set up to be used properly.
  80.  * Whereas analytical propagators are configured only thanks to their various
  81.  * constructors and can be used immediately after construction, numerical propagators
  82.  * configuration involve setting several parameters between construction time
  83.  * and propagation time.</p>
  84.  * <p>The configuration parameters that can be set are:</p>
  85.  * <ul>
  86.  *   <li>the initial spacecraft state ({@link #setInitialState(SpacecraftState)})</li>
  87.  *   <li>the central attraction coefficient ({@link #setMu(double)})</li>
  88.  *   <li>the various force models ({@link #addForceModel(ForceModel)},
  89.  *   {@link #removeForceModels()})</li>
  90.  *   <li>the {@link OrbitType type} of orbital parameters to be used for propagation
  91.  *   ({@link #setOrbitType(OrbitType)}),</li>
  92.  *   <li>the {@link PositionAngle type} of position angle to be used in orbital parameters
  93.  *   to be used for propagation where it is relevant ({@link
  94.  *   #setPositionAngleType(PositionAngle)}),</li>
  95.  *   <li>whether {@link MatricesHarvester state transition matrices and Jacobians matrices}
  96.  *   should be propagated along with orbital state ({@link
  97.  *   #setupMatricesComputation(String, RealMatrix, DoubleArrayDictionary)}),</li>
  98.  *   <li>whether {@link org.orekit.propagation.integration.AdditionalDerivativesProvider additional derivatives}
  99.  *   should be propagated along with orbital state ({@link
  100.  *   #addAdditionalDerivativesProvider(AdditionalDerivativesProvider)}),</li>
  101.  *   <li>the discrete events that should be triggered during propagation
  102.  *   ({@link #addEventDetector(EventDetector)},
  103.  *   {@link #clearEventsDetectors()})</li>
  104.  *   <li>the binding logic with the rest of the application ({@link #getMultiplexer()})</li>
  105.  * </ul>
  106.  * <p>From these configuration parameters, only the initial state is mandatory. The default
  107.  * propagation settings are in {@link OrbitType#EQUINOCTIAL equinoctial} parameters with
  108.  * {@link PositionAngle#TRUE true} longitude argument. If the central attraction coefficient
  109.  * is not explicitly specified, the one used to define the initial orbit will be used.
  110.  * However, specifying only the initial state and perhaps the central attraction coefficient
  111.  * would mean the propagator would use only Keplerian forces. In this case, the simpler {@link
  112.  * org.orekit.propagation.analytical.KeplerianPropagator KeplerianPropagator} class would
  113.  * perhaps be more effective.</p>
  114.  * <p>The underlying numerical integrator set up in the constructor may also have its own
  115.  * configuration parameters. Typical configuration parameters for adaptive stepsize integrators
  116.  * are the min, max and perhaps start step size as well as the absolute and/or relative errors
  117.  * thresholds.</p>
  118.  * <p>The state that is seen by the integrator is a simple seven elements double array.
  119.  * The six first elements are either:
  120.  * <ul>
  121.  *   <li>the {@link org.orekit.orbits.EquinoctialOrbit equinoctial orbit parameters} (a, e<sub>x</sub>,
  122.  *   e<sub>y</sub>, h<sub>x</sub>, h<sub>y</sub>, λ<sub>M</sub> or λ<sub>E</sub>
  123.  *   or λ<sub>v</sub>) in meters and radians,</li>
  124.  *   <li>the {@link org.orekit.orbits.KeplerianOrbit Keplerian orbit parameters} (a, e, i, ω, Ω,
  125.  *   M or E or v) in meters and radians,</li>
  126.  *   <li>the {@link org.orekit.orbits.CircularOrbit circular orbit parameters} (a, e<sub>x</sub>, e<sub>y</sub>, i,
  127.  *   Ω, α<sub>M</sub> or α<sub>E</sub> or α<sub>v</sub>) in meters
  128.  *   and radians,</li>
  129.  *   <li>the {@link org.orekit.orbits.CartesianOrbit Cartesian orbit parameters} (x, y, z, v<sub>x</sub>,
  130.  *   v<sub>y</sub>, v<sub>z</sub>) in meters and meters per seconds.
  131.  * </ul>
  132.  * <p> The last element is the mass in kilograms and changes only during thrusters firings
  133.  *
  134.  * <p>The following code snippet shows a typical setting for Low Earth Orbit propagation in
  135.  * equinoctial parameters and true longitude argument:</p>
  136.  * <pre>
  137.  * final double dP       = 0.001;
  138.  * final double minStep  = 0.001;
  139.  * final double maxStep  = 500;
  140.  * final double initStep = 60;
  141.  * final double[][] tolerance = NumericalPropagator.tolerances(dP, orbit, OrbitType.EQUINOCTIAL);
  142.  * AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, tolerance[0], tolerance[1]);
  143.  * integrator.setInitialStepSize(initStep);
  144.  * propagator = new NumericalPropagator(integrator);
  145.  * </pre>
  146.  * <p>By default, at the end of the propagation, the propagator resets the initial state to the final state,
  147.  * thus allowing a new propagation to be started from there without recomputing the part already performed.
  148.  * This behaviour can be chenged by calling {@link #setResetAtEnd(boolean)}.
  149.  * </p>
  150.  * <p>Beware the same instance cannot be used simultaneously by different threads, the class is <em>not</em>
  151.  * thread-safe.</p>
  152.  *
  153.  * @see SpacecraftState
  154.  * @see ForceModel
  155.  * @see org.orekit.propagation.sampling.OrekitStepHandler
  156.  * @see org.orekit.propagation.sampling.OrekitFixedStepHandler
  157.  * @see org.orekit.propagation.integration.IntegratedEphemeris
  158.  * @see TimeDerivativesEquations
  159.  *
  160.  * @author Mathieu Rom&eacute;ro
  161.  * @author Luc Maisonobe
  162.  * @author Guylaine Prat
  163.  * @author Fabien Maussion
  164.  * @author V&eacute;ronique Pommier-Maurussane
  165.  */
  166. public class NumericalPropagator extends AbstractIntegratedPropagator {

  168.     /** Space dimension. */
  169.     private static final int SPACE_DIMENSION = 3;

  171.     /** State dimension. */
  172.     private static final int STATE_DIMENSION = 2 * SPACE_DIMENSION;

  174.     /** Threshold for matrix solving. */
  175.     private static final double THRESHOLD = Precision.SAFE_MIN;

  177.     /** Force models used during the extrapolation of the orbit. */
  178.     private final List<ForceModel> forceModels;

  180.     /** boolean to ignore or not the creation of a NewtonianAttraction. */
  181.     private boolean ignoreCentralAttraction;

  183.     /** Create a new instance of NumericalPropagator, based on orbit definition mu.
  184.      * After creation, the instance is empty, i.e. the attitude provider is set to an
  185.      * unspecified default law and there are no perturbing forces at all.
  186.      * This means that if {@link #addForceModel addForceModel} is not
  187.      * called after creation, the integrated orbit will follow a Keplerian
  188.      * evolution only. The defaults are {@link OrbitType#EQUINOCTIAL}
  189.      * for {@link #setOrbitType(OrbitType) propagation
  190.      * orbit type} and {@link PositionAngle#TRUE} for {@link
  191.      * #setPositionAngleType(PositionAngle) position angle type}.
  192.      *
  193.      * <p>This constructor uses the {@link DataContext#getDefault() default data context}.
  194.      *
  195.      * @param integrator numerical integrator to use for propagation.
  196.      * @see #NumericalPropagator(ODEIntegrator, AttitudeProvider)
  197.      */
  198.     @DefaultDataContext
  199.     public NumericalPropagator(final ODEIntegrator integrator) {
  200.         this(integrator,
  201.                 Propagator.getDefaultLaw(DataContext.getDefault().getFrames()));
  202.     }

  204.     /** Create a new instance of NumericalPropagator, based on orbit definition mu.
  205.      * After creation, the instance is empty, i.e. the attitude provider is set to an
  206.      * unspecified default law and there are no perturbing forces at all.
  207.      * This means that if {@link #addForceModel addForceModel} is not
  208.      * called after creation, the integrated orbit will follow a Keplerian
  209.      * evolution only. The defaults are {@link OrbitType#EQUINOCTIAL}
  210.      * for {@link #setOrbitType(OrbitType) propagation
  211.      * orbit type} and {@link PositionAngle#TRUE} for {@link
  212.      * #setPositionAngleType(PositionAngle) position angle type}.
  213.      * @param integrator numerical integrator to use for propagation.
  214.      * @param attitudeProvider the attitude law.
  215.      * @since 10.1
  216.      */
  217.     public NumericalPropagator(final ODEIntegrator integrator,
  218.                                final AttitudeProvider attitudeProvider) {
  219.         super(integrator, PropagationType.OSCULATING);
  220.         forceModels             = new ArrayList<ForceModel>();
  221.         ignoreCentralAttraction = false;
  222.         initMapper();
  223.         setAttitudeProvider(attitudeProvider);
  224.         clearStepHandlers();
  225.         setOrbitType(OrbitType.EQUINOCTIAL);
  226.         setPositionAngleType(PositionAngle.TRUE);
  227.     }

  229.     /** Set the flag to ignore or not the creation of a {@link NewtonianAttraction}.
  230.      * @param ignoreCentralAttraction if true, {@link NewtonianAttraction} is <em>not</em>
  231.      * added automatically if missing
  232.      */
  233.     public void setIgnoreCentralAttraction(final boolean ignoreCentralAttraction) {
  234.         this.ignoreCentralAttraction = ignoreCentralAttraction;
  235.     }

  237.      /** Set the central attraction coefficient μ.
  238.       * <p>
  239.       * Setting the central attraction coefficient is
  240.       * equivalent to {@link #addForceModel(ForceModel) add}
  241.       * a {@link NewtonianAttraction} force model.
  242.       * </p>
  243.      * @param mu central attraction coefficient (m³/s²)
  244.      * @see #addForceModel(ForceModel)
  245.      * @see #getAllForceModels()
  246.      */
  247.     public void setMu(final double mu) {
  248.         if (ignoreCentralAttraction) {
  249.             superSetMu(mu);
  250.         } else {
  251.             addForceModel(new NewtonianAttraction(mu));
  252.         }
  253.     }

  255.     /** Set the central attraction coefficient μ only in upper class.
  256.      * @param mu central attraction coefficient (m³/s²)
  257.      */
  258.     private void superSetMu(final double mu) {
  259.         super.setMu(mu);
  260.     }

  262.     /** Check if Newtonian attraction force model is available.
  263.      * <p>
  264.      * Newtonian attraction is always the last force model in the list.
  265.      * </p>
  266.      * @return true if Newtonian attraction force model is available
  267.      */
  268.     private boolean hasNewtonianAttraction() {
  269.         final int last = forceModels.size() - 1;
  270.         return last >= 0 && forceModels.get(last) instanceof NewtonianAttraction;
  271.     }

  273.     /** Add a force model.
  274.      * <p>If this method is not called at all, the integrated orbit will follow
  275.      * a Keplerian evolution only.</p>
  276.      * @param model {@link ForceModel} to add (it can be either a perturbing force
  277.      * model or an instance of {@link NewtonianAttraction})
  278.      * @see #removeForceModels()
  279.      * @see #setMu(double)
  280.      */
  281.     public void addForceModel(final ForceModel model) {

  283.         if (model instanceof NewtonianAttraction) {
  284.             // we want to add the central attraction force model

  286.             try {
  287.                 // ensure we are notified of any mu change
  288.                 model.getParametersDrivers().get(0).addObserver(new ParameterObserver() {
  289.                     /** {@inheritDoc} */
  290.                     @Override
  291.                     public void valueChanged(final double previousValue, final ParameterDriver driver) {
  292.                         superSetMu(driver.getValue());
  293.                     }
  294.                 });
  295.             } catch (OrekitException oe) {
  296.                 // this should never happen
  297.                 throw new OrekitInternalError(oe);
  298.             }

  300.             if (hasNewtonianAttraction()) {
  301.                 // there is already a central attraction model, replace it
  302.                 forceModels.set(forceModels.size() - 1, model);
  303.             } else {
  304.                 // there are no central attraction model yet, add it at the end of the list
  305.                 forceModels.add(model);
  306.             }
  307.         } else {
  308.             // we want to add a perturbing force model
  309.             if (hasNewtonianAttraction()) {
  310.                 // insert the new force model before Newtonian attraction,
  311.                 // which should always be the last one in the list
  312.                 forceModels.add(forceModels.size() - 1, model);
  313.             } else {
  314.                 // we only have perturbing force models up to now, just append at the end of the list
  315.                 forceModels.add(model);
  316.             }
  317.         }

  319.     }

  321.     /** Remove all force models (except central attraction).
  322.      * <p>Once all perturbing forces have been removed (and as long as no new force
  323.      * model is added), the integrated orbit will follow a Keplerian evolution
  324.      * only.</p>
  325.      * @see #addForceModel(ForceModel)
  326.      */
  327.     public void removeForceModels() {
  328.         final int last = forceModels.size() - 1;
  329.         if (hasNewtonianAttraction()) {
  330.             // preserve the Newtonian attraction model at the end
  331.             final ForceModel newton = forceModels.get(last);
  332.             forceModels.clear();
  333.             forceModels.add(newton);
  334.         } else {
  335.             forceModels.clear();
  336.         }
  337.     }

  339.     /** Get all the force models, perturbing forces and Newtonian attraction included.
  340.      * @return list of perturbing force models, with Newtonian attraction being the
  341.      * last one
  342.      * @see #addForceModel(ForceModel)
  343.      * @see #setMu(double)
  344.      */
  345.     public List<ForceModel> getAllForceModels() {
  346.         return Collections.unmodifiableList(forceModels);
  347.     }

  349.     /** Set propagation orbit type.
  350.      * @param orbitType orbit type to use for propagation, null for
  351.      * propagating using {@link org.orekit.utils.AbsolutePVCoordinates} rather than {@link Orbit}
  352.      */
  353.     public void setOrbitType(final OrbitType orbitType) {
  354.         super.setOrbitType(orbitType);
  355.     }

  357.     /** Get propagation parameter type.
  358.      * @return orbit type used for propagation, null for
  359.      * propagating using {@link org.orekit.utils.AbsolutePVCoordinates} rather than {@link Orbit}
  360.      */
  361.     public OrbitType getOrbitType() {
  362.         return super.getOrbitType();
  363.     }

  365.     /** Set position angle type.
  366.      * <p>
  367.      * The position parameter type is meaningful only if {@link
  368.      * #getOrbitType() propagation orbit type}
  369.      * support it. As an example, it is not meaningful for propagation
  370.      * in {@link OrbitType#CARTESIAN Cartesian} parameters.
  371.      * </p>
  372.      * @param positionAngleType angle type to use for propagation
  373.      */
  374.     public void setPositionAngleType(final PositionAngle positionAngleType) {
  375.         super.setPositionAngleType(positionAngleType);
  376.     }

  378.     /** Get propagation parameter type.
  379.      * @return angle type to use for propagation
  380.      */
  381.     public PositionAngle getPositionAngleType() {
  382.         return super.getPositionAngleType();
  383.     }

  385.     /** Set the initial state.
  386.      * @param initialState initial state
  387.      */
  388.     public void setInitialState(final SpacecraftState initialState) {
  389.         resetInitialState(initialState);
  390.     }

  392.     /** {@inheritDoc} */
  393.     public void resetInitialState(final SpacecraftState state) {
  394.         super.resetInitialState(state);
  395.         if (!hasNewtonianAttraction()) {
  396.             // use the state to define central attraction
  397.             setMu(state.getMu());
  398.         }
  399.         setStartDate(state.getDate());
  400.     }

  402.     /** Get the names of the parameters in the matrix returned by {@link MatricesHarvester#getParametersJacobian}.
  403.      * @return names of the parameters (i.e. columns) of the Jacobian matrix
  404.      */
  405.     List<String> getJacobiansColumnsNames() {
  406.         final List<String> columnsNames = new ArrayList<>();
  407.         for (final ForceModel forceModel : getAllForceModels()) {
  408.             for (final ParameterDriver driver : forceModel.getParametersDrivers()) {
  409.                 if (driver.isSelected() && !columnsNames.contains(driver.getName())) {
  410.                     columnsNames.add(driver.getName());
  411.                 }
  412.             }
  413.         }
  414.         Collections.sort(columnsNames);
  415.         return columnsNames;
  416.     }

  418.     /** {@inheritDoc} */
  419.     @Override
  420.     protected AbstractMatricesHarvester createHarvester(final String stmName, final RealMatrix initialStm,
  421.                                                         final DoubleArrayDictionary initialJacobianColumns) {
  422.         return new NumericalPropagationHarvester(this, stmName, initialStm, initialJacobianColumns);
  423.     }

  425.     /** {@inheritDoc} */
  426.     @Override
  427.     protected void setUpStmAndJacobianGenerators() {

  429.         final AbstractMatricesHarvester harvester = getHarvester();
  430.         if (harvester != null) {

  432.             // set up the additional equations and additional state providers
  433.             final StateTransitionMatrixGenerator stmGenerator = setUpStmGenerator();
  434.             final List<String> triggersDates = setUpTriggerDatesJacobiansColumns(stmGenerator.getName());
  435.             setUpRegularParametersJacobiansColumns(stmGenerator, triggersDates);

  437.             // as we are now starting the propagation, everything is configured
  438.             // we can freeze the names in the harvester
  439.             harvester.freezeColumnsNames();

  441.         }

  443.     }

  445.     /** Set up the State Transition Matrix Generator.
  446.      * @return State Transition Matrix Generator
  447.      * @since 11.1
  448.      */
  449.     private StateTransitionMatrixGenerator setUpStmGenerator() {

  451.         final AbstractMatricesHarvester harvester = getHarvester();

  453.         // add the STM generator corresponding to the current settings, and setup state accordingly
  454.         StateTransitionMatrixGenerator stmGenerator = null;
  455.         for (final AdditionalDerivativesProvider equations : getAdditionalDerivativesProviders()) {
  456.             if (equations instanceof StateTransitionMatrixGenerator &&
  457.                 equations.getName().equals(harvester.getStmName())) {
  458.                 // the STM generator has already been set up in a previous propagation
  459.                 stmGenerator = (StateTransitionMatrixGenerator) equations;
  460.                 break;
  461.             }
  462.         }
  463.         if (stmGenerator == null) {
  464.             // this is the first time we need the STM generate, create it
  465.             stmGenerator = new StateTransitionMatrixGenerator(harvester.getStmName(), getAllForceModels(), getAttitudeProvider());
  466.             addAdditionalDerivativesProvider(stmGenerator);
  467.         }

  469.         if (!getInitialIntegrationState().hasAdditionalState(harvester.getStmName())) {
  470.             // add the initial State Transition Matrix if it is not already there
  471.             // (perhaps due to a previous propagation)
  472.             setInitialState(stmGenerator.setInitialStateTransitionMatrix(getInitialState(),
  473.                                                                          harvester.getInitialStateTransitionMatrix(),
  474.                                                                          getOrbitType(),
  475.                                                                          getPositionAngleType()));
  476.         }

  478.         return stmGenerator;

  480.     }

  482.     /** Set up the Jacobians columns generator dedicated to trigger dates.
  483.      * @param stmName name of the State Transition Matrix state
  484.      * @return names of the columns corresponding to trigger dates
  485.      * @since 11.1
  486.      */
  487.     private List<String> setUpTriggerDatesJacobiansColumns(final String stmName) {

  489.         final List<String> names = new ArrayList<>();
  490.         for (final ForceModel forceModel : getAllForceModels()) {
  491.             if (forceModel instanceof Maneuver) {
  492.                 final Maneuver         maneuver         = (Maneuver) forceModel;
  493.                 final ManeuverTriggers maneuverTriggers = maneuver.getManeuverTriggers();
  494.                 if (maneuverTriggers instanceof AbstractManeuverTriggers) {

  496.                     // FIXME: when issue https://gitlab.orekit.org/orekit/orekit/-/issues/854 is solved
  497.                     // the previous if statement and the following cast should be removed as the following
  498.                     // code should really be done for all ManeuverTriggers and not only AbstractManeuverTriggers
  499.                     final AbstractManeuverTriggers amt = (AbstractManeuverTriggers) maneuverTriggers;

  501.                     amt.getEventsDetectors().
  502.                         filter(d -> d instanceof ParameterDrivenDateIntervalDetector).
  503.                         map (d -> (ParameterDrivenDateIntervalDetector) d).
  504.                         forEach(d -> {
  505.                             if (d.getStartDriver().isSelected() || d.getMedianDriver().isSelected() || d.getDurationDriver().isSelected()) {
  506.                                 final TriggerDate start =
  507.                                                 manageTriggerDate(stmName, maneuver, amt, d.getStartDriver().getName(), true,  d.getThreshold());
  508.                                 names.add(start.getName());
  509.                             }
  510.                             if (d.getStopDriver().isSelected() || d.getMedianDriver().isSelected() || d.getDurationDriver().isSelected()) {
  511.                                 final TriggerDate stop =
  512.                                                 manageTriggerDate(stmName, maneuver, amt, d.getStopDriver().getName(),  false, d.getThreshold());
  513.                                 names.add(stop.getName());
  514.                             }
  515.                             if (d.getMedianDriver().isSelected()) {
  516.                                 final MedianDate median =
  517.                                                 manageMedianDate(d.getStartDriver().getName(), d.getStopDriver().getName(), d.getMedianDriver().getName());
  518.                                 names.add(median.getName());
  519.                             }
  520.                             if (d.getDurationDriver().isSelected()) {
  521.                                 final Duration duration =
  522.                                                 manageManeuverDuration(d.getStartDriver().getName(), d.getStopDriver().getName(), d.getDurationDriver().getName());
  523.                                 names.add(duration.getName());
  524.                             }
  525.                         });

  527.                 }
  528.             }
  529.         }

  531.         return names;

  533.     }

  535.     /** Manage a maneuver trigger date.
  536.      * @param stmName name of the State Transition Matrix state
  537.      * @param maneuver maneuver force model
  538.      * @param amt trigger to which the driver is bound
  539.      * @param driverName name of the date driver
  540.      * @param start if true, the driver is a maneuver start
  541.      * @param threshold event detector threshold
  542.      * @return generator for the date driver
  543.      * @since 11.1
  544.      */
  545.     private TriggerDate manageTriggerDate(final String stmName,
  546.                                           final Maneuver maneuver,
  547.                                           final AbstractManeuverTriggers amt,
  548.                                           final String driverName,
  549.                                           final boolean start,
  550.                                           final double threshold) {

  552.         TriggerDate triggerGenerator = null;

  554.         // check if we already have set up the provider
  555.         for (final AdditionalStateProvider provider : getAdditionalStateProviders()) {
  556.             if (provider instanceof TriggerDate &&
  557.                 provider.getName().equals(driverName)) {
  558.                 // the Jacobian column generator has already been set up in a previous propagation
  559.                 triggerGenerator = (TriggerDate) provider;
  560.                 break;
  561.             }
  562.         }

  564.         if (triggerGenerator == null) {
  565.             // this is the first time we need the Jacobian column generator, create it
  566.             triggerGenerator = new TriggerDate(stmName, driverName, start, maneuver, threshold);
  567.             amt.addResetter(triggerGenerator);
  568.             addAdditionalDerivativesProvider(triggerGenerator.getMassDepletionDelay());
  569.             addAdditionalStateProvider(triggerGenerator);
  570.         }

  572.         if (!getInitialIntegrationState().hasAdditionalState(driverName)) {
  573.             // add the initial Jacobian column if it is not already there
  574.             // (perhaps due to a previous propagation)
  575.             setInitialColumn(triggerGenerator.getMassDepletionDelay().getName(), new double[6]);
  576.             setInitialColumn(driverName, getHarvester().getInitialJacobianColumn(driverName));
  577.         }

  579.         return triggerGenerator;

  581.     }

  583.     /** Manage a maneuver median date.
  584.      * @param startName name of the start driver
  585.      * @param stopName name of the stop driver
  586.      * @param medianName name of the median driver
  587.      * @return generator for the median driver
  588.      * @since 11.1
  589.      */
  590.     private MedianDate manageMedianDate(final String startName, final String stopName, final String medianName) {

  592.         MedianDate medianGenerator = null;

  594.         // check if we already have set up the provider
  595.         for (final AdditionalStateProvider provider : getAdditionalStateProviders()) {
  596.             if (provider instanceof MedianDate &&
  597.                 provider.getName().equals(medianName)) {
  598.                 // the Jacobian column generator has already been set up in a previous propagation
  599.                 medianGenerator = (MedianDate) provider;
  600.                 break;
  601.             }
  602.         }

  604.         if (medianGenerator == null) {
  605.             // this is the first time we need the Jacobian column generator, create it
  606.             medianGenerator = new MedianDate(startName, stopName, medianName);
  607.             addAdditionalStateProvider(medianGenerator);
  608.         }

  610.         if (!getInitialIntegrationState().hasAdditionalState(medianName)) {
  611.             // add the initial Jacobian column if it is not already there
  612.             // (perhaps due to a previous propagation)
  613.             setInitialColumn(medianName, getHarvester().getInitialJacobianColumn(medianName));
  614.         }

  616.         return medianGenerator;

  618.     }

  620.     /** Manage a maneuver duration.
  621.      * @param startName name of the start driver
  622.      * @param stopName name of the stop driver
  623.      * @param durationName name of the duration driver
  624.      * @return generator for the median driver
  625.      * @since 11.1
  626.      */
  627.     private Duration manageManeuverDuration(final String startName, final String stopName, final String durationName) {

  629.         Duration durationGenerator = null;

  631.         // check if we already have set up the provider
  632.         for (final AdditionalStateProvider provider : getAdditionalStateProviders()) {
  633.             if (provider instanceof Duration &&
  634.                 provider.getName().equals(durationName)) {
  635.                 // the Jacobian column generator has already been set up in a previous propagation
  636.                 durationGenerator = (Duration) provider;
  637.                 break;
  638.             }
  639.         }

  641.         if (durationGenerator == null) {
  642.             // this is the first time we need the Jacobian column generator, create it
  643.             durationGenerator = new Duration(startName, stopName, durationName);
  644.             addAdditionalStateProvider(durationGenerator);
  645.         }

  647.         if (!getInitialIntegrationState().hasAdditionalState(durationName)) {
  648.             // add the initial Jacobian column if it is not already there
  649.             // (perhaps due to a previous propagation)
  650.             setInitialColumn(durationName, getHarvester().getInitialJacobianColumn(durationName));
  651.         }

  653.         return durationGenerator;

  655.     }

  657.     /** Set up the Jacobians columns generator for regular parameters.
  658.      * @param stmGenerator generator for the State Transition Matrix
  659.      * @param triggerDates names of the columns already managed as trigger dates
  660.      * @since 11.1
  661.      */
  662.     private void setUpRegularParametersJacobiansColumns(final StateTransitionMatrixGenerator stmGenerator,
  663.                                                         final List<String> triggerDates) {

  665.         // first pass: gather all parameters (excluding trigger dates), binding similar names together
  666.         final ParameterDriversList selected = new ParameterDriversList();
  667.         for (final ForceModel forceModel : getAllForceModels()) {
  668.             for (final ParameterDriver driver : forceModel.getParametersDrivers()) {
  669.                 if (!triggerDates.contains(driver.getName())) {
  670.                     selected.add(driver);
  671.                 }
  672.             }
  673.         }

  675.         // second pass: now that shared parameter names are bound together,
  676.         // their selections status have been synchronized, we can filter them
  677.         selected.filter(true);

  679.         // third pass: sort parameters lexicographically
  680.         selected.sort();

  682.         // add the Jacobians column generators corresponding to parameters, and setup state accordingly
  683.         for (final DelegatingDriver driver : selected.getDrivers()) {

  685.             IntegrableJacobianColumnGenerator generator = null;

  687.             // check if we already have set up the providers
  688.             for (final AdditionalDerivativesProvider provider : getAdditionalDerivativesProviders()) {
  689.                 if (provider instanceof IntegrableJacobianColumnGenerator &&
  690.                     provider.getName().equals(driver.getName())) {
  691.                     // the Jacobian column generator has already been set up in a previous propagation
  692.                     generator = (IntegrableJacobianColumnGenerator) provider;
  693.                     break;
  694.                 }
  695.             }

  697.             if (generator == null) {
  698.                 // this is the first time we need the Jacobian column generator, create it
  699.                 generator = new IntegrableJacobianColumnGenerator(stmGenerator, driver.getName());
  700.                 addAdditionalDerivativesProvider(generator);
  701.             }

  703.             if (!getInitialIntegrationState().hasAdditionalState(driver.getName())) {
  704.                 // add the initial Jacobian column if it is not already there
  705.                 // (perhaps due to a previous propagation)
  706.                 setInitialColumn(driver.getName(), getHarvester().getInitialJacobianColumn(driver.getName()));
  707.             }

  709.         }

  711.     }

  713.     /** Add the initial value of the column to the initial state.
  714.      * <p>
  715.      * The initial state must already contain the Cartesian State Transition Matrix.
  716.      * </p>
  717.      * @param columnName name of the column
  718.      * @param dYdQ column of the Jacobian ∂Y/∂qₘ with respect to propagation type,
  719.      * if null (which is the most frequent case), assumed to be 0
  720.      * @since 11.1
  721.      */
  722.     private void setInitialColumn(final String columnName, final double[] dYdQ) {

  724.         final SpacecraftState state = getInitialState();

  726.         if (dYdQ.length != STATE_DIMENSION) {
  727.             throw new OrekitException(LocalizedCoreFormats.DIMENSIONS_MISMATCH,
  728.                                       dYdQ.length, STATE_DIMENSION);
  729.         }

  731.         // convert to Cartesian Jacobian
  732.         final double[][] dYdC = new double[STATE_DIMENSION][STATE_DIMENSION];
  733.         getOrbitType().convertType(state.getOrbit()).getJacobianWrtCartesian(getPositionAngleType(), dYdC);
  734.         final double[] column = new QRDecomposition(MatrixUtils.createRealMatrix(dYdC), THRESHOLD).
  735.                         getSolver().
  736.                         solve(MatrixUtils.createRealVector(dYdQ)).
  737.                         toArray();

  739.         // set additional state
  740.         setInitialState(state.addAdditionalState(columnName, column));

  742.     }

  744.     /** {@inheritDoc} */
  745.     @Override
  746.     public TimeStampedPVCoordinates getPVCoordinates(final AbsoluteDate date, final Frame frame) {
  747.         return propagate(date).getPVCoordinates(frame);
  748.     }

  750.     /** {@inheritDoc} */
  751.     @Override
  752.     protected StateMapper createMapper(final AbsoluteDate referenceDate, final double mu,
  753.                                        final OrbitType orbitType, final PositionAngle positionAngleType,
  754.                                        final AttitudeProvider attitudeProvider, final Frame frame) {
  755.         return new OsculatingMapper(referenceDate, mu, orbitType, positionAngleType, attitudeProvider, frame);
  756.     }

  758.     /** Internal mapper using directly osculating parameters. */
  759.     private static class OsculatingMapper extends StateMapper {

  761.         /** Simple constructor.
  762.          * <p>
  763.          * The position parameter type is meaningful only if {@link
  764.          * #getOrbitType() propagation orbit type}
  765.          * support it. As an example, it is not meaningful for propagation
  766.          * in {@link OrbitType#CARTESIAN Cartesian} parameters.
  767.          * </p>
  768.          * @param referenceDate reference date
  769.          * @param mu central attraction coefficient (m³/s²)
  770.          * @param orbitType orbit type to use for mapping (can be null for {@link AbsolutePVCoordinates})
  771.          * @param positionAngleType angle type to use for propagation
  772.          * @param attitudeProvider attitude provider
  773.          * @param frame inertial frame
  774.          */
  775.         OsculatingMapper(final AbsoluteDate referenceDate, final double mu,
  776.                          final OrbitType orbitType, final PositionAngle positionAngleType,
  777.                          final AttitudeProvider attitudeProvider, final Frame frame) {
  778.             super(referenceDate, mu, orbitType, positionAngleType, attitudeProvider, frame);
  779.         }

  781.         /** {@inheritDoc} */
  782.         public SpacecraftState mapArrayToState(final AbsoluteDate date, final double[] y, final double[] yDot,
  783.                                                final PropagationType type) {
  784.             // the parameter type is ignored for the Numerical Propagator

  786.             final double mass = y[6];
  787.             if (mass <= 0.0) {
  788.                 throw new OrekitException(OrekitMessages.SPACECRAFT_MASS_BECOMES_NEGATIVE, mass);
  789.             }

  791.             if (getOrbitType() == null) {
  792.                 // propagation uses absolute position-velocity-acceleration
  793.                 final Vector3D p = new Vector3D(y[0],    y[1],    y[2]);
  794.                 final Vector3D v = new Vector3D(y[3],    y[4],    y[5]);
  795.                 final Vector3D a;
  796.                 final AbsolutePVCoordinates absPva;
  797.                 if (yDot == null) {
  798.                     absPva = new AbsolutePVCoordinates(getFrame(), new TimeStampedPVCoordinates(date, p, v));
  799.                 } else {
  800.                     a = new Vector3D(yDot[3], yDot[4], yDot[5]);
  801.                     absPva = new AbsolutePVCoordinates(getFrame(), new TimeStampedPVCoordinates(date, p, v, a));
  802.                 }

  804.                 final Attitude attitude = getAttitudeProvider().getAttitude(absPva, date, getFrame());
  805.                 return new SpacecraftState(absPva, attitude, mass);
  806.             } else {
  807.                 // propagation uses regular orbits
  808.                 final Orbit orbit       = getOrbitType().mapArrayToOrbit(y, yDot, getPositionAngleType(), date, getMu(), getFrame());
  809.                 final Attitude attitude = getAttitudeProvider().getAttitude(orbit, date, getFrame());

  811.                 return new SpacecraftState(orbit, attitude, mass);
  812.             }

  814.         }

  816.         /** {@inheritDoc} */
  817.         public void mapStateToArray(final SpacecraftState state, final double[] y, final double[] yDot) {
  818.             if (getOrbitType() == null) {
  819.                 // propagation uses absolute position-velocity-acceleration
  820.                 final Vector3D p = state.getAbsPVA().getPosition();
  821.                 final Vector3D v = state.getAbsPVA().getVelocity();
  822.                 y[0] = p.getX();
  823.                 y[1] = p.getY();
  824.                 y[2] = p.getZ();
  825.                 y[3] = v.getX();
  826.                 y[4] = v.getY();
  827.                 y[5] = v.getZ();
  828.                 y[6] = state.getMass();
  829.             }
  830.             else {
  831.                 getOrbitType().mapOrbitToArray(state.getOrbit(), getPositionAngleType(), y, yDot);
  832.                 y[6] = state.getMass();
  833.             }
  834.         }

  836.     }

  838.     /** {@inheritDoc} */
  839.     protected MainStateEquations getMainStateEquations(final ODEIntegrator integrator) {
  840.         return new Main(integrator);
  841.     }

  843.     /** Internal class for osculating parameters integration. */
  844.     private class Main implements MainStateEquations, TimeDerivativesEquations {

  846.         /** Derivatives array. */
  847.         private final double[] yDot;

  849.         /** Current state. */
  850.         private SpacecraftState currentState;

  852.         /** Jacobian of the orbital parameters with respect to the Cartesian parameters. */
  853.         private double[][] jacobian;

  855.         /** Simple constructor.
  856.          * @param integrator numerical integrator to use for propagation.
  857.          */
  858.         Main(final ODEIntegrator integrator) {

  860.             this.yDot     = new double[7];
  861.             this.jacobian = new double[6][6];

  863.             for (final ForceModel forceModel : forceModels) {
  864.                 forceModel.getEventsDetectors().forEach(detector -> setUpEventDetector(integrator, detector));
  865.             }

  867.             if (getOrbitType() == null) {
  868.                 // propagation uses absolute position-velocity-acceleration
  869.                 // we can set Jacobian once and for all
  870.                 for (int i = 0; i < jacobian.length; ++i) {
  871.                     Arrays.fill(jacobian[i], 0.0);
  872.                     jacobian[i][i] = 1.0;
  873.                 }
  874.             }

  876.         }

  878.         /** {@inheritDoc} */
  879.         @Override
  880.         public void init(final SpacecraftState initialState, final AbsoluteDate target) {
  881.             forceModels.forEach(fm -> fm.init(initialState, target));
  882.         }

  884.         /** {@inheritDoc} */
  885.         @Override
  886.         public double[] computeDerivatives(final SpacecraftState state) {

  888.             currentState = state;
  889.             Arrays.fill(yDot, 0.0);
  890.             if (getOrbitType() != null) {
  891.                 // propagation uses regular orbits
  892.                 currentState.getOrbit().getJacobianWrtCartesian(getPositionAngleType(), jacobian);
  893.             }

  895.             // compute the contributions of all perturbing forces,
  896.             // using the Kepler contribution at the end since
  897.             // NewtonianAttraction is always the last instance in the list
  898.             for (final ForceModel forceModel : forceModels) {
  899.                 forceModel.addContribution(state, this);
  900.             }

  902.             if (getOrbitType() == null) {
  903.                 // position derivative is velocity, and was not added above in the force models
  904.                 // (it is added when orbit type is non-null because NewtonianAttraction considers it)
  905.                 final Vector3D velocity = currentState.getPVCoordinates().getVelocity();
  906.                 yDot[0] += velocity.getX();
  907.                 yDot[1] += velocity.getY();
  908.                 yDot[2] += velocity.getZ();
  909.             }


  912.             return yDot.clone();

  914.         }

  916.         /** {@inheritDoc} */
  917.         @Override
  918.         public void addKeplerContribution(final double mu) {
  919.             if (getOrbitType() == null) {

  921.                 // if mu is neither 0 nor NaN, we want to include Newtonian acceleration
  922.                 if (mu > 0) {
  923.                     // velocity derivative is Newtonian acceleration
  924.                     final Vector3D position = currentState.getPVCoordinates().getPosition();
  925.                     final double r2         = position.getNormSq();
  926.                     final double coeff      = -mu / (r2 * FastMath.sqrt(r2));
  927.                     yDot[3] += coeff * position.getX();
  928.                     yDot[4] += coeff * position.getY();
  929.                     yDot[5] += coeff * position.getZ();
  930.                 }

  932.             } else {
  933.                 // propagation uses regular orbits
  934.                 currentState.getOrbit().addKeplerContribution(getPositionAngleType(), mu, yDot);
  935.             }
  936.         }

  938.         /** {@inheritDoc} */
  939.         public void addNonKeplerianAcceleration(final Vector3D gamma) {
  940.             for (int i = 0; i < 6; ++i) {
  941.                 final double[] jRow = jacobian[i];
  942.                 yDot[i] += jRow[3] * gamma.getX() + jRow[4] * gamma.getY() + jRow[5] * gamma.getZ();
  943.             }
  944.         }

  946.         /** {@inheritDoc} */
  947.         @Override
  948.         public void addMassDerivative(final double q) {
  949.             if (q > 0) {
  950.                 throw new OrekitIllegalArgumentException(OrekitMessages.POSITIVE_FLOW_RATE, q);
  951.             }
  952.             yDot[6] += q;
  953.         }

  955.     }

  957.     /** Estimate tolerance vectors for integrators when propagating in absolute position-velocity-acceleration.
  958.      * @param dP user specified position error
  959.      * @param absPva reference absolute position-velocity-acceleration
  960.      * @return a two rows array, row 0 being the absolute tolerance error and row 1
  961.      * being the relative tolerance error
  962.      * @see NumericalPropagator#tolerances(double, Orbit, OrbitType)
  963.      */
  964.     public static double[][] tolerances(final double dP, final AbsolutePVCoordinates absPva) {

  966.         final double relative = dP / absPva.getPosition().getNorm();
  967.         final double dV = relative * absPva.getVelocity().getNorm();

  969.         final double[] absTol = new double[7];
  970.         final double[] relTol = new double[7];

  972.         absTol[0] = dP;
  973.         absTol[1] = dP;
  974.         absTol[2] = dP;
  975.         absTol[3] = dV;
  976.         absTol[4] = dV;
  977.         absTol[5] = dV;

  979.         // we set the mass tolerance arbitrarily to 1.0e-6 kg, as mass evolves linearly
  980.         // with trust, this often has no influence at all on propagation
  981.         absTol[6] = 1.0e-6;

  983.         Arrays.fill(relTol, relative);

  985.         return new double[][] {
  986.             absTol, relTol
  987.         };

  989.     }

  991.     /** Estimate tolerance vectors for integrators when propagating in orbits.
  992.      * <p>
  993.      * The errors are estimated from partial derivatives properties of orbits,
  994.      * starting from a scalar position error specified by the user.
  995.      * Considering the energy conservation equation V = sqrt(mu (2/r - 1/a)),
  996.      * we get at constant energy (i.e. on a Keplerian trajectory):
  997.      * <pre>
  998.      * V r² |dV| = mu |dr|
  999.      * </pre>
  1000.      * <p> So we deduce a scalar velocity error consistent with the position error.
  1001.      * From here, we apply orbits Jacobians matrices to get consistent errors
  1002.      * on orbital parameters.
  1003.      *
  1004.      * <p>
  1005.      * The tolerances are only <em>orders of magnitude</em>, and integrator tolerances
  1006.      * are only local estimates, not global ones. So some care must be taken when using
  1007.      * these tolerances. Setting 1mm as a position error does NOT mean the tolerances
  1008.      * will guarantee a 1mm error position after several orbits integration.
  1009.      * </p>
  1010.      * @param dP user specified position error
  1011.      * @param orbit reference orbit
  1012.      * @param type propagation type for the meaning of the tolerance vectors elements
  1013.      * (it may be different from {@code orbit.getType()})
  1014.      * @return a two rows array, row 0 being the absolute tolerance error and row 1
  1015.      * being the relative tolerance error
  1016.      */
  1017.     public static double[][] tolerances(final double dP, final Orbit orbit, final OrbitType type) {

  1019.         // estimate the scalar velocity error
  1020.         final PVCoordinates pv = orbit.getPVCoordinates();
  1021.         final double r2 = pv.getPosition().getNormSq();
  1022.         final double v  = pv.getVelocity().getNorm();
  1023.         final double dV = orbit.getMu() * dP / (v * r2);

  1025.         return tolerances(dP, dV, orbit, type);

  1027.     }

  1029.     /** Estimate tolerance vectors for integrators when propagating in orbits.
  1030.      * <p>
  1031.      * The errors are estimated from partial derivatives properties of orbits,
  1032.      * starting from scalar position and velocity errors specified by the user.
  1033.      * <p>
  1034.      * The tolerances are only <em>orders of magnitude</em>, and integrator tolerances
  1035.      * are only local estimates, not global ones. So some care must be taken when using
  1036.      * these tolerances. Setting 1mm as a position error does NOT mean the tolerances
  1037.      * will guarantee a 1mm error position after several orbits integration.
  1038.      * </p>
  1039.      * @param dP user specified position error
  1040.      * @param dV user specified velocity error
  1041.      * @param orbit reference orbit
  1042.      * @param type propagation type for the meaning of the tolerance vectors elements
  1043.      * (it may be different from {@code orbit.getType()})
  1044.      * @return a two rows array, row 0 being the absolute tolerance error and row 1
  1045.      * being the relative tolerance error
  1046.      * @since 10.3
  1047.      */
  1048.     public static double[][] tolerances(final double dP, final double dV,
  1049.                                         final Orbit orbit, final OrbitType type) {

  1051.         final double[] absTol = new double[7];
  1052.         final double[] relTol = new double[7];

  1054.         // we set the mass tolerance arbitrarily to 1.0e-6 kg, as mass evolves linearly
  1055.         // with trust, this often has no influence at all on propagation
  1056.         absTol[6] = 1.0e-6;

  1058.         if (type == OrbitType.CARTESIAN) {
  1059.             absTol[0] = dP;
  1060.             absTol[1] = dP;
  1061.             absTol[2] = dP;
  1062.             absTol[3] = dV;
  1063.             absTol[4] = dV;
  1064.             absTol[5] = dV;
  1065.         } else {

  1067.             // convert the orbit to the desired type
  1068.             final double[][] jacobian = new double[6][6];
  1069.             final Orbit converted = type.convertType(orbit);
  1070.             converted.getJacobianWrtCartesian(PositionAngle.TRUE, jacobian);

  1072.             for (int i = 0; i < 6; ++i) {
  1073.                 final double[] row = jacobian[i];
  1074.                 absTol[i] = FastMath.abs(row[0]) * dP +
  1075.                             FastMath.abs(row[1]) * dP +
  1076.                             FastMath.abs(row[2]) * dP +
  1077.                             FastMath.abs(row[3]) * dV +
  1078.                             FastMath.abs(row[4]) * dV +
  1079.                             FastMath.abs(row[5]) * dV;
  1080.                 if (Double.isNaN(absTol[i])) {
  1081.                     throw new OrekitException(OrekitMessages.SINGULAR_JACOBIAN_FOR_ORBIT_TYPE, type);
  1082.                 }
  1083.             }

  1085.         }

  1087.         Arrays.fill(relTol, dP / FastMath.sqrt(orbit.getPVCoordinates().getPosition().getNormSq()));

  1089.         return new double[][] {
  1090.             absTol, relTol
  1091.         };

  1093.     }

  1095.     /** {@inheritDoc} */
  1096.     @Override
  1097.     protected void beforeIntegration(final SpacecraftState initialState, final AbsoluteDate tEnd) {

  1099.         if (!getFrame().isPseudoInertial()) {

  1101.             // inspect all force models to find InertialForces
  1102.             for (ForceModel force : forceModels) {
  1103.                 if (force instanceof InertialForces) {
  1104.                     return;
  1105.                 }
  1106.             }

  1108.             // throw exception if no inertial forces found
  1109.             throw new OrekitException(OrekitMessages.INERTIAL_FORCE_MODEL_MISSING, getFrame().getName());

  1111.         }

  1113.     }

  1115. }
Created with JaCoCo 0.8.8.202204050719