ThirdBodyAttractionEpoch.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.forces.gravity;

  18. import org.hipparchus.analysis.differentiation.Gradient;
  19. import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
  20. import org.hipparchus.geometry.euclidean.threed.Vector3D;
  21. import org.orekit.bodies.CelestialBody;
  22. import org.orekit.propagation.SpacecraftState;

  23. /** Third body attraction force model.
  24.  * This class is a copy of {@link ThirdBodyAttraction} class.
  25.  * The computation of derivatives of the acceleration
  26.  * w.r.t. the Epoch has been added.
  27.  *
  28.  * @author Fabien Maussion
  29.  * @author Véronique Pommier-Maurussane
  30.  * @since 10.2
  31.  */
  32. public class ThirdBodyAttractionEpoch extends ThirdBodyAttraction {

  33.     /** The body to consider. */
  34.     private final CelestialBody body;

  35.     /** Simple constructor.
  36.      * @param body the third body to consider
  37.      * (ex: {@link org.orekit.bodies.CelestialBodyFactory#getSun()} or
  38.      * {@link org.orekit.bodies.CelestialBodyFactory#getMoon()})
  39.      */
  40.     public ThirdBodyAttractionEpoch(final CelestialBody body) {
  41.         super(body);
  42.         this.body = body;
  43.     }

  44.     /** Compute acceleration.
  45.      * @param s current state information: date, kinematics, attitude
  46.      * @param parameters values of the force model parameters
  47.      * @return acceleration in same frame as state
  48.      */
  49.     private FieldVector3D<Gradient> accelerationToEpoch(final SpacecraftState s, final double[] parameters) {

  50.         final double gm = parameters[0];

  51.         // compute bodies separation vectors and squared norm
  52.         final Vector3D centralToBody = body.getPVCoordinates(s.getDate(), s.getFrame()).getPosition();

  53.         // Spacecraft Position
  54.         final double rx = centralToBody.getX();
  55.         final double ry = centralToBody.getY();
  56.         final double rz = centralToBody.getZ();

  57.         final int freeParameters = 3;
  58.         final Gradient fpx = Gradient.variable(freeParameters, 0, rx);
  59.         final Gradient fpy = Gradient.variable(freeParameters, 1, ry);
  60.         final Gradient fpz = Gradient.variable(freeParameters, 2, rz);

  61.         final FieldVector3D<Gradient> centralToBodyFV = new FieldVector3D<>(new Gradient[] {fpx, fpy, fpz});


  62.         final Gradient                r2Central = centralToBodyFV.getNormSq();
  63.         final FieldVector3D<Gradient> satToBody = centralToBodyFV.subtract(s.getPVCoordinates().getPosition());
  64.         final Gradient                r2Sat     = satToBody.getNormSq();

  65.         return new FieldVector3D<>(gm, satToBody.scalarMultiply(r2Sat.multiply(r2Sat.sqrt()).reciprocal()),
  66.                                   -gm, centralToBodyFV.scalarMultiply(r2Central.multiply(r2Central.sqrt()).reciprocal()));
  67.     }

  68.     /** Compute derivatives of the state w.r.t epoch.
  69.      * @param s current state information: date, kinematics, attitude
  70.      * @param parameters values of the force model parameters
  71.      * @return derivatives
  72.      */
  73.     public double[] getDerivativesToEpoch(final SpacecraftState s, final double[] parameters) {

  74.         final FieldVector3D<Gradient> acc = accelerationToEpoch(s, parameters);
  75.         final Vector3D centralToBodyVelocity = body.getPVCoordinates(s.getDate(), s.getFrame()).getVelocity();

  76.         final double[] dAccxdR1i = acc.getX().getGradient();
  77.         final double[] dAccydR1i = acc.getY().getGradient();
  78.         final double[] dAcczdR1i = acc.getZ().getGradient();
  79.         final double[] v = centralToBodyVelocity.toArray();

  80.         return new double[] {
  81.             dAccxdR1i[0] * v[0] + dAccxdR1i[1] * v[1] + dAccxdR1i[2] * v[2],
  82.             dAccydR1i[0] * v[0] + dAccydR1i[1] * v[1] + dAccydR1i[2] * v[2],
  83.             dAcczdR1i[0] * v[0] + dAcczdR1i[1] * v[1] + dAcczdR1i[2] * v[2]
  84.         };

  85.     }

  86. }