Examples with BrentSolver org.apache.commons.math3.analysis.solvers.BrentSolver used on opensource projects

Example 1 with BrentSolver

use of org.apache.commons.math3.analysis.solvers.BrentSolver in project gatk-protected by broadinstitute.

the class RobustBrentSolver method doSolve.

@Override
protected double doSolve() throws TooManyEvaluationsException, NoBracketingException {
    final double min = getMin();
    final double max = getMax();
    final double[] xSearchGrid = createHybridSearchGrid(min, max, numBisections, depth);
    final double[] fSearchGrid = Arrays.stream(xSearchGrid).map(this::computeObjectiveValue).toArray();
    /* find bracketing intervals on the search grid */
    final List<Bracket> bracketsList = detectBrackets(xSearchGrid, fSearchGrid);
    if (bracketsList.isEmpty()) {
        throw new NoBracketingException(min, max, fSearchGrid[0], fSearchGrid[fSearchGrid.length - 1]);
    }
    final BrentSolver solver = new BrentSolver(getRelativeAccuracy(), getAbsoluteAccuracy(), getFunctionValueAccuracy());
    final List<Double> roots = bracketsList.stream().map(b -> solver.solve(getMaxEvaluations(), this::computeObjectiveValue, b.min, b.max, 0.5 * (b.min + b.max))).collect(Collectors.toList());
    if (roots.size() == 1 || meritFunc == null) {
        return roots.get(0);
    }
    final double[] merits = roots.stream().mapToDouble(meritFunc::value).toArray();
    final int bestRootIndex = IntStream.range(0, roots.size()).boxed().max((i, j) -> (int) (merits[i] - merits[j])).get();
    return roots.get(bestRootIndex);
}

Example 2 with BrentSolver

use of org.apache.commons.math3.analysis.solvers.BrentSolver in project gatk-protected by broadinstitute.

the class RobustBrentSolverUnitTest method simpleTest.

/**
     * Test on a 4th degree polynomial with 4 real roots at x = 0, 1, 2, 3. This objective function is positive for
     * large enough positive and negative values of its arguments. Therefore, the simple Brent solver complains that
     * the search interval does not bracket a root. The robust Brent solver, however, subdivides the given search
     * interval and finds a bracketing sub-interval.
     *
     * The "best" root according to the given merit function (set to the anti-derivative of the objective function)
     * is in fact the one at x = 0. We require the robust solver to output x = 0, and the simple solver to fail.
     */
@Test
public void simpleTest() {
    final UnivariateFunction objFunc = x -> 30 * x * (x - 1) * (x - 2) * (x - 3);
    final UnivariateFunction meritFunc = x -> 6 * FastMath.pow(x, 5) - 45 * FastMath.pow(x, 4) + 110 * FastMath.pow(x, 3) - 90 * FastMath.pow(x, 2);
    final RobustBrentSolver solverRobust = new RobustBrentSolver(DEF_REL_ACC, DEF_REL_ACC, DEF_F_ACC, meritFunc, 4, 1);
    final BrentSolver solverSimple = new BrentSolver(DEF_REL_ACC, DEF_REL_ACC, DEF_F_ACC);
    final double xRobust = solverRobust.solve(100, objFunc, -1, 4);
    Assert.assertEquals(xRobust, 0, DEF_ABS_ACC);
    boolean simpleSolverFails = false;
    try {
        /* this will fail */
        solverSimple.solve(100, objFunc, -1, 4);
    } catch (final NoBracketingException ex) {
        simpleSolverFails = true;
    }
    Assert.assertTrue(simpleSolverFails);
}

Example 3 with BrentSolver

use of org.apache.commons.math3.analysis.solvers.BrentSolver in project gatk by broadinstitute.

the class RobustBrentSolver method doSolve.

@Override
protected double doSolve() throws TooManyEvaluationsException, NoBracketingException {
    final double min = getMin();
    final double max = getMax();
    final double[] xSearchGrid = createHybridSearchGrid(min, max, numBisections, depth);
    final double[] fSearchGrid = Arrays.stream(xSearchGrid).map(this::computeObjectiveValue).toArray();
    /* find bracketing intervals on the search grid */
    final List<Bracket> bracketsList = detectBrackets(xSearchGrid, fSearchGrid);
    if (bracketsList.isEmpty()) {
        throw new NoBracketingException(min, max, fSearchGrid[0], fSearchGrid[fSearchGrid.length - 1]);
    }
    final BrentSolver solver = new BrentSolver(getRelativeAccuracy(), getAbsoluteAccuracy(), getFunctionValueAccuracy());
    final List<Double> roots = bracketsList.stream().map(b -> solver.solve(getMaxEvaluations(), this::computeObjectiveValue, b.min, b.max, 0.5 * (b.min + b.max))).collect(Collectors.toList());
    if (roots.size() == 1 || meritFunc == null) {
        return roots.get(0);
    }
    final double[] merits = roots.stream().mapToDouble(meritFunc::value).toArray();
    final int bestRootIndex = IntStream.range(0, roots.size()).boxed().max((i, j) -> (int) (merits[i] - merits[j])).get();
    return roots.get(bestRootIndex);
}

Example 4 with BrentSolver

use of org.apache.commons.math3.analysis.solvers.BrentSolver in project gatk by broadinstitute.

the class RobustBrentSolverUnitTest method simpleTest.

/**
     * Test on a 4th degree polynomial with 4 real roots at x = 0, 1, 2, 3. This objective function is positive for
     * large enough positive and negative values of its arguments. Therefore, the simple Brent solver complains that
     * the search interval does not bracket a root. The robust Brent solver, however, subdivides the given search
     * interval and finds a bracketing sub-interval.
     *
     * The "best" root according to the given merit function (set to the anti-derivative of the objective function)
     * is in fact the one at x = 0. We require the robust solver to output x = 0, and the simple solver to fail.
     */
@Test
public void simpleTest() {
    final UnivariateFunction objFunc = x -> 30 * x * (x - 1) * (x - 2) * (x - 3);
    final UnivariateFunction meritFunc = x -> 6 * FastMath.pow(x, 5) - 45 * FastMath.pow(x, 4) + 110 * FastMath.pow(x, 3) - 90 * FastMath.pow(x, 2);
    final RobustBrentSolver solverRobust = new RobustBrentSolver(DEF_REL_ACC, DEF_REL_ACC, DEF_F_ACC, meritFunc, 4, 1);
    final BrentSolver solverSimple = new BrentSolver(DEF_REL_ACC, DEF_REL_ACC, DEF_F_ACC);
    final double xRobust = solverRobust.solve(100, objFunc, -1, 4);
    Assert.assertEquals(xRobust, 0, DEF_ABS_ACC);
    boolean simpleSolverFails = false;
    try {
        /* this will fail */
        solverSimple.solve(100, objFunc, -1, 4);
    } catch (final NoBracketingException ex) {
        simpleSolverFails = true;
    }
    Assert.assertTrue(simpleSolverFails);
}