org.apache.commons.math.optimization.univariate.BrentOptimizer Maven / Gradle / Ivy
Go to download
Show more of this group Show more artifacts with this name
Show all versions of commons-math Show documentation
Show all versions of commons-math Show documentation
The Math project is a library of lightweight, self-contained mathematics and statistics components addressing the most common practical problems not immediately available in the Java programming language or commons-lang.
/*
* 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.math.optimization.univariate;
import org.apache.commons.math.FunctionEvaluationException;
import org.apache.commons.math.MaxIterationsExceededException;
import org.apache.commons.math.analysis.UnivariateRealFunction;
import org.apache.commons.math.optimization.GoalType;
/**
* Implements Richard Brent's algorithm (from his book "Algorithms for
* Minimization without Derivatives", p. 79) for finding minima of real
* univariate functions.
*
* @version $Revision: 799857 $ $Date: 2009-08-01 09:07:12 -0400 (Sat, 01 Aug 2009) $
* @since 2.0
*/
public class BrentOptimizer extends AbstractUnivariateRealOptimizer {
/**
* Golden section.
*/
private static final double c = 0.5 * (3 - Math.sqrt(5));
/**
* Construct a solver.
*/
public BrentOptimizer() {
super(100, 1E-10);
}
/** {@inheritDoc} */
public double optimize(final UnivariateRealFunction f, final GoalType goalType,
final double min, final double max, final double startValue)
throws MaxIterationsExceededException, FunctionEvaluationException {
return optimize(f, goalType, min, max);
}
/** {@inheritDoc} */
public double optimize(final UnivariateRealFunction f, final GoalType goalType,
final double min, final double max)
throws MaxIterationsExceededException, FunctionEvaluationException {
clearResult();
return localMin(f, goalType, min, max, relativeAccuracy, absoluteAccuracy);
}
/**
* Find the minimum of the function {@code f} within the interval {@code (a, b)}.
*
* If the function {@code f} is defined on the interval {@code (a, b)}, then
* this method finds an approximation {@code x} to the point at which {@code f}
* attains its minimum.
* {@code t} and {@code eps} define a tolerance {@code tol = eps |x| + t} and
* {@code f} is never evaluated at two points closer together than {@code tol}.
* {@code eps} should be no smaller than 2 macheps and preferable not
* much less than sqrt(macheps), where macheps is the relative
* machine precision. {@code t} should be positive.
* @param f the function to solve
* @param goalType type of optimization goal: either {@link GoalType#MAXIMIZE}
* or {@link GoalType#MINIMIZE}
* @param a Lower bound of the interval
* @param b Higher bound of the interval
* @param eps Relative accuracy
* @param t Absolute accuracy
* @return the point at which the function is minimal.
* @throws MaxIterationsExceededException if the maximum iteration count
* is exceeded.
* @throws FunctionEvaluationException if an error occurs evaluating
* the function.
*/
private double localMin(final UnivariateRealFunction f, final GoalType goalType,
double a, double b, final double eps, final double t)
throws MaxIterationsExceededException, FunctionEvaluationException {
double x = a + c * (b - a);
double v = x;
double w = x;
double e = 0;
double fx = computeObjectiveValue(f, x);
if (goalType == GoalType.MAXIMIZE) {
fx = -fx;
}
double fv = fx;
double fw = fx;
int count = 0;
while (count < maximalIterationCount) {
double m = 0.5 * (a + b);
double tol = eps * Math.abs(x) + t;
double t2 = 2 * tol;
// Check stopping criterion.
if (Math.abs(x - m) > t2 - 0.5 * (b - a)) {
double p = 0;
double q = 0;
double r = 0;
double d = 0;
double u = 0;
if (Math.abs(e) > tol) { // Fit parabola.
r = (x - w) * (fx - fv);
q = (x - v) * (fx - fw);
p = (x - v) * q - (x - w) * r;
q = 2 * (q - r);
if (q > 0) {
p = -p;
} else {
q = -q;
}
r = e;
e = d;
}
if (Math.abs(p) < Math.abs(0.5 * q * r) &&
(p < q * (a - x)) && (p < q * (b - x))) { // Parabolic interpolation step.
d = p / q;
u = x + d;
// f must not be evaluated too close to a or b.
if (((u - a) < t2) || ((b - u) < t2)) {
d = (x < m) ? tol : -tol;
}
} else { // Golden section step.
e = ((x < m) ? b : a) - x;
d = c * e;
}
// f must not be evaluated too close to a or b.
u = x + ((Math.abs(d) > tol) ? d : ((d > 0) ? tol : -tol));
double fu = computeObjectiveValue(f, u);
if (goalType == GoalType.MAXIMIZE) {
fu = -fu;
}
// Update a, b, v, w and x.
if (fu <= fx) {
if (u < x) {
b = x;
} else {
a = x;
}
v = w;
fv = fw;
w = x;
fw = fx;
x = u;
fx = fu;
} else {
if (u < x) {
a = u;
} else {
b = u;
}
if ((fu <= fw) || (w == x)) {
v = w;
fv = fw;
w = u;
fw = fu;
} else if ((fu <= fv) || (v == x) || (v == w)) {
v = u;
fv = fu;
}
}
} else { // termination
setResult(x, (goalType == GoalType.MAXIMIZE) ? -fx : fx, count);
return x;
}
++count;
}
throw new MaxIterationsExceededException(maximalIterationCount);
}
}
© 2015 - 2024 Weber Informatics LLC | Privacy Policy