You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by tn...@apache.org on 2015/02/16 23:40:09 UTC
[39/82] [partial] [math] Update for next development iteration:
commons-math4
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/BracketingNthOrderBrentSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/BracketingNthOrderBrentSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/BracketingNthOrderBrentSolver.java
deleted file mode 100644
index 4981592..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/BracketingNthOrderBrentSolver.java
+++ /dev/null
@@ -1,412 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-package org.apache.commons.math3.analysis.solvers;
-
-
-import org.apache.commons.math3.analysis.UnivariateFunction;
-import org.apache.commons.math3.exception.MathInternalError;
-import org.apache.commons.math3.exception.NoBracketingException;
-import org.apache.commons.math3.exception.NumberIsTooSmallException;
-import org.apache.commons.math3.exception.NumberIsTooLargeException;
-import org.apache.commons.math3.exception.TooManyEvaluationsException;
-import org.apache.commons.math3.util.FastMath;
-import org.apache.commons.math3.util.Precision;
-
-/**
- * This class implements a modification of the <a
- * href="http://mathworld.wolfram.com/BrentsMethod.html"> Brent algorithm</a>.
- * <p>
- * The changes with respect to the original Brent algorithm are:
- * <ul>
- * <li>the returned value is chosen in the current interval according
- * to user specified {@link AllowedSolution},</li>
- * <li>the maximal order for the invert polynomial root search is
- * user-specified instead of being invert quadratic only</li>
- * </ul>
- * </p>
- * The given interval must bracket the root.
- *
- */
-public class BracketingNthOrderBrentSolver
- extends AbstractUnivariateSolver
- implements BracketedUnivariateSolver<UnivariateFunction> {
-
- /** Default absolute accuracy. */
- private static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
-
- /** Default maximal order. */
- private static final int DEFAULT_MAXIMAL_ORDER = 5;
-
- /** Maximal aging triggering an attempt to balance the bracketing interval. */
- private static final int MAXIMAL_AGING = 2;
-
- /** Reduction factor for attempts to balance the bracketing interval. */
- private static final double REDUCTION_FACTOR = 1.0 / 16.0;
-
- /** Maximal order. */
- private final int maximalOrder;
-
- /** The kinds of solutions that the algorithm may accept. */
- private AllowedSolution allowed;
-
- /**
- * Construct a solver with default accuracy and maximal order (1e-6 and 5 respectively)
- */
- public BracketingNthOrderBrentSolver() {
- this(DEFAULT_ABSOLUTE_ACCURACY, DEFAULT_MAXIMAL_ORDER);
- }
-
- /**
- * Construct a solver.
- *
- * @param absoluteAccuracy Absolute accuracy.
- * @param maximalOrder maximal order.
- * @exception NumberIsTooSmallException if maximal order is lower than 2
- */
- public BracketingNthOrderBrentSolver(final double absoluteAccuracy,
- final int maximalOrder)
- throws NumberIsTooSmallException {
- super(absoluteAccuracy);
- if (maximalOrder < 2) {
- throw new NumberIsTooSmallException(maximalOrder, 2, true);
- }
- this.maximalOrder = maximalOrder;
- this.allowed = AllowedSolution.ANY_SIDE;
- }
-
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- * @param maximalOrder maximal order.
- * @exception NumberIsTooSmallException if maximal order is lower than 2
- */
- public BracketingNthOrderBrentSolver(final double relativeAccuracy,
- final double absoluteAccuracy,
- final int maximalOrder)
- throws NumberIsTooSmallException {
- super(relativeAccuracy, absoluteAccuracy);
- if (maximalOrder < 2) {
- throw new NumberIsTooSmallException(maximalOrder, 2, true);
- }
- this.maximalOrder = maximalOrder;
- this.allowed = AllowedSolution.ANY_SIDE;
- }
-
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- * @param functionValueAccuracy Function value accuracy.
- * @param maximalOrder maximal order.
- * @exception NumberIsTooSmallException if maximal order is lower than 2
- */
- public BracketingNthOrderBrentSolver(final double relativeAccuracy,
- final double absoluteAccuracy,
- final double functionValueAccuracy,
- final int maximalOrder)
- throws NumberIsTooSmallException {
- super(relativeAccuracy, absoluteAccuracy, functionValueAccuracy);
- if (maximalOrder < 2) {
- throw new NumberIsTooSmallException(maximalOrder, 2, true);
- }
- this.maximalOrder = maximalOrder;
- this.allowed = AllowedSolution.ANY_SIDE;
- }
-
- /** Get the maximal order.
- * @return maximal order
- */
- public int getMaximalOrder() {
- return maximalOrder;
- }
-
- /**
- * {@inheritDoc}
- */
- @Override
- protected double doSolve()
- throws TooManyEvaluationsException,
- NumberIsTooLargeException,
- NoBracketingException {
- // prepare arrays with the first points
- final double[] x = new double[maximalOrder + 1];
- final double[] y = new double[maximalOrder + 1];
- x[0] = getMin();
- x[1] = getStartValue();
- x[2] = getMax();
- verifySequence(x[0], x[1], x[2]);
-
- // evaluate initial guess
- y[1] = computeObjectiveValue(x[1]);
- if (Precision.equals(y[1], 0.0, 1)) {
- // return the initial guess if it is a perfect root.
- return x[1];
- }
-
- // evaluate first endpoint
- y[0] = computeObjectiveValue(x[0]);
- if (Precision.equals(y[0], 0.0, 1)) {
- // return the first endpoint if it is a perfect root.
- return x[0];
- }
-
- int nbPoints;
- int signChangeIndex;
- if (y[0] * y[1] < 0) {
-
- // reduce interval if it brackets the root
- nbPoints = 2;
- signChangeIndex = 1;
-
- } else {
-
- // evaluate second endpoint
- y[2] = computeObjectiveValue(x[2]);
- if (Precision.equals(y[2], 0.0, 1)) {
- // return the second endpoint if it is a perfect root.
- return x[2];
- }
-
- if (y[1] * y[2] < 0) {
- // use all computed point as a start sampling array for solving
- nbPoints = 3;
- signChangeIndex = 2;
- } else {
- throw new NoBracketingException(x[0], x[2], y[0], y[2]);
- }
-
- }
-
- // prepare a work array for inverse polynomial interpolation
- final double[] tmpX = new double[x.length];
-
- // current tightest bracketing of the root
- double xA = x[signChangeIndex - 1];
- double yA = y[signChangeIndex - 1];
- double absYA = FastMath.abs(yA);
- int agingA = 0;
- double xB = x[signChangeIndex];
- double yB = y[signChangeIndex];
- double absYB = FastMath.abs(yB);
- int agingB = 0;
-
- // search loop
- while (true) {
-
- // check convergence of bracketing interval
- final double xTol = getAbsoluteAccuracy() +
- getRelativeAccuracy() * FastMath.max(FastMath.abs(xA), FastMath.abs(xB));
- if (((xB - xA) <= xTol) || (FastMath.max(absYA, absYB) < getFunctionValueAccuracy())) {
- switch (allowed) {
- case ANY_SIDE :
- return absYA < absYB ? xA : xB;
- case LEFT_SIDE :
- return xA;
- case RIGHT_SIDE :
- return xB;
- case BELOW_SIDE :
- return (yA <= 0) ? xA : xB;
- case ABOVE_SIDE :
- return (yA < 0) ? xB : xA;
- default :
- // this should never happen
- throw new MathInternalError();
- }
- }
-
- // target for the next evaluation point
- double targetY;
- if (agingA >= MAXIMAL_AGING) {
- // we keep updating the high bracket, try to compensate this
- final int p = agingA - MAXIMAL_AGING;
- final double weightA = (1 << p) - 1;
- final double weightB = p + 1;
- targetY = (weightA * yA - weightB * REDUCTION_FACTOR * yB) / (weightA + weightB);
- } else if (agingB >= MAXIMAL_AGING) {
- // we keep updating the low bracket, try to compensate this
- final int p = agingB - MAXIMAL_AGING;
- final double weightA = p + 1;
- final double weightB = (1 << p) - 1;
- targetY = (weightB * yB - weightA * REDUCTION_FACTOR * yA) / (weightA + weightB);
- } else {
- // bracketing is balanced, try to find the root itself
- targetY = 0;
- }
-
- // make a few attempts to guess a root,
- double nextX;
- int start = 0;
- int end = nbPoints;
- do {
-
- // guess a value for current target, using inverse polynomial interpolation
- System.arraycopy(x, start, tmpX, start, end - start);
- nextX = guessX(targetY, tmpX, y, start, end);
-
- if (!((nextX > xA) && (nextX < xB))) {
- // the guessed root is not strictly inside of the tightest bracketing interval
-
- // the guessed root is either not strictly inside the interval or it
- // is a NaN (which occurs when some sampling points share the same y)
- // we try again with a lower interpolation order
- if (signChangeIndex - start >= end - signChangeIndex) {
- // we have more points before the sign change, drop the lowest point
- ++start;
- } else {
- // we have more points after sign change, drop the highest point
- --end;
- }
-
- // we need to do one more attempt
- nextX = Double.NaN;
-
- }
-
- } while (Double.isNaN(nextX) && (end - start > 1));
-
- if (Double.isNaN(nextX)) {
- // fall back to bisection
- nextX = xA + 0.5 * (xB - xA);
- start = signChangeIndex - 1;
- end = signChangeIndex;
- }
-
- // evaluate the function at the guessed root
- final double nextY = computeObjectiveValue(nextX);
- if (Precision.equals(nextY, 0.0, 1)) {
- // we have found an exact root, since it is not an approximation
- // we don't need to bother about the allowed solutions setting
- return nextX;
- }
-
- if ((nbPoints > 2) && (end - start != nbPoints)) {
-
- // we have been forced to ignore some points to keep bracketing,
- // they are probably too far from the root, drop them from now on
- nbPoints = end - start;
- System.arraycopy(x, start, x, 0, nbPoints);
- System.arraycopy(y, start, y, 0, nbPoints);
- signChangeIndex -= start;
-
- } else if (nbPoints == x.length) {
-
- // we have to drop one point in order to insert the new one
- nbPoints--;
-
- // keep the tightest bracketing interval as centered as possible
- if (signChangeIndex >= (x.length + 1) / 2) {
- // we drop the lowest point, we have to shift the arrays and the index
- System.arraycopy(x, 1, x, 0, nbPoints);
- System.arraycopy(y, 1, y, 0, nbPoints);
- --signChangeIndex;
- }
-
- }
-
- // insert the last computed point
- //(by construction, we know it lies inside the tightest bracketing interval)
- System.arraycopy(x, signChangeIndex, x, signChangeIndex + 1, nbPoints - signChangeIndex);
- x[signChangeIndex] = nextX;
- System.arraycopy(y, signChangeIndex, y, signChangeIndex + 1, nbPoints - signChangeIndex);
- y[signChangeIndex] = nextY;
- ++nbPoints;
-
- // update the bracketing interval
- if (nextY * yA <= 0) {
- // the sign change occurs before the inserted point
- xB = nextX;
- yB = nextY;
- absYB = FastMath.abs(yB);
- ++agingA;
- agingB = 0;
- } else {
- // the sign change occurs after the inserted point
- xA = nextX;
- yA = nextY;
- absYA = FastMath.abs(yA);
- agingA = 0;
- ++agingB;
-
- // update the sign change index
- signChangeIndex++;
-
- }
-
- }
-
- }
-
- /** Guess an x value by n<sup>th</sup> order inverse polynomial interpolation.
- * <p>
- * The x value is guessed by evaluating polynomial Q(y) at y = targetY, where Q
- * is built such that for all considered points (x<sub>i</sub>, y<sub>i</sub>),
- * Q(y<sub>i</sub>) = x<sub>i</sub>.
- * </p>
- * @param targetY target value for y
- * @param x reference points abscissas for interpolation,
- * note that this array <em>is</em> modified during computation
- * @param y reference points ordinates for interpolation
- * @param start start index of the points to consider (inclusive)
- * @param end end index of the points to consider (exclusive)
- * @return guessed root (will be a NaN if two points share the same y)
- */
- private double guessX(final double targetY, final double[] x, final double[] y,
- final int start, final int end) {
-
- // compute Q Newton coefficients by divided differences
- for (int i = start; i < end - 1; ++i) {
- final int delta = i + 1 - start;
- for (int j = end - 1; j > i; --j) {
- x[j] = (x[j] - x[j-1]) / (y[j] - y[j - delta]);
- }
- }
-
- // evaluate Q(targetY)
- double x0 = 0;
- for (int j = end - 1; j >= start; --j) {
- x0 = x[j] + x0 * (targetY - y[j]);
- }
-
- return x0;
-
- }
-
- /** {@inheritDoc} */
- public double solve(int maxEval, UnivariateFunction f, double min,
- double max, AllowedSolution allowedSolution)
- throws TooManyEvaluationsException,
- NumberIsTooLargeException,
- NoBracketingException {
- this.allowed = allowedSolution;
- return super.solve(maxEval, f, min, max);
- }
-
- /** {@inheritDoc} */
- public double solve(int maxEval, UnivariateFunction f, double min,
- double max, double startValue,
- AllowedSolution allowedSolution)
- throws TooManyEvaluationsException,
- NumberIsTooLargeException,
- NoBracketingException {
- this.allowed = allowedSolution;
- return super.solve(maxEval, f, min, max, startValue);
- }
-
-}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/BrentSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/BrentSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/BrentSolver.java
deleted file mode 100644
index 0cc8750..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/BrentSolver.java
+++ /dev/null
@@ -1,244 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-package org.apache.commons.math3.analysis.solvers;
-
-
-import org.apache.commons.math3.exception.NoBracketingException;
-import org.apache.commons.math3.exception.TooManyEvaluationsException;
-import org.apache.commons.math3.exception.NumberIsTooLargeException;
-import org.apache.commons.math3.util.FastMath;
-import org.apache.commons.math3.util.Precision;
-
-/**
- * This class implements the <a href="http://mathworld.wolfram.com/BrentsMethod.html">
- * Brent algorithm</a> for finding zeros of real univariate functions.
- * The function should be continuous but not necessarily smooth.
- * The {@code solve} method returns a zero {@code x} of the function {@code f}
- * in the given interval {@code [a, b]} to within a tolerance
- * {@code 2 eps abs(x) + t} where {@code eps} is the relative accuracy and
- * {@code t} is the absolute accuracy.
- * The given interval must bracket the root.
- * <p>
- * The reference implementation is given in chapter 4 of
- * <blockquote>
- * <b>Algorithms for Minimization Without Derivatives</b><br>
- * <em>Richard P. Brent</em><br>
- * Dover, 2002<br>
- * </blockquote>
- * </p>
- *
- * @see BaseAbstractUnivariateSolver
- */
-public class BrentSolver extends AbstractUnivariateSolver {
-
- /** Default absolute accuracy. */
- private static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
-
- /**
- * Construct a solver with default absolute accuracy (1e-6).
- */
- public BrentSolver() {
- this(DEFAULT_ABSOLUTE_ACCURACY);
- }
- /**
- * Construct a solver.
- *
- * @param absoluteAccuracy Absolute accuracy.
- */
- public BrentSolver(double absoluteAccuracy) {
- super(absoluteAccuracy);
- }
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- */
- public BrentSolver(double relativeAccuracy,
- double absoluteAccuracy) {
- super(relativeAccuracy, absoluteAccuracy);
- }
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- * @param functionValueAccuracy Function value accuracy.
- *
- * @see BaseAbstractUnivariateSolver#BaseAbstractUnivariateSolver(double,double,double)
- */
- public BrentSolver(double relativeAccuracy,
- double absoluteAccuracy,
- double functionValueAccuracy) {
- super(relativeAccuracy, absoluteAccuracy, functionValueAccuracy);
- }
-
- /**
- * {@inheritDoc}
- */
- @Override
- protected double doSolve()
- throws NoBracketingException,
- TooManyEvaluationsException,
- NumberIsTooLargeException {
- double min = getMin();
- double max = getMax();
- final double initial = getStartValue();
- final double functionValueAccuracy = getFunctionValueAccuracy();
-
- verifySequence(min, initial, max);
-
- // Return the initial guess if it is good enough.
- double yInitial = computeObjectiveValue(initial);
- if (FastMath.abs(yInitial) <= functionValueAccuracy) {
- return initial;
- }
-
- // Return the first endpoint if it is good enough.
- double yMin = computeObjectiveValue(min);
- if (FastMath.abs(yMin) <= functionValueAccuracy) {
- return min;
- }
-
- // Reduce interval if min and initial bracket the root.
- if (yInitial * yMin < 0) {
- return brent(min, initial, yMin, yInitial);
- }
-
- // Return the second endpoint if it is good enough.
- double yMax = computeObjectiveValue(max);
- if (FastMath.abs(yMax) <= functionValueAccuracy) {
- return max;
- }
-
- // Reduce interval if initial and max bracket the root.
- if (yInitial * yMax < 0) {
- return brent(initial, max, yInitial, yMax);
- }
-
- throw new NoBracketingException(min, max, yMin, yMax);
- }
-
- /**
- * Search for a zero inside the provided interval.
- * This implementation is based on the algorithm described at page 58 of
- * the book
- * <blockquote>
- * <b>Algorithms for Minimization Without Derivatives</b>
- * <it>Richard P. Brent</it>
- * Dover 0-486-41998-3
- * </blockquote>
- *
- * @param lo Lower bound of the search interval.
- * @param hi Higher bound of the search interval.
- * @param fLo Function value at the lower bound of the search interval.
- * @param fHi Function value at the higher bound of the search interval.
- * @return the value where the function is zero.
- */
- private double brent(double lo, double hi,
- double fLo, double fHi) {
- double a = lo;
- double fa = fLo;
- double b = hi;
- double fb = fHi;
- double c = a;
- double fc = fa;
- double d = b - a;
- double e = d;
-
- final double t = getAbsoluteAccuracy();
- final double eps = getRelativeAccuracy();
-
- while (true) {
- if (FastMath.abs(fc) < FastMath.abs(fb)) {
- a = b;
- b = c;
- c = a;
- fa = fb;
- fb = fc;
- fc = fa;
- }
-
- final double tol = 2 * eps * FastMath.abs(b) + t;
- final double m = 0.5 * (c - b);
-
- if (FastMath.abs(m) <= tol ||
- Precision.equals(fb, 0)) {
- return b;
- }
- if (FastMath.abs(e) < tol ||
- FastMath.abs(fa) <= FastMath.abs(fb)) {
- // Force bisection.
- d = m;
- e = d;
- } else {
- double s = fb / fa;
- double p;
- double q;
- // The equality test (a == c) is intentional,
- // it is part of the original Brent's method and
- // it should NOT be replaced by proximity test.
- if (a == c) {
- // Linear interpolation.
- p = 2 * m * s;
- q = 1 - s;
- } else {
- // Inverse quadratic interpolation.
- q = fa / fc;
- final double r = fb / fc;
- p = s * (2 * m * q * (q - r) - (b - a) * (r - 1));
- q = (q - 1) * (r - 1) * (s - 1);
- }
- if (p > 0) {
- q = -q;
- } else {
- p = -p;
- }
- s = e;
- e = d;
- if (p >= 1.5 * m * q - FastMath.abs(tol * q) ||
- p >= FastMath.abs(0.5 * s * q)) {
- // Inverse quadratic interpolation gives a value
- // in the wrong direction, or progress is slow.
- // Fall back to bisection.
- d = m;
- e = d;
- } else {
- d = p / q;
- }
- }
- a = b;
- fa = fb;
-
- if (FastMath.abs(d) > tol) {
- b += d;
- } else if (m > 0) {
- b += tol;
- } else {
- b -= tol;
- }
- fb = computeObjectiveValue(b);
- if ((fb > 0 && fc > 0) ||
- (fb <= 0 && fc <= 0)) {
- c = a;
- fc = fa;
- d = b - a;
- e = d;
- }
- }
- }
-}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/DifferentiableUnivariateSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/DifferentiableUnivariateSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/DifferentiableUnivariateSolver.java
deleted file mode 100644
index b9ae158..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/DifferentiableUnivariateSolver.java
+++ /dev/null
@@ -1,30 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-package org.apache.commons.math3.analysis.solvers;
-
-import org.apache.commons.math3.analysis.DifferentiableUnivariateFunction;
-
-
-/**
- * Interface for (univariate real) rootfinding algorithms.
- * Implementations will search for only one zero in the given interval.
- *
- * @deprecated as of 3.1, replaced by {@link UnivariateDifferentiableSolver}
- */
-@Deprecated
-public interface DifferentiableUnivariateSolver
- extends BaseUnivariateSolver<DifferentiableUnivariateFunction> {}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/IllinoisSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/IllinoisSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/IllinoisSolver.java
deleted file mode 100644
index bd3bc71..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/IllinoisSolver.java
+++ /dev/null
@@ -1,82 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-
-package org.apache.commons.math3.analysis.solvers;
-
-
-/**
- * Implements the <em>Illinois</em> method for root-finding (approximating
- * a zero of a univariate real function). It is a modified
- * {@link RegulaFalsiSolver <em>Regula Falsi</em>} method.
- *
- * <p>Like the <em>Regula Falsi</em> method, convergence is guaranteed by
- * maintaining a bracketed solution. The <em>Illinois</em> method however,
- * should converge much faster than the original <em>Regula Falsi</em>
- * method. Furthermore, this implementation of the <em>Illinois</em> method
- * should not suffer from the same implementation issues as the <em>Regula
- * Falsi</em> method, which may fail to convergence in certain cases.</p>
- *
- * <p>The <em>Illinois</em> method assumes that the function is continuous,
- * but not necessarily smooth.</p>
- *
- * <p>Implementation based on the following article: M. Dowell and P. Jarratt,
- * <em>A modified regula falsi method for computing the root of an
- * equation</em>, BIT Numerical Mathematics, volume 11, number 2,
- * pages 168-174, Springer, 1971.</p>
- *
- * @since 3.0
- */
-public class IllinoisSolver extends BaseSecantSolver {
-
- /** Construct a solver with default accuracy (1e-6). */
- public IllinoisSolver() {
- super(DEFAULT_ABSOLUTE_ACCURACY, Method.ILLINOIS);
- }
-
- /**
- * Construct a solver.
- *
- * @param absoluteAccuracy Absolute accuracy.
- */
- public IllinoisSolver(final double absoluteAccuracy) {
- super(absoluteAccuracy, Method.ILLINOIS);
- }
-
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- */
- public IllinoisSolver(final double relativeAccuracy,
- final double absoluteAccuracy) {
- super(relativeAccuracy, absoluteAccuracy, Method.ILLINOIS);
- }
-
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- * @param functionValueAccuracy Maximum function value error.
- */
- public IllinoisSolver(final double relativeAccuracy,
- final double absoluteAccuracy,
- final double functionValueAccuracy) {
- super(relativeAccuracy, absoluteAccuracy, functionValueAccuracy, Method.PEGASUS);
- }
-}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/LaguerreSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/LaguerreSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/LaguerreSolver.java
deleted file mode 100644
index c127b42..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/LaguerreSolver.java
+++ /dev/null
@@ -1,390 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-package org.apache.commons.math3.analysis.solvers;
-
-import org.apache.commons.math3.complex.Complex;
-import org.apache.commons.math3.complex.ComplexUtils;
-import org.apache.commons.math3.analysis.polynomials.PolynomialFunction;
-import org.apache.commons.math3.exception.NoBracketingException;
-import org.apache.commons.math3.exception.NullArgumentException;
-import org.apache.commons.math3.exception.NoDataException;
-import org.apache.commons.math3.exception.TooManyEvaluationsException;
-import org.apache.commons.math3.exception.NumberIsTooLargeException;
-import org.apache.commons.math3.exception.util.LocalizedFormats;
-import org.apache.commons.math3.util.FastMath;
-
-/**
- * Implements the <a href="http://mathworld.wolfram.com/LaguerresMethod.html">
- * Laguerre's Method</a> for root finding of real coefficient polynomials.
- * For reference, see
- * <blockquote>
- * <b>A First Course in Numerical Analysis</b><br>
- * ISBN 048641454X, chapter 8.<br>
- * </blockquote>
- * Laguerre's method is global in the sense that it can start with any initial
- * approximation and be able to solve all roots from that point.
- * The algorithm requires a bracketing condition.
- *
- * @since 1.2
- */
-public class LaguerreSolver extends AbstractPolynomialSolver {
- /** Default absolute accuracy. */
- private static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
- /** Complex solver. */
- private final ComplexSolver complexSolver = new ComplexSolver();
-
- /**
- * Construct a solver with default accuracy (1e-6).
- */
- public LaguerreSolver() {
- this(DEFAULT_ABSOLUTE_ACCURACY);
- }
- /**
- * Construct a solver.
- *
- * @param absoluteAccuracy Absolute accuracy.
- */
- public LaguerreSolver(double absoluteAccuracy) {
- super(absoluteAccuracy);
- }
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- */
- public LaguerreSolver(double relativeAccuracy,
- double absoluteAccuracy) {
- super(relativeAccuracy, absoluteAccuracy);
- }
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- * @param functionValueAccuracy Function value accuracy.
- */
- public LaguerreSolver(double relativeAccuracy,
- double absoluteAccuracy,
- double functionValueAccuracy) {
- super(relativeAccuracy, absoluteAccuracy, functionValueAccuracy);
- }
-
- /**
- * {@inheritDoc}
- */
- @Override
- public double doSolve()
- throws TooManyEvaluationsException,
- NumberIsTooLargeException,
- NoBracketingException {
- final double min = getMin();
- final double max = getMax();
- final double initial = getStartValue();
- final double functionValueAccuracy = getFunctionValueAccuracy();
-
- verifySequence(min, initial, max);
-
- // Return the initial guess if it is good enough.
- final double yInitial = computeObjectiveValue(initial);
- if (FastMath.abs(yInitial) <= functionValueAccuracy) {
- return initial;
- }
-
- // Return the first endpoint if it is good enough.
- final double yMin = computeObjectiveValue(min);
- if (FastMath.abs(yMin) <= functionValueAccuracy) {
- return min;
- }
-
- // Reduce interval if min and initial bracket the root.
- if (yInitial * yMin < 0) {
- return laguerre(min, initial, yMin, yInitial);
- }
-
- // Return the second endpoint if it is good enough.
- final double yMax = computeObjectiveValue(max);
- if (FastMath.abs(yMax) <= functionValueAccuracy) {
- return max;
- }
-
- // Reduce interval if initial and max bracket the root.
- if (yInitial * yMax < 0) {
- return laguerre(initial, max, yInitial, yMax);
- }
-
- throw new NoBracketingException(min, max, yMin, yMax);
- }
-
- /**
- * Find a real root in the given interval.
- *
- * Despite the bracketing condition, the root returned by
- * {@link LaguerreSolver.ComplexSolver#solve(Complex[],Complex)} may
- * not be a real zero inside {@code [min, max]}.
- * For example, <code>p(x) = x<sup>3</sup> + 1,</code>
- * with {@code min = -2}, {@code max = 2}, {@code initial = 0}.
- * When it occurs, this code calls
- * {@link LaguerreSolver.ComplexSolver#solveAll(Complex[],Complex)}
- * in order to obtain all roots and picks up one real root.
- *
- * @param lo Lower bound of the search interval.
- * @param hi Higher bound of the search interval.
- * @param fLo Function value at the lower bound of the search interval.
- * @param fHi Function value at the higher bound of the search interval.
- * @return the point at which the function value is zero.
- * @deprecated This method should not be part of the public API: It will
- * be made private in version 4.0.
- */
- @Deprecated
- public double laguerre(double lo, double hi,
- double fLo, double fHi) {
- final Complex c[] = ComplexUtils.convertToComplex(getCoefficients());
-
- final Complex initial = new Complex(0.5 * (lo + hi), 0);
- final Complex z = complexSolver.solve(c, initial);
- if (complexSolver.isRoot(lo, hi, z)) {
- return z.getReal();
- } else {
- double r = Double.NaN;
- // Solve all roots and select the one we are seeking.
- Complex[] root = complexSolver.solveAll(c, initial);
- for (int i = 0; i < root.length; i++) {
- if (complexSolver.isRoot(lo, hi, root[i])) {
- r = root[i].getReal();
- break;
- }
- }
- return r;
- }
- }
-
- /**
- * Find all complex roots for the polynomial with the given
- * coefficients, starting from the given initial value.
- * <br/>
- * Note: This method is not part of the API of {@link BaseUnivariateSolver}.
- *
- * @param coefficients Polynomial coefficients.
- * @param initial Start value.
- * @return the point at which the function value is zero.
- * @throws org.apache.commons.math3.exception.TooManyEvaluationsException
- * if the maximum number of evaluations is exceeded.
- * @throws NullArgumentException if the {@code coefficients} is
- * {@code null}.
- * @throws NoDataException if the {@code coefficients} array is empty.
- * @since 3.1
- */
- public Complex[] solveAllComplex(double[] coefficients,
- double initial)
- throws NullArgumentException,
- NoDataException,
- TooManyEvaluationsException {
- setup(Integer.MAX_VALUE,
- new PolynomialFunction(coefficients),
- Double.NEGATIVE_INFINITY,
- Double.POSITIVE_INFINITY,
- initial);
- return complexSolver.solveAll(ComplexUtils.convertToComplex(coefficients),
- new Complex(initial, 0d));
- }
-
- /**
- * Find a complex root for the polynomial with the given coefficients,
- * starting from the given initial value.
- * <br/>
- * Note: This method is not part of the API of {@link BaseUnivariateSolver}.
- *
- * @param coefficients Polynomial coefficients.
- * @param initial Start value.
- * @return the point at which the function value is zero.
- * @throws org.apache.commons.math3.exception.TooManyEvaluationsException
- * if the maximum number of evaluations is exceeded.
- * @throws NullArgumentException if the {@code coefficients} is
- * {@code null}.
- * @throws NoDataException if the {@code coefficients} array is empty.
- * @since 3.1
- */
- public Complex solveComplex(double[] coefficients,
- double initial)
- throws NullArgumentException,
- NoDataException,
- TooManyEvaluationsException {
- setup(Integer.MAX_VALUE,
- new PolynomialFunction(coefficients),
- Double.NEGATIVE_INFINITY,
- Double.POSITIVE_INFINITY,
- initial);
- return complexSolver.solve(ComplexUtils.convertToComplex(coefficients),
- new Complex(initial, 0d));
- }
-
- /**
- * Class for searching all (complex) roots.
- */
- private class ComplexSolver {
- /**
- * Check whether the given complex root is actually a real zero
- * in the given interval, within the solver tolerance level.
- *
- * @param min Lower bound for the interval.
- * @param max Upper bound for the interval.
- * @param z Complex root.
- * @return {@code true} if z is a real zero.
- */
- public boolean isRoot(double min, double max, Complex z) {
- if (isSequence(min, z.getReal(), max)) {
- double tolerance = FastMath.max(getRelativeAccuracy() * z.abs(), getAbsoluteAccuracy());
- return (FastMath.abs(z.getImaginary()) <= tolerance) ||
- (z.abs() <= getFunctionValueAccuracy());
- }
- return false;
- }
-
- /**
- * Find all complex roots for the polynomial with the given
- * coefficients, starting from the given initial value.
- *
- * @param coefficients Polynomial coefficients.
- * @param initial Start value.
- * @return the point at which the function value is zero.
- * @throws org.apache.commons.math3.exception.TooManyEvaluationsException
- * if the maximum number of evaluations is exceeded.
- * @throws NullArgumentException if the {@code coefficients} is
- * {@code null}.
- * @throws NoDataException if the {@code coefficients} array is empty.
- */
- public Complex[] solveAll(Complex coefficients[], Complex initial)
- throws NullArgumentException,
- NoDataException,
- TooManyEvaluationsException {
- if (coefficients == null) {
- throw new NullArgumentException();
- }
- final int n = coefficients.length - 1;
- if (n == 0) {
- throw new NoDataException(LocalizedFormats.POLYNOMIAL);
- }
- // Coefficients for deflated polynomial.
- final Complex c[] = new Complex[n + 1];
- for (int i = 0; i <= n; i++) {
- c[i] = coefficients[i];
- }
-
- // Solve individual roots successively.
- final Complex root[] = new Complex[n];
- for (int i = 0; i < n; i++) {
- final Complex subarray[] = new Complex[n - i + 1];
- System.arraycopy(c, 0, subarray, 0, subarray.length);
- root[i] = solve(subarray, initial);
- // Polynomial deflation using synthetic division.
- Complex newc = c[n - i];
- Complex oldc = null;
- for (int j = n - i - 1; j >= 0; j--) {
- oldc = c[j];
- c[j] = newc;
- newc = oldc.add(newc.multiply(root[i]));
- }
- }
-
- return root;
- }
-
- /**
- * Find a complex root for the polynomial with the given coefficients,
- * starting from the given initial value.
- *
- * @param coefficients Polynomial coefficients.
- * @param initial Start value.
- * @return the point at which the function value is zero.
- * @throws org.apache.commons.math3.exception.TooManyEvaluationsException
- * if the maximum number of evaluations is exceeded.
- * @throws NullArgumentException if the {@code coefficients} is
- * {@code null}.
- * @throws NoDataException if the {@code coefficients} array is empty.
- */
- public Complex solve(Complex coefficients[], Complex initial)
- throws NullArgumentException,
- NoDataException,
- TooManyEvaluationsException {
- if (coefficients == null) {
- throw new NullArgumentException();
- }
-
- final int n = coefficients.length - 1;
- if (n == 0) {
- throw new NoDataException(LocalizedFormats.POLYNOMIAL);
- }
-
- final double absoluteAccuracy = getAbsoluteAccuracy();
- final double relativeAccuracy = getRelativeAccuracy();
- final double functionValueAccuracy = getFunctionValueAccuracy();
-
- final Complex nC = new Complex(n, 0);
- final Complex n1C = new Complex(n - 1, 0);
-
- Complex z = initial;
- Complex oldz = new Complex(Double.POSITIVE_INFINITY,
- Double.POSITIVE_INFINITY);
- while (true) {
- // Compute pv (polynomial value), dv (derivative value), and
- // d2v (second derivative value) simultaneously.
- Complex pv = coefficients[n];
- Complex dv = Complex.ZERO;
- Complex d2v = Complex.ZERO;
- for (int j = n-1; j >= 0; j--) {
- d2v = dv.add(z.multiply(d2v));
- dv = pv.add(z.multiply(dv));
- pv = coefficients[j].add(z.multiply(pv));
- }
- d2v = d2v.multiply(new Complex(2.0, 0.0));
-
- // Check for convergence.
- final double tolerance = FastMath.max(relativeAccuracy * z.abs(),
- absoluteAccuracy);
- if ((z.subtract(oldz)).abs() <= tolerance) {
- return z;
- }
- if (pv.abs() <= functionValueAccuracy) {
- return z;
- }
-
- // Now pv != 0, calculate the new approximation.
- final Complex G = dv.divide(pv);
- final Complex G2 = G.multiply(G);
- final Complex H = G2.subtract(d2v.divide(pv));
- final Complex delta = n1C.multiply((nC.multiply(H)).subtract(G2));
- // Choose a denominator larger in magnitude.
- final Complex deltaSqrt = delta.sqrt();
- final Complex dplus = G.add(deltaSqrt);
- final Complex dminus = G.subtract(deltaSqrt);
- final Complex denominator = dplus.abs() > dminus.abs() ? dplus : dminus;
- // Perturb z if denominator is zero, for instance,
- // p(x) = x^3 + 1, z = 0.
- if (denominator.equals(new Complex(0.0, 0.0))) {
- z = z.add(new Complex(absoluteAccuracy, absoluteAccuracy));
- oldz = new Complex(Double.POSITIVE_INFINITY,
- Double.POSITIVE_INFINITY);
- } else {
- oldz = z;
- z = z.subtract(nC.divide(denominator));
- }
- incrementEvaluationCount();
- }
- }
- }
-}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/MullerSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/MullerSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/MullerSolver.java
deleted file mode 100644
index 06a7b6b..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/MullerSolver.java
+++ /dev/null
@@ -1,202 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-package org.apache.commons.math3.analysis.solvers;
-
-import org.apache.commons.math3.util.FastMath;
-import org.apache.commons.math3.exception.NumberIsTooLargeException;
-import org.apache.commons.math3.exception.NoBracketingException;
-import org.apache.commons.math3.exception.TooManyEvaluationsException;
-
-/**
- * This class implements the <a href="http://mathworld.wolfram.com/MullersMethod.html">
- * Muller's Method</a> for root finding of real univariate functions. For
- * reference, see <b>Elementary Numerical Analysis</b>, ISBN 0070124477,
- * chapter 3.
- * <p>
- * Muller's method applies to both real and complex functions, but here we
- * restrict ourselves to real functions.
- * This class differs from {@link MullerSolver} in the way it avoids complex
- * operations.</p>
- * Muller's original method would have function evaluation at complex point.
- * Since our f(x) is real, we have to find ways to avoid that. Bracketing
- * condition is one way to go: by requiring bracketing in every iteration,
- * the newly computed approximation is guaranteed to be real.</p>
- * <p>
- * Normally Muller's method converges quadratically in the vicinity of a
- * zero, however it may be very slow in regions far away from zeros. For
- * example, f(x) = exp(x) - 1, min = -50, max = 100. In such case we use
- * bisection as a safety backup if it performs very poorly.</p>
- * <p>
- * The formulas here use divided differences directly.</p>
- *
- * @since 1.2
- * @see MullerSolver2
- */
-public class MullerSolver extends AbstractUnivariateSolver {
-
- /** Default absolute accuracy. */
- private static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
-
- /**
- * Construct a solver with default accuracy (1e-6).
- */
- public MullerSolver() {
- this(DEFAULT_ABSOLUTE_ACCURACY);
- }
- /**
- * Construct a solver.
- *
- * @param absoluteAccuracy Absolute accuracy.
- */
- public MullerSolver(double absoluteAccuracy) {
- super(absoluteAccuracy);
- }
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- */
- public MullerSolver(double relativeAccuracy,
- double absoluteAccuracy) {
- super(relativeAccuracy, absoluteAccuracy);
- }
-
- /**
- * {@inheritDoc}
- */
- @Override
- protected double doSolve()
- throws TooManyEvaluationsException,
- NumberIsTooLargeException,
- NoBracketingException {
- final double min = getMin();
- final double max = getMax();
- final double initial = getStartValue();
-
- final double functionValueAccuracy = getFunctionValueAccuracy();
-
- verifySequence(min, initial, max);
-
- // check for zeros before verifying bracketing
- final double fMin = computeObjectiveValue(min);
- if (FastMath.abs(fMin) < functionValueAccuracy) {
- return min;
- }
- final double fMax = computeObjectiveValue(max);
- if (FastMath.abs(fMax) < functionValueAccuracy) {
- return max;
- }
- final double fInitial = computeObjectiveValue(initial);
- if (FastMath.abs(fInitial) < functionValueAccuracy) {
- return initial;
- }
-
- verifyBracketing(min, max);
-
- if (isBracketing(min, initial)) {
- return solve(min, initial, fMin, fInitial);
- } else {
- return solve(initial, max, fInitial, fMax);
- }
- }
-
- /**
- * Find a real root in the given interval.
- *
- * @param min Lower bound for the interval.
- * @param max Upper bound for the interval.
- * @param fMin function value at the lower bound.
- * @param fMax function value at the upper bound.
- * @return the point at which the function value is zero.
- * @throws TooManyEvaluationsException if the allowed number of calls to
- * the function to be solved has been exhausted.
- */
- 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();
-
- // [x0, x2] is the bracketing interval in each iteration
- // x1 is the last approximation and an interpolation point in (x0, x2)
- // x is the new root approximation and new x1 for next round
- // d01, d12, d012 are divided differences
-
- 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) {
- // Muller's method employs quadratic interpolation through
- // x0, x1, x2 and x is the zero of the interpolating parabola.
- // Due to bracketing condition, this parabola must have two
- // real roots and we choose one in [x0, x2] to be x.
- 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));
- // xplus and xminus are two roots of parabola and at least
- // one of them should lie in (x0, x2)
- final double x = isSequence(x0, xplus, x2) ? xplus : xminus;
- final double y = computeObjectiveValue(x);
-
- // check for convergence
- final double tolerance = FastMath.max(relativeAccuracy * FastMath.abs(x), absoluteAccuracy);
- if (FastMath.abs(x - oldx) <= tolerance ||
- FastMath.abs(y) <= functionValueAccuracy) {
- return x;
- }
-
- // Bisect if convergence is too slow. Bisection would waste
- // our calculation of x, hopefully it won't happen often.
- // the real number equality test x == x1 is intentional and
- // completes the proximity tests above it
- boolean bisect = (x < x1 && (x1 - x0) > 0.95 * (x2 - x0)) ||
- (x > x1 && (x2 - x1) > 0.95 * (x2 - x0)) ||
- (x == x1);
- // prepare the new bracketing interval for next iteration
- 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;
- }
- }
- }
-}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/MullerSolver2.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/MullerSolver2.java b/src/main/java/org/apache/commons/math3/analysis/solvers/MullerSolver2.java
deleted file mode 100644
index 2a79c11..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/MullerSolver2.java
+++ /dev/null
@@ -1,168 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-package org.apache.commons.math3.analysis.solvers;
-
-import org.apache.commons.math3.exception.NoBracketingException;
-import org.apache.commons.math3.exception.NumberIsTooLargeException;
-import org.apache.commons.math3.exception.TooManyEvaluationsException;
-import org.apache.commons.math3.util.FastMath;
-
-/**
- * This class implements the <a href="http://mathworld.wolfram.com/MullersMethod.html">
- * Muller's Method</a> for root finding of real univariate functions. For
- * reference, see <b>Elementary Numerical Analysis</b>, ISBN 0070124477,
- * chapter 3.
- * <p>
- * Muller's method applies to both real and complex functions, but here we
- * restrict ourselves to real functions.
- * This class differs from {@link MullerSolver} in the way it avoids complex
- * operations.</p>
- * Except for the initial [min, max], it does not require bracketing
- * condition, e.g. f(x0), f(x1), f(x2) can have the same sign. If complex
- * number arises in the computation, we simply use its modulus as real
- * approximation.</p>
- * <p>
- * Because the interval may not be bracketing, bisection alternative is
- * not applicable here. However in practice our treatment usually works
- * well, especially near real zeroes where the imaginary part of complex
- * approximation is often negligible.</p>
- * <p>
- * The formulas here do not use divided differences directly.</p>
- *
- * @since 1.2
- * @see MullerSolver
- */
-public class MullerSolver2 extends AbstractUnivariateSolver {
-
- /** Default absolute accuracy. */
- private static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
-
- /**
- * Construct a solver with default accuracy (1e-6).
- */
- public MullerSolver2() {
- this(DEFAULT_ABSOLUTE_ACCURACY);
- }
- /**
- * Construct a solver.
- *
- * @param absoluteAccuracy Absolute accuracy.
- */
- public MullerSolver2(double absoluteAccuracy) {
- super(absoluteAccuracy);
- }
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- */
- public MullerSolver2(double relativeAccuracy,
- double absoluteAccuracy) {
- super(relativeAccuracy, absoluteAccuracy);
- }
-
- /**
- * {@inheritDoc}
- */
- @Override
- protected double doSolve()
- throws TooManyEvaluationsException,
- NumberIsTooLargeException,
- NoBracketingException {
- final double min = getMin();
- final double max = getMax();
-
- verifyInterval(min, max);
-
- final double relativeAccuracy = getRelativeAccuracy();
- final double absoluteAccuracy = getAbsoluteAccuracy();
- final double functionValueAccuracy = getFunctionValueAccuracy();
-
- // x2 is the last root approximation
- // x is the new approximation and new x2 for next round
- // x0 < x1 < x2 does not hold here
-
- double x0 = min;
- double y0 = computeObjectiveValue(x0);
- if (FastMath.abs(y0) < functionValueAccuracy) {
- return x0;
- }
- double x1 = max;
- double y1 = computeObjectiveValue(x1);
- if (FastMath.abs(y1) < functionValueAccuracy) {
- return x1;
- }
-
- if(y0 * y1 > 0) {
- throw new NoBracketingException(x0, x1, y0, y1);
- }
-
- double x2 = 0.5 * (x0 + x1);
- double y2 = computeObjectiveValue(x2);
-
- double oldx = Double.POSITIVE_INFINITY;
- while (true) {
- // quadratic interpolation through x0, x1, x2
- final double q = (x2 - x1) / (x1 - x0);
- final double a = q * (y2 - (1 + q) * y1 + q * y0);
- final double b = (2 * q + 1) * y2 - (1 + q) * (1 + q) * y1 + q * q * y0;
- final double c = (1 + q) * y2;
- final double delta = b * b - 4 * a * c;
- double x;
- final double denominator;
- if (delta >= 0.0) {
- // choose a denominator larger in magnitude
- double dplus = b + FastMath.sqrt(delta);
- double dminus = b - FastMath.sqrt(delta);
- denominator = FastMath.abs(dplus) > FastMath.abs(dminus) ? dplus : dminus;
- } else {
- // take the modulus of (B +/- FastMath.sqrt(delta))
- denominator = FastMath.sqrt(b * b - delta);
- }
- if (denominator != 0) {
- x = x2 - 2.0 * c * (x2 - x1) / denominator;
- // perturb x if it exactly coincides with x1 or x2
- // the equality tests here are intentional
- while (x == x1 || x == x2) {
- x += absoluteAccuracy;
- }
- } else {
- // extremely rare case, get a random number to skip it
- x = min + FastMath.random() * (max - min);
- oldx = Double.POSITIVE_INFINITY;
- }
- final double y = computeObjectiveValue(x);
-
- // check for convergence
- final double tolerance = FastMath.max(relativeAccuracy * FastMath.abs(x), absoluteAccuracy);
- if (FastMath.abs(x - oldx) <= tolerance ||
- FastMath.abs(y) <= functionValueAccuracy) {
- return x;
- }
-
- // prepare the next iteration
- x0 = x1;
- y0 = y1;
- x1 = x2;
- y1 = y2;
- x2 = x;
- y2 = y;
- oldx = x;
- }
- }
-}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/NewtonRaphsonSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/NewtonRaphsonSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/NewtonRaphsonSolver.java
deleted file mode 100644
index 4cf2688..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/NewtonRaphsonSolver.java
+++ /dev/null
@@ -1,92 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-
-package org.apache.commons.math3.analysis.solvers;
-
-import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
-import org.apache.commons.math3.analysis.differentiation.UnivariateDifferentiableFunction;
-import org.apache.commons.math3.util.FastMath;
-import org.apache.commons.math3.exception.TooManyEvaluationsException;
-
-/**
- * Implements <a href="http://mathworld.wolfram.com/NewtonsMethod.html">
- * Newton's Method</a> for finding zeros of real univariate differentiable
- * functions.
- *
- * @since 3.1
- */
-public class NewtonRaphsonSolver extends AbstractUnivariateDifferentiableSolver {
- /** Default absolute accuracy. */
- private static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
-
- /**
- * Construct a solver.
- */
- public NewtonRaphsonSolver() {
- this(DEFAULT_ABSOLUTE_ACCURACY);
- }
- /**
- * Construct a solver.
- *
- * @param absoluteAccuracy Absolute accuracy.
- */
- public NewtonRaphsonSolver(double absoluteAccuracy) {
- super(absoluteAccuracy);
- }
-
- /**
- * Find a zero near the midpoint of {@code min} and {@code max}.
- *
- * @param f Function to solve.
- * @param min Lower bound for the interval.
- * @param max Upper bound for the interval.
- * @param maxEval Maximum number of evaluations.
- * @return the value where the function is zero.
- * @throws org.apache.commons.math3.exception.TooManyEvaluationsException
- * if the maximum evaluation count is exceeded.
- * @throws org.apache.commons.math3.exception.NumberIsTooLargeException
- * if {@code min >= max}.
- */
- @Override
- public double solve(int maxEval, final UnivariateDifferentiableFunction f,
- final double min, final double max)
- throws TooManyEvaluationsException {
- return super.solve(maxEval, f, UnivariateSolverUtils.midpoint(min, max));
- }
-
- /**
- * {@inheritDoc}
- */
- @Override
- protected double doSolve()
- throws TooManyEvaluationsException {
- final double startValue = getStartValue();
- final double absoluteAccuracy = getAbsoluteAccuracy();
-
- double x0 = startValue;
- double x1;
- while (true) {
- final DerivativeStructure y0 = computeObjectiveValueAndDerivative(x0);
- x1 = x0 - (y0.getValue() / y0.getPartialDerivative(1));
- if (FastMath.abs(x1 - x0) <= absoluteAccuracy) {
- return x1;
- }
-
- x0 = x1;
- }
- }
-}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/NewtonSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/NewtonSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/NewtonSolver.java
deleted file mode 100644
index 3ba7bf2..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/NewtonSolver.java
+++ /dev/null
@@ -1,92 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-
-package org.apache.commons.math3.analysis.solvers;
-
-import org.apache.commons.math3.analysis.DifferentiableUnivariateFunction;
-import org.apache.commons.math3.util.FastMath;
-import org.apache.commons.math3.exception.TooManyEvaluationsException;
-
-/**
- * Implements <a href="http://mathworld.wolfram.com/NewtonsMethod.html">
- * Newton's Method</a> for finding zeros of real univariate functions.
- * <p>
- * The function should be continuous but not necessarily smooth.</p>
- *
- * @deprecated as of 3.1, replaced by {@link NewtonRaphsonSolver}
- */
-@Deprecated
-public class NewtonSolver extends AbstractDifferentiableUnivariateSolver {
- /** Default absolute accuracy. */
- private static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
-
- /**
- * Construct a solver.
- */
- public NewtonSolver() {
- this(DEFAULT_ABSOLUTE_ACCURACY);
- }
- /**
- * Construct a solver.
- *
- * @param absoluteAccuracy Absolute accuracy.
- */
- public NewtonSolver(double absoluteAccuracy) {
- super(absoluteAccuracy);
- }
-
- /**
- * Find a zero near the midpoint of {@code min} and {@code max}.
- *
- * @param f Function to solve.
- * @param min Lower bound for the interval.
- * @param max Upper bound for the interval.
- * @param maxEval Maximum number of evaluations.
- * @return the value where the function is zero.
- * @throws org.apache.commons.math3.exception.TooManyEvaluationsException
- * if the maximum evaluation count is exceeded.
- * @throws org.apache.commons.math3.exception.NumberIsTooLargeException
- * if {@code min >= max}.
- */
- @Override
- public double solve(int maxEval, final DifferentiableUnivariateFunction f,
- final double min, final double max)
- throws TooManyEvaluationsException {
- return super.solve(maxEval, f, UnivariateSolverUtils.midpoint(min, max));
- }
-
- /**
- * {@inheritDoc}
- */
- @Override
- protected double doSolve()
- throws TooManyEvaluationsException {
- final double startValue = getStartValue();
- final double absoluteAccuracy = getAbsoluteAccuracy();
-
- double x0 = startValue;
- double x1;
- while (true) {
- x1 = x0 - (computeObjectiveValue(x0) / computeDerivativeObjectiveValue(x0));
- if (FastMath.abs(x1 - x0) <= absoluteAccuracy) {
- return x1;
- }
-
- x0 = x1;
- }
- }
-}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/PegasusSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/PegasusSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/PegasusSolver.java
deleted file mode 100644
index 0d80895..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/PegasusSolver.java
+++ /dev/null
@@ -1,84 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-
-package org.apache.commons.math3.analysis.solvers;
-
-/**
- * Implements the <em>Pegasus</em> method for root-finding (approximating
- * a zero of a univariate real function). It is a modified
- * {@link RegulaFalsiSolver <em>Regula Falsi</em>} method.
- *
- * <p>Like the <em>Regula Falsi</em> method, convergence is guaranteed by
- * maintaining a bracketed solution. The <em>Pegasus</em> method however,
- * should converge much faster than the original <em>Regula Falsi</em>
- * method. Furthermore, this implementation of the <em>Pegasus</em> method
- * should not suffer from the same implementation issues as the <em>Regula
- * Falsi</em> method, which may fail to convergence in certain cases. Also,
- * the <em>Pegasus</em> method should converge faster than the
- * {@link IllinoisSolver <em>Illinois</em>} method, another <em>Regula
- * Falsi</em>-based method.</p>
- *
- * <p>The <em>Pegasus</em> method assumes that the function is continuous,
- * but not necessarily smooth.</p>
- *
- * <p>Implementation based on the following article: M. Dowell and P. Jarratt,
- * <em>The "Pegasus" method for computing the root of an equation</em>,
- * BIT Numerical Mathematics, volume 12, number 4, pages 503-508, Springer,
- * 1972.</p>
- *
- * @since 3.0
- */
-public class PegasusSolver extends BaseSecantSolver {
-
- /** Construct a solver with default accuracy (1e-6). */
- public PegasusSolver() {
- super(DEFAULT_ABSOLUTE_ACCURACY, Method.PEGASUS);
- }
-
- /**
- * Construct a solver.
- *
- * @param absoluteAccuracy Absolute accuracy.
- */
- public PegasusSolver(final double absoluteAccuracy) {
- super(absoluteAccuracy, Method.PEGASUS);
- }
-
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- */
- public PegasusSolver(final double relativeAccuracy,
- final double absoluteAccuracy) {
- super(relativeAccuracy, absoluteAccuracy, Method.PEGASUS);
- }
-
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- * @param functionValueAccuracy Maximum function value error.
- */
- public PegasusSolver(final double relativeAccuracy,
- final double absoluteAccuracy,
- final double functionValueAccuracy) {
- super(relativeAccuracy, absoluteAccuracy, functionValueAccuracy, Method.PEGASUS);
- }
-}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/PolynomialSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/PolynomialSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/PolynomialSolver.java
deleted file mode 100644
index c21f076..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/PolynomialSolver.java
+++ /dev/null
@@ -1,28 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-package org.apache.commons.math3.analysis.solvers;
-
-import org.apache.commons.math3.analysis.polynomials.PolynomialFunction;
-
-/**
- * Interface for (polynomial) root-finding algorithms.
- * Implementations will search for only one zero in the given interval.
- *
- * @since 3.0
- */
-public interface PolynomialSolver
- extends BaseUnivariateSolver<PolynomialFunction> {}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/RegulaFalsiSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/RegulaFalsiSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/RegulaFalsiSolver.java
deleted file mode 100644
index cfb7055..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/RegulaFalsiSolver.java
+++ /dev/null
@@ -1,94 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-
-package org.apache.commons.math3.analysis.solvers;
-
-/**
- * Implements the <em>Regula Falsi</em> or <em>False position</em> method for
- * root-finding (approximating a zero of a univariate real function). It is a
- * modified {@link SecantSolver <em>Secant</em>} method.
- *
- * <p>The <em>Regula Falsi</em> method is included for completeness, for
- * testing purposes, for educational purposes, for comparison to other
- * algorithms, etc. It is however <strong>not</strong> intended to be used
- * for actual problems, as one of the bounds often remains fixed, resulting
- * in very slow convergence. Instead, one of the well-known modified
- * <em>Regula Falsi</em> algorithms can be used ({@link IllinoisSolver
- * <em>Illinois</em>} or {@link PegasusSolver <em>Pegasus</em>}). These two
- * algorithms solve the fundamental issues of the original <em>Regula
- * Falsi</em> algorithm, and greatly out-performs it for most, if not all,
- * (practical) functions.
- *
- * <p>Unlike the <em>Secant</em> method, the <em>Regula Falsi</em> guarantees
- * convergence, by maintaining a bracketed solution. Note however, that due to
- * the finite/limited precision of Java's {@link Double double} type, which is
- * used in this implementation, the algorithm may get stuck in a situation
- * where it no longer makes any progress. Such cases are detected and result
- * in a {@code ConvergenceException} exception being thrown. In other words,
- * the algorithm theoretically guarantees convergence, but the implementation
- * does not.</p>
- *
- * <p>The <em>Regula Falsi</em> method assumes that the function is continuous,
- * but not necessarily smooth.</p>
- *
- * <p>Implementation based on the following article: M. Dowell and P. Jarratt,
- * <em>A modified regula falsi method for computing the root of an
- * equation</em>, BIT Numerical Mathematics, volume 11, number 2,
- * pages 168-174, Springer, 1971.</p>
- *
- * @since 3.0
- */
-public class RegulaFalsiSolver extends BaseSecantSolver {
-
- /** Construct a solver with default accuracy (1e-6). */
- public RegulaFalsiSolver() {
- super(DEFAULT_ABSOLUTE_ACCURACY, Method.REGULA_FALSI);
- }
-
- /**
- * Construct a solver.
- *
- * @param absoluteAccuracy Absolute accuracy.
- */
- public RegulaFalsiSolver(final double absoluteAccuracy) {
- super(absoluteAccuracy, Method.REGULA_FALSI);
- }
-
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- */
- public RegulaFalsiSolver(final double relativeAccuracy,
- final double absoluteAccuracy) {
- super(relativeAccuracy, absoluteAccuracy, Method.REGULA_FALSI);
- }
-
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- * @param functionValueAccuracy Maximum function value error.
- */
- public RegulaFalsiSolver(final double relativeAccuracy,
- final double absoluteAccuracy,
- final double functionValueAccuracy) {
- super(relativeAccuracy, absoluteAccuracy, functionValueAccuracy, Method.REGULA_FALSI);
- }
-}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/analysis/solvers/RiddersSolver.java
----------------------------------------------------------------------
diff --git a/src/main/java/org/apache/commons/math3/analysis/solvers/RiddersSolver.java b/src/main/java/org/apache/commons/math3/analysis/solvers/RiddersSolver.java
deleted file mode 100644
index d83f595..0000000
--- a/src/main/java/org/apache/commons/math3/analysis/solvers/RiddersSolver.java
+++ /dev/null
@@ -1,142 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-package org.apache.commons.math3.analysis.solvers;
-
-import org.apache.commons.math3.util.FastMath;
-import org.apache.commons.math3.exception.NoBracketingException;
-import org.apache.commons.math3.exception.TooManyEvaluationsException;
-
-/**
- * Implements the <a href="http://mathworld.wolfram.com/RiddersMethod.html">
- * Ridders' Method</a> for root finding of real univariate functions. For
- * reference, see C. Ridders, <i>A new algorithm for computing a single root
- * of a real continuous function </i>, IEEE Transactions on Circuits and
- * Systems, 26 (1979), 979 - 980.
- * <p>
- * The function should be continuous but not necessarily smooth.</p>
- *
- * @since 1.2
- */
-public class RiddersSolver extends AbstractUnivariateSolver {
- /** Default absolute accuracy. */
- private static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
-
- /**
- * Construct a solver with default accuracy (1e-6).
- */
- public RiddersSolver() {
- this(DEFAULT_ABSOLUTE_ACCURACY);
- }
- /**
- * Construct a solver.
- *
- * @param absoluteAccuracy Absolute accuracy.
- */
- public RiddersSolver(double absoluteAccuracy) {
- super(absoluteAccuracy);
- }
- /**
- * Construct a solver.
- *
- * @param relativeAccuracy Relative accuracy.
- * @param absoluteAccuracy Absolute accuracy.
- */
- public RiddersSolver(double relativeAccuracy,
- double absoluteAccuracy) {
- super(relativeAccuracy, absoluteAccuracy);
- }
-
- /**
- * {@inheritDoc}
- */
- @Override
- protected double doSolve()
- throws TooManyEvaluationsException,
- NoBracketingException {
- double min = getMin();
- double max = getMax();
- // [x1, x2] is the bracketing interval in each iteration
- // x3 is the midpoint of [x1, x2]
- // x is the new root approximation and an endpoint of the new interval
- double x1 = min;
- double y1 = computeObjectiveValue(x1);
- double x2 = max;
- double y2 = computeObjectiveValue(x2);
-
- // check for zeros before verifying bracketing
- if (y1 == 0) {
- return min;
- }
- if (y2 == 0) {
- return max;
- }
- verifyBracketing(min, max);
-
- final double absoluteAccuracy = getAbsoluteAccuracy();
- final double functionValueAccuracy = getFunctionValueAccuracy();
- final double relativeAccuracy = getRelativeAccuracy();
-
- double oldx = Double.POSITIVE_INFINITY;
- while (true) {
- // calculate the new root approximation
- final double x3 = 0.5 * (x1 + x2);
- final double y3 = computeObjectiveValue(x3);
- if (FastMath.abs(y3) <= functionValueAccuracy) {
- return x3;
- }
- final double delta = 1 - (y1 * y2) / (y3 * y3); // delta > 1 due to bracketing
- final double correction = (FastMath.signum(y2) * FastMath.signum(y3)) *
- (x3 - x1) / FastMath.sqrt(delta);
- final double x = x3 - correction; // correction != 0
- final double y = computeObjectiveValue(x);
-
- // check for convergence
- final double tolerance = FastMath.max(relativeAccuracy * FastMath.abs(x), absoluteAccuracy);
- if (FastMath.abs(x - oldx) <= tolerance) {
- return x;
- }
- if (FastMath.abs(y) <= functionValueAccuracy) {
- return x;
- }
-
- // prepare the new interval for next iteration
- // Ridders' method guarantees x1 < x < x2
- if (correction > 0.0) { // x1 < x < x3
- if (FastMath.signum(y1) + FastMath.signum(y) == 0.0) {
- x2 = x;
- y2 = y;
- } else {
- x1 = x;
- x2 = x3;
- y1 = y;
- y2 = y3;
- }
- } else { // x3 < x < x2
- if (FastMath.signum(y2) + FastMath.signum(y) == 0.0) {
- x1 = x;
- y1 = y;
- } else {
- x1 = x3;
- x2 = x;
- y1 = y3;
- y2 = y;
- }
- }
- oldx = x;
- }
- }
-}