Search in sources :

Example 86 with DerivativeStructure

use of org.hipparchus.analysis.differentiation.DerivativeStructure in project Orekit by CS-SI.

the class AssociatedLegendreFunction method value.

public DerivativeStructure value(DerivativeStructure t) {
    DerivativeStructure y1 = t.getField().getOne().multiply(polynomial[polynomial.length - 1].toDouble());
    for (int j = polynomial.length - 2; j >= 0; j--) {
        y1 = y1.multiply(t).add(polynomial[j].toDouble());
    }
    DerivativeStructure oneMinusT2 = t.getField().getOne().subtract(t.multiply(t));
    DerivativeStructure y2 = oneMinusT2.pow(m).sqrt();
    return y1.multiply(y2).multiply(normalization.toDouble());
}
Also used : DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure)

Example 87 with DerivativeStructure

use of org.hipparchus.analysis.differentiation.DerivativeStructure in project Orekit by CS-SI.

the class HolmesFeatherstoneAttractionModelTest method RealFieldExpectErrorTest.

/**
 *Same test as the previous one but not adding the ForceModel to the NumericalPropagator
 *    it is a test to validate the previous test.
 *    (to test if the ForceModel it's actually
 *    doing something in the Propagator and the FieldPropagator)
 */
@Test
public void RealFieldExpectErrorTest() throws OrekitException {
    DSFactory factory = new DSFactory(6, 0);
    DerivativeStructure a_0 = factory.variable(0, 7201009.7124401);
    DerivativeStructure e_0 = factory.variable(1, 1e-3);
    DerivativeStructure i_0 = factory.variable(2, 98.7 * FastMath.PI / 180);
    DerivativeStructure R_0 = factory.variable(3, 15.0 * 22.5 * FastMath.PI / 180);
    DerivativeStructure O_0 = factory.variable(4, 93.0 * FastMath.PI / 180);
    DerivativeStructure n_0 = factory.variable(5, 0.1);
    Field<DerivativeStructure> field = a_0.getField();
    DerivativeStructure zero = field.getZero();
    FieldAbsoluteDate<DerivativeStructure> J2000 = new FieldAbsoluteDate<>(field);
    Frame EME = FramesFactory.getEME2000();
    FieldKeplerianOrbit<DerivativeStructure> FKO = new FieldKeplerianOrbit<>(a_0, e_0, i_0, R_0, O_0, n_0, PositionAngle.MEAN, EME, J2000, Constants.EIGEN5C_EARTH_MU);
    FieldSpacecraftState<DerivativeStructure> initialState = new FieldSpacecraftState<>(FKO);
    SpacecraftState iSR = initialState.toSpacecraftState();
    OrbitType type = OrbitType.EQUINOCTIAL;
    double[][] tolerance = NumericalPropagator.tolerances(10.0, FKO.toOrbit(), type);
    AdaptiveStepsizeFieldIntegrator<DerivativeStructure> integrator = new DormandPrince853FieldIntegrator<>(field, 0.001, 200, tolerance[0], tolerance[1]);
    integrator.setInitialStepSize(zero.add(60));
    AdaptiveStepsizeIntegrator RIntegrator = new DormandPrince853Integrator(0.001, 200, tolerance[0], tolerance[1]);
    RIntegrator.setInitialStepSize(60);
    FieldNumericalPropagator<DerivativeStructure> FNP = new FieldNumericalPropagator<>(field, integrator);
    FNP.setOrbitType(type);
    FNP.setInitialState(initialState);
    NumericalPropagator NP = new NumericalPropagator(RIntegrator);
    NP.setOrbitType(type);
    NP.setInitialState(iSR);
    double[][] c = new double[3][1];
    c[0][0] = 0.0;
    c[2][0] = normalizedC20;
    double[][] s = new double[3][1];
    NormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getNormalizedProvider(6378136.460, mu, TideSystem.UNKNOWN, c, s);
    HolmesFeatherstoneAttractionModel forceModel = new HolmesFeatherstoneAttractionModel(itrf, provider);
    // FNP.addForceModel(forceModel);
    NP.addForceModel(forceModel);
    FieldAbsoluteDate<DerivativeStructure> target = J2000.shiftedBy(100.);
    FieldSpacecraftState<DerivativeStructure> finalState_DS = FNP.propagate(target);
    SpacecraftState finalState_R = NP.propagate(target.toAbsoluteDate());
    FieldPVCoordinates<DerivativeStructure> finPVC_DS = finalState_DS.getPVCoordinates();
    PVCoordinates finPVC_R = finalState_R.getPVCoordinates();
    Assert.assertFalse(FastMath.abs(finPVC_DS.toPVCoordinates().getPosition().getX() - finPVC_R.getPosition().getX()) < FastMath.abs(finPVC_R.getPosition().getX()) * 1e-11);
    Assert.assertFalse(FastMath.abs(finPVC_DS.toPVCoordinates().getPosition().getY() - finPVC_R.getPosition().getY()) < FastMath.abs(finPVC_R.getPosition().getY()) * 1e-11);
    Assert.assertFalse(FastMath.abs(finPVC_DS.toPVCoordinates().getPosition().getZ() - finPVC_R.getPosition().getZ()) < FastMath.abs(finPVC_R.getPosition().getZ()) * 1e-11);
}
Also used : DormandPrince853FieldIntegrator(org.hipparchus.ode.nonstiff.DormandPrince853FieldIntegrator) Frame(org.orekit.frames.Frame) FieldSpacecraftState(org.orekit.propagation.FieldSpacecraftState) AdaptiveStepsizeIntegrator(org.hipparchus.ode.nonstiff.AdaptiveStepsizeIntegrator) DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) DSFactory(org.hipparchus.analysis.differentiation.DSFactory) PVCoordinates(org.orekit.utils.PVCoordinates) FieldPVCoordinates(org.orekit.utils.FieldPVCoordinates) FieldKeplerianOrbit(org.orekit.orbits.FieldKeplerianOrbit) SpacecraftState(org.orekit.propagation.SpacecraftState) FieldSpacecraftState(org.orekit.propagation.FieldSpacecraftState) FieldNumericalPropagator(org.orekit.propagation.numerical.FieldNumericalPropagator) NumericalPropagator(org.orekit.propagation.numerical.NumericalPropagator) FieldNumericalPropagator(org.orekit.propagation.numerical.FieldNumericalPropagator) OrbitType(org.orekit.orbits.OrbitType) DormandPrince853Integrator(org.hipparchus.ode.nonstiff.DormandPrince853Integrator) NormalizedSphericalHarmonicsProvider(org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) AbstractLegacyForceModelTest(org.orekit.forces.AbstractLegacyForceModelTest) Test(org.junit.Test)

