Search in sources :

Example 6 with Field

use of org.hipparchus.Field in project Orekit by CS-SI.

the class AbstractForceModelTest method checkStateJacobianVsFiniteDifferences.

protected void checkStateJacobianVsFiniteDifferences(final SpacecraftState state0, final ForceModel forceModel, final AttitudeProvider provider, final double dP, final double checkTolerance, final boolean print) throws OrekitException {
    double[][] finiteDifferencesJacobian = Differentiation.differentiate(state -> forceModel.acceleration(state, forceModel.getParameters()).toArray(), 3, provider, OrbitType.CARTESIAN, PositionAngle.MEAN, dP, 5).value(state0);
    DSFactory factory = new DSFactory(6, 1);
    Field<DerivativeStructure> field = factory.getDerivativeField();
    final FieldAbsoluteDate<DerivativeStructure> fDate = new FieldAbsoluteDate<>(field, state0.getDate());
    final Vector3D p = state0.getPVCoordinates().getPosition();
    final Vector3D v = state0.getPVCoordinates().getVelocity();
    final Vector3D a = state0.getPVCoordinates().getAcceleration();
    final TimeStampedFieldPVCoordinates<DerivativeStructure> fPVA = new TimeStampedFieldPVCoordinates<>(fDate, new FieldVector3D<>(factory.variable(0, p.getX()), factory.variable(1, p.getY()), factory.variable(2, p.getZ())), new FieldVector3D<>(factory.variable(3, v.getX()), factory.variable(4, v.getY()), factory.variable(5, v.getZ())), new FieldVector3D<>(factory.constant(a.getX()), factory.constant(a.getY()), factory.constant(a.getZ())));
    final TimeStampedFieldAngularCoordinates<DerivativeStructure> fAC = new TimeStampedFieldAngularCoordinates<>(fDate, new FieldRotation<>(field, state0.getAttitude().getRotation()), new FieldVector3D<>(field, state0.getAttitude().getSpin()), new FieldVector3D<>(field, state0.getAttitude().getRotationAcceleration()));
    final FieldSpacecraftState<DerivativeStructure> fState = new FieldSpacecraftState<>(new FieldCartesianOrbit<>(fPVA, state0.getFrame(), state0.getMu()), new FieldAttitude<>(state0.getFrame(), fAC), field.getZero().add(state0.getMass()));
    FieldVector3D<DerivativeStructure> dsJacobian = forceModel.acceleration(fState, forceModel.getParameters(fState.getDate().getField()));
    Vector3D dFdPXRef = new Vector3D(finiteDifferencesJacobian[0][0], finiteDifferencesJacobian[1][0], finiteDifferencesJacobian[2][0]);
    Vector3D dFdPXRes = new Vector3D(dsJacobian.getX().getPartialDerivative(1, 0, 0, 0, 0, 0), dsJacobian.getY().getPartialDerivative(1, 0, 0, 0, 0, 0), dsJacobian.getZ().getPartialDerivative(1, 0, 0, 0, 0, 0));
    Vector3D dFdPYRef = new Vector3D(finiteDifferencesJacobian[0][1], finiteDifferencesJacobian[1][1], finiteDifferencesJacobian[2][1]);
    Vector3D dFdPYRes = new Vector3D(dsJacobian.getX().getPartialDerivative(0, 1, 0, 0, 0, 0), dsJacobian.getY().getPartialDerivative(0, 1, 0, 0, 0, 0), dsJacobian.getZ().getPartialDerivative(0, 1, 0, 0, 0, 0));
    Vector3D dFdPZRef = new Vector3D(finiteDifferencesJacobian[0][2], finiteDifferencesJacobian[1][2], finiteDifferencesJacobian[2][2]);
    Vector3D dFdPZRes = new Vector3D(dsJacobian.getX().getPartialDerivative(0, 0, 1, 0, 0, 0), dsJacobian.getY().getPartialDerivative(0, 0, 1, 0, 0, 0), dsJacobian.getZ().getPartialDerivative(0, 0, 1, 0, 0, 0));
    Vector3D dFdVXRef = new Vector3D(finiteDifferencesJacobian[0][3], finiteDifferencesJacobian[1][3], finiteDifferencesJacobian[2][3]);
    Vector3D dFdVXRes = new Vector3D(dsJacobian.getX().getPartialDerivative(0, 0, 0, 1, 0, 0), dsJacobian.getY().getPartialDerivative(0, 0, 0, 1, 0, 0), dsJacobian.getZ().getPartialDerivative(0, 0, 0, 1, 0, 0));
    Vector3D dFdVYRef = new Vector3D(finiteDifferencesJacobian[0][4], finiteDifferencesJacobian[1][4], finiteDifferencesJacobian[2][4]);
    Vector3D dFdVYRes = new Vector3D(dsJacobian.getX().getPartialDerivative(0, 0, 0, 0, 1, 0), dsJacobian.getY().getPartialDerivative(0, 0, 0, 0, 1, 0), dsJacobian.getZ().getPartialDerivative(0, 0, 0, 0, 1, 0));
    Vector3D dFdVZRef = new Vector3D(finiteDifferencesJacobian[0][5], finiteDifferencesJacobian[1][5], finiteDifferencesJacobian[2][5]);
    Vector3D dFdVZRes = new Vector3D(dsJacobian.getX().getPartialDerivative(0, 0, 0, 0, 0, 1), dsJacobian.getY().getPartialDerivative(0, 0, 0, 0, 0, 1), dsJacobian.getZ().getPartialDerivative(0, 0, 0, 0, 0, 1));
    if (print) {
        System.out.println("dF/dPX ref: " + dFdPXRef.getX() + " " + dFdPXRef.getY() + " " + dFdPXRef.getZ());
        System.out.println("dF/dPX res: " + dFdPXRes.getX() + " " + dFdPXRes.getY() + " " + dFdPXRes.getZ());
        System.out.println("dF/dPY ref: " + dFdPYRef.getX() + " " + dFdPYRef.getY() + " " + dFdPYRef.getZ());
        System.out.println("dF/dPY res: " + dFdPYRes.getX() + " " + dFdPYRes.getY() + " " + dFdPYRes.getZ());
        System.out.println("dF/dPZ ref: " + dFdPZRef.getX() + " " + dFdPZRef.getY() + " " + dFdPZRef.getZ());
        System.out.println("dF/dPZ res: " + dFdPZRes.getX() + " " + dFdPZRes.getY() + " " + dFdPZRes.getZ());
        System.out.println("dF/dPX ref norm: " + dFdPXRef.getNorm() + ", abs error: " + Vector3D.distance(dFdPXRef, dFdPXRes) + ", rel error: " + (Vector3D.distance(dFdPXRef, dFdPXRes) / dFdPXRef.getNorm()));
        System.out.println("dF/dPY ref norm: " + dFdPYRef.getNorm() + ", abs error: " + Vector3D.distance(dFdPYRef, dFdPYRes) + ", rel error: " + (Vector3D.distance(dFdPYRef, dFdPYRes) / dFdPYRef.getNorm()));
        System.out.println("dF/dPZ ref norm: " + dFdPZRef.getNorm() + ", abs error: " + Vector3D.distance(dFdPZRef, dFdPZRes) + ", rel error: " + (Vector3D.distance(dFdPZRef, dFdPZRes) / dFdPZRef.getNorm()));
        System.out.println("dF/dVX ref norm: " + dFdVXRef.getNorm() + ", abs error: " + Vector3D.distance(dFdVXRef, dFdVXRes) + ", rel error: " + (Vector3D.distance(dFdVXRef, dFdVXRes) / dFdVXRef.getNorm()));
        System.out.println("dF/dVY ref norm: " + dFdVYRef.getNorm() + ", abs error: " + Vector3D.distance(dFdVYRef, dFdVYRes) + ", rel error: " + (Vector3D.distance(dFdVYRef, dFdVYRes) / dFdVYRef.getNorm()));
        System.out.println("dF/dVZ ref norm: " + dFdVZRef.getNorm() + ", abs error: " + Vector3D.distance(dFdVZRef, dFdVZRes) + ", rel error: " + (Vector3D.distance(dFdVZRef, dFdVZRes) / dFdVZRef.getNorm()));
    }
    checkdFdP(dFdPXRef, dFdPXRes, checkTolerance);
    checkdFdP(dFdPYRef, dFdPYRes, checkTolerance);
    checkdFdP(dFdPZRef, dFdPZRes, checkTolerance);
    checkdFdP(dFdVXRef, dFdVXRes, checkTolerance);
    checkdFdP(dFdVYRef, dFdVYRes, checkTolerance);
    checkdFdP(dFdVZRef, dFdVZRes, checkTolerance);
}
Also used : ParameterDriver(org.orekit.utils.ParameterDriver) OrekitStepHandler(org.orekit.propagation.sampling.OrekitStepHandler) FieldRotation(org.hipparchus.geometry.euclidean.threed.FieldRotation) Differentiation(org.orekit.utils.Differentiation) FieldAttitude(org.orekit.attitudes.FieldAttitude) AttitudeProvider(org.orekit.attitudes.AttitudeProvider) FieldSpacecraftState(org.orekit.propagation.FieldSpacecraftState) DSFactory(org.hipparchus.analysis.differentiation.DSFactory) Precision(org.hipparchus.util.Precision) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) SpacecraftState(org.orekit.propagation.SpacecraftState) OrbitType(org.orekit.orbits.OrbitType) DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) PositionAngle(org.orekit.orbits.PositionAngle) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) OrekitStepInterpolator(org.orekit.propagation.sampling.OrekitStepInterpolator) JacobiansMapper(org.orekit.propagation.numerical.JacobiansMapper) TimeStampedFieldAngularCoordinates(org.orekit.utils.TimeStampedFieldAngularCoordinates) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) PartialDerivativesEquations(org.orekit.propagation.numerical.PartialDerivativesEquations) FieldCartesianOrbit(org.orekit.orbits.FieldCartesianOrbit) Field(org.hipparchus.Field) OrekitException(org.orekit.errors.OrekitException) TimeStampedFieldPVCoordinates(org.orekit.utils.TimeStampedFieldPVCoordinates) NumericalPropagator(org.orekit.propagation.numerical.NumericalPropagator) Assert(org.junit.Assert) AbsoluteDate(org.orekit.time.AbsoluteDate) FieldSpacecraftState(org.orekit.propagation.FieldSpacecraftState) DerivativeStructure(org.hipparchus.analysis.differentiation.DerivativeStructure) DSFactory(org.hipparchus.analysis.differentiation.DSFactory) TimeStampedFieldAngularCoordinates(org.orekit.utils.TimeStampedFieldAngularCoordinates) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) TimeStampedFieldPVCoordinates(org.orekit.utils.TimeStampedFieldPVCoordinates)

