diff options
Diffstat (limited to 'src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/PowellOptimizer.java')
-rw-r--r-- | src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/PowellOptimizer.java | 299 |
1 files changed, 299 insertions, 0 deletions
diff --git a/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/PowellOptimizer.java b/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/PowellOptimizer.java new file mode 100644 index 0000000..afa8426 --- /dev/null +++ b/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/PowellOptimizer.java @@ -0,0 +1,299 @@ +/* + * Licensed to the Apache Software Foundation (ASF) under one or more + * contributor license agreements. See the NOTICE file distributed with + * this work for additional information regarding copyright ownership. + * The ASF licenses this file to You under the Apache License, Version 2.0 + * (the "License"); you may not use this file except in compliance with + * the License. You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + */ +package org.apache.commons.math3.optim.nonlinear.scalar.noderiv; + +import org.apache.commons.math3.util.FastMath; +import org.apache.commons.math3.util.MathArrays; +import org.apache.commons.math3.exception.NumberIsTooSmallException; +import org.apache.commons.math3.exception.NotStrictlyPositiveException; +import org.apache.commons.math3.exception.MathUnsupportedOperationException; +import org.apache.commons.math3.exception.util.LocalizedFormats; +import org.apache.commons.math3.optim.nonlinear.scalar.GoalType; +import org.apache.commons.math3.optim.PointValuePair; +import org.apache.commons.math3.optim.ConvergenceChecker; +import org.apache.commons.math3.optim.nonlinear.scalar.MultivariateOptimizer; +import org.apache.commons.math3.optim.nonlinear.scalar.LineSearch; +import org.apache.commons.math3.optim.univariate.UnivariatePointValuePair; + +/** + * Powell's 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 default stopping criterion is based on the differences of the + * function value between two successive iterations. It is however possible + * to define a custom convergence checker that might terminate the algorithm + * earlier. + * <br/> + * Line search is performed by the {@link LineSearch} class. + * <br/> + * Constraints are not supported: the call to + * {@link #optimize(OptimizationData[]) optimize} will throw + * {@link MathUnsupportedOperationException} if bounds are passed to it. + * In order to impose simple constraints, the objective function must be + * wrapped in an adapter like + * {@link org.apache.commons.math3.optim.nonlinear.scalar.MultivariateFunctionMappingAdapter + * MultivariateFunctionMappingAdapter} or + * {@link org.apache.commons.math3.optim.nonlinear.scalar.MultivariateFunctionPenaltyAdapter + * MultivariateFunctionPenaltyAdapter}. + * + * @since 2.2 + */ +public class PowellOptimizer + extends MultivariateOptimizer { + /** + * Minimum relative tolerance. + */ + private static final double MIN_RELATIVE_TOLERANCE = 2 * FastMath.ulp(1d); + /** + * Relative threshold. + */ + private final double relativeThreshold; + /** + * Absolute threshold. + */ + private final double absoluteThreshold; + /** + * Line search. + */ + private final LineSearch line; + + /** + * This constructor allows to specify a user-defined convergence checker, + * in addition to the parameters that control the default convergence + * checking procedure. + * <br/> + * The internal line search tolerances are set to the square-root of their + * corresponding value in the multivariate optimizer. + * + * @param rel Relative threshold. + * @param abs Absolute threshold. + * @param checker Convergence checker. + * @throws NotStrictlyPositiveException if {@code abs <= 0}. + * @throws NumberIsTooSmallException if {@code rel < 2 * Math.ulp(1d)}. + */ + public PowellOptimizer(double rel, + double abs, + ConvergenceChecker<PointValuePair> checker) { + this(rel, abs, FastMath.sqrt(rel), FastMath.sqrt(abs), checker); + } + + /** + * This constructor allows to specify a user-defined convergence checker, + * in addition to the parameters that control the default convergence + * checking procedure and the line search tolerances. + * + * @param rel Relative threshold for this optimizer. + * @param abs Absolute threshold for this optimizer. + * @param lineRel Relative threshold for the internal line search optimizer. + * @param lineAbs Absolute threshold for the internal line search optimizer. + * @param checker Convergence checker. + * @throws NotStrictlyPositiveException if {@code abs <= 0}. + * @throws NumberIsTooSmallException if {@code rel < 2 * Math.ulp(1d)}. + */ + public PowellOptimizer(double rel, + double abs, + double lineRel, + double lineAbs, + ConvergenceChecker<PointValuePair> checker) { + super(checker); + + if (rel < MIN_RELATIVE_TOLERANCE) { + throw new NumberIsTooSmallException(rel, MIN_RELATIVE_TOLERANCE, true); + } + if (abs <= 0) { + throw new NotStrictlyPositiveException(abs); + } + relativeThreshold = rel; + absoluteThreshold = abs; + + // Create the line search optimizer. + line = new LineSearch(this, + lineRel, + lineAbs, + 1d); + } + + /** + * The parameters control the default convergence checking procedure. + * <br/> + * The internal line search tolerances are set to the square-root of their + * corresponding value in the multivariate optimizer. + * + * @param rel Relative threshold. + * @param abs Absolute threshold. + * @throws NotStrictlyPositiveException if {@code abs <= 0}. + * @throws NumberIsTooSmallException if {@code rel < 2 * Math.ulp(1d)}. + */ + public PowellOptimizer(double rel, + double abs) { + this(rel, abs, null); + } + + /** + * Builds an instance with the default convergence checking procedure. + * + * @param rel Relative threshold. + * @param abs Absolute threshold. + * @param lineRel Relative threshold for the internal line search optimizer. + * @param lineAbs Absolute threshold for the internal line search optimizer. + * @throws NotStrictlyPositiveException if {@code abs <= 0}. + * @throws NumberIsTooSmallException if {@code rel < 2 * Math.ulp(1d)}. + */ + public PowellOptimizer(double rel, + double abs, + double lineRel, + double lineAbs) { + this(rel, abs, lineRel, lineAbs, null); + } + + /** {@inheritDoc} */ + @Override + protected PointValuePair doOptimize() { + checkParameters(); + + final GoalType goal = getGoalType(); + final double[] guess = getStartPoint(); + final int n = guess.length; + + final double[][] direc = new double[n][n]; + for (int i = 0; i < n; i++) { + direc[i][i] = 1; + } + + final ConvergenceChecker<PointValuePair> checker + = getConvergenceChecker(); + + double[] x = guess; + double fVal = computeObjectiveValue(x); + double[] x1 = x.clone(); + while (true) { + incrementIterationCount(); + + double fX = fVal; + double fX2 = 0; + double delta = 0; + int bigInd = 0; + double alphaMin = 0; + + for (int i = 0; i < n; i++) { + final double[] d = MathArrays.copyOf(direc[i]); + + fX2 = fVal; + + final UnivariatePointValuePair optimum = line.search(x, d); + fVal = optimum.getValue(); + alphaMin = optimum.getPoint(); + final double[][] result = newPointAndDirection(x, d, alphaMin); + x = result[0]; + + if ((fX2 - fVal) > delta) { + delta = fX2 - fVal; + bigInd = i; + } + } + + // Default convergence check. + boolean stop = 2 * (fX - fVal) <= + (relativeThreshold * (FastMath.abs(fX) + FastMath.abs(fVal)) + + absoluteThreshold); + + final PointValuePair previous = new PointValuePair(x1, fX); + final PointValuePair current = new PointValuePair(x, fVal); + if (!stop && checker != null) { // User-defined stopping criteria. + stop = checker.converged(getIterations(), previous, current); + } + if (stop) { + if (goal == GoalType.MINIMIZE) { + return (fVal < fX) ? current : previous; + } else { + return (fVal > fX) ? current : previous; + } + } + + final double[] d = new double[n]; + final double[] x2 = new double[n]; + for (int i = 0; i < n; i++) { + d[i] = x[i] - x1[i]; + x2[i] = 2 * x[i] - x1[i]; + } + + x1 = x.clone(); + fX2 = computeObjectiveValue(x2); + + if (fX > fX2) { + double t = 2 * (fX + fX2 - 2 * fVal); + double temp = fX - fVal - delta; + t *= temp * temp; + temp = fX - fX2; + t -= delta * temp * temp; + + if (t < 0.0) { + final UnivariatePointValuePair optimum = line.search(x, d); + fVal = optimum.getValue(); + alphaMin = optimum.getPoint(); + final double[][] result = newPointAndDirection(x, d, alphaMin); + x = result[0]; + + final int lastInd = n - 1; + direc[bigInd] = direc[lastInd]; + direc[lastInd] = result[1]; + } + } + } + } + + /** + * Compute a new point (in the original space) and a new direction + * vector, resulting from the line search. + * + * @param p Point used in the line search. + * @param d Direction used in the line search. + * @param optimum Optimum found by the line search. + * @return a 2-element array containing the new point (at index 0) and + * the new direction (at index 1). + */ + private double[][] newPointAndDirection(double[] p, + double[] d, + double optimum) { + final int n = p.length; + final double[] nP = new double[n]; + final double[] nD = new double[n]; + for (int i = 0; i < n; i++) { + nD[i] = d[i] * optimum; + nP[i] = p[i] + nD[i]; + } + + final double[][] result = new double[2][]; + result[0] = nP; + result[1] = nD; + + return result; + } + + /** + * @throws MathUnsupportedOperationException if bounds were passed to the + * {@link #optimize(OptimizationData[]) optimize} method. + */ + private void checkParameters() { + if (getLowerBound() != null || + getUpperBound() != null) { + throw new MathUnsupportedOperationException(LocalizedFormats.CONSTRAINT); + } + } +} |