Search in sources :

Example 6 with ODEIntegrator

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

the class NumericalPropagatorTest method createPropagator.

private static synchronized NumericalPropagator createPropagator(SpacecraftState spacecraftState, OrbitType orbitType, PositionAngle angleType) throws OrekitException {
    final double minStep = 0.001;
    final double maxStep = 120.0;
    final double positionTolerance = 0.1;
    final int degree = 20;
    final int order = 20;
    final double spacecraftArea = 1.0;
    final double spacecraftDragCoefficient = 2.0;
    final double spacecraftReflectionCoefficient = 2.0;
    // propagator main configuration
    final double[][] tol = NumericalPropagator.tolerances(positionTolerance, spacecraftState.getOrbit(), orbitType);
    final ODEIntegrator integrator = new DormandPrince853Integrator(minStep, maxStep, tol[0], tol[1]);
    final NumericalPropagator np = new NumericalPropagator(integrator);
    np.setOrbitType(orbitType);
    np.setPositionAngleType(angleType);
    np.setInitialState(spacecraftState);
    // Earth gravity field
    final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
    final NormalizedSphericalHarmonicsProvider harmonicsGravityProvider = GravityFieldFactory.getNormalizedProvider(degree, order);
    np.addForceModel(new HolmesFeatherstoneAttractionModel(earth.getBodyFrame(), harmonicsGravityProvider));
    // Sun and Moon attraction
    np.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getSun()));
    np.addForceModel(new ThirdBodyAttraction(CelestialBodyFactory.getMoon()));
    // atmospheric drag
    MarshallSolarActivityFutureEstimation msafe = new MarshallSolarActivityFutureEstimation("Jan2000F10-edited-data\\.txt", MarshallSolarActivityFutureEstimation.StrengthLevel.AVERAGE);
    DataProvidersManager.getInstance().feed(msafe.getSupportedNames(), msafe);
    DTM2000 atmosphere = new DTM2000(msafe, CelestialBodyFactory.getSun(), earth);
    np.addForceModel(new DragForce(atmosphere, new IsotropicDrag(spacecraftArea, spacecraftDragCoefficient)));
    // solar radiation pressure
    np.addForceModel(new SolarRadiationPressure(CelestialBodyFactory.getSun(), earth.getEquatorialRadius(), new IsotropicRadiationSingleCoefficient(spacecraftArea, spacecraftReflectionCoefficient)));
    return np;
}
Also used : OneAxisEllipsoid(org.orekit.bodies.OneAxisEllipsoid) IsotropicDrag(org.orekit.forces.drag.IsotropicDrag) DTM2000(org.orekit.forces.drag.atmosphere.DTM2000) SolarRadiationPressure(org.orekit.forces.radiation.SolarRadiationPressure) MarshallSolarActivityFutureEstimation(org.orekit.forces.drag.atmosphere.data.MarshallSolarActivityFutureEstimation) ThirdBodyAttraction(org.orekit.forces.gravity.ThirdBodyAttraction) ODEIntegrator(org.hipparchus.ode.ODEIntegrator) DragForce(org.orekit.forces.drag.DragForce) DormandPrince853Integrator(org.hipparchus.ode.nonstiff.DormandPrince853Integrator) NormalizedSphericalHarmonicsProvider(org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider) HolmesFeatherstoneAttractionModel(org.orekit.forces.gravity.HolmesFeatherstoneAttractionModel) IsotropicRadiationSingleCoefficient(org.orekit.forces.radiation.IsotropicRadiationSingleCoefficient)

Example 7 with ODEIntegrator

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

the class AngularCoordinatesTest method testShiftWithAcceleration.

