HansenZonalLinear.java

  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.propagation.semianalytical.dsst.utilities.hansen;

  19. import org.hipparchus.analysis.polynomials.PolynomialFunction;
  20. import org.hipparchus.util.FastMath;

  22. /**
  23.  * Hansen coefficients K(t,n,s) for t=0 and n < 0.
  24.  * <p>
  25.  *Implements Collins 4-242 or echivalently, Danielson 2.7.3-(6) for Hansen Coefficients and
  26.  * Collins 4-245 or Danielson 3.1-(7) for derivatives. The recursions are transformed into
  27.  * composition of linear transformations to obtain the associated polynomials
  28.  * for coefficients and their derivatives - see Petre's paper
  29.  *
  30.  * @author Petre Bazavan
  31.  * @author Lucian Barbulescu
  32.  */
  33. public class HansenZonalLinear {

  35.     /** The number of coefficients that will be computed with a set of roots. */
  36.     private static final int SLICE = 10;

  38.     /**
  39.      * The first vector of polynomials associated to Hansen coefficients and
  40.      * derivatives.
  41.      */
  42.     private final PolynomialFunction[][] mpvec;

  44.     /** The second vector of polynomials associated only to derivatives. */
  45.     private final PolynomialFunction[][] mpvecDeriv;

  47.     /** The Hansen coefficients used as roots. */
  48.     private final double[][] hansenRoot;

  50.     /** The derivatives of the Hansen coefficients used as roots. */
  51.     private final double[][] hansenDerivRoot;

  53.     /** The s coefficient. */
  54.     private final int s;

  56.     /**
  57.      * The offset used to identify the polynomial that corresponds to a negative
  58.      * value of n in the internal array that starts at 0.
  59.      */
  60.     private final int offset;

  62.     /** The number of slices needed to compute the coefficients. */
  63.     private final int numSlices;

  65.     /** 2<sup>s</sup>. */
  66.     private final double twots;

  68.     /** 2*s+1. */
  69.     private final int twosp1;

  71.     /** 2*s. */
  72.     private final int twos;

  74.     /** (2*s+1) / 2<sup>s</sup>. */
  75.     private final double twosp1otwots;

  77.     /**
  78.      * Constructor.
  79.      *
  80.      * @param nMax the maximum (absolute) value of n coefficient
  81.      * @param s s coefficient
  82.      */
  83.     public HansenZonalLinear(final int nMax, final int s) {

  85.         //Initialize fields
  86.         final int Nmin = -nMax - 1;
  87.         final int N0 = -(s + 2);
  88.         this.offset = nMax + 1;
  89.         this.s = s;
  90.         this.twots = FastMath.pow(2., s);
  91.         this.twos = 2 * s;
  92.         this.twosp1 = this.twos + 1;
  93.         this.twosp1otwots = (double) this.twosp1 / this.twots;

  95.         // prepare structures for stored data
  96.         final int size = nMax - s - 1;
  97.         mpvec      = new PolynomialFunction[size][];
  98.         mpvecDeriv = new PolynomialFunction[size][];

  100.         this.numSlices  = FastMath.max((int) FastMath.ceil(((double) size) / SLICE), 1);
  101.         hansenRoot      = new double[numSlices][2];
  102.         hansenDerivRoot = new double[numSlices][2];

  104.         // Prepare the data base of associated polynomials
  105.         HansenUtilities.generateZonalPolynomials(N0, Nmin, offset, SLICE, s,
  106.                                                  mpvec, mpvecDeriv);

  108.     }

  110.     /**
  111.      * Compute the roots for the Hansen coefficients and their derivatives.
  112.      *
  113.      * @param chi 1 / sqrt(1 - e²)
  114.      */
  115.     public void computeInitValues(final double chi) {
  116.         // compute the values for n=s and n=s+1
  117.         // See Danielson 2.7.3-(6a,b)
  118.         hansenRoot[0][0] = 0;
  119.         hansenRoot[0][1] = FastMath.pow(chi, this.twosp1) / this.twots;
  120.         hansenDerivRoot[0][0] = 0;
  121.         hansenDerivRoot[0][1] = this.twosp1otwots * FastMath.pow(chi, this.twos);

  123.         final int st = -s - 1;
  124.         for (int i = 1; i < numSlices; i++) {
  125.             for (int j = 0; j < 2; j++) {
  126.                 // Get the required polynomials
  127.                 final PolynomialFunction[] mv = mpvec[st - (i * SLICE) - j + offset];
  128.                 final PolynomialFunction[] sv = mpvecDeriv[st - (i * SLICE) - j + offset];

  130.                 //Compute the root derivatives
  131.                 hansenDerivRoot[i][j] = mv[1].value(chi) * hansenDerivRoot[i - 1][1] +
  132.                                         mv[0].value(chi) * hansenDerivRoot[i - 1][0] +
  133.                                         (sv[1].value(chi) * hansenRoot[i - 1][1] +
  134.                                          sv[0].value(chi) * hansenRoot[i - 1][0]
  135.                                         ) / chi;
  136.                 hansenRoot[i][j] =     mv[1].value(chi) * hansenRoot[i - 1][1] +
  137.                                        mv[0].value(chi) * hansenRoot[i - 1][0];

  139.             }

  141.         }
  142.     }

  144.     /**
  145.      * Get the K₀<sup>-n-1,s</sup> coefficient value.
  146.      *
  147.      * <p> The s value is given in the class constructor
  148.      *
  149.      * @param mnm1 (-n-1) coefficient
  150.      * @param chi The value of χ
  151.      * @return K₀<sup>-n-1,s</sup>
  152.      */
  153.     public double getValue(final int mnm1, final double chi) {

  155.         //Compute n
  156.         final int n = -mnm1 - 1;

  158.         //Compute the potential slice
  159.         int sliceNo = (n - s) / SLICE;
  160.         if (sliceNo < numSlices) {
  161.             //Compute the index within the slice
  162.             final int indexInSlice = (n - s) % SLICE;

  164.             //Check if a root must be returned
  165.             if (indexInSlice <= 1) {
  166.                 return hansenRoot[sliceNo][indexInSlice];
  167.             }
  168.         } else {
  169.             //the value was a potential root for a slice, but that slice was not required
  170.             //Decrease the slice number
  171.             sliceNo--;
  172.         }

  174.         // Danielson 2.7.3-(6c)/Collins 4-242 and Petre's paper
  175.         final PolynomialFunction[] v = mpvec[mnm1 + offset];
  176.         double ret = v[1].value(chi) * hansenRoot[sliceNo][1];
  177.         if (hansenRoot[sliceNo][0] != 0) {
  178.             ret += v[0].value(chi) * hansenRoot[sliceNo][0];
  179.         }
  180.         return  ret;
  181.     }

  183.     /**
  184.      * Get the dK₀<sup>-n-1,s</sup> / d&Chi; coefficient derivative.
  185.      *
  186.      * <p> The s value is given in the class constructor.
  187.      *
  188.      * @param mnm1 (-n-1) coefficient
  189.      * @param chi The value of χ
  190.      * @return dK₀<sup>-n-1,s</sup> / d&Chi;
  191.      */
  192.     public double getDerivative(final int mnm1, final double chi) {

  194.         //Compute n
  195.         final int n = -mnm1 - 1;

  197.         //Compute the potential slice
  198.         int sliceNo = (n - s) / SLICE;
  199.         if (sliceNo < numSlices) {
  200.             //Compute the index within the slice
  201.             final int indexInSlice = (n - s) % SLICE;

  203.             //Check if a root must be returned
  204.             if (indexInSlice <= 1) {
  205.                 return hansenDerivRoot[sliceNo][indexInSlice];
  206.             }
  207.         } else {
  208.             //the value was a potential root for a slice, but that slice was not required
  209.             //Decrease the slice number
  210.             sliceNo--;
  211.         }

  213.         // Danielson 3.1-(7c) and Petre's paper
  214.         final PolynomialFunction[] v = mpvec[mnm1 + offset];
  215.         double ret = v[1].value(chi) * hansenDerivRoot[sliceNo][1];
  216.         if (hansenDerivRoot[sliceNo][0] != 0) {
  217.             ret += v[0].value(chi) * hansenDerivRoot[sliceNo][0];
  218.         }

  220.         // Danielson 2.7.3-(6b)
  221.         final PolynomialFunction[] v1 = mpvecDeriv[mnm1 + offset];
  222.         double hret = v1[1].value(chi) * hansenRoot[sliceNo][1];
  223.         if (hansenRoot[sliceNo][0] != 0) {
  224.             hret += v1[0].value(chi) * hansenRoot[sliceNo][0];
  225.         }
  226.         ret += hret / chi;

  228.         return ret;

  230.     }

  232. }
Created with JaCoCo 0.8.12.202403310830