Example 7 with Field

use of org.hipparchus.Field in project Orekit by CS-SI.

the class FieldCartesianOrbitTest method doTestNonKeplerianDerivatives.

private <T extends RealFieldElement<T>> void doTestNonKeplerianDerivatives(Field<T> field) throws OrekitException {
    final FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field, "2003-05-01T00:00:20.000", TimeScalesFactory.getUTC());
    final FieldVector3D<T> position = new FieldVector3D<>(field.getZero().add(6896874.444705), field.getZero().add(1956581.072644), field.getZero().add(-147476.245054));
    final FieldVector3D<T> velocity = new FieldVector3D<>(field.getZero().add(166.816407662), field.getZero().add(-1106.783301861), field.getZero().add(-7372.745712770));
    final FieldVector3D<T> acceleration = new FieldVector3D<>(field.getZero().add(-7.466182457944), field.getZero().add(-2.118153357345), field.getZero().add(0.160004048437));
    final TimeStampedFieldPVCoordinates<T> pv = new TimeStampedFieldPVCoordinates<>(date, position, velocity, acceleration);
    final Frame frame = FramesFactory.getEME2000();
    final double mu = Constants.EIGEN5C_EARTH_MU;
    final FieldCartesianOrbit<T> orbit = new FieldCartesianOrbit<>(pv, frame, mu);
    Assert.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getA()), orbit.getADot().getReal(), 4.3e-8);
    Assert.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getEquinoctialEx()), orbit.getEquinoctialExDot().getReal(), 2.1e-15);
    Assert.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getEquinoctialEy()), orbit.getEquinoctialEyDot().getReal(), 5.3e-16);
    Assert.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getHx()), orbit.getHxDot().getReal(), 4.4e-15);
    Assert.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getHy()), orbit.getHyDot().getReal(), 8.0e-16);
    Assert.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getLv()), orbit.getLvDot().getReal(), 1.2e-15);
    Assert.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getLE()), orbit.getLEDot().getReal(), 7.8e-16);
    Assert.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getLM()), orbit.getLMDot().getReal(), 8.8e-16);
    Assert.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getE()), orbit.getEDot().getReal(), 7.0e-16);
    Assert.assertEquals(differentiate(pv, frame, mu, shifted -> shifted.getI()), orbit.getIDot().getReal(), 5.7e-16);
}
Also used : Frame(org.orekit.frames.Frame) FieldPVCoordinates(org.orekit.utils.FieldPVCoordinates) Function(java.util.function.Function) FieldMatrixPreservingVisitor(org.hipparchus.linear.FieldMatrixPreservingVisitor) FieldEcksteinHechlerPropagator(org.orekit.propagation.analytical.FieldEcksteinHechlerPropagator) ArrayList(java.util.ArrayList) DSFactory(org.hipparchus.analysis.differentiation.DSFactory) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) MathArrays(org.hipparchus.util.MathArrays) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) FastMath(org.hipparchus.util.FastMath) Method(java.lang.reflect.Method) Utils(org.orekit.Utils) Before(org.junit.Before) Constants(org.orekit.utils.Constants) UnivariateFunction(org.hipparchus.analysis.UnivariateFunction) OrekitMatchers.relativelyCloseTo(org.orekit.OrekitMatchers.relativelyCloseTo) FramesFactory(org.orekit.frames.FramesFactory) Test(org.junit.Test) InvocationTargetException(java.lang.reflect.InvocationTargetException) Field(org.hipparchus.Field) MathUtils(org.hipparchus.util.MathUtils) List(java.util.List) UnivariateDifferentiableFunction(org.hipparchus.analysis.differentiation.UnivariateDifferentiableFunction) FiniteDifferencesDifferentiator(org.hipparchus.analysis.differentiation.FiniteDifferencesDifferentiator) OrekitException(org.orekit.errors.OrekitException) TimeStampedFieldPVCoordinates(org.orekit.utils.TimeStampedFieldPVCoordinates) RealFieldElement(org.hipparchus.RealFieldElement) TimeScalesFactory(org.orekit.time.TimeScalesFactory) MatrixUtils(org.hipparchus.linear.MatrixUtils) Decimal64Field(org.hipparchus.util.Decimal64Field) Transform(org.orekit.frames.Transform) Assert(org.junit.Assert) Frame(org.orekit.frames.Frame) FieldAbsoluteDate(org.orekit.time.FieldAbsoluteDate) FieldVector3D(org.hipparchus.geometry.euclidean.threed.FieldVector3D) TimeStampedFieldPVCoordinates(org.orekit.utils.TimeStampedFieldPVCoordinates)