@Test
public void testShiftWithAcceleration() throws OrekitException {
    double rate = 2 * FastMath.PI / (12 * 60);
    double acc = 0.001;
    double dt = 1.0;
    int n = 2000;
    final AngularCoordinates quadratic = new AngularCoordinates(Rotation.IDENTITY, new Vector3D(rate, Vector3D.PLUS_K), new Vector3D(acc, Vector3D.PLUS_J));
    final AngularCoordinates linear = new AngularCoordinates(quadratic.getRotation(), quadratic.getRotationRate(), Vector3D.ZERO);
    final OrdinaryDifferentialEquation ode = new OrdinaryDifferentialEquation() {

        public int getDimension() {
            return 4;
        }

        public double[] computeDerivatives(final double t, final double[] q) {
            final double omegaX = quadratic.getRotationRate().getX() + t * quadratic.getRotationAcceleration().getX();
            final double omegaY = quadratic.getRotationRate().getY() + t * quadratic.getRotationAcceleration().getY();
            final double omegaZ = quadratic.getRotationRate().getZ() + t * quadratic.getRotationAcceleration().getZ();
            return new double[] { 0.5 * MathArrays.linearCombination(-q[1], omegaX, -q[2], omegaY, -q[3], omegaZ), 0.5 * MathArrays.linearCombination(q[0], omegaX, -q[3], omegaY, q[2], omegaZ), 0.5 * MathArrays.linearCombination(q[3], omegaX, q[0], omegaY, -q[1], omegaZ), 0.5 * MathArrays.linearCombination(-q[2], omegaX, q[1], omegaY, q[0], omegaZ) };
        }
    };
    ODEIntegrator integrator = new DormandPrince853Integrator(1.0e-6, 1.0, 1.0e-12, 1.0e-12);
    integrator.addStepHandler(new StepNormalizer(dt / n, new ODEFixedStepHandler() {

        public void handleStep(ODEStateAndDerivative s, boolean isLast) {
            final double t = s.getTime();
            final double[] y = s.getPrimaryState();
            Rotation reference = new Rotation(y[0], y[1], y[2], y[3], true);
            // the error in shiftedBy taking acceleration into account is cubic
            double expectedCubicError = 1.4544e-6 * t * t * t;
            Assert.assertEquals(expectedCubicError, Rotation.distance(reference, quadratic.shiftedBy(t).getRotation()), 0.0001 * expectedCubicError);
            // the error in shiftedBy not taking acceleration into account is quadratic
            double expectedQuadraticError = 5.0e-4 * t * t;
            Assert.assertEquals(expectedQuadraticError, Rotation.distance(reference, linear.shiftedBy(t).getRotation()), 0.00001 * expectedQuadraticError);
        }
    }));
    double[] y = new double[] { quadratic.getRotation().getQ0(), quadratic.getRotation().getQ1(), quadratic.getRotation().getQ2(), quadratic.getRotation().getQ3() };
    integrator.integrate(ode, new ODEState(0, y), dt);
}
Also used : ODEFixedStepHandler(org.hipparchus.ode.sampling.ODEFixedStepHandler) ODEState(org.hipparchus.ode.ODEState) Rotation(org.hipparchus.geometry.euclidean.threed.Rotation) Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) OrdinaryDifferentialEquation(org.hipparchus.ode.OrdinaryDifferentialEquation) ODEIntegrator(org.hipparchus.ode.ODEIntegrator) DormandPrince853Integrator(org.hipparchus.ode.nonstiff.DormandPrince853Integrator) StepNormalizer(org.hipparchus.ode.sampling.StepNormalizer) ODEStateAndDerivative(org.hipparchus.ode.ODEStateAndDerivative) Test(org.junit.Test)

Aggregations

ODEIntegrator (org.hipparchus.ode.ODEIntegrator)7 DormandPrince853Integrator (org.hipparchus.ode.nonstiff.DormandPrince853Integrator)4 ArrayList (java.util.ArrayList)3 Rotation (org.hipparchus.geometry.euclidean.threed.Rotation)2 Vector3D (org.hipparchus.geometry.euclidean.threed.Vector3D)2 ODEState (org.hipparchus.ode.ODEState)2 ODEStateAndDerivative (org.hipparchus.ode.ODEStateAndDerivative)2 OrdinaryDifferentialEquation (org.hipparchus.ode.OrdinaryDifferentialEquation)2 ODEFixedStepHandler (org.hipparchus.ode.sampling.ODEFixedStepHandler)2 ODEStepHandler (org.hipparchus.ode.sampling.ODEStepHandler)2 StepNormalizer (org.hipparchus.ode.sampling.StepNormalizer)2 DragForce (org.orekit.forces.drag.DragForce)2 IsotropicDrag (org.orekit.forces.drag.IsotropicDrag)2 DTM2000 (org.orekit.forces.drag.atmosphere.DTM2000)2 HolmesFeatherstoneAttractionModel (org.orekit.forces.gravity.HolmesFeatherstoneAttractionModel)2 NumericalPropagator (org.orekit.propagation.numerical.NumericalPropagator)2 ClassicalRungeKuttaIntegrator (org.hipparchus.ode.nonstiff.ClassicalRungeKuttaIntegrator)1 Test (org.junit.Test)1 OneAxisEllipsoid (org.orekit.bodies.OneAxisEllipsoid)1 ForceModel (org.orekit.forces.ForceModel)1