FieldGHmsjPolynomials.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;

  18. import org.hipparchus.CalculusFieldElement;
  19. import org.hipparchus.Field;
  20. import org.hipparchus.util.FastMath;

  21. /** Compute the G<sub>ms</sub><sup>j</sup> and the H<sub>ms</sub><sup>j</sup>
  22.  *  polynomials in the equinoctial elements h, k and the direction cosines α and β
  23.  *  and their partial derivatives with respect to k, h, α and β.
  24.  *  <p>
  25.  *  The expressions used are equations 2.7.5-(1)(2) from the Danielson paper.
  26.  *  </p>
  27.  *  @author Romain Di Costanzo
  28.  *  @author Bryan Cazabonne (field translation)
  29.  * @param <T> type of the field elements
  30.  */
  31. public class FieldGHmsjPolynomials <T extends CalculusFieldElement<T>> {
  32.     /** C<sub>j</sub>(k, h), S<sub>j</sub>(k, h) coefficient.
  33.      * (k, h) are the (x, y) component of the eccentricity vector in equinoctial elements
  34.      */
  35.     private final FieldCjSjCoefficient<T> cjsjKH;

  36.     /** C<sub>j</sub>(α, β), S<sub>j</sub>(α, β) coefficient.
  37.      * (α, β) are the direction cosines
  38.      */
  39.     private final FieldCjSjCoefficient<T> cjsjAB;

  40.     /** Is the orbit represented as a retrograde orbit.
  41.      *  I = -1 if yes, +1 otherwise.
  42.      */
  43.     private int                   I;

  44.     /** Create a set of G<sub>ms</sub><sup>j</sup> and H<sub>ms</sub><sup>j</sup> polynomials.
  45.      *  @param k X component of the eccentricity vector
  46.      *  @param h Y component of the eccentricity vector
  47.      *  @param alpha direction cosine α
  48.      *  @param beta direction cosine β
  49.      *  @param retroFactor -1 if the orbit is represented as retrograde, +1 otherwise
  50.      *  @param field field element
  51.      **/
  52.     public FieldGHmsjPolynomials(final T k, final T h,
  53.                             final T alpha, final T beta,
  54.                             final int retroFactor,
  55.                             final Field<T> field) {
  56.         this.cjsjKH = new FieldCjSjCoefficient<>(k, h, field);
  57.         this.cjsjAB = new FieldCjSjCoefficient<>(alpha, beta, field);
  58.         this.I = retroFactor;
  59.     }

  60.     /** Get the G<sub>ms</sub><sup>j</sup> coefficient.
  61.      * @param m m subscript
  62.      * @param s s subscript
  63.      * @param j order
  64.      * @return the G<sub>ms</sub><sup>j</sup>
  65.      */
  66.     public T getGmsj(final int m, final int s, final int j) {
  67.         final int sMj = FastMath.abs(s - j);
  68.         final T gms;
  69.         if (FastMath.abs(s) <= m) {
  70.             final int mMis = m - I * s;
  71.             gms = cjsjKH.getCj(sMj).multiply(cjsjAB.getCj(mMis)).
  72.                   subtract(cjsjKH.getSj(sMj).multiply(cjsjAB.getSj(mMis)).multiply(sgn(s - j)).multiply(I));
  73.         } else {
  74.             final int sMim = FastMath.abs(s - I * m);
  75.             gms = cjsjKH.getCj(sMj).multiply(cjsjAB.getCj(sMim)).
  76.                   add(cjsjKH.getSj(sMj).multiply(cjsjAB.getSj(sMim)).multiply(sgn(s - j)).multiply(sgn(s - m)));
  77.         }
  78.         return gms;
  79.     }

  80.     /** Get the H<sub>ms</sub><sup>j</sup> coefficient.
  81.      * @param m m subscript
  82.      * @param s s subscript
  83.      * @param j order
  84.      * @return the H<sub>ms</sub><sup>j</sup>
  85.      */
  86.     public T getHmsj(final int m, final int s, final int j) {
  87.         final int sMj = FastMath.abs(s - j);
  88.         final T hms;
  89.         if (FastMath.abs(s) <= m) {
  90.             final int mMis = m - I * s;
  91.             hms = cjsjKH.getCj(sMj).multiply(cjsjAB.getSj(mMis)).multiply(I).
  92.                             add(cjsjKH.getSj(sMj).multiply(cjsjAB.getCj(mMis)).multiply(sgn(s - j)));
  93.         } else {
  94.             final int sMim = FastMath.abs(s - I * m);
  95.             hms = cjsjKH.getCj(sMj).multiply(cjsjAB.getSj(sMim)).multiply(-sgn(s - m)).
  96.                   add(cjsjKH.getSj(sMj).multiply(cjsjAB.getCj(sMim)).multiply(sgn(s - j)));
  97.         }
  98.         return hms;
  99.     }

  100.     /** Get the dG<sub>ms</sub><sup>j</sup> / d<sub>k</sub> coefficient.
  101.      * @param m m subscript
  102.      * @param s s subscript
  103.      * @param j order
  104.      * @return dG<sub>ms</sub><sup>j</sup> / d<sub>k</sub>
  105.      */
  106.     public T getdGmsdk(final int m, final int s, final int j) {
  107.         final int sMj = FastMath.abs(s - j);
  108.         final T dGmsdk;
  109.         if (FastMath.abs(s) <= m) {
  110.             final int mMis = m - I * s;
  111.             dGmsdk = cjsjKH.getDcjDk(sMj).multiply(cjsjAB.getCj(mMis)).
  112.                      subtract(cjsjKH.getDsjDk(sMj).multiply(cjsjAB.getSj(mMis)).multiply(I).multiply(sgn(s - j)));
  113.         } else {
  114.             final int sMim = FastMath.abs(s - I * m);
  115.             dGmsdk = cjsjKH.getDcjDk(sMj).multiply(cjsjAB.getCj(sMim)).
  116.                      add(cjsjKH.getDsjDk(sMj).multiply(cjsjAB.getSj(sMim)).multiply(sgn(s - m)).multiply(sgn(s - j)));
  117.         }
  118.         return dGmsdk;
  119.     }

  120.     /** Get the dG<sub>ms</sub><sup>j</sup> / d<sub>h</sub> coefficient.
  121.      * @param m m subscript
  122.      * @param s s subscript
  123.      * @param j order
  124.      * @return dG<sub>ms</sub><sup>j</sup> / d<sub>h</sub>
  125.      */
  126.     public T getdGmsdh(final int m, final int s, final int j) {
  127.         final int sMj = FastMath.abs(s - j);
  128.         final T dGmsdh;
  129.         if (FastMath.abs(s) <= m) {
  130.             final int mMis = m - I * s;
  131.             dGmsdh = cjsjKH.getDcjDh(sMj).multiply(cjsjAB.getCj(mMis)).
  132.                      subtract(cjsjKH.getDsjDh(sMj).multiply(cjsjAB.getSj(mMis)).multiply(I).multiply(sgn(s - j)));
  133.         } else {
  134.             final int sMim = FastMath.abs(s - I * m);
  135.             dGmsdh = cjsjKH.getDcjDh(sMj).multiply(cjsjAB.getCj(sMim)).
  136.                      add(cjsjKH.getDsjDh(sMj).multiply(cjsjAB.getSj(sMim)).multiply(sgn(s - m)).multiply(sgn(s - j)));
  137.         }
  138.         return dGmsdh;
  139.     }

  140.     /** Get the dG<sub>ms</sub><sup>j</sup> / d<sub>α</sub> coefficient.
  141.      * @param m m subscript
  142.      * @param s s subscript
  143.      * @param j order
  144.      * @return dG<sub>ms</sub><sup>j</sup> / d<sub>α</sub>
  145.      */
  146.     public T getdGmsdAlpha(final int m, final int s, final int j) {
  147.         final int sMj  = FastMath.abs(s - j);
  148.         final T dGmsdAl;
  149.         if (FastMath.abs(s) <= m) {
  150.             final int mMis = m - I * s;
  151.             dGmsdAl = cjsjKH.getCj(sMj).multiply(cjsjAB.getDcjDk(mMis)).
  152.                       subtract(cjsjKH.getSj(sMj).multiply(cjsjAB.getDsjDk(mMis)).multiply(I).multiply(sgn(s - j)));
  153.         } else {
  154.             final int sMim = FastMath.abs(s - I * m);
  155.             dGmsdAl = cjsjKH.getCj(sMj).multiply(cjsjAB.getDcjDk(sMim)).
  156.                       add(cjsjKH.getSj(sMj).multiply(cjsjAB.getDsjDk(sMim)).multiply(sgn(s - j)).multiply(sgn(s - m)));
  157.         }
  158.         return dGmsdAl;
  159.     }

  160.     /** Get the dG<sub>ms</sub><sup>j</sup> / d<sub>β</sub> coefficient.
  161.      * @param m m subscript
  162.      * @param s s subscript
  163.      * @param j order
  164.      * @return dG<sub>ms</sub><sup>j</sup> / d<sub>β</sub>
  165.      */
  166.     public T getdGmsdBeta(final int m, final int s, final int j) {
  167.         final int sMj = FastMath.abs(s - j);
  168.         final T dGmsdBe;
  169.         if (FastMath.abs(s) <= m) {
  170.             final int mMis = m - I * s;
  171.             dGmsdBe = cjsjKH.getCj(sMj).multiply(cjsjAB.getDcjDh(mMis)).
  172.                       subtract(cjsjKH.getSj(sMj).multiply(cjsjAB.getDsjDh(mMis)).multiply(I).multiply(sgn(s - j)));
  173.         } else {
  174.             final int sMim = FastMath.abs(s - I * m);
  175.             dGmsdBe = cjsjKH.getCj(sMj).multiply(cjsjAB.getDcjDh(sMim)).
  176.                       add(cjsjKH.getSj(sMj).multiply(cjsjAB.getDsjDh(sMim)).multiply(sgn(s - j)).multiply(sgn(s - m)));
  177.         }
  178.         return dGmsdBe;
  179.     }

  180.     /** Get the dH<sub>ms</sub><sup>j</sup> / d<sub>k</sub> coefficient.
  181.      * @param m m subscript
  182.      * @param s s subscript
  183.      * @param j order
  184.      * @return dH<sub>ms</sub><sup>j</sup> / d<sub>k</sub>
  185.      */
  186.     public T getdHmsdk(final int m, final int s, final int j) {
  187.         final int sMj = FastMath.abs(s - j);
  188.         final T dHmsdk;
  189.         if (FastMath.abs(s) <= m) {
  190.             final int mMis = m - I * s;
  191.             dHmsdk = cjsjKH.getDcjDk(sMj).multiply(cjsjAB.getSj(mMis)).multiply(I).
  192.                      add(cjsjKH.getDsjDk(sMj).multiply(cjsjAB.getCj(mMis)).multiply(sgn(s - j)));
  193.         } else {
  194.             final int sMim = FastMath.abs(s - I * m);
  195.             dHmsdk = cjsjKH.getDcjDk(sMj).multiply(cjsjAB.getSj(sMim)).multiply(-sgn(s - m)).
  196.                      add(cjsjKH.getDsjDk(sMj).multiply(cjsjAB.getCj(sMim)).multiply(sgn(s - j)));
  197.         }
  198.         return dHmsdk;
  199.     }

  200.     /** Get the dH<sub>ms</sub><sup>j</sup> / d<sub>h</sub> coefficient.
  201.      * @param m m subscript
  202.      * @param s s subscript
  203.      * @param j order
  204.      * @return dH<sub>ms</sub><sup>j</sup> / d<sub>h</sub>
  205.      */
  206.     public T getdHmsdh(final int m,  final int s, final int j) {
  207.         final int sMj = FastMath.abs(s - j);
  208.         final T dHmsdh;
  209.         if (FastMath.abs(s) <= m) {
  210.             final int mMis = m - I * s;
  211.             dHmsdh = cjsjKH.getDcjDh(sMj).multiply(cjsjAB.getSj(mMis)).multiply(I).
  212.                      add(cjsjKH.getDsjDh(sMj).multiply(cjsjAB.getCj(mMis)).multiply(sgn(s - j)));
  213.         } else {
  214.             final int sMim = FastMath.abs(s - I * m);
  215.             dHmsdh = cjsjKH.getDcjDh(sMj).multiply(cjsjAB.getSj(sMim)).multiply(-sgn(s - m)).
  216.                      add(cjsjKH.getDsjDh(sMj).multiply(cjsjAB.getCj(sMim)).multiply(sgn(s - j)));
  217.         }
  218.         return dHmsdh;
  219.     }

  220.     /** Get the dH<sub>ms</sub><sup>j</sup> / d<sub>α</sub> coefficient.
  221.      * @param m m subscript
  222.      * @param s s subscript
  223.      * @param j order
  224.      * @return dH<sub>ms</sub><sup>j</sup> / d<sub>α</sub>
  225.      */
  226.     public T getdHmsdAlpha(final int m, final int s, final int j) {
  227.         final int sMj  = FastMath.abs(s - j);
  228.         final T dHmsdAl;
  229.         if (FastMath.abs(s) <= m) {
  230.             final int mMis = m - I * s;
  231.             dHmsdAl = cjsjKH.getCj(sMj).multiply(cjsjAB.getDsjDk(mMis)).multiply(I).
  232.                       add(cjsjKH.getSj(sMj).multiply(cjsjAB.getDcjDk(mMis)).multiply(sgn(s - j)));
  233.         } else {
  234.             final int sMim = FastMath.abs(s - I * m);
  235.             dHmsdAl = cjsjKH.getCj(sMj).multiply(cjsjAB.getDsjDk(sMim)).multiply(-sgn(s - m)).
  236.                       add(cjsjKH.getSj(sMj).multiply(cjsjAB.getDcjDk(sMim)).multiply(sgn(s - j)));
  237.         }
  238.         return dHmsdAl;
  239.     }

  240.     /** Get the dH<sub>ms</sub><sup>j</sup> / d<sub>β</sub> coefficient.
  241.      * @param m m subscript
  242.      * @param s s subscript
  243.      * @param j order
  244.      * @return dH<sub>ms</sub><sup>j</sup> / d<sub>β</sub>
  245.      */
  246.     public T getdHmsdBeta(final int m, final int s, final int j) {
  247.         final int sMj = FastMath.abs(s - j);
  248.         final T dHmsdBe;
  249.         if (FastMath.abs(s) <= m) {
  250.             final int mMis = m - I * s;
  251.             dHmsdBe = cjsjKH.getCj(sMj).multiply(cjsjAB.getDsjDh(mMis)).multiply(I).
  252.                       add(cjsjKH.getSj(sMj).multiply(cjsjAB.getDcjDh(mMis)).multiply(sgn(s - j)));
  253.         } else {
  254.             final int sMim = FastMath.abs(s - I * m);
  255.             dHmsdBe = cjsjKH.getCj(sMj).multiply(cjsjAB.getDsjDh(sMim)).multiply(-sgn(s - m)).
  256.                       add(cjsjKH.getSj(sMj).multiply(cjsjAB.getDcjDh(sMim)).multiply(sgn(s - j)));
  257.         }
  258.         return dHmsdBe;
  259.     }

  260.     /** Get the sign of an integer.
  261.      *  @param i number on which evaluation is done
  262.      *  @return -1 or +1 depending on sign of i
  263.      */
  264.     private int sgn(final int i) {
  265.         return (i < 0) ? -1 : 1;
  266.     }
  267. }