Example usage for org.apache.commons.math3.analysis.solvers BracketingNthOrderBrentSolver solve

Introduction

In this page you can find the example usage for org.apache.commons.math3.analysis.solvers BracketingNthOrderBrentSolver solve.

Prototype

public double solve(int maxEval, UnivariateFunction f, double min, double max, AllowedSolution allowedSolution)
        throws TooManyEvaluationsException, NumberIsTooLargeException, NoBracketingException 

Source Link

Usage

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 a  va 2 s  . c om
 * @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();
    }

}

---