You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by er...@apache.org on 2010/08/30 15:06:24 UTC
svn commit: r990792 [3/5] - in /commons/proper/math/trunk/src:
main/java/org/apache/commons/math/
main/java/org/apache/commons/math/optimization/
main/java/org/apache/commons/math/optimization/direct/
main/java/org/apache/commons/math/optimization/fitt...
Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/GaussNewtonOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/GaussNewtonOptimizer.java?rev=990792&r1=990791&r2=990792&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/GaussNewtonOptimizer.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/GaussNewtonOptimizer.java Mon Aug 30 13:06:22 2010
@@ -25,7 +25,7 @@ import org.apache.commons.math.linear.In
import org.apache.commons.math.linear.LUDecompositionImpl;
import org.apache.commons.math.linear.QRDecompositionImpl;
import org.apache.commons.math.linear.RealMatrix;
-import org.apache.commons.math.optimization.OptimizationException;
+import org.apache.commons.math.exception.ConvergenceException;
import org.apache.commons.math.optimization.VectorialPointValuePair;
/**
@@ -43,17 +43,17 @@ import org.apache.commons.math.optimizat
*/
public class GaussNewtonOptimizer extends AbstractLeastSquaresOptimizer {
-
/** Indicator for using LU decomposition. */
private final boolean useLU;
- /** Simple constructor with default settings.
- * <p>The convergence check is set to a {@link
- * org.apache.commons.math.optimization.SimpleVectorialValueChecker}
- * and the maximal number of evaluation is set to
- * {@link AbstractLeastSquaresOptimizer#DEFAULT_MAX_ITERATIONS}.
- * @param useLU if true, the normal equations will be solved using LU
- * decomposition, otherwise they will be solved using QR decomposition
+ /**
+ * Simple constructor with default settings.
+ * The convergence check is set to a {@link
+ * org.apache.commons.math.optimization.SimpleVectorialValueChecker}.
+ *
+ * @param useLU if {@code true}, the normal equations will be solved
+ * using LU decomposition, otherwise they will be solved using QR
+ * decomposition.
*/
public GaussNewtonOptimizer(final boolean useLU) {
this.useLU = useLU;
@@ -62,13 +62,13 @@ public class GaussNewtonOptimizer extend
/** {@inheritDoc} */
@Override
public VectorialPointValuePair doOptimize()
- throws FunctionEvaluationException, OptimizationException, IllegalArgumentException {
+ throws FunctionEvaluationException {
// iterate until convergence is reached
VectorialPointValuePair current = null;
+ int iter = 0;
for (boolean converged = false; !converged;) {
-
- incrementIterationsCounter();
+ ++iter;
// evaluate the objective function and its jacobian
VectorialPointValuePair previous = current;
@@ -76,6 +76,9 @@ public class GaussNewtonOptimizer extend
updateJacobian();
current = new VectorialPointValuePair(point, objective);
+ final double[] targetValues = getTargetRef();
+ final double[] residualsWeights = getWeightRef();
+
// build the linear problem
final double[] b = new double[cols];
final double[][] a = new double[cols][cols];
@@ -99,11 +102,9 @@ public class GaussNewtonOptimizer extend
ak[l] += wgk * grad[l];
}
}
-
}
try {
-
// solve the linearized least squares problem
RealMatrix mA = new BlockRealMatrix(a);
DecompositionSolver solver = useLU ?
@@ -115,21 +116,16 @@ public class GaussNewtonOptimizer extend
for (int i = 0; i < cols; ++i) {
point[i] += dX[i];
}
-
- } catch(InvalidMatrixException e) {
- throw new OptimizationException(LocalizedFormats.UNABLE_TO_SOLVE_SINGULAR_PROBLEM);
+ } catch (InvalidMatrixException e) {
+ throw new ConvergenceException(LocalizedFormats.UNABLE_TO_SOLVE_SINGULAR_PROBLEM);
}
// check convergence
if (previous != null) {
- converged = checker.converged(getIterations(), previous, current);
+ converged = getConvergenceChecker().converged(iter, previous, current);
}
-
}
-
// we have converged
return current;
-
}
-
}
Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizer.java?rev=990792&r1=990791&r2=990792&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizer.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizer.java Mon Aug 30 13:06:22 2010
@@ -19,9 +19,10 @@ package org.apache.commons.math.optimiza
import java.util.Arrays;
import org.apache.commons.math.FunctionEvaluationException;
+import org.apache.commons.math.exception.ConvergenceException;
import org.apache.commons.math.exception.util.LocalizedFormats;
-import org.apache.commons.math.optimization.OptimizationException;
import org.apache.commons.math.optimization.VectorialPointValuePair;
+import org.apache.commons.math.optimization.ConvergenceChecker;
import org.apache.commons.math.util.MathUtils;
import org.apache.commons.math.util.FastMath;
@@ -150,9 +151,8 @@ public class LevenbergMarquardtOptimizer
* Build an optimizer for least squares problems.
* <p>The default values for the algorithm settings are:
* <ul>
- * <li>{@link #setConvergenceChecker(VectorialConvergenceChecker) vectorial convergence checker}: null</li>
+ * <li>{@link #setConvergenceChecker(ConvergenceChecker) vectorial convergence checker}: null</li>
* <li>{@link #setInitialStepBoundFactor(double) initial step bound factor}: 100.0</li>
- * <li>{@link #setMaxIterations(int) maximal iterations}: 1000</li>
* <li>{@link #setCostRelativeTolerance(double) cost relative tolerance}: 1.0e-10</li>
* <li>{@link #setParRelativeTolerance(double) parameters relative tolerance}: 1.0e-10</li>
* <li>{@link #setOrthoTolerance(double) orthogonality tolerance}: 1.0e-10</li>
@@ -165,10 +165,6 @@ public class LevenbergMarquardtOptimizer
* and {@link #setParRelativeTolerance parameters relative tolerance} settings.
*/
public LevenbergMarquardtOptimizer() {
-
- // set up the superclass with a default max cost evaluations setting
- setMaxIterations(1000);
-
// default values for the tuning parameters
setConvergenceChecker(null);
setInitialStepBoundFactor(100.0);
@@ -176,7 +172,6 @@ public class LevenbergMarquardtOptimizer
setParRelativeTolerance(1.0e-10);
setOrthoTolerance(1.0e-10);
setQRRankingThreshold(MathUtils.SAFE_MIN);
-
}
/**
@@ -239,8 +234,7 @@ public class LevenbergMarquardtOptimizer
/** {@inheritDoc} */
@Override
- protected VectorialPointValuePair doOptimize()
- throws FunctionEvaluationException, OptimizationException, IllegalArgumentException {
+ protected VectorialPointValuePair doOptimize() throws FunctionEvaluationException {
// arrays shared with the other private methods
solvedCols = FastMath.min(rows, cols);
@@ -269,11 +263,14 @@ public class LevenbergMarquardtOptimizer
lmPar = 0;
boolean firstIteration = true;
VectorialPointValuePair current = new VectorialPointValuePair(point, objective);
+ int iter = 0;
+ final ConvergenceChecker<VectorialPointValuePair> checker = getConvergenceChecker();
while (true) {
+ ++iter;
+
for (int i=0;i<rows;i++) {
qtf[i]=weightedResiduals[i];
}
- incrementIterationsCounter();
// compute the Q.R. decomposition of the jacobian matrix
VectorialPointValuePair previous = current;
@@ -290,7 +287,6 @@ public class LevenbergMarquardtOptimizer
}
if (firstIteration) {
-
// scale the point according to the norms of the columns
// of the initial jacobian
xNorm = 0;
@@ -307,7 +303,6 @@ public class LevenbergMarquardtOptimizer
// initialize the step bound delta
delta = (xNorm == 0) ? initialStepBoundFactor : (initialStepBoundFactor * xNorm);
-
}
// check orthogonality between function vector and jacobian columns
@@ -425,17 +420,17 @@ public class LevenbergMarquardtOptimizer
xNorm = 0;
for (int k = 0; k < cols; ++k) {
double xK = diag[k] * point[k];
- xNorm += xK * xK;
+ xNorm += xK * xK;
}
xNorm = FastMath.sqrt(xNorm);
current = new VectorialPointValuePair(point, objective);
// tests for convergence.
if (checker != null) {
- // we use the vectorial convergence checker
- if (checker.converged(getIterations(), previous, current)) {
- return current;
- }
+ // we use the vectorial convergence checker
+ if (checker.converged(iter, previous, current)) {
+ return current;
+ }
}
} else {
// failed iteration, reset the previous values
@@ -462,20 +457,17 @@ public class LevenbergMarquardtOptimizer
// tests for termination and stringent tolerances
// (2.2204e-16 is the machine epsilon for IEEE754)
if ((FastMath.abs(actRed) <= 2.2204e-16) && (preRed <= 2.2204e-16) && (ratio <= 2.0)) {
- throw new OptimizationException(LocalizedFormats.TOO_SMALL_COST_RELATIVE_TOLERANCE,
+ throw new ConvergenceException(LocalizedFormats.TOO_SMALL_COST_RELATIVE_TOLERANCE,
costRelativeTolerance);
} else if (delta <= 2.2204e-16 * xNorm) {
- throw new OptimizationException(LocalizedFormats.TOO_SMALL_PARAMETERS_RELATIVE_TOLERANCE,
+ throw new ConvergenceException(LocalizedFormats.TOO_SMALL_PARAMETERS_RELATIVE_TOLERANCE,
parRelativeTolerance);
} else if (maxCosine <= 2.2204e-16) {
- throw new OptimizationException(LocalizedFormats.TOO_SMALL_ORTHOGONALITY_TOLERANCE,
+ throw new ConvergenceException(LocalizedFormats.TOO_SMALL_ORTHOGONALITY_TOLERANCE,
orthoTolerance);
}
-
}
-
}
-
}
/**
@@ -733,7 +725,6 @@ public class LevenbergMarquardtOptimizer
lmDiag[i] = -sin * rik + cos * lmDiag[i];
weightedResidualJacobian[i][pk] = temp2;
}
-
}
}
@@ -741,7 +732,6 @@ public class LevenbergMarquardtOptimizer
// the corresponding diagonal element of R
lmDiag[j] = weightedResidualJacobian[j][permutation[j]];
weightedResidualJacobian[j][permutation[j]] = lmDir[j];
-
}
// solve the triangular system for z, if the system is
@@ -770,7 +760,6 @@ public class LevenbergMarquardtOptimizer
for (int j = 0; j < lmDir.length; ++j) {
lmDir[permutation[j]] = work[j];
}
-
}
/**
@@ -793,9 +782,9 @@ public class LevenbergMarquardtOptimizer
* are performed in non-increasing columns norms order thanks to columns
* pivoting. The diagonal elements of the R matrix are therefore also in
* non-increasing absolute values order.</p>
- * @exception OptimizationException if the decomposition cannot be performed
+ * @exception ConvergenceException if the decomposition cannot be performed
*/
- private void qrDecomposition() throws OptimizationException {
+ private void qrDecomposition() throws ConvergenceException {
// initializations
for (int k = 0; k < cols; ++k) {
@@ -821,7 +810,7 @@ public class LevenbergMarquardtOptimizer
norm2 += aki * aki;
}
if (Double.isInfinite(norm2) || Double.isNaN(norm2)) {
- throw new OptimizationException(LocalizedFormats.UNABLE_TO_PERFORM_QR_DECOMPOSITION_ON_JACOBIAN,
+ throw new ConvergenceException(LocalizedFormats.UNABLE_TO_PERFORM_QR_DECOMPOSITION_ON_JACOBIAN,
rows, cols);
}
if (norm2 > ak2) {
@@ -858,11 +847,8 @@ public class LevenbergMarquardtOptimizer
weightedResidualJacobian[j][permutation[k + dk]] -= gamma * weightedResidualJacobian[j][pk];
}
}
-
}
-
rank = solvedCols;
-
}
/**
@@ -883,5 +869,4 @@ public class LevenbergMarquardtOptimizer
}
}
}
-
}
Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/NonLinearConjugateGradientOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/NonLinearConjugateGradientOptimizer.java?rev=990792&r1=990791&r2=990792&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/NonLinearConjugateGradientOptimizer.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/NonLinearConjugateGradientOptimizer.java Mon Aug 30 13:06:22 2010
@@ -17,14 +17,14 @@
package org.apache.commons.math.optimization.general;
-import org.apache.commons.math.ConvergenceException;
import org.apache.commons.math.FunctionEvaluationException;
+import org.apache.commons.math.exception.MathIllegalStateException;
+import org.apache.commons.math.exception.ConvergenceException;
import org.apache.commons.math.analysis.UnivariateRealFunction;
import org.apache.commons.math.analysis.solvers.BrentSolver;
import org.apache.commons.math.analysis.solvers.UnivariateRealSolver;
import org.apache.commons.math.exception.util.LocalizedFormats;
import org.apache.commons.math.optimization.GoalType;
-import org.apache.commons.math.optimization.OptimizationException;
import org.apache.commons.math.optimization.RealPointValuePair;
import org.apache.commons.math.util.FastMath;
@@ -54,11 +54,11 @@ public class NonLinearConjugateGradientO
/** Current point. */
private double[] point;
- /** Simple constructor with default settings.
- * <p>The convergence check is set to a {@link
- * org.apache.commons.math.optimization.SimpleVectorialValueChecker}
- * and the maximal number of iterations is set to
- * {@link AbstractScalarDifferentiableOptimizer#DEFAULT_MAX_ITERATIONS}.
+ /**
+ * Simple constructor with default settings.
+ * The convergence check is set to a {@link
+ * org.apache.commons.math.optimization.SimpleVectorialValueChecker}.
+ *
* @param updateFormula formula to use for updating the β parameter,
* must be one of {@link ConjugateGradientFormula#FLETCHER_REEVES} or {@link
* ConjugateGradientFormula#POLAK_RIBIERE}
@@ -110,104 +110,104 @@ public class NonLinearConjugateGradientO
/** {@inheritDoc} */
@Override
protected RealPointValuePair doOptimize()
- throws FunctionEvaluationException, OptimizationException, IllegalArgumentException {
- try {
- // initialization
- if (preconditioner == null) {
- preconditioner = new IdentityPreconditioner();
- }
- if (solver == null) {
- solver = new BrentSolver();
- }
- point = getStartPoint();
- final GoalType goal = getGoalType();
- final int n = point.length;
- double[] r = computeObjectiveGradient(point);
- if (goal == GoalType.MINIMIZE) {
- for (int i = 0; i < n; ++i) {
- r[i] = -r[i];
- }
- }
-
- // initial search direction
- double[] steepestDescent = preconditioner.precondition(point, r);
- double[] searchDirection = steepestDescent.clone();
-
- double delta = 0;
+ throws FunctionEvaluationException {
+ // Initialization.
+ if (preconditioner == null) {
+ preconditioner = new IdentityPreconditioner();
+ }
+ if (solver == null) {
+ solver = new BrentSolver();
+ }
+ point = getStartPoint();
+ final GoalType goal = getGoalType();
+ final int n = point.length;
+ double[] r = computeObjectiveGradient(point);
+ if (goal == GoalType.MINIMIZE) {
for (int i = 0; i < n; ++i) {
- delta += r[i] * searchDirection[i];
+ r[i] = -r[i];
}
+ }
- RealPointValuePair current = null;
- while (true) {
+ // Initial search direction.
+ double[] steepestDescent = preconditioner.precondition(point, r);
+ double[] searchDirection = steepestDescent.clone();
+
+ double delta = 0;
+ for (int i = 0; i < n; ++i) {
+ delta += r[i] * searchDirection[i];
+ }
- final double objective = computeObjectiveValue(point);
- RealPointValuePair previous = current;
- current = new RealPointValuePair(point, objective);
- if (previous != null) {
- if (getConvergenceChecker().converged(getIterations(), previous, current)) {
- // we have found an optimum
- return current;
- }
+ RealPointValuePair current = null;
+ int iter = 0;
+ while (true) {
+ ++iter;
+
+ final double objective = computeObjectiveValue(point);
+ RealPointValuePair previous = current;
+ current = new RealPointValuePair(point, objective);
+ if (previous != null) {
+ if (getConvergenceChecker().converged(iter, previous, current)) {
+ // We have found an optimum.
+ return current;
}
+ }
- incrementIterationsCounter();
-
- double dTd = 0;
- for (final double di : searchDirection) {
- dTd += di * di;
- }
+ double dTd = 0;
+ for (final double di : searchDirection) {
+ dTd += di * di;
+ }
- // find the optimal step in the search direction
- final UnivariateRealFunction lsf = new LineSearchFunction(searchDirection);
+ // Find the optimal step in the search direction.
+ final UnivariateRealFunction lsf = new LineSearchFunction(searchDirection);
+ try {
final double step = solver.solve(lsf, 0, findUpperBound(lsf, 0, initialStep));
- // validate new point
+ // Validate new point.
for (int i = 0; i < point.length; ++i) {
point[i] += step * searchDirection[i];
}
- r = computeObjectiveGradient(point);
- if (goal == GoalType.MINIMIZE) {
- for (int i = 0; i < n; ++i) {
- r[i] = -r[i];
- }
- }
+ } catch (org.apache.commons.math.ConvergenceException e) {
+ throw new ConvergenceException(); // XXX ugly workaround.
+ }
- // compute beta
- final double deltaOld = delta;
- final double[] newSteepestDescent = preconditioner.precondition(point, r);
- delta = 0;
+ r = computeObjectiveGradient(point);
+ if (goal == GoalType.MINIMIZE) {
for (int i = 0; i < n; ++i) {
- delta += r[i] * newSteepestDescent[i];
+ r[i] = -r[i];
}
+ }
- final double beta;
- if (updateFormula == ConjugateGradientFormula.FLETCHER_REEVES) {
- beta = delta / deltaOld;
- } else {
- double deltaMid = 0;
- for (int i = 0; i < r.length; ++i) {
- deltaMid += r[i] * steepestDescent[i];
- }
- beta = (delta - deltaMid) / deltaOld;
- }
- steepestDescent = newSteepestDescent;
+ // Compute beta.
+ final double deltaOld = delta;
+ final double[] newSteepestDescent = preconditioner.precondition(point, r);
+ delta = 0;
+ for (int i = 0; i < n; ++i) {
+ delta += r[i] * newSteepestDescent[i];
+ }
- // compute conjugate search direction
- if ((getIterations() % n == 0) || (beta < 0)) {
- // break conjugation: reset search direction
- searchDirection = steepestDescent.clone();
- } else {
- // compute new conjugate search direction
- for (int i = 0; i < n; ++i) {
- searchDirection[i] = steepestDescent[i] + beta * searchDirection[i];
- }
+ final double beta;
+ if (updateFormula == ConjugateGradientFormula.FLETCHER_REEVES) {
+ beta = delta / deltaOld;
+ } else {
+ double deltaMid = 0;
+ for (int i = 0; i < r.length; ++i) {
+ deltaMid += r[i] * steepestDescent[i];
}
-
+ beta = (delta - deltaMid) / deltaOld;
}
+ steepestDescent = newSteepestDescent;
- } catch (ConvergenceException ce) {
- throw new OptimizationException(ce);
+ // Compute conjugate search direction.
+ if (iter % n == 0 ||
+ beta < 0) {
+ // Break conjugation: reset search direction.
+ searchDirection = steepestDescent.clone();
+ } else {
+ // Compute new conjugate search direction.
+ for (int i = 0; i < n; ++i) {
+ searchDirection[i] = steepestDescent[i] + beta * searchDirection[i];
+ }
+ }
}
}
@@ -218,11 +218,12 @@ public class NonLinearConjugateGradientO
* @param h initial step to try
* @return b such that f(a) and f(b) have opposite signs
* @exception FunctionEvaluationException if the function cannot be computed
- * @exception OptimizationException if no bracket can be found
+ * @exception MathIllegalStateException if no bracket can be found
+ * @deprecated in 2.2 (must be replaced with "BracketFinder").
*/
private double findUpperBound(final UnivariateRealFunction f,
final double a, final double h)
- throws FunctionEvaluationException, OptimizationException {
+ throws FunctionEvaluationException {
final double yA = f.value(a);
double yB = yA;
for (double step = h; step < Double.MAX_VALUE; step *= FastMath.max(2, yA / yB)) {
@@ -232,7 +233,7 @@ public class NonLinearConjugateGradientO
return b;
}
}
- throw new OptimizationException(LocalizedFormats.UNABLE_TO_BRACKET_OPTIMUM_IN_LINE_SEARCH);
+ throw new MathIllegalStateException(LocalizedFormats.UNABLE_TO_BRACKET_OPTIMUM_IN_LINE_SEARCH);
}
/** Default identity preconditioner. */
Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/PowellOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/PowellOptimizer.java?rev=990792&r1=990791&r2=990792&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/PowellOptimizer.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/general/PowellOptimizer.java Mon Aug 30 13:06:22 2010
@@ -20,20 +20,24 @@ package org.apache.commons.math.optimiza
import java.util.Arrays;
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;
-import org.apache.commons.math.optimization.OptimizationException;
import org.apache.commons.math.optimization.RealPointValuePair;
+import org.apache.commons.math.optimization.ConvergenceChecker;
import org.apache.commons.math.optimization.univariate.AbstractUnivariateRealOptimizer;
import org.apache.commons.math.optimization.univariate.BracketFinder;
import org.apache.commons.math.optimization.univariate.BrentOptimizer;
+import org.apache.commons.math.optimization.univariate.UnivariateRealPointValuePair;
/**
* Powell algorithm.
* This code is translated and adapted from the Python version of this
* algorithm (as implemented in module {@code optimize.py} v0.5 of
* <em>SciPy</em>).
+ * <br/>
+ * The user is responsible for calling {@link
+ * #setConvergenceChecker(ConvergenceChecker) ConvergenceChecker}
+ * prior to using the optimizer.
*
* @version $Revision$ $Date$
* @since 2.2
@@ -41,56 +45,44 @@ import org.apache.commons.math.optimizat
public class PowellOptimizer
extends AbstractScalarOptimizer {
/**
- * Default relative tolerance for line search ({@value}).
- */
- public static final double DEFAULT_LS_RELATIVE_TOLERANCE = 1e-7;
- /**
- * Default absolute tolerance for line search ({@value}).
- */
- public static final double DEFAULT_LS_ABSOLUTE_TOLERANCE = 1e-11;
- /**
* Line search.
*/
- private final LineSearch line;
+ private LineSearch line = new LineSearch();
/**
- * Constructor with default line search tolerances (see the
- * {@link #PowellOptimizer(double,double) other constructor}).
- */
- public PowellOptimizer() {
- this(DEFAULT_LS_RELATIVE_TOLERANCE,
- DEFAULT_LS_ABSOLUTE_TOLERANCE);
- }
-
- /**
- * Constructor with default absolute line search tolerances (see
- * the {@link #PowellOptimizer(double,double) other constructor}).
+ * Set the convergence checker.
+ * It also indirectly sets the line search tolerances to the square-root
+ * of the correponding tolerances in the checker.
*
- * @param lsRelativeTolerance Relative error tolerance for
- * the line search algorithm ({@link BrentOptimizer}).
+ * @param checker Convergence checker.
*/
- public PowellOptimizer(double lsRelativeTolerance) {
- this(lsRelativeTolerance,
- DEFAULT_LS_ABSOLUTE_TOLERANCE);
+ public void setConvergenceChecker(ConvergenceChecker<RealPointValuePair> checker) {
+ super.setConvergenceChecker(checker);
+
+ // Line search tolerances can be much lower than the tolerances
+ // required for the optimizer itself.
+ final double minTol = 1e-4;
+ final double rel = Math.min(Math.sqrt(checker.getRelativeThreshold()), minTol);
+ final double abs = Math.min(Math.sqrt(checker.getAbsoluteThreshold()), minTol);
+ line.setConvergenceChecker(new BrentOptimizer.BrentConvergenceChecker(rel, abs));
}
- /**
- * @param lsRelativeTolerance Relative error tolerance for
- * the line search algorithm ({@link BrentOptimizer}).
- * @param lsAbsoluteTolerance Relative error tolerance for
- * the line search algorithm ({@link BrentOptimizer}).
- */
- public PowellOptimizer(double lsRelativeTolerance,
- double lsAbsoluteTolerance) {
- line = new LineSearch(lsRelativeTolerance,
- lsAbsoluteTolerance);
+ /** {@inheritDoc} */
+ @Override
+ public void setMaxEvaluations(int maxEvaluations) {
+ super.setMaxEvaluations(maxEvaluations);
+
+ // We must allow at least as many iterations to the underlying line
+ // search optimizer. Because the line search inner class will call
+ // "computeObjectiveValue" in this class, we ensure that this class
+ // will be the first to eventually throw "TooManyEvaluationsException".
+ line.setMaxEvaluations(maxEvaluations);
}
/** {@inheritDoc} */
@Override
protected RealPointValuePair doOptimize()
- throws FunctionEvaluationException,
- OptimizationException {
+ throws FunctionEvaluationException {
final GoalType goal = getGoalType();
final double[] guess = getStartPoint();
final int n = guess.length;
@@ -103,8 +95,9 @@ public class PowellOptimizer
double[] x = guess;
double fVal = computeObjectiveValue(x);
double[] x1 = x.clone();
+ int iter = 0;
while (true) {
- incrementIterationsCounter();
+ ++iter;
double fX = fVal;
double fX2 = 0;
@@ -117,9 +110,9 @@ public class PowellOptimizer
fX2 = fVal;
- line.search(x, d);
- fVal = line.getValueAtOptimum();
- alphaMin = line.getOptimum();
+ final UnivariateRealPointValuePair optimum = line.search(x, d);
+ fVal = optimum.getValue();
+ alphaMin = optimum.getPoint();
final double[][] result = newPointAndDirection(x, d, alphaMin);
x = result[0];
@@ -131,7 +124,7 @@ public class PowellOptimizer
final RealPointValuePair previous = new RealPointValuePair(x1, fX);
final RealPointValuePair current = new RealPointValuePair(x, fVal);
- if (getConvergenceChecker().converged(getIterations(), previous, current)) {
+ if (getConvergenceChecker().converged(iter, previous, current)) {
if (goal == GoalType.MINIMIZE) {
return (fVal < fX) ? current : previous;
} else {
@@ -157,9 +150,9 @@ public class PowellOptimizer
t -= delta * temp * temp;
if (t < 0.0) {
- line.search(x, d);
- fVal = line.getValueAtOptimum();
- alphaMin = line.getOptimum();
+ final UnivariateRealPointValuePair optimum = line.search(x, d);
+ fVal = optimum.getValue();
+ alphaMin = optimum.getPoint();
final double[][] result = newPointAndDirection(x, d, alphaMin);
x = result[0];
@@ -200,11 +193,7 @@ public class PowellOptimizer
* Class for finding the minimum of the objective function along a given
* direction.
*/
- private class LineSearch {
- /**
- * Optimizer.
- */
- private final AbstractUnivariateRealOptimizer optim = new BrentOptimizer();
+ private class LineSearch extends BrentOptimizer {
/**
* Automatic bracketing.
*/
@@ -212,78 +201,41 @@ public class PowellOptimizer
/**
* Value of the optimum.
*/
- private double optimum = Double.NaN;
- /**
- * Value of the objective function at the optimum.
- */
- private double valueAtOptimum = Double.NaN;
-
- /**
- * @param relativeTolerance Relative tolerance.
- * @param absoluteTolerance Absolute tolerance.
- */
- public LineSearch(double relativeTolerance,
- double absoluteTolerance) {
- optim.setRelativeAccuracy(relativeTolerance);
- optim.setAbsoluteAccuracy(absoluteTolerance);
- }
+ private UnivariateRealPointValuePair optimum;
/**
* Find the minimum of the function {@code f(p + alpha * d)}.
*
* @param p Starting point.
* @param d Search direction.
- * @throws OptimizationException if function cannot be evaluated at some test point
- * or algorithm fails to converge
- */
- public void search(final double[] p,
- final double[] d)
- throws OptimizationException {
-
- // Reset.
- optimum = Double.NaN;
- valueAtOptimum = Double.NaN;
-
- try {
- final int n = p.length;
- final UnivariateRealFunction f = new UnivariateRealFunction() {
- public double value(double alpha)
- throws FunctionEvaluationException {
-
- final double[] x = new double[n];
- for (int i = 0; i < n; i++) {
- x[i] = p[i] + alpha * d[i];
- }
- final double obj = computeObjectiveValue(x);
- return obj;
- }
- };
-
- final GoalType goal = getGoalType();
- bracket.search(f, goal, 0, 1);
- optimum = optim.optimize(f, goal,
- bracket.getLo(),
- bracket.getHi(),
- bracket.getMid());
- valueAtOptimum = optim.getFunctionValue();
- } catch (FunctionEvaluationException e) {
- throw new OptimizationException(e);
- } catch (MaxIterationsExceededException e) {
- throw new OptimizationException(e);
- }
- }
-
- /**
* @return the optimum.
- */
- public double getOptimum() {
- return optimum;
- }
- /**
- * @return the value of the function at the optimum.
- */
- public double getValueAtOptimum() {
- return valueAtOptimum;
+ * @throws FunctionEvaluationException if the function evaluation
+ * fails.
+ * @throws TooManyEvaluationsException if the number of evaluations is
+ * exceeded.
+ */
+ public UnivariateRealPointValuePair search(final double[] p,
+ final double[] d)
+ throws FunctionEvaluationException {
+
+ final int n = p.length;
+ final UnivariateRealFunction f = new UnivariateRealFunction() {
+ public double value(double alpha)
+ throws FunctionEvaluationException {
+
+ final double[] x = new double[n];
+ for (int i = 0; i < n; i++) {
+ x[i] = p[i] + alpha * d[i];
+ }
+ final double obj = PowellOptimizer.this.computeObjectiveValue(x);
+ return obj;
+ }
+ };
+
+ final GoalType goal = PowellOptimizer.this.getGoalType();
+ bracket.search(f, goal, 0, 1);
+ return optimize(f, goal, bracket.getLo(), bracket.getHi(),
+ bracket.getMid());
}
}
}
Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/package.html
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/package.html?rev=990792&r1=990791&r2=990792&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/package.html (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/package.html Mon Aug 30 13:06:22 2010
@@ -33,7 +33,7 @@ by changing its input variables set unti
interfaces defining the common behavior of optimizers, one for each supported type of objective
function:
<ul>
- <li>{@link org.apache.commons.math.optimization.UnivariateRealOptimizer
+ <li>{@link org.apache.commons.math.optimization.univariate.UnivariateRealOptimizer
UnivariateRealOptimizer} for {@link org.apache.commons.math.analysis.UnivariateRealFunction
univariate real functions}</li>
<li>{@link org.apache.commons.math.optimization.MultivariateRealOptimizer
Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/AbstractUnivariateRealOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/AbstractUnivariateRealOptimizer.java?rev=990792&r1=990791&r2=990792&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/AbstractUnivariateRealOptimizer.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/AbstractUnivariateRealOptimizer.java Mon Aug 30 13:06:22 2010
@@ -17,14 +17,13 @@
package org.apache.commons.math.optimization.univariate;
-import org.apache.commons.math.ConvergingAlgorithmImpl;
import org.apache.commons.math.FunctionEvaluationException;
-import org.apache.commons.math.MaxEvaluationsExceededException;
-import org.apache.commons.math.MaxIterationsExceededException;
+import org.apache.commons.math.util.Incrementor;
+import org.apache.commons.math.exception.MaxCountExceededException;
+import org.apache.commons.math.exception.TooManyEvaluationsException;
import org.apache.commons.math.analysis.UnivariateRealFunction;
-import org.apache.commons.math.exception.NoDataException;
import org.apache.commons.math.optimization.GoalType;
-import org.apache.commons.math.optimization.UnivariateRealOptimizer;
+import org.apache.commons.math.optimization.ConvergenceChecker;
/**
* Provide a default implementation for several functions useful to generic
@@ -34,19 +33,13 @@ import org.apache.commons.math.optimizat
* @since 2.0
*/
public abstract class AbstractUnivariateRealOptimizer
- extends ConvergingAlgorithmImpl implements UnivariateRealOptimizer {
- /** Indicates where a root has been computed. */
- private boolean resultComputed;
- /** The last computed root. */
- private double result;
- /** Value of the function at the last computed result. */
- private double functionValue;
- /** Maximal number of evaluations allowed. */
- private int maxEvaluations;
- /** Number of evaluations already performed. */
- private int evaluations;
+ implements UnivariateRealOptimizer {
+ /** Convergence checker. */
+ private ConvergenceChecker<UnivariateRealPointValuePair> checker;
+ /** Evaluations counter. */
+ private final Incrementor evaluations = new Incrementor();
/** Optimization type */
- private GoalType optimizationGoal;
+ private GoalType goal;
/** Lower end of search interval. */
private double searchMin;
/** Higher end of search interval. */
@@ -56,115 +49,35 @@ public abstract class AbstractUnivariate
/** Function to optimize. */
private UnivariateRealFunction function;
- /**
- * Construct a solver with given iteration count and accuracy.
- * FunctionEvaluationExceptionFunctionEvaluationException
- * @param defaultAbsoluteAccuracy maximum absolute error
- * @param defaultMaximalIterationCount maximum number of iterations
- * @throws IllegalArgumentException if f is null or the
- * defaultAbsoluteAccuracy is not valid
- * @deprecated in 2.2. Please use the "setter" methods to assign meaningful
- * values to the maximum numbers of iterations and evaluations, and to the
- * absolute and relative accuracy thresholds.
- */
- protected AbstractUnivariateRealOptimizer(final int defaultMaximalIterationCount,
- final double defaultAbsoluteAccuracy) {
- super(defaultMaximalIterationCount, defaultAbsoluteAccuracy);
- resultComputed = false;
- setMaxEvaluations(Integer.MAX_VALUE);
- }
-
- /**
- * Default constructor.
- * To be removed once the single non-default one has been removed.
- */
- protected AbstractUnivariateRealOptimizer() {}
-
- /**
- * Check whether a result has been computed.
- * @throws NoDataException if no result has been computed
- * @deprecated in 2.2 (no alternative).
- */
- protected void checkResultComputed() {
- if (!resultComputed) {
- throw new NoDataException();
- }
- }
-
- /** {@inheritDoc} */
- public double getResult() {
- if (!resultComputed) {
- throw new NoDataException();
- }
- return result;
- }
-
- /** {@inheritDoc} */
- public double getFunctionValue() {
- if (functionValue == Double.NaN) {
- final double opt = getResult();
- try {
- functionValue = function.value(opt);
- } catch (FunctionEvaluationException e) {
- throw new RuntimeException(e);
- }
- }
- return functionValue;
- }
-
- /**
- * Convenience function for implementations.
- *
- * @param x the result to set
- * @param fx the result to set
- * @param iterationCount the iteration count to set
- * @deprecated in 2.2 (no alternative).
- */
- protected final void setResult(final double x, final double fx,
- final int iterationCount) {
- this.result = x;
- this.functionValue = fx;
- this.iterationCount = iterationCount;
- this.resultComputed = true;
- }
-
- /**
- * Convenience function for implementations.
- * @deprecated in 2.2 (no alternative).
- */
- protected final void clearResult() {
- this.resultComputed = false;
- }
-
/** {@inheritDoc} */
public void setMaxEvaluations(int maxEvaluations) {
- this.maxEvaluations = maxEvaluations;
+ evaluations.setMaximalCount(maxEvaluations);
}
/** {@inheritDoc} */
public int getMaxEvaluations() {
- return maxEvaluations;
+ return evaluations.getMaximalCount();
}
/** {@inheritDoc} */
public int getEvaluations() {
- return evaluations;
+ return evaluations.getCount();
}
/**
* @return the optimization type.
*/
public GoalType getGoalType() {
- return optimizationGoal;
+ return goal;
}
/**
- * @return the lower of the search interval.
+ * @return the lower end of the search interval.
*/
public double getMin() {
return searchMin;
}
/**
- * @return the higher of the search interval.
+ * @return the higher end of the search interval.
*/
public double getMax() {
return searchMax;
@@ -178,91 +91,74 @@ public abstract class AbstractUnivariate
/**
* Compute the objective function value.
- * @param f objective function
- * @param point point at which the objective function must be evaluated
- * @return objective function value at specified point
- * @exception FunctionEvaluationException if the function cannot be evaluated
- * or the maximal number of iterations is exceeded
- * @deprecated in 2.2. Use this {@link #computeObjectiveValue(double)
- * replacement} instead.
- */
- protected double computeObjectiveValue(final UnivariateRealFunction f,
- final double point)
- throws FunctionEvaluationException {
- if (++evaluations > maxEvaluations) {
- throw new FunctionEvaluationException(new MaxEvaluationsExceededException(maxEvaluations),
- point);
- }
- return f.value(point);
- }
-
- /**
- * Compute the objective function value.
*
* @param point Point at which the objective function must be evaluated.
* @return the objective function value at specified point.
- * @exception FunctionEvaluationException if the function cannot be evaluated
- * or the maximal number of iterations is exceeded.
+ * @throws FunctionEvaluationException if the function cannot be
+ * evaluated.
+ * @throws TooManyEvaluationsException if the maximal number of evaluations
+ * is exceeded.
*/
protected double computeObjectiveValue(double point)
throws FunctionEvaluationException {
- if (++evaluations > maxEvaluations) {
- resultComputed = false;
- throw new FunctionEvaluationException(new MaxEvaluationsExceededException(maxEvaluations),
- point);
+ try {
+ evaluations.incrementCount();
+ } catch (MaxCountExceededException e) {
+ throw new TooManyEvaluationsException(e.getMax());
}
return function.value(point);
}
/** {@inheritDoc} */
- public double optimize(UnivariateRealFunction f, GoalType goal,
- double min, double max, double startValue)
- throws MaxIterationsExceededException, FunctionEvaluationException {
- // Initialize.
- this.searchMin = min;
- this.searchMax = max;
- this.searchStart = startValue;
- this.optimizationGoal = goal;
- this.function = f;
-
+ public UnivariateRealPointValuePair optimize(UnivariateRealFunction f,
+ GoalType goalType,
+ double min, double max,
+ double startValue)
+ throws FunctionEvaluationException {
// Reset.
- functionValue = Double.NaN;
- evaluations = 0;
- resetIterationsCounter();
+ searchMin = min;
+ searchMax = max;
+ searchStart = startValue;
+ goal = goalType;
+ function = f;
+ evaluations.resetCount();
// Perform computation.
- result = doOptimize();
- resultComputed = true;
+ return doOptimize();
+ }
- return result;
+ /** {@inheritDoc} */
+ public UnivariateRealPointValuePair optimize(UnivariateRealFunction f,
+ GoalType goal,
+ double min, double max)
+ throws FunctionEvaluationException {
+ return optimize(f, goal, min, max, min + 0.5 * (max - min));
}
/**
- * Set the value at the optimum.
- *
- * @param functionValue Value of the objective function at the optimum.
+ * {@inheritDoc}
*/
- protected void setFunctionValue(double functionValue) {
- this.functionValue = functionValue;
+ public void setConvergenceChecker(ConvergenceChecker<UnivariateRealPointValuePair> checker) {
+ this.checker = checker;
}
- /** {@inheritDoc} */
- public double optimize(UnivariateRealFunction f, GoalType goal,
- double min, double max)
- throws MaxIterationsExceededException, FunctionEvaluationException {
- return optimize(f, goal, min, max, min + 0.5 * (max - min));
+ /**
+ * {@inheritDoc}
+ */
+ public ConvergenceChecker<UnivariateRealPointValuePair> getConvergenceChecker() {
+ return checker;
}
/**
* Method for implementing actual optimization algorithms in derived
* classes.
*
- * @return the optimum.
- * @throws MaxIterationsExceededException if the maximum iteration count
+ * @return the optimum and its corresponding function value.
+ * @throws TooManyEvaluationsException if the maximal number of evaluations
* is exceeded.
* @throws FunctionEvaluationException if an error occurs evaluating
* the function.
*/
- protected abstract double doOptimize()
- throws MaxIterationsExceededException, FunctionEvaluationException;
+ protected abstract UnivariateRealPointValuePair doOptimize()
+ throws FunctionEvaluationException;
}
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BaseUnivariateRealOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BaseUnivariateRealOptimizer.java?rev=990792&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BaseUnivariateRealOptimizer.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BaseUnivariateRealOptimizer.java Mon Aug 30 13:06:22 2010
@@ -0,0 +1,88 @@
+/*
+ * 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.analysis.UnivariateRealFunction;
+import org.apache.commons.math.optimization.BaseOptimizer;
+import org.apache.commons.math.optimization.GoalType;
+
+/**
+ * This interface is mainly intended to enforce the internal coherence of
+ * Commons-Math. Users of the API are advised to base their code on
+ * the following interfaces:
+ * <ul>
+ * <li>{@link org.apache.commons.math.optimization.univariate.UnivariateRealOptimizer}</li>
+ * </ul>
+ *
+ * @param <FUNC> Type of the objective function to be optimized.
+ *
+ * @version $Revision$ $Date$
+ * @since 3.0
+ */
+public interface BaseUnivariateRealOptimizer<FUNC extends UnivariateRealFunction>
+ extends BaseOptimizer<UnivariateRealPointValuePair> {
+ /**
+ * Find an optimum in the given interval.
+ *
+ * An optimizer may require that the interval brackets a single optimum.
+ *
+ * @param f Function to optimize.
+ * @param goalType Type of optimization goal: either
+ * {@link GoalType#MAXIMIZE} or {@link GoalType#MINIMIZE}.
+ * @param min Lower bound for the interval.
+ * @param max Upper bound for the interval.
+ * @return a (point, value) pair where the function is optimum.
+ * @throws {@link org.apache.commons.math.exception.TooManyEvaluationsException}
+ * if the maximum evaluation count is exceeded.
+ * @throws {@link org.apache.commons.math.exception.ConvergenceException}
+ * if the optimizer detects a convergence problem.
+ * @throws FunctionEvaluationException if an error occurs evaluating the
+ * function.
+ * @throws IllegalArgumentException if {@code min > max} or the endpoints
+ * do not satisfy the requirements specified by the optimizer.
+ */
+ UnivariateRealPointValuePair optimize(FUNC f, GoalType goalType,
+ double min, double max)
+ throws FunctionEvaluationException;
+
+ /**
+ * Find an optimum in the given interval, start at startValue.
+ * An optimizer may require that the interval brackets a single optimum.
+ *
+ * @param f Function to optimize.
+ * @param goalType Type of optimization goal: either
+ * {@link GoalType#MAXIMIZE} or {@link GoalType#MINIMIZE}.
+ * @param min Lower bound for the interval.
+ * @param max Upper bound for the interval.
+ * @param startValue Start value to use.
+ * @return a (point, value) pair where the function is optimum.
+ * @throws {@link org.apache.commons.math.exception.TooManyEvaluationsException}
+ * if the maximum evaluation count is exceeded.
+ * @throws {@link org.apache.commons.math.exception.ConvergenceException}
+ * if the optimizer detects a convergence problem.
+ * @throws FunctionEvaluationException if an error occurs evaluating the
+ * function.
+ * @throws IllegalArgumentException if {@code min > max} or the endpoints
+ * do not satisfy the requirements specified by the optimizer.
+ */
+ UnivariateRealPointValuePair optimize(FUNC f, GoalType goalType,
+ double min, double max,
+ double startValue)
+ throws FunctionEvaluationException;
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BaseUnivariateRealOptimizer.java
------------------------------------------------------------------------------
svn:eol-style = native
Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BracketFinder.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BracketFinder.java?rev=990792&r1=990791&r2=990792&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BracketFinder.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BracketFinder.java Mon Aug 30 13:06:22 2010
@@ -16,9 +16,11 @@
*/
package org.apache.commons.math.optimization.univariate;
+import org.apache.commons.math.util.Incrementor;
import org.apache.commons.math.exception.NotStrictlyPositiveException;
+import org.apache.commons.math.exception.TooManyEvaluationsException;
+import org.apache.commons.math.exception.MaxCountExceededException;
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;
@@ -26,6 +28,7 @@ import org.apache.commons.math.optimizat
* Provide an interval that brackets a local optimum of a function.
* This code is based on a Python implementation (from <em>SciPy</em>,
* module {@code optimize.py} v0.5).
+ *
* @version $Revision$ $Date$
* @since 2.2
*/
@@ -41,17 +44,9 @@ public class BracketFinder {
*/
private final double growLimit;
/**
- * Maximum number of iterations.
+ * Counter for function evaluations.
*/
- private final int maxIterations;
- /**
- * Number of iterations.
- */
- private int iterations;
- /**
- * Number of function evaluations.
- */
- private int evaluations;
+ private final Incrementor evaluations = new Incrementor();
/**
* Lower bound of the bracket.
*/
@@ -89,20 +84,20 @@ public class BracketFinder {
* Create a bracketing interval finder.
*
* @param growLimit Expanding factor.
- * @param maxIterations Maximum number of iterations allowed for finding
+ * @param maxEvaluations Maximum number of evaluations allowed for finding
* a bracketing interval.
*/
public BracketFinder(double growLimit,
- int maxIterations) {
+ int maxEvaluations) {
if (growLimit <= 0) {
throw new NotStrictlyPositiveException(growLimit);
}
- if (maxIterations <= 0) {
- throw new NotStrictlyPositiveException(maxIterations);
+ if (maxEvaluations <= 0) {
+ throw new NotStrictlyPositiveException(maxEvaluations);
}
this.growLimit = growLimit;
- this.maxIterations = maxIterations;
+ evaluations.setMaximalCount(maxEvaluations);
}
/**
@@ -112,7 +107,7 @@ public class BracketFinder {
* @param goal {@link GoalType Goal type}.
* @param xA Initial point.
* @param xB Initial point.
- * @throws MaxIterationsExceededException if the maximum iteration count
+ * @throws TooManyEvaluationsException if the maximum number of evaluations
* is exceeded.
* @throws FunctionEvaluationException if an error occurs evaluating
* the function.
@@ -121,9 +116,8 @@ public class BracketFinder {
GoalType goal,
double xA,
double xB)
- throws MaxIterationsExceededException,
- FunctionEvaluationException {
- reset();
+ throws FunctionEvaluationException {
+ evaluations.resetCount();
final boolean isMinim = goal == GoalType.MINIMIZE;
double fA = eval(func, xA);
@@ -131,6 +125,7 @@ public class BracketFinder {
if (isMinim ?
fA < fB :
fA > fB) {
+
double tmp = xA;
xA = xB;
xB = tmp;
@@ -144,10 +139,6 @@ public class BracketFinder {
double fC = eval(func, xC);
while (isMinim ? fC < fB : fC > fB) {
- if (++iterations > maxIterations) {
- throw new MaxIterationsExceededException(maxIterations);
- }
-
double tmp1 = (xB - xA) * (fB - fC);
double tmp2 = (xB - xC) * (fB - fA);
@@ -187,7 +178,7 @@ public class BracketFinder {
fW > fC) {
xB = xC;
xC = w;
- w = xC + GOLD * (xC -xB);
+ w = xC + GOLD * (xC - xB);
fB = fC;
fC =fW;
fW = eval(func, w);
@@ -198,37 +189,48 @@ public class BracketFinder {
}
xA = xB;
- xB = xC;
- xC = w;
fA = fB;
+ xB = xC;
fB = fC;
+ xC = w;
fC = fW;
}
lo = xA;
- mid = xB;
- hi = xC;
fLo = fA;
+ mid = xB;
fMid = fB;
+ hi = xC;
fHi = fC;
+
+ if (lo > hi) {
+ double tmp = lo;
+ lo = hi;
+ hi = tmp;
+
+ tmp = fLo;
+ fLo = fHi;
+ fHi = tmp;
+ }
}
/**
- * @return the number of iterations.
+ * @return the number of evalutations.
*/
- public int getIterations() {
- return iterations;
+ public int getMaxEvaluations() {
+ return evaluations.getMaximalCount();
}
+
/**
* @return the number of evalutations.
*/
public int getEvaluations() {
- return evaluations;
+ return evaluations.getCount();
}
/**
* @return the lower bound of the bracket.
- * @see #getFLow()
+ * @see #getFLo()
*/
public double getLo() {
return lo;
@@ -238,7 +240,7 @@ public class BracketFinder {
* Get function value at {@link #getLo()}.
* @return function value at {@link #getLo()}
*/
- public double getFLow() {
+ public double getFLo() {
return fLo;
}
@@ -278,21 +280,18 @@ public class BracketFinder {
* @param f Function.
* @param x Argument.
* @return {@code f(x)}
- * @throws FunctionEvaluationException if function cannot be evaluated at x
+ * @throws FunctionEvaluationException if function cannot be evaluated.
+ * @throws TooManyEvaluationsException if the maximal number of evaluations is
+ * exceeded.
*/
private double eval(UnivariateRealFunction f,
double x)
throws FunctionEvaluationException {
-
- ++evaluations;
+ try {
+ evaluations.incrementCount();
+ } catch (MaxCountExceededException e) {
+ throw new TooManyEvaluationsException(e.getMax());
+ }
return f.value(x);
}
-
- /**
- * Reset internal state.
- */
- private void reset() {
- iterations = 0;
- evaluations = 0;
- }
}
Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BrentOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BrentOptimizer.java?rev=990792&r1=990791&r2=990792&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BrentOptimizer.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/univariate/BrentOptimizer.java Mon Aug 30 13:06:22 2010
@@ -17,16 +17,28 @@
package org.apache.commons.math.optimization.univariate;
import org.apache.commons.math.FunctionEvaluationException;
-import org.apache.commons.math.MaxIterationsExceededException;
+import org.apache.commons.math.util.MathUtils;
+import org.apache.commons.math.util.FastMath;
+import org.apache.commons.math.exception.DimensionMismatchException;
+import org.apache.commons.math.exception.NumberIsTooSmallException;
import org.apache.commons.math.exception.NotStrictlyPositiveException;
+import org.apache.commons.math.exception.MathUnsupportedOperationException;
+import org.apache.commons.math.optimization.ConvergenceChecker;
+import org.apache.commons.math.optimization.AbstractConvergenceChecker;
import org.apache.commons.math.optimization.GoalType;
-import org.apache.commons.math.util.FastMath;
/**
* Implements Richard Brent's algorithm (from his book "Algorithms for
* Minimization without Derivatives", p. 79) for finding minima of real
* univariate functions. This implementation is an adaptation partly
* based on the Python code from SciPy (module "optimize.py" v0.5).
+ * If the function is defined on some interval {@code (lo, hi)}, then
+ * this method finds an approximation {@code x} to the point at which
+ * the function attains its minimum.
+ * <br/>
+ * The user is responsible for calling {@link
+ * #setConvergenceChecker(ConvergenceChecker) ConvergenceChecker}
+ * prior to using the optimizer.
*
* @version $Revision$ $Date$
* @since 2.0
@@ -38,57 +50,105 @@ public class BrentOptimizer extends Abst
private static final double GOLDEN_SECTION = 0.5 * (3 - FastMath.sqrt(5));
/**
- * Construct a solver.
+ * Convergence checker that implements the original stopping criterion
+ * of Brent's algorithm.
+ * {@code abs} and {@code rel} define a tolerance
+ * {@code tol = rel |x| + abs}. {@code rel} should be no smaller than
+ * <em>2 macheps</em> and preferably not much less than <em>sqrt(macheps)</em>,
+ * where <em>macheps</em> is the relative machine precision. {@code abs} must
+ * be positive.
+ *
+ * @since 3.0
*/
- public BrentOptimizer() {
- setMaxEvaluations(1000);
- setMaximalIterationCount(100);
- setAbsoluteAccuracy(1e-11);
- setRelativeAccuracy(1e-9);
- }
+ public static class BrentConvergenceChecker
+ extends AbstractConvergenceChecker<UnivariateRealPointValuePair> {
+ /**
+ * Minimum relative tolerance.
+ */
+ private static final double MIN_RELATIVE_TOLERANCE = 2 * FastMath.ulp(1d);
+
+ /**
+ * Build an instance with specified thresholds.
+ *
+ * @param rel Relative tolerance threshold
+ * @param abs Absolute tolerance threshold
+ */
+ public BrentConvergenceChecker(final double rel,
+ final double abs) {
+ super(rel, abs);
+
+ if (rel < MIN_RELATIVE_TOLERANCE) {
+ throw new NumberIsTooSmallException(rel, MIN_RELATIVE_TOLERANCE, true);
+ }
+ if (abs <= 0) {
+ throw new NotStrictlyPositiveException(abs);
+ }
+ }
- /** {@inheritDoc} */
- protected double doOptimize()
- throws MaxIterationsExceededException, FunctionEvaluationException {
- return localMin(getGoalType() == GoalType.MINIMIZE,
- getMin(), getStartValue(), getMax(),
- getRelativeAccuracy(), getAbsoluteAccuracy());
+ /**
+ * Convergence criterion.
+ *
+ * @param iteration Current iteration.
+ * @param points Points used for checking the stopping criterion. The list
+ * must contain 3 points (in the following order):
+ * <ul>
+ * <li>the lower end of the current interval</li>
+ * <li>the current best point</li>
+ * <li>the higher end of the current interval</li>
+ * </ul>
+ * @return {@code true} if the stopping criterion is satisfied.
+ * @throws DimensionMismatchException if the length of the {@code points}
+ * list is not equal to 3.
+ */
+ public boolean converged(final int iteration,
+ final UnivariateRealPointValuePair ... points) {
+ if (points.length != 3) {
+ throw new DimensionMismatchException(points.length, 3);
+ }
+
+ final double a = points[0].getPoint();
+ final double x = points[1].getPoint();
+ final double b = points[2].getPoint();
+
+ final double tol1 = getRelativeThreshold() * FastMath.abs(x) + getAbsoluteThreshold();
+ final double tol2 = 2 * tol1;
+
+ final double m = 0.5 * (a + b);
+ return FastMath.abs(x - m) <= tol2 - 0.5 * (b - a);
+ }
}
/**
- * Find the minimum of the function within the interval {@code (lo, hi)}.
+ * Set the convergence checker.
+ * Since this algorithm requires a specific checker, this method will throw
+ * an {@code UnsupportedOperationexception} if the argument type is not
+ * {@link BrentConvergenceChecker}.
*
- * If the function is defined on the interval {@code (lo, hi)}, then
- * this method finds an approximation {@code x} to the point at which
- * the function attains its minimum.<br/>
- * {@code t} and {@code eps} define a tolerance {@code tol = eps |x| + t}
- * and the function is never evaluated at two points closer together than
- * {@code tol}. {@code eps} should be no smaller than <em>2 macheps</em> and
- * preferable not much less than <em>sqrt(macheps)</em>, where
- * <em>macheps</em> is the relative machine precision. {@code t} should be
- * positive.
- * @param isMinim {@code true} when minimizing the function.
- * @param lo Lower bound of the interval.
- * @param mid Point inside the interval {@code [lo, hi]}.
- * @param hi Higher bound of the interval.
- * @param eps Relative accuracy.
- * @param t Absolute accuracy.
- * @return the optimum point.
- * @throws MaxIterationsExceededException if the maximum iteration count
- * is exceeded.
- * @throws FunctionEvaluationException if an error occurs evaluating
- * the function.
+ * @throws MathUnsupportedOperationexception if the checker is not an
+ * instance of {@link BrentConvergenceChecker}.
*/
- private double localMin(boolean isMinim,
- double lo, double mid, double hi,
- double eps, double t)
- throws MaxIterationsExceededException, FunctionEvaluationException {
- if (eps <= 0) {
- throw new NotStrictlyPositiveException(eps);
- }
- if (t <= 0) {
- throw new NotStrictlyPositiveException(t);
+ @Override
+ public void setConvergenceChecker(ConvergenceChecker<UnivariateRealPointValuePair> checker) {
+ if (checker instanceof BrentConvergenceChecker) {
+ super.setConvergenceChecker(checker);
+ } else {
+ throw new MathUnsupportedOperationException();
}
+ }
+
+ /** {@inheritDoc} */
+ protected UnivariateRealPointValuePair doOptimize()
+ throws FunctionEvaluationException {
+ final boolean isMinim = (getGoalType() == GoalType.MINIMIZE);
+ final double lo = getMin();
+ final double mid = getStartValue();
+ final double hi = getMax();
+
+ final ConvergenceChecker<UnivariateRealPointValuePair> checker
+ = getConvergenceChecker();
+ final double eps = checker.getRelativeThreshold();
+ final double t = checker.getAbsoluteThreshold();
+
double a;
double b;
if (lo < hi) {
@@ -111,13 +171,19 @@ public class BrentOptimizer extends Abst
double fv = fx;
double fw = fx;
+ int iter = 0;
while (true) {
double m = 0.5 * (a + b);
final double tol1 = eps * FastMath.abs(x) + t;
final double tol2 = 2 * tol1;
// Check stopping criterion.
- if (FastMath.abs(x - m) > tol2 - 0.5 * (b - a)) {
+ // This test will work only if the "checker" is an instance of
+ // "BrentOptimizer.BrentConvergenceChecker".
+ if (!getConvergenceChecker().converged(iter,
+ new UnivariateRealPointValuePair(a, Double.NaN),
+ new UnivariateRealPointValuePair(x, Double.NaN),
+ new UnivariateRealPointValuePair(b, Double.NaN))) {
double p = 0;
double q = 0;
double r = 0;
@@ -217,11 +283,10 @@ public class BrentOptimizer extends Abst
fv = fu;
}
}
- } else { // termination
- setFunctionValue(isMinim ? fx : -fx);
- return x;
+ } else { // Termination.
+ return new UnivariateRealPointValuePair(x, (isMinim ? fx : -fx));
}
- incrementIterationsCounter();
+ ++iter;
}
}
}