List of usage examples for org.apache.commons.math3.analysis UnivariateFunction UnivariateFunction
UnivariateFunction
From source file:org.orekit.orbits.KeplerianParametersTest.java
@Test public void testKeplerianDerivatives() { final KeplerianOrbit o = new KeplerianOrbit(new PVCoordinates(new Vector3D(-4947831., -3765382., -3708221.), new Vector3D(-2079., 5291., -7842.)), FramesFactory.getEME2000(), date, 3.9860047e14); final Vector3D p = o.getPVCoordinates().getPosition(); final Vector3D v = o.getPVCoordinates().getVelocity(); final Vector3D a = o.getPVCoordinates().getAcceleration(); // check that despite we did not provide acceleration, it got recomputed Assert.assertEquals(7.605422, a.getNorm(), 1.0e-6); FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(8, 0.1); // check velocity is the derivative of position double vx = differentiator.differentiate(new UnivariateFunction() { public double value(double dt) { return o.shiftedBy(dt).getPVCoordinates().getPosition().getX(); }//w w w . j a v a 2 s.co m }).value(new DerivativeStructure(1, 1, 0, 0.0)).getPartialDerivative(1); Assert.assertEquals(o.getPVCoordinates().getVelocity().getX(), vx, 3.0e-12 * v.getNorm()); double vy = differentiator.differentiate(new UnivariateFunction() { public double value(double dt) { return o.shiftedBy(dt).getPVCoordinates().getPosition().getY(); } }).value(new DerivativeStructure(1, 1, 0, 0.0)).getPartialDerivative(1); Assert.assertEquals(o.getPVCoordinates().getVelocity().getY(), vy, 3.0e-12 * v.getNorm()); double vz = differentiator.differentiate(new UnivariateFunction() { public double value(double dt) { return o.shiftedBy(dt).getPVCoordinates().getPosition().getZ(); } }).value(new DerivativeStructure(1, 1, 0, 0.0)).getPartialDerivative(1); Assert.assertEquals(o.getPVCoordinates().getVelocity().getZ(), vz, 3.0e-12 * v.getNorm()); // check acceleration is the derivative of velocity double ax = differentiator.differentiate(new UnivariateFunction() { public double value(double dt) { return o.shiftedBy(dt).getPVCoordinates().getVelocity().getX(); } }).value(new DerivativeStructure(1, 1, 0, 0.0)).getPartialDerivative(1); Assert.assertEquals(o.getPVCoordinates().getAcceleration().getX(), ax, 3.0e-12 * a.getNorm()); double ay = differentiator.differentiate(new UnivariateFunction() { public double value(double dt) { return o.shiftedBy(dt).getPVCoordinates().getVelocity().getY(); } }).value(new DerivativeStructure(1, 1, 0, 0.0)).getPartialDerivative(1); Assert.assertEquals(o.getPVCoordinates().getAcceleration().getY(), ay, 3.0e-12 * a.getNorm()); double az = differentiator.differentiate(new UnivariateFunction() { public double value(double dt) { return o.shiftedBy(dt).getPVCoordinates().getVelocity().getZ(); } }).value(new DerivativeStructure(1, 1, 0, 0.0)).getPartialDerivative(1); Assert.assertEquals(o.getPVCoordinates().getAcceleration().getZ(), az, 3.0e-12 * a.getNorm()); // check jerk is the derivative of acceleration final double r2 = p.getNormSq(); final double r = FastMath.sqrt(r2); Vector3D keplerianJerk = new Vector3D(-3 * Vector3D.dotProduct(p, v) / r2, a, -a.getNorm() / r, v); double jx = differentiator.differentiate(new UnivariateFunction() { public double value(double dt) { return o.shiftedBy(dt).getPVCoordinates().getAcceleration().getX(); } }).value(new DerivativeStructure(1, 1, 0, 0.0)).getPartialDerivative(1); Assert.assertEquals(keplerianJerk.getX(), jx, 3.0e-12 * keplerianJerk.getNorm()); double jy = differentiator.differentiate(new UnivariateFunction() { public double value(double dt) { return o.shiftedBy(dt).getPVCoordinates().getAcceleration().getY(); } }).value(new DerivativeStructure(1, 1, 0, 0.0)).getPartialDerivative(1); Assert.assertEquals(keplerianJerk.getY(), jy, 3.0e-12 * keplerianJerk.getNorm()); double jz = differentiator.differentiate(new UnivariateFunction() { public double value(double dt) { return o.shiftedBy(dt).getPVCoordinates().getAcceleration().getZ(); } }).value(new DerivativeStructure(1, 1, 0, 0.0)).getPartialDerivative(1); Assert.assertEquals(keplerianJerk.getZ(), jz, 3.0e-12 * keplerianJerk.getNorm()); }
From source file:org.orekit.propagation.events.EventState.java
/** Evaluate the impact of the proposed step on the event detector. * @param interpolator step interpolator for the proposed step * @return true if the event detector triggers an event before * the end of the proposed step (this implies the step should be * rejected)/*from w w w . j av a 2s . com*/ * @exception OrekitException if the switching function * cannot be evaluated * @exception TooManyEvaluationsException if an event cannot be located * @exception NoBracketingException if bracketing cannot be performed */ public boolean evaluateStep(final OrekitStepInterpolator interpolator) throws OrekitException, TooManyEvaluationsException, NoBracketingException { try { final double convergence = detector.getThreshold(); final int maxIterationcount = detector.getMaxIterationCount(); if (forward ^ interpolator.isForward()) { forward = !forward; pendingEvent = false; pendingEventTime = null; previousEventTime = null; } final AbsoluteDate t1 = interpolator.getCurrentDate(); final double dt = t1.durationFrom(t0); if (FastMath.abs(dt) < convergence) { // we cannot do anything on such a small step, don't trigger any events return false; } final int n = FastMath.max(1, (int) FastMath.ceil(FastMath.abs(dt) / detector.getMaxCheckInterval())); final double h = dt / n; final UnivariateFunction f = new UnivariateFunction() { public double value(final double t) throws LocalWrapperException { try { interpolator.setInterpolatedDate(t0.shiftedBy(t)); return g(interpolator.getInterpolatedState()); } catch (OrekitException oe) { throw new LocalWrapperException(oe); } } }; final BracketingNthOrderBrentSolver solver = new BracketingNthOrderBrentSolver(convergence, 5); AbsoluteDate ta = t0; double ga = g0; for (int i = 0; i < n; ++i) { // evaluate detector value at the end of the substep final AbsoluteDate tb = t0.shiftedBy((i + 1) * h); interpolator.setInterpolatedDate(tb); final double gb = g(interpolator.getInterpolatedState()); // check events occurrence if (g0Positive ^ (gb >= 0)) { // there is a sign change: an event is expected during this step // variation direction, with respect to the integration direction increasing = gb >= ga; // find the event time making sure we select a solution just at or past the exact root final double dtA = ta.durationFrom(t0); final double dtB = tb.durationFrom(t0); final double dtRoot = forward ? solver.solve(maxIterationcount, f, dtA, dtB, AllowedSolution.RIGHT_SIDE) : solver.solve(maxIterationcount, f, dtB, dtA, AllowedSolution.LEFT_SIDE); final AbsoluteDate root = t0.shiftedBy(dtRoot); if ((previousEventTime != null) && (FastMath.abs(root.durationFrom(ta)) <= convergence) && (FastMath.abs(root.durationFrom(previousEventTime)) <= convergence)) { // we have either found nothing or found (again ?) a past event, // retry the substep excluding this value, and taking care to have the // required sign in case the g function is noisy around its zero and // crosses the axis several times do { ta = forward ? ta.shiftedBy(convergence) : ta.shiftedBy(-convergence); ga = f.value(ta.durationFrom(t0)); } while ((g0Positive ^ (ga >= 0)) && (forward ^ (ta.compareTo(tb) >= 0))); if (forward ^ (ta.compareTo(tb) >= 0)) { // we were able to skip this spurious root --i; } else { // we can't avoid this root before the end of the step, // we have to handle it despite it is close to the former one // maybe we have two very close roots pendingEventTime = root; pendingEvent = true; return true; } } else if ((previousEventTime == null) || (FastMath.abs(previousEventTime.durationFrom(root)) > convergence)) { pendingEventTime = root; pendingEvent = true; return true; } else { // no sign change: there is no event for now ta = tb; ga = gb; } } else { // no sign change: there is no event for now ta = tb; ga = gb; } } // no event during the whole step pendingEvent = false; pendingEventTime = null; return false; } catch (LocalWrapperException lwe) { throw lwe.getWrappedException(); } }
From source file:org.orekit.utils.FieldPVCoordinatesTest.java
@Test public void testNormalize() { RandomGenerator generator = new Well19937a(0x7ede9376e4e1ab5al); FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(5, 1.0e-3); for (int i = 0; i < 200; ++i) { final FieldPVCoordinates<DerivativeStructure> pv = randomPVCoordinates(generator, 1e6, 1e3, 1.0); DerivativeStructure x = differentiator.differentiate(new UnivariateFunction() { public double value(double t) { return pv.shiftedBy(t).getPosition().normalize().getX().getValue(); }/* ww w .j a v a 2 s .c om*/ }).value(new DerivativeStructure(1, 2, 0, 0.0)); DerivativeStructure y = differentiator.differentiate(new UnivariateFunction() { public double value(double t) { return pv.shiftedBy(t).getPosition().normalize().getY().getValue(); } }).value(new DerivativeStructure(1, 2, 0, 0.0)); DerivativeStructure z = differentiator.differentiate(new UnivariateFunction() { public double value(double t) { return pv.shiftedBy(t).getPosition().normalize().getZ().getValue(); } }).value(new DerivativeStructure(1, 2, 0, 0.0)); FieldPVCoordinates<DerivativeStructure> normalized = pv.normalize(); Assert.assertEquals(x.getValue(), normalized.getPosition().getX().getValue(), 1.0e-16); Assert.assertEquals(y.getValue(), normalized.getPosition().getY().getValue(), 1.0e-16); Assert.assertEquals(z.getValue(), normalized.getPosition().getZ().getValue(), 1.0e-16); Assert.assertEquals(x.getPartialDerivative(1), normalized.getVelocity().getX().getValue(), 3.0e-13); Assert.assertEquals(y.getPartialDerivative(1), normalized.getVelocity().getY().getValue(), 3.0e-13); Assert.assertEquals(z.getPartialDerivative(1), normalized.getVelocity().getZ().getValue(), 3.0e-13); Assert.assertEquals(x.getPartialDerivative(2), normalized.getAcceleration().getX().getValue(), 6.0e-10); Assert.assertEquals(y.getPartialDerivative(2), normalized.getAcceleration().getY().getValue(), 6.0e-10); Assert.assertEquals(z.getPartialDerivative(2), normalized.getAcceleration().getZ().getValue(), 6.0e-10); } }
From source file:org.orekit.utils.IERSConventionsTest.java
private void checkDerivative(final TimeFunction<DerivativeStructure> function, final AbsoluteDate date, final double span, final double sampleStep, final double h, final double tolerance) { UnivariateDifferentiableFunction differentiated = new FiniteDifferencesDifferentiator(4, h) .differentiate(new UnivariateFunction() { @Override/* w w w . j a va2 s .co m*/ public double value(final double dt) { return function.value(date.shiftedBy(dt)).getValue(); } }); for (double dt = 0; dt < span; dt += sampleStep) { DerivativeStructure yRef = differentiated.value(new DerivativeStructure(1, 1, 0, dt)); DerivativeStructure y = function.value(date.shiftedBy(dt)); Assert.assertEquals(yRef.getPartialDerivative(1), y.getPartialDerivative(1), tolerance); } }
From source file:org.orekit.utils.PVCoordinatesTest.java
@Test public void testNormalize() { RandomGenerator generator = new Well19937a(0xb2011ffd25412067l); FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(5, 1.0e-3); for (int i = 0; i < 200; ++i) { final PVCoordinates pv = randomPVCoordinates(generator, 1e6, 1e3, 1.0); DerivativeStructure x = differentiator.differentiate(new UnivariateFunction() { public double value(double t) { return pv.shiftedBy(t).getPosition().normalize().getX(); }/* w w w. j av a 2s .c om*/ }).value(new DerivativeStructure(1, 2, 0, 0.0)); DerivativeStructure y = differentiator.differentiate(new UnivariateFunction() { public double value(double t) { return pv.shiftedBy(t).getPosition().normalize().getY(); } }).value(new DerivativeStructure(1, 2, 0, 0.0)); DerivativeStructure z = differentiator.differentiate(new UnivariateFunction() { public double value(double t) { return pv.shiftedBy(t).getPosition().normalize().getZ(); } }).value(new DerivativeStructure(1, 2, 0, 0.0)); PVCoordinates normalized = pv.normalize(); Assert.assertEquals(x.getValue(), normalized.getPosition().getX(), 1.0e-16); Assert.assertEquals(y.getValue(), normalized.getPosition().getY(), 1.0e-16); Assert.assertEquals(z.getValue(), normalized.getPosition().getZ(), 1.0e-16); Assert.assertEquals(x.getPartialDerivative(1), normalized.getVelocity().getX(), 3.0e-13); Assert.assertEquals(y.getPartialDerivative(1), normalized.getVelocity().getY(), 3.0e-13); Assert.assertEquals(z.getPartialDerivative(1), normalized.getVelocity().getZ(), 3.0e-13); Assert.assertEquals(x.getPartialDerivative(2), normalized.getAcceleration().getX(), 6.0e-10); Assert.assertEquals(y.getPartialDerivative(2), normalized.getAcceleration().getY(), 6.0e-10); Assert.assertEquals(z.getPartialDerivative(2), normalized.getAcceleration().getZ(), 6.0e-10); } }
From source file:org.orekit.utils.PVCoordinatesTest.java
@Test public void testCrossProduct() { RandomGenerator generator = new Well19937a(0x85c592b3be733d23l); FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(5, 1.0e-3); for (int i = 0; i < 200; ++i) { final PVCoordinates pv1 = randomPVCoordinates(generator, 1.0, 1.0, 1.0); final PVCoordinates pv2 = randomPVCoordinates(generator, 1.0, 1.0, 1.0); DerivativeStructure x = differentiator.differentiate(new UnivariateFunction() { public double value(double t) { return Vector3D.crossProduct(pv1.shiftedBy(t).getPosition(), pv2.shiftedBy(t).getPosition()) .getX();/*from w w w . ja v a 2 s . co m*/ } }).value(new DerivativeStructure(1, 2, 0, 0.0)); DerivativeStructure y = differentiator.differentiate(new UnivariateFunction() { public double value(double t) { return Vector3D.crossProduct(pv1.shiftedBy(t).getPosition(), pv2.shiftedBy(t).getPosition()) .getY(); } }).value(new DerivativeStructure(1, 2, 0, 0.0)); DerivativeStructure z = differentiator.differentiate(new UnivariateFunction() { public double value(double t) { return Vector3D.crossProduct(pv1.shiftedBy(t).getPosition(), pv2.shiftedBy(t).getPosition()) .getZ(); } }).value(new DerivativeStructure(1, 2, 0, 0.0)); PVCoordinates product = PVCoordinates.crossProduct(pv1, pv2); Assert.assertEquals(x.getValue(), product.getPosition().getX(), 1.0e-16); Assert.assertEquals(y.getValue(), product.getPosition().getY(), 1.0e-16); Assert.assertEquals(z.getValue(), product.getPosition().getZ(), 1.0e-16); Assert.assertEquals(x.getPartialDerivative(1), product.getVelocity().getX(), 9.0e-10); Assert.assertEquals(y.getPartialDerivative(1), product.getVelocity().getY(), 9.0e-10); Assert.assertEquals(z.getPartialDerivative(1), product.getVelocity().getZ(), 9.0e-10); Assert.assertEquals(x.getPartialDerivative(2), product.getAcceleration().getX(), 3.0e-9); Assert.assertEquals(y.getPartialDerivative(2), product.getAcceleration().getY(), 3.0e-9); Assert.assertEquals(z.getPartialDerivative(2), product.getAcceleration().getZ(), 3.0e-9); } }
From source file:org.orekit.utils.SecularAndHarmonicTest.java
private SpacecraftState findLatitudeCrossing(final double latitude, final AbsoluteDate guessDate, final AbsoluteDate endDate, final double shift, final double maxShift, final Propagator propagator) throws OrekitException, NoBracketingException { // function evaluating to 0 at latitude crossings final UnivariateFunction latitudeFunction = new UnivariateFunction() { /** {@inheritDoc} */ public double value(double x) { try { final SpacecraftState state = propagator.propagate(guessDate.shiftedBy(x)); final Vector3D position = state.getPVCoordinates(earth.getBodyFrame()).getPosition(); final GeodeticPoint point = earth.transform(position, earth.getBodyFrame(), state.getDate()); return point.getLatitude() - latitude; } catch (OrekitException oe) { throw new RuntimeException(oe); }// w w w. j av a 2s . co m } }; // try to bracket the encounter double span; if (guessDate.shiftedBy(shift).compareTo(endDate) > 0) { // Take a 1e-3 security margin span = endDate.durationFrom(guessDate) - 1e-3; } else { span = shift; } while (!UnivariateSolverUtils.isBracketing(latitudeFunction, -span, span)) { if (2 * span > maxShift) { // let the Apache Commons Math exception be thrown UnivariateSolverUtils.verifyBracketing(latitudeFunction, -span, span); } else if (guessDate.shiftedBy(2 * span).compareTo(endDate) > 0) { // Out of range : return null; } // expand the search interval span *= 2; } // find the encounter in the bracketed interval final BaseUnivariateSolver<UnivariateFunction> solver = new BracketingNthOrderBrentSolver(0.1, 5); final double dt = solver.solve(1000, latitudeFunction, -span, span); return propagator.propagate(guessDate.shiftedBy(dt)); }