BracketedSingleOptimizer.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.optimization;
import com.irurueta.numerical.EvaluationException;
import com.irurueta.numerical.InvalidBracketRangeException;
import com.irurueta.numerical.LockedException;
import com.irurueta.numerical.NotAvailableException;
import com.irurueta.numerical.NotReadyException;
import com.irurueta.numerical.SingleDimensionFunctionEvaluatorListener;
/**
* This class searches for brackets of values containing a minimum in a single
* dimension function.
* A bracket is a set of points: "a" a minimum evaluation point,
* "b" a middle evaluation point and "c" a maximum evaluation where a <= b
* <= c, and where f(b) <= f(a) and f(b) <= f(c).
* This class uses a downhill algorithm that is better suited to continuous
* functions. Other functions might not obtain reliable results when using this
* algorithm to obtain a bracket of points.
* Some subclasses of this class will implement algorithms to refine the
* solution obtained in a bracket in order to find an accurate estimation of a
* minimum.
* Some algorithms might not need to previously compute a bracket and will
* simply search for a minimum in all the range of possible values, whereas
* other algorithms will require first the computation of a bracket.
* In either case, computing a bracket prior estimating a minimum will always
* ensure that a more reliable solution will be found.
* Besides, bracket computation is required when a function contains several
* minima and search of an accurate minimum estimation is desired to be
* restricted to a certain range of values.
*/
public abstract class BracketedSingleOptimizer extends SingleOptimizer {
/**
* The default ratio by which intervals are magnified and.
*/
public static final double GOLD = 1.618034;
/**
* The maximum magnification allowed for a parabolic-fit step.
*/
public static final double GLIMIT = 100.0;
/**
* Small value representing machine precision.
*/
public static final double TINY = 1e-20;
/**
* Default minimum evaluation point where the bracket is supposed to start
* By default, if no bracket is computed, the whole range of values is used
* for minimum estimation.
*/
public static final double DEFAULT_MIN_EVAL_POINT = -Double.MAX_VALUE;
/**
* Default middle evaluation point where the bracket is supposed to start
* By default, if no bracket is computed, the whole range of values is used
* for minimum estimation.
*/
public static final double DEFAULT_MIDDLE_EVAL_POINT = 0.0;
/**
* Default maximum evaluation point where the bracket is supposed to start
* By default, if no bracket is computed, the whole range of values is used
* for minimum estimation.
*/
public static final double DEFAULT_MAX_EVAL_POINT = Double.MAX_VALUE;
/**
* Minimum evaluation point inside the bracket.
*/
protected double ax;
/**
* Middle evaluation point inside the bracket.
*/
protected double bx;
/**
* Maximum evaluation point inside the bracket.
*/
protected double cx;
/**
* Boolean indicating whether a bracket has been provided or computed.
*/
private boolean bracketAvailable;
/**
* Function evaluation value at minimum evaluation point inside the bracket.
*/
private double fa;
/**
* Function evaluation value at middle evaluation point inside the bracket.
*/
private double fb;
/**
* Function evaluation value at maximum evaluation point inside the bracket.
*/
private double fc;
/**
* Boolean indicating whether function evaluation at bracket limits and
* middle point are available or not.
*/
private boolean bracketEvaluationAvailable;
/**
* Constructor. Creates an instance with provided bracket of values.
*
* @param minEvalPoint Minimum bracket evaluation point.
* @param middleEvalPoint Middle bracket evaluation point.
* @param maxEvalPoint Maximum bracket evaluation point.
* @throws InvalidBracketRangeException Raised if the following condition is
* not met: minEvalPoint <= middleEvalPoint <= maxEvalPoint.
*/
protected BracketedSingleOptimizer(
final double minEvalPoint, final double middleEvalPoint, final double maxEvalPoint)
throws InvalidBracketRangeException {
internalSetBracket(minEvalPoint, middleEvalPoint, maxEvalPoint);
}
/**
* Empty Constructor. Creates an instance using default bracket values.
*/
protected BracketedSingleOptimizer() {
ax = DEFAULT_MIN_EVAL_POINT;
bx = DEFAULT_MIDDLE_EVAL_POINT;
cx = DEFAULT_MAX_EVAL_POINT;
bracketAvailable = true;
}
/**
* Constructor. Creates an instance with provided bracket of values and a
* listener to get single dimension function evaluations.
*
* @param listener Listener to evaluate a function.
* @param minEvalPoint Minimum bracket evaluation point.
* @param middleEvalPoint Middle bracket evaluation point.
* @param maxEvalPoint Maximum bracket evaluation point.
* @throws InvalidBracketRangeException Raised if the following condition is
* not met: minEvalPoint <= middleEvalPoint <= maxEvalPoint.
*/
protected BracketedSingleOptimizer(
final SingleDimensionFunctionEvaluatorListener listener, final double minEvalPoint,
final double middleEvalPoint, final double maxEvalPoint) throws InvalidBracketRangeException {
super(listener);
internalSetBracket(minEvalPoint, middleEvalPoint, maxEvalPoint);
}
/**
* Sets a bracket of values to later search for a minimum. A local minimum
* will only be search within the minimum and maximum evaluation points of
* a given bracket.
* If bracket is not provided, it can also be computed from a default or
* coarse set of points in order to obtain a more refined bracket so that
* a minimum search can be estimated more precisely.
*
* @param minEvalPoint Minimum bracket evaluation point.
* @param middleEvalPoint Middle bracket evaluation point.
* @param maxEvalPoint Maximum bracket evaluation point.
* @throws InvalidBracketRangeException Raised if the following condition is
* not met: minEvalPoint <= middleEvalPoint <= maxEvalPoint.
* @throws LockedException Raised if this instance is locked. This instance
* will be locked while doing some operations. Attempting to change any
* parameter while being locked will raise this exception.
*/
public void setBracket(
final double minEvalPoint, final double middleEvalPoint, final double maxEvalPoint) throws LockedException,
InvalidBracketRangeException {
if (isLocked()) {
throw new LockedException();
}
internalSetBracket(minEvalPoint, middleEvalPoint, maxEvalPoint);
}
/**
* Returns boolean indicating whether a bracket has been provided or
* computed and is available for retrieval.
*
* @return true if a bracket has been provided, false otherwise.
*/
public boolean isBracketAvailable() {
return bracketAvailable;
}
/**
* Returns minimum evaluation point where the bracket starts
*
* @return Minimum evaluation point.
* @throws NotAvailableException Raised if not provided or computed.
*/
public double getMinEvaluationPoint() throws NotAvailableException {
if (!isBracketAvailable()) {
throw new NotAvailableException();
}
return ax;
}
/**
* Returns middle evaluation point within the bracket.
*
* @return Middle evaluation point.
* @throws NotAvailableException Raised if not provided or computed.
*/
public double getMiddleEvaluationPoint() throws NotAvailableException {
if (!isBracketAvailable()) {
throw new NotAvailableException();
}
return bx;
}
/**
* Returns maximum evaluation point whether the bracket finishes.
*
* @return Maximum evaluation point.
* @throws NotAvailableException Raised if not provided or computed.
*/
public double getMaxEvaluationPoint() throws NotAvailableException {
if (!isBracketAvailable()) {
throw new NotAvailableException();
}
return cx;
}
/**
* Returns single dimension function evaluation at provided or computed
* minimum evaluation point where the bracket starts.
*
* @return Function evaluation at bracket's minimum evaluation point.
* @throws NotAvailableException Raised if bracket evaluations are not
* available.
*/
public double getEvaluationAtMin() throws NotAvailableException {
if (!areBracketEvaluationsAvailable()) {
throw new NotAvailableException();
}
return fa;
}
/**
* Returns single dimension function evaluation at provided or computed
* middle evaluation point within the bracket.
*
* @return Function evaluation at bracket's middle evaluation point.
* @throws NotAvailableException Raised if bracket evaluations are not
* available.
*/
public double getEvaluationAtMiddle() throws NotAvailableException {
if (!areBracketEvaluationsAvailable()) {
throw new NotAvailableException();
}
return fb;
}
/**
* Returns single dimension function evaluation at provided or computed
* maximum evaluation point where the bracket finishes.
*
* @return Function evaluation at bracket's maximum evaluation point.
* @throws NotAvailableException Raised if bracket evaluations are not
* available.
*/
public double getEvaluationAtMax() throws NotAvailableException {
if (!areBracketEvaluationsAvailable()) {
throw new NotAvailableException();
}
return fc;
}
/**
* Computes a bracket of values using provided values as a starting point.
* Given a function f, and given distinct initial points ax and bx, this
* routine searches in the downhill direction (defined by the function as
* evaluated at the initial points) and returns.
* ax (minimum evaluation point), bx (middle evaluation point), cx (maximum
* evaluation point) that bracket a minimum of the function. Also returned
* are the function values at the three points fa, fb, and fc, which are the
* function evaluations at minimum, middle and maximum bracket points.
*
* @param minEvalPoint Initial minimum evaluation point of bracket.
* @param middleEvalPoint Initial middle evaluation point of bracket.
* @throws LockedException Raised if this instance is locked. This instance
* will be locked while doing some operations. Attempting to change any
* parameter while being locked will raise this exception.
* @throws NotReadyException Raised if this instance is not ready because a
* listener has not yet been provided.
* @throws InvalidBracketRangeException Raised if minEvalPoint <
* middleEvalPoint.
* @throws OptimizationException Raised if a bracket couldn't be found .
* because convergence was not achieved or function evaluation failed.
*/
@SuppressWarnings("Duplicates")
public void computeBracket(final double minEvalPoint, final double middleEvalPoint) throws LockedException,
NotReadyException, InvalidBracketRangeException, OptimizationException {
if (isLocked()) {
throw new LockedException();
}
if (!isReady()) {
throw new NotReadyException();
}
if (minEvalPoint > middleEvalPoint) {
throw new InvalidBracketRangeException();
}
locked = true;
final var a = new double[1];
final var b = new double[1];
final var c = new double[1];
try {
ax = minEvalPoint;
bx = middleEvalPoint;
double fu;
fa = listener.evaluate(ax);
fb = listener.evaluate(bx);
//switch roles of a and b so that we can go downhill in the
//direction from a to b
if (fb > fa) {
a[0] = ax;
b[0] = bx;
swap(a, b);
ax = a[0];
bx = b[0];
a[0] = fa;
b[0] = fb;
swap(a, b);
a[0] = fa;
b[0] = fb;
}
//First guess for c
cx = bx + GOLD * (bx - ax);
fc = listener.evaluate(cx);
//Keep returning here until we bracket.
while (fb > fc) {
//Compute u by parabolic extrapolation from a, b, c. TINY is
//used to prevent any possible division by zero.
final var r = (bx - ax) * (fb - fc);
final var q = (bx - cx) * (fb - fa);
var u = bx - ((bx - cx) * q - (bx - ax) * r) /
(2.0 * sign(Math.max(Math.abs(q - r), TINY), q - r));
final var ulim = bx + GLIMIT * (cx - bx);
//We won't go farther than this. Test various possibilities:
if ((bx - u) * (u - cx) > 0.0) {
//Parabolic u is between b and c: try it.
fu = listener.evaluate(u);
if (fu < fc) {
//Got a minimum between b and c.
ax = bx;
bx = u;
fa = fb;
fb = fu;
break;
} else if (fu > fb) {
//Got a minimum between a and u
cx = u;
fc = fu;
break;
}
//Parabolic fit was no use. Use default magnification.
u = cx + GOLD * (cx - bx);
fu = listener.evaluate(u);
} else if ((cx - u) * (u - ulim) > 0.0) {
//Parabolic fit is between c and its allowed limit
fu = listener.evaluate(u);
if (fu < fc) {
a[0] = bx;
b[0] = cx;
c[0] = u;
shift3(a, b, c, u + GOLD * (u - cx));
bx = a[0];
cx = b[0];
u = c[0];
a[0] = fb;
b[0] = fc;
c[0] = fu;
shift3(a, b, c, listener.evaluate(u));
fb = a[0];
fc = b[0];
fu = c[0];
}
} else if ((u - ulim) * (ulim - cx) >= 0.0) {
//Limit parabolic u to maximum allowed value.
u = ulim;
fu = listener.evaluate(u);
} else {
//Reject parabolic u, use default magnification
u = cx + GOLD * (cx - bx);
fu = listener.evaluate(u);
}
//Eliminate the oldest point and continue
a[0] = ax;
b[0] = bx;
c[0] = cx;
shift3(a, b, c, u);
ax = a[0];
bx = b[0];
cx = c[0];
a[0] = fa;
b[0] = fb;
c[0] = fc;
shift3(a, b, c, fu);
fa = a[0];
fb = b[0];
fc = c[0];
}
} catch (final EvaluationException e) {
throw new OptimizationException(e);
} finally {
locked = false;
}
bracketAvailable = true;
bracketEvaluationAvailable = true;
}
/**
* Computes a bracket of values using provided value as a starting point,
* and assuming that bracket finishes at Double.MAX_VALUE.
* Given a function f, and given distinct initial points ax and bx = 0.0,
* this routine searches in the downhill direction (defined by the function
* as evaluated at the initial points) and returns
* ax (minimum evaluation point), bx (middle evaluation point), cx (maximum
* evaluation point) that bracket a minimum of the function. Also returned
* are the function values at the three points fa, fb, and fc, which are the
* function evaluations at minimum, middle and maximum bracket points.
*
* @param minEvalPoint Initial minimum evaluation point of bracket.
* @throws LockedException Raised if this instance is locked. This instance
* will be locked while doing some operations. Attempting to change any
* parameter while being locked will raise this exception.
* @throws NotReadyException Raised if this instance is not ready because a
* listener has not yet been provided.
* @throws InvalidBracketRangeException Raised if minEvalPoint < 0.0.
* @throws OptimizationException Raised if a bracket couldn't be found
* because convergence was not achieved or function evaluation failed.
*/
public void computeBracket(final double minEvalPoint) throws LockedException, NotReadyException,
OptimizationException, InvalidBracketRangeException {
computeBracket(minEvalPoint, DEFAULT_MIDDLE_EVAL_POINT);
}
/**
* Computes a bracket of values using the whole range of possible values as
* an initial guess.
* Given a function f, and given distinct initial points ax =
* -Double.MAX_VALUE and bx = 0.0, this
* routine searches in the downhill direction (defined by the function as
* evaluated at the initial points) and returns
* ax (minimum evaluation point), bx (middle evaluation point), cx (maximum
* evaluation point) that bracket a minimum of the function. Also returned
* are the function values at the three points fa, fb, and fc, which are the
* function evaluations at minimum, middle and maximum bracket points
*
* @throws LockedException Raised if this instance is locked. This instance
* will be locked while doing some operations. Attempting to change any
* parameter while being locked will raise this exception.
* @throws NotReadyException Raised if this instance is not ready because a
* listener has not yet been provided.
* @throws OptimizationException Raised if a bracket couldn't be found
* because convergence was not achieved or function evaluation failed.
*/
public void computeBracket() throws LockedException, NotReadyException, OptimizationException {
try {
computeBracket(DEFAULT_MIN_EVAL_POINT, DEFAULT_MIDDLE_EVAL_POINT);
} catch (InvalidBracketRangeException ignore) {
//never happens
}
}
/**
* Computes function evaluations at provided or estimated bracket locations.
* After calling this method bracket evaluations will be available.
*
* @throws LockedException Raised if this instance is locked. This instance
* will be locked while doing some operations. Attempting to change any
* parameter while being locked will raise this exception.
* @throws NotReadyException Raised if this instance is not ready because a
* listener has not yet been provided.
* @throws OptimizationException Raised if function evaluation failed.
*/
public void evaluateBracket() throws LockedException, NotReadyException, OptimizationException {
if (isLocked()) {
throw new LockedException();
}
if (!isReady()) {
throw new NotReadyException();
}
locked = true;
try {
fa = listener.evaluate(ax);
fb = listener.evaluate(bx);
fc = listener.evaluate(cx);
} catch (final EvaluationException e) {
throw new OptimizationException(e);
} finally {
locked = false;
}
bracketEvaluationAvailable = true;
}
/**
* Returns boolean indicating whether bracket evaluations are available for
* retrieval.
*
* @return True if bracket evaluations are available, false otherwise.
*/
public boolean areBracketEvaluationsAvailable() {
return bracketEvaluationAvailable;
}
/**
* Internal method to determine whether a and b have the same sign.
*
* @param a Value to be compared.
* @param b Value to be compared.
* @return Returns "a" if "a" and "b" have the same sign or "-a" otherwise.
*/
protected double sign(final double a, final double b) {
if (b >= 0.0) {
return a >= 0.0 ? a : -a;
} else {
return a >= 0.0 ? -a : a;
}
}
/**
* Pushes "b" value into "a", and "c" value into "b". "a" and "b" are in/out parameters.
* Results will be available at a[0] and b[0] after executing this method.
*
* @param a a value to be lost.
* @param b a value to be shifted into "a".
* @param c a value to be shifted into "b".
*/
protected void shift2(final double[] a, final double[] b, final double c) {
a[0] = b[0];
b[0] = c;
}
/**
* Pushes "b" value into "a", and "c" value into "b" and "d" value into "c". "a", "b" and "c"
* are in/out parameters.
* Results will be available at a[0], b[0] and c[0] after executing this
* method.
*
* @param a a value to be lost.
* @param b a value to be shifted into "a".
* @param c a value to be shifted into "b".
* @param d a value to be shifted into "c".
*/
protected void shift3(final double[] a, final double[] b, final double[] c, final double d) {
a[0] = b[0];
b[0] = c[0];
c[0] = d;
}
/**
* Moves d, e and f into a[0], b[0] and c[0]. Previously existing values
* into a, b, c will be lost after executing this method.
*
* @param a a value to be set.
* @param b a value to be set.
* @param c a value to be set.
* @param d a value to be copied.
* @param e a value to be copied.
* @param f a value to be copied.
*/
protected void mov3(final double[] a, final double[] b, final double[] c, final double d, final double e,
final double f) {
a[0] = d;
b[0] = e;
c[0] = f;
}
/**
* Internal method to swap two values. Value inside a[0] will be swapped
* with value provided in b[0].
*
* @param a Value to be swapped.
* @param b Value to be swapped.
*/
private void swap(final double[] a, final double[] b) {
final var tmp = a[0];
a[0] = b[0];
b[0] = tmp;
}
/**
* Internal method to set a bracket of values. This method does not check
* whether this instance is locked.
*
* @param minEvalPoint Minimum bracket evaluation point.
* @param middleEvalPoint Middle bracket evaluation point.
* @param maxEvalPoint Maximum bracket evaluation point.
* @throws InvalidBracketRangeException Raised if the following condition is
* not met: minEvalPoint <= middleEvalPoint <= maxEvalPoint.
*/
private void internalSetBracket(final double minEvalPoint, final double middleEvalPoint,
final double maxEvalPoint) throws InvalidBracketRangeException {
if ((minEvalPoint > middleEvalPoint) || (middleEvalPoint > maxEvalPoint)) {
//which also means || (minEvalPoint > maxEvalPoint))
throw new InvalidBracketRangeException();
}
ax = minEvalPoint;
bx = middleEvalPoint;
cx = maxEvalPoint;
bracketAvailable = true;
}
}