BisectionSingleRootEstimator.java
/*
* Copyright (C) 2012 Alberto Irurueta Carro (alberto@irurueta.com)
*
* Licensed 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 com.irurueta.numerical.roots;
import com.irurueta.numerical.EvaluationException;
import com.irurueta.numerical.InvalidBracketRangeException;
import com.irurueta.numerical.LockedException;
import com.irurueta.numerical.NotReadyException;
import com.irurueta.numerical.SingleDimensionFunctionEvaluatorListener;
/**
* This class searches for a single REAL root on a single dimension function
* (i.e. f(x) ). The root to be found must be inside a given or computed bracket
* of values.
* <p>
* This class uses the bisection method, which is one of the slowest but most
* reliable existing methods (i.e. the method doubles the accuracy on each
* iteration, so in 64 iterations the highest accuracy provided by a double
* value can be obtained).
* <p>
* In comparison to other methods, the Bisection one can find roots even in
* situations where isn't a zero crossing. In other words, bracket estimation
* might fail for this class, but even then a root might be obtained.
* <p>
* In order to find a root around a given range of values, a bracket of values
* can be provided. Otherwise, a bracket can be computed by attempting to find
* a zero-crossing while expanding an initial range of values.
* <p>
* Bracket computation estimates a larger range of values, which can later be
* refined in order to estimate a given root.
* <p>
* This class is based on the implementation of Numerical Recipes 3rd ed.
* Section 9.1.1 page 447.
*/
public class BisectionSingleRootEstimator extends BracketedSingleRootEstimator {
/**
* Constant defining maximum number of iterations to estimate a root.
*/
public static final int JMAX = 50;
/**
* Constant defining default tolerance to find a root.
*/
public static final double DEFAULT_TOLERANCE = 1e-6;
/**
* Minimum allowed tolerance that can be set.
*/
public static final double MIN_TOLERANCE = 0.0;
/**
* Tolerance to find a root. Whenever the variation of the estimated root is
* smaller than the provided tolerance, then the algorithm is assumed to be
* converged.
*/
private double tolerance;
/**
* Empty constructor.
*/
public BisectionSingleRootEstimator() {
super();
tolerance = DEFAULT_TOLERANCE;
}
/**
* Constructor.
*
* @param listener Listener to evaluate a single dimension function f(x)
* to find its roots.
* @param minEvalPoint Smallest value inside the bracket of values where the
* root will be searched.
* @param maxEvalPoint Largest value inside the bracket of values where the
* root will be searched.
* @param tolerance Tolerance to find a root. During the estimation of a
* root, if the variation of the estimated root is smaller than the provided
* tolerance, then the algorithm is assumed to be converged, and the root
* is ensured to have an accuracy that equals tolerance.
* @throws InvalidBracketRangeException Raised if minEvalPoint <
* maxEvalPoint.
* @throws IllegalArgumentException Raised if provided tolerance is negative.
*/
public BisectionSingleRootEstimator(
final SingleDimensionFunctionEvaluatorListener listener, final double minEvalPoint,
final double maxEvalPoint, final double tolerance) throws InvalidBracketRangeException {
super(listener, minEvalPoint, maxEvalPoint);
internalSetTolerance(tolerance);
}
/**
* Returns tolerance to find a root. Whenever the variation of the estimated
* root is smaller than returned tolerance, then the algorithm is assumed to
* be converged, and the estimated root is ensured to have an accuracy that
* equals the returned tolerance.
*
* @return Tolerance to find a root.
*/
public double getTolerance() {
return tolerance;
}
/**
* Sets tolerance to find a root. Whenever the variation of the estimated
* root is smaller than provided tolerance, then the algorithm is assumed to
* be converged, and the estimated root is ensured to have an accuracy that
* equals provided tolerance.
*
* @param tolerance Tolerance to find a root.
* @throws LockedException Raised if this instance is locked while doing
* some computations.
* @throws IllegalArgumentException Raised if provided tolerance is negative.
*/
public void setTolerance(final double tolerance) throws LockedException {
if (isLocked()) {
throw new LockedException();
}
internalSetTolerance(tolerance);
}
/**
* Internal method to set tolerance to find a root. This method does not
* check whether this instance is locked.
* Whenever the variation of the estimated root is smaller than provided
* tolerance, then the algorithm is assumed to be converged, and the
* estimated root is ensured to have an accuracy that equals provided
* tolerance.
*
* @param tolerance Tolerance to find a root.
* @throws IllegalArgumentException Raised if provided tolerance is negative.
*/
private void internalSetTolerance(final double tolerance) {
if (tolerance < MIN_TOLERANCE) {
throw new IllegalArgumentException();
}
this.tolerance = tolerance;
}
/**
* Estimates a single root of the provided single dimension function
* contained within a given bracket of values.
*
* @throws LockedException Exception raised if this instance is already
* locked.
* @throws NotReadyException Exception raised if not enough parameters have
* been provided in order to start the estimation.
* @throws RootEstimationException Raised if the root estimation failed for
* some other reason (usually inability to evaluate the function,
* numerical instability or convergence problems, or no roots are found).
*/
@Override
public void estimate() throws LockedException, NotReadyException, RootEstimationException {
if (isLocked()) {
throw new LockedException();
}
if (!isReady()) {
throw new NotReadyException();
}
locked = true;
rootAvailable = false;
try {
final var x1 = minEvalPoint;
final var x2 = maxEvalPoint;
final var xacc = tolerance;
double dx;
double xmid;
double rtb;
final var f = listener.evaluate(x1);
var fmid = listener.evaluate(x2);
if (f * fmid >= 0.0) {
// check that bracket contains a sign change in function
locked = false;
throw new RootEstimationException();
}
if (f < 0.0) {
dx = x2 - x1;
rtb = x1;
} else {
dx = x1 - x2;
rtb = x2;
}
for (var j = 0; j < JMAX; j++) {
dx *= 0.5;
xmid = rtb + dx;
fmid = listener.evaluate(xmid);
if (fmid <= 0.0) {
rtb = xmid;
}
if (Math.abs(dx) < xacc || fmid == 0.0) {
// result obtained
root = rtb;
rootAvailable = true;
locked = false;
return;
}
}
} catch (final EvaluationException e) {
throw new RootEstimationException(e);
} finally {
locked = false;
}
// too many iterations and error exceeds desired tolerance
throw new RootEstimationException();
}
/**
* Returns boolean indicating whether this instance is ready to compute a
* root.
* This instance is ready as soon as a listener is provided and a bracket
* is provided or computed.
*
* @return True if instance is ready, false otherwise
*/
@Override
public boolean isReady() {
return isListenerAvailable() && isBracketAvailable();
}
}