JacobiansMapper.java

  1. /* Copyright 2002-2020 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;

  18. import org.hipparchus.linear.Array2DRowRealMatrix;
  19. import org.hipparchus.linear.DecompositionSolver;
  20. import org.hipparchus.linear.QRDecomposition;
  21. import org.hipparchus.linear.RealMatrix;
  22. import org.orekit.orbits.Orbit;
  23. import org.orekit.orbits.OrbitType;
  24. import org.orekit.orbits.PositionAngle;
  25. import org.orekit.propagation.SpacecraftState;
  26. import org.orekit.propagation.integration.AbstractJacobiansMapper;
  27. import org.orekit.utils.ParameterDriversList;

  28. /** Mapper between two-dimensional Jacobian matrices and one-dimensional {@link
  29.  * SpacecraftState#getAdditionalState(String) additional state arrays}.
  30.  * <p>
  31.  * This class does not hold the states by itself. Instances of this class are guaranteed
  32.  * to be immutable.
  33.  * </p>
  34.  * @author Luc Maisonobe
  35.  * @see org.orekit.propagation.numerical.PartialDerivativesEquations
  36.  * @see org.orekit.propagation.numerical.NumericalPropagator
  37.  * @see SpacecraftState#getAdditionalState(String)
  38.  * @see org.orekit.propagation.AbstractPropagator
  39.  */
  40. public class JacobiansMapper extends AbstractJacobiansMapper {

  41.     /** State dimension, fixed to 6.
  42.      * @since 9.0
  43.      */
  44.     public static final int STATE_DIMENSION = 6;

  45.     /** Selected parameters for Jacobian computation. */
  46.     private final ParameterDriversList parameters;

  47.     /** Name. */
  48.     private String name;

  49.     /** Orbit type. */
  50.     private final OrbitType orbitType;

  51.     /** Position angle type. */
  52.     private final PositionAngle angleType;

  53.     /** Simple constructor.
  54.      * @param name name of the Jacobians
  55.      * @param parameters selected parameters for Jacobian computation
  56.      * @param orbitType orbit type
  57.      * @param angleType position angle type
  58.      */
  59.     JacobiansMapper(final String name, final ParameterDriversList parameters,
  60.                     final OrbitType orbitType, final PositionAngle angleType) {
  61.         super(name, parameters);
  62.         this.orbitType  = orbitType;
  63.         this.angleType  = angleType;
  64.         this.parameters = parameters;
  65.         this.name = name;
  66.     }

  67.     /** {@inheritDoc} */
  68.     protected double[][] getConversionJacobian(final SpacecraftState state) {

  69.         final double[][] dYdC = new double[STATE_DIMENSION][STATE_DIMENSION];

  70.         // make sure the state is in the desired orbit type
  71.         final Orbit orbit = orbitType.convertType(state.getOrbit());

  72.         // compute the Jacobian, taking the position angle type into account
  73.         orbit.getJacobianWrtCartesian(angleType, dYdC);

  74.         return dYdC;

  75.     }

  76.     /** {@inheritDoc}
  77.      * <p>
  78.      * This method converts the Jacobians to Cartesian parameters and put the converted data
  79.      * in the one-dimensional {@code p} array.
  80.      * </p>
  81.      */
  82.     public void setInitialJacobians(final SpacecraftState state, final double[][] dY1dY0,
  83.                              final double[][] dY1dP, final double[] p) {

  84.         // set up a converter
  85.         final RealMatrix dY1dC1 = new Array2DRowRealMatrix(getConversionJacobian(state), false);
  86.         final DecompositionSolver solver = new QRDecomposition(dY1dC1).getSolver();

  87.         // convert the provided state Jacobian
  88.         final RealMatrix dC1dY0 = solver.solve(new Array2DRowRealMatrix(dY1dY0, false));

  89.         // map the converted state Jacobian to one-dimensional array
  90.         int index = 0;
  91.         for (int i = 0; i < STATE_DIMENSION; ++i) {
  92.             for (int j = 0; j < STATE_DIMENSION; ++j) {
  93.                 p[index++] = dC1dY0.getEntry(i, j);
  94.             }
  95.         }

  96.         if (parameters.getNbParams() != 0) {
  97.             // convert the provided state Jacobian
  98.             final RealMatrix dC1dP = solver.solve(new Array2DRowRealMatrix(dY1dP, false));

  99.             // map the converted parameters Jacobian to one-dimensional array
  100.             for (int i = 0; i < STATE_DIMENSION; ++i) {
  101.                 for (int j = 0; j < parameters.getNbParams(); ++j) {
  102.                     p[index++] = dC1dP.getEntry(i, j);
  103.                 }
  104.             }
  105.         }

  106.     }

  107.     /** {@inheritDoc} */
  108.     public void getStateJacobian(final SpacecraftState state,  final double[][] dYdY0) {

  109.         // get the conversion Jacobian
  110.         final double[][] dYdC = getConversionJacobian(state);

  111.         // extract the additional state
  112.         final double[] p = state.getAdditionalState(name);

  113.         // compute dYdY0 = dYdC * dCdY0, without allocating new arrays
  114.         for (int i = 0; i < STATE_DIMENSION; i++) {
  115.             final double[] rowC = dYdC[i];
  116.             final double[] rowD = dYdY0[i];
  117.             for (int j = 0; j < STATE_DIMENSION; ++j) {
  118.                 double sum = 0;
  119.                 int pIndex = j;
  120.                 for (int k = 0; k < STATE_DIMENSION; ++k) {
  121.                     sum += rowC[k] * p[pIndex];
  122.                     pIndex += STATE_DIMENSION;
  123.                 }
  124.                 rowD[j] = sum;
  125.             }
  126.         }

  127.     }

  128.     /** {@inheritDoc} */
  129.     public void getParametersJacobian(final SpacecraftState state, final double[][] dYdP) {

  130.         if (parameters.getNbParams() != 0) {

  131.             // get the conversion Jacobian
  132.             final double[][] dYdC = getConversionJacobian(state);

  133.             // extract the additional state
  134.             final double[] p = state.getAdditionalState(name);

  135.             // compute dYdP = dYdC * dCdP, without allocating new arrays
  136.             for (int i = 0; i < STATE_DIMENSION; i++) {
  137.                 final double[] rowC = dYdC[i];
  138.                 final double[] rowD = dYdP[i];
  139.                 for (int j = 0; j < parameters.getNbParams(); ++j) {
  140.                     double sum = 0;
  141.                     int pIndex = j + STATE_DIMENSION * STATE_DIMENSION;
  142.                     for (int k = 0; k < STATE_DIMENSION; ++k) {
  143.                         sum += rowC[k] * p[pIndex];
  144.                         pIndex += parameters.getNbParams();
  145.                     }
  146.                     rowD[j] = sum;
  147.                 }
  148.             }

  149.         }

  150.     }

  151.     /** {@inheritDoc} */
  152.     @Override
  153.     public int getAdditionalStateDimension() {
  154.         return STATE_DIMENSION * (STATE_DIMENSION + parameters.getNbParams());
  155.     }

  156. }