Example 88 with DerivativeStructure

use of org.hipparchus.analysis.differentiation.DerivativeStructure in project Orekit by CS-SI.

the class HolmesFeatherstoneAttractionModelTest method accelerationDerivatives.

@Override
protected FieldVector3D<DerivativeStructure> accelerationDerivatives(final ForceModel forceModel, final AbsoluteDate date, final Frame frame, final FieldVector3D<DerivativeStructure> position, final FieldVector3D<DerivativeStructure> velocity, final FieldRotation<DerivativeStructure> rotation, final DerivativeStructure mass) throws OrekitException {
    try {
        java.lang.reflect.Field bodyFrameField = HolmesFeatherstoneAttractionModel.class.getDeclaredField("bodyFrame");
        bodyFrameField.setAccessible(true);
        Frame bodyFrame = (Frame) bodyFrameField.get(forceModel);
        // get the position in body frame
        final Transform fromBodyFrame = bodyFrame.getTransformTo(frame, date);
        final Transform toBodyFrame = fromBodyFrame.getInverse();
        final Vector3D positionBody = toBodyFrame.transformPosition(position.toVector3D());
        // compute gradient and Hessian
        final GradientHessian gh = gradientHessian((HolmesFeatherstoneAttractionModel) forceModel, date, positionBody);
        // gradient of the non-central part of the gravity field
        final double[] gInertial = fromBodyFrame.transformVector(new Vector3D(gh.getGradient())).toArray();
        // Hessian of the non-central part of the gravity field
        final RealMatrix hBody = new Array2DRowRealMatrix(gh.getHessian(), false);
        final RealMatrix rot = new Array2DRowRealMatrix(toBodyFrame.getRotation().getMatrix());
        final RealMatrix hInertial = rot.transpose().multiply(hBody).multiply(rot);
        // distribute all partial derivatives in a compact acceleration vector
        final double[] derivatives = new double[1 + mass.getFreeParameters()];
        final DerivativeStructure[] accDer = new DerivativeStructure[3];
        for (int i = 0; i < 3; ++i) {
            // first element is value of acceleration (i.e. gradient of field)
            derivatives[0] = gInertial[i];
            // next three elements are one row of the Jacobian of acceleration (i.e. Hessian of field)
            derivatives[1] = hInertial.getEntry(i, 0);
            derivatives[2] = hInertial.getEntry(i, 1);
            derivatives[3] = hInertial.getEntry(i, 2);
            // next elements (three or four depending on mass being used or not) are left as 0
            accDer[i] = mass.getFactory().build(derivatives);
        }
        return new FieldVector3D<>(accDer);
    } catch (IllegalArgumentException | IllegalAccessException | NoSuchFieldException | SecurityException e) {
        return null;
    }
}
Also used : Frame(org.orekit.frames.Frame) DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) RealMatrix(org.hipparchus.linear.RealMatrix) Array2DRowRealMatrix(org.hipparchus.linear.Array2DRowRealMatrix) Array2DRowRealMatrix(org.hipparchus.linear.Array2DRowRealMatrix) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) Transform(org.orekit.frames.Transform)

Example 89 with DerivativeStructure

use of org.hipparchus.analysis.differentiation.DerivativeStructure in project Orekit by CS-SI.

the class RelativityTest method accelerationDerivatives.

