MullerSolver extends AbstractUnivariateSolver { private double solve(double min, double max, double fMin, double fMax) throws TooManyEvaluationsException { final double relativeAccuracy = getRelativeAccuracy(); final double absoluteAccuracy = getAbsoluteAccuracy(); final double functionValueAccuracy = getFunctionValueAccuracy(); double x0 = min; double y0 = fMin; double x2 = max; double y2 = fMax; double x1 = 0.5 * (x0 + x2); double y1 = computeObjectiveValue(x1); double oldx = Double.POSITIVE_INFINITY; while (true) { final double d01 = (y1 - y0) / (x1 - x0); final double d12 = (y2 - y1) / (x2 - x1); final double d012 = (d12 - d01) / (x2 - x0); final double c1 = d01 + (x1 - x0) * d012; final double delta = c1 * c1 - 4 * y1 * d012; final double xplus = x1 + (-2.0 * y1) / (c1 + FastMath.sqrt(delta)); final double xminus = x1 + (-2.0 * y1) / (c1 - FastMath.sqrt(delta)); final double x = isSequence(x0, xplus, x2) ? xplus : xminus; final double y = computeObjectiveValue(x); final double tolerance = FastMath.max(relativeAccuracy * FastMath.abs(x), absoluteAccuracy); if (FastMath.abs(x - oldx) <= tolerance || FastMath.abs(y) <= functionValueAccuracy) { return x; } boolean bisect = (x < x1 && (x1 - x0) > 0.95 * (x2 - x0)) || (x > x1 && (x2 - x1) > 0.95 * (x2 - x0)) || (x == x1); if (!bisect) { x0 = x < x1 ? x0 : x1; y0 = x < x1 ? y0 : y1; x2 = x > x1 ? x2 : x1; y2 = x > x1 ? y2 : y1; x1 = x; y1 = y; oldx = x; } else { double xm = 0.5 * (x0 + x2); double ym = computeObjectiveValue(xm); if (FastMath.signum(y0) + FastMath.signum(ym) == 0.0) { x2 = xm; y2 = ym; } else { x0 = xm; y0 = ym; } x1 = 0.5 * (x0 + x2); y1 = computeObjectiveValue(x1); oldx = Double.POSITIVE_INFINITY; } } } MullerSolver();  MullerSolver(double absoluteAccuracy);  MullerSolver(double relativeAccuracy,
                        double absoluteAccuracy);   }
@Test public void testSinFunction() { UnivariateFunction f = new Sin(); UnivariateSolver solver = new MullerSolver(); double min, max, expected, result, tolerance; min = 3.0; max = 4.0; expected = FastMath.PI; tolerance = FastMath.max(solver.getAbsoluteAccuracy(), FastMath.abs(expected * solver.getRelativeAccuracy())); result = solver.solve(100, f, min, max); Assert.assertEquals(expected, result, tolerance); min = -1.0; max = 1.5; expected = 0.0; tolerance = FastMath.max(solver.getAbsoluteAccuracy(), FastMath.abs(expected * solver.getRelativeAccuracy())); result = solver.solve(100, f, min, max); Assert.assertEquals(expected, result, tolerance); }
@Test public void testQuinticFunction() { UnivariateFunction f = new QuinticFunction(); UnivariateSolver solver = new MullerSolver(); double min, max, expected, result, tolerance; min = -0.4; max = 0.2; expected = 0.0; tolerance = FastMath.max(solver.getAbsoluteAccuracy(), FastMath.abs(expected * solver.getRelativeAccuracy())); result = solver.solve(100, f, min, max); Assert.assertEquals(expected, result, tolerance); min = 0.75; max = 1.5; expected = 1.0; tolerance = FastMath.max(solver.getAbsoluteAccuracy(), FastMath.abs(expected * solver.getRelativeAccuracy())); result = solver.solve(100, f, min, max); Assert.assertEquals(expected, result, tolerance); min = -0.9; max = -0.2; expected = -0.5; tolerance = FastMath.max(solver.getAbsoluteAccuracy(), FastMath.abs(expected * solver.getRelativeAccuracy())); result = solver.solve(100, f, min, max); Assert.assertEquals(expected, result, tolerance); }
@Test public void testExpm1Function() { UnivariateFunction f = new Expm1(); UnivariateSolver solver = new MullerSolver(); double min, max, expected, result, tolerance; min = -1.0; max = 2.0; expected = 0.0; tolerance = FastMath.max(solver.getAbsoluteAccuracy(), FastMath.abs(expected * solver.getRelativeAccuracy())); result = solver.solve(100, f, min, max); Assert.assertEquals(expected, result, tolerance); min = -20.0; max = 10.0; expected = 0.0; tolerance = FastMath.max(solver.getAbsoluteAccuracy(), FastMath.abs(expected * solver.getRelativeAccuracy())); result = solver.solve(100, f, min, max); Assert.assertEquals(expected, result, tolerance); min = -50.0; max = 100.0; expected = 0.0; tolerance = FastMath.max(solver.getAbsoluteAccuracy(), FastMath.abs(expected * solver.getRelativeAccuracy())); result = solver.solve(100, f, min, max); Assert.assertEquals(expected, result, tolerance); }
@Test public void testParameters() { UnivariateFunction f = new Sin(); UnivariateSolver solver = new MullerSolver(); try { double root = solver.solve(100, f, 1, -1); System.out.println("root=" + root); Assert.fail("Expecting NumberIsTooLargeException - bad interval"); } catch (NumberIsTooLargeException ex) { } try { solver.solve(100, f, 2, 3); Assert.fail("Expecting NoBracketingException - no bracketing"); } catch (NoBracketingException ex) { } }