Aggregations

Field (org.hipparchus.Field)7 FieldVector3D (org.hipparchus.geometry.euclidean.threed.FieldVector3D)7 OrekitException (org.orekit.errors.OrekitException)7 FieldAbsoluteDate (org.orekit.time.FieldAbsoluteDate)7 RealFieldElement (org.hipparchus.RealFieldElement)6 UnivariateFunction (org.hipparchus.analysis.UnivariateFunction)6 DSFactory (org.hipparchus.analysis.differentiation.DSFactory)6 FastMath (org.hipparchus.util.FastMath)6 Assert (org.junit.Assert)6 Frame (org.orekit.frames.Frame)6 Transform (org.orekit.frames.Transform)6 TimeStampedFieldPVCoordinates (org.orekit.utils.TimeStampedFieldPVCoordinates)6 ArrayList (java.util.ArrayList)5 List (java.util.List)5 Function (java.util.function.Function)5 FiniteDifferencesDifferentiator (org.hipparchus.analysis.differentiation.FiniteDifferencesDifferentiator)5 UnivariateDifferentiableFunction (org.hipparchus.analysis.differentiation.UnivariateDifferentiableFunction)5 FieldMatrixPreservingVisitor (org.hipparchus.linear.FieldMatrixPreservingVisitor)5 MatrixUtils (org.hipparchus.linear.MatrixUtils)5 Decimal64Field (org.hipparchus.util.Decimal64Field)5