@Override
protected FieldVector3D<DerivativeStructure> accelerationDerivatives(final ForceModel forceModel, final AbsoluteDate date, final Frame frame, final FieldVector3D<DerivativeStructure> position, final FieldVector3D<DerivativeStructure> velocity, final FieldRotation<DerivativeStructure> rotation, final DerivativeStructure mass) throws OrekitException {
    try {
        double gm = forceModel.getParameterDriver(NewtonianAttraction.CENTRAL_ATTRACTION_COEFFICIENT).getValue();
        // radius
        final DerivativeStructure r2 = position.getNormSq();
        final DerivativeStructure r = r2.sqrt();
        // speed squared
        final DerivativeStructure s2 = velocity.getNormSq();
        final double c2 = Constants.SPEED_OF_LIGHT * Constants.SPEED_OF_LIGHT;
        // eq. 3.146
        return new FieldVector3D<>(r.reciprocal().multiply(4 * gm).subtract(s2), position, position.dotProduct(velocity).multiply(4), velocity).scalarMultiply(r2.multiply(r).multiply(c2).reciprocal().multiply(gm));
    } catch (IllegalArgumentException | SecurityException e) {
        return null;
    }
}
Also used : DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure)

Example 90 with DerivativeStructure

use of org.hipparchus.analysis.differentiation.DerivativeStructure in project Orekit by CS-SI.

the class RelativityTest method testAcceleration.

/**
 * check the acceleration from relativity
 *
 * @throws OrekitException on error
 */
@Test
public void testAcceleration() throws OrekitException {
    double gm = Constants.EIGEN5C_EARTH_MU;
    Relativity relativity = new Relativity(gm);
    Assert.assertFalse(relativity.dependsOnPositionOnly());
    final Vector3D p = new Vector3D(3777828.75000531, -5543949.549783845, 2563117.448578311);
    final Vector3D v = new Vector3D(489.0060271721, -2849.9328929417, -6866.4671013153);
    SpacecraftState s = new SpacecraftState(new CartesianOrbit(new PVCoordinates(p, v), frame, date, gm));
    // action
    Vector3D acceleration = relativity.acceleration(s, relativity.getParameters());
    // verify
    // force is ~1e-8 so this give ~3 sig figs.
    double tol = 2e-11;
    Vector3D circularApproximation = p.normalize().scalarMultiply(gm / p.getNormSq() * 3 * v.getNormSq() / (c * c));
    Assert.assertEquals(0, acceleration.subtract(circularApproximation).getNorm(), tol);
    // check derivatives
    FieldSpacecraftState<DerivativeStructure> sDS = toDS(s, new LofOffset(s.getFrame(), LOFType.VVLH));
    final Vector3D actualDerivatives = relativity.acceleration(sDS, relativity.getParameters(sDS.getDate().getField())).toVector3D();
    Assert.assertEquals(0, actualDerivatives.subtract(circularApproximation).getNorm(), tol);
}
Also used : FieldSpacecraftState(org.orekit.propagation.FieldSpacecraftState) SpacecraftState(org.orekit.propagation.SpacecraftState) CartesianOrbit(org.orekit.orbits.CartesianOrbit) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) FieldPVCoordinates(org.orekit.utils.FieldPVCoordinates) PVCoordinates(org.orekit.utils.PVCoordinates) LofOffset(org.orekit.attitudes.LofOffset) AbstractLegacyForceModelTest(org.orekit.forces.AbstractLegacyForceModelTest) Test(org.junit.Test)

Aggregations

DerivativeStructure (org.hipparchus.analysis.differentiation.DerivativeStructure)140 Test (org.junit.Test)69 DSFactory (org.hipparchus.analysis.differentiation.DSFactory)63 FieldVector3D (org.hipparchus.geometry.euclidean.threed.FieldVector3D)42 FieldAbsoluteDate (org.orekit.time.FieldAbsoluteDate)40 Vector3D (org.hipparchus.geometry.euclidean.threed.Vector3D)33 SpacecraftState (org.orekit.propagation.SpacecraftState)30 AbsoluteDate (org.orekit.time.AbsoluteDate)25 RandomGenerator (org.hipparchus.random.RandomGenerator)22 Frame (org.orekit.frames.Frame)22 PVCoordinates (org.orekit.utils.PVCoordinates)21 FieldSpacecraftState (org.orekit.propagation.FieldSpacecraftState)20 FieldPVCoordinates (org.orekit.utils.FieldPVCoordinates)18 OrekitException (org.orekit.errors.OrekitException)16 FiniteDifferencesDifferentiator (org.hipparchus.analysis.differentiation.FiniteDifferencesDifferentiator)15 AbstractLegacyForceModelTest (org.orekit.forces.AbstractLegacyForceModelTest)15 OrbitType (org.orekit.orbits.OrbitType)15 ParameterDriver (org.orekit.utils.ParameterDriver)15 FieldKeplerianOrbit (org.orekit.orbits.FieldKeplerianOrbit)14 FieldNumericalPropagator (org.orekit.propagation.numerical.FieldNumericalPropagator)14