directsearchoptimizer.java

Apache的common math数学软件包

/*
 * 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;

import java.util.Arrays;
import java.util.Comparator;

import org.apache.commons.math.ConvergenceException;
import org.apache.commons.math.DimensionMismatchException;
import org.apache.commons.math.linear.RealMatrix;
import org.apache.commons.math.random.CorrelatedRandomVectorGenerator;
import org.apache.commons.math.random.JDKRandomGenerator;
import org.apache.commons.math.random.NotPositiveDefiniteMatrixException;
import org.apache.commons.math.random.RandomGenerator;
import org.apache.commons.math.random.RandomVectorGenerator;
import org.apache.commons.math.random.UncorrelatedRandomVectorGenerator;
import org.apache.commons.math.random.UniformRandomGenerator;
import org.apache.commons.math.stat.descriptive.moment.VectorialCovariance;
import org.apache.commons.math.stat.descriptive.moment.VectorialMean;

/** 
 * This class implements simplex-based direct search optimization
 * algorithms.
 *
 * <p>Direct search methods only use cost function values, they don't
 * need derivatives and don't either try to compute approximation of
 * the derivatives. According to a 1996 paper by Margaret H. Wright
 * (<a href="http://cm.bell-labs.com/cm/cs/doc/96/4-02.ps.gz">Direct
 * Search Methods: Once Scorned, Now Respectable</a>), they are used
 * when either the computation of the derivative is impossible (noisy
 * functions, unpredictable dicontinuities) or difficult (complexity,
 * computation cost). In the first cases, rather than an optimum, a
 * <em>not too bad</em> point is desired. In the latter cases, an
 * optimum is desired but cannot be reasonably found. In all cases
 * direct search methods can be useful.</p>
 *
 * <p>Simplex-based direct search methods are based on comparison of
 * the cost function values at the vertices of a simplex (which is a
 * set of n+1 points in dimension n) that is updated by the algorithms
 * steps.</p>
 *
 * <p>Minimization can be attempted either in single-start or in
 * multi-start mode. Multi-start is a traditional way to try to avoid
 * being trapped in a local minimum and miss the global minimum of a
 * function. It can also be used to verify the convergence of an
 * algorithm. The various multi-start-enabled <code>minimize</code>
 * methods return the best minimum found after all starts, and the
 * {@link #getMinima getMinima} method can be used to retrieve all
 * minima from all starts (including the one already provided by the
 * {@link #minimize(CostFunction, int, ConvergenceChecker, double[],
 * double[]) minimize} method).</p>
 *
 * <p>This class is the base class performing the boilerplate simplex
 * initialization and handling. The simplex update by itself is
 * performed by the derived classes according to the implemented
 * algorithms.</p>
 *
 * @see CostFunction
 * @see NelderMead
 * @see MultiDirectional
 * @version $Revision: 628000 $ $Date: 2008-02-15 03:31:48 -0700 (Fri, 15 Feb 2008) $
 * @since 1.2
 */
public abstract class DirectSearchOptimizer {

    /** Simple constructor.
     */
    protected DirectSearchOptimizer() {
    }

    /** Minimizes a cost function.
     * <p>The initial simplex is built from two vertices that are
     * considered to represent two opposite vertices of a box parallel
     * to the canonical axes of the space. The simplex is the subset of
     * vertices encountered while going from vertexA to vertexB
     * traveling along the box edges only. This can be seen as a scaled
     * regular simplex using the projected separation between the given
     * points as the scaling factor along each coordinate axis.</p>
     * <p>The optimization is performed in single-start mode.</p>
     * @param f cost function
     * @param maxEvaluations maximal number of function calls for each
     * start (note that the number will be checked <em>after</em>
     * complete simplices have been evaluated, this means that in some
     * cases this number will be exceeded by a few units, depending on
     * the dimension of the problem)
     * @param checker object to use to check for convergence
     * @param vertexA first vertex
     * @param vertexB last vertex
     * @return the point/cost pairs giving the minimal cost
     * @exception CostException if the cost function throws one during
     * the search
     * @exception ConvergenceException if none of the starts did
     * converge (it is not thrown if at least one start did converge)
     */
    public PointCostPair minimize(CostFunction f, int maxEvaluations,
                                  ConvergenceChecker checker,
                                  double[] vertexA, double[] vertexB)
    throws CostException, ConvergenceException {

        // set up optimizer
        buildSimplex(vertexA, vertexB);
        setSingleStart();

        // compute minimum
        return minimize(f, maxEvaluations, checker);

    }

    /** Minimizes a cost function.
     * <p>The initial simplex is built from two vertices that are
     * considered to represent two opposite vertices of a box parallel
     * to the canonical axes of the space. The simplex is the subset of
     * vertices encountered while going from vertexA to vertexB
     * traveling along the box edges only. This can be seen as a scaled
     * regular simplex using the projected separation between the given
     * points as the scaling factor along each coordinate axis.</p>
     * <p>The optimization is performed in multi-start mode.</p>
     * @param f cost function
     * @param maxEvaluations maximal number of function calls for each
     * start (note that the number will be checked <em>after</em>
     * complete simplices have been evaluated, this means that in some
     * cases this number will be exceeded by a few units, depending on
     * the dimension of the problem)
     * @param checker object to use to check for convergence
     * @param vertexA first vertex
     * @param vertexB last vertex
     * @param starts number of starts to perform (including the
     * first one), multi-start is disabled if value is less than or
     * equal to 1
     * @param seed seed for the random vector generator
     * @return the point/cost pairs giving the minimal cost
     * @exception CostException if the cost function throws one during
     * the search
     * @exception ConvergenceException if none of the starts did
     * converge (it is not thrown if at least one start did converge)
     */
    public PointCostPair minimize(CostFunction f, int maxEvaluations,
                                  ConvergenceChecker checker,
                                  double[] vertexA, double[] vertexB,
                                  int starts, long seed)
    throws CostException, ConvergenceException {

        // set up the simplex traveling around the box
        buildSimplex(vertexA, vertexB);

        // we consider the simplex could have been produced by a generator
        // having its mean value at the center of the box, the standard
        // deviation along each axe being the corresponding half size
        double[] mean              = new double[vertexA.length];
        double[] standardDeviation = new double[vertexA.length];
        for (int i = 0; i < vertexA.length; ++i) {
            mean[i]              = 0.5 * (vertexA[i] + vertexB[i]);
            standardDeviation[i] = 0.5 * Math.abs(vertexA[i] - vertexB[i]);
        }

        RandomGenerator rg = new JDKRandomGenerator();
        rg.setSeed(seed);
        UniformRandomGenerator urg = new UniformRandomGenerator(rg);
        RandomVectorGenerator rvg =
            new UncorrelatedRandomVectorGenerator(mean, standardDeviation, urg);
        setMultiStart(starts, rvg);

        // compute minimum
        return minimize(f, maxEvaluations, checker);

    }

    /** Minimizes a cost function.
     * <p>The simplex is built from all its vertices.</p>
     * <p>The optimization is performed in single-start mode.</p>
     * @param f cost function
     * @param maxEvaluations maximal number of function calls for each
     * start (note that the number will be checked <em>after</em>
     * complete simplices have been evaluated, this means that in some
     * cases this number will be exceeded by a few units, depending on
     * the dimension of the problem)
     * @param checker object to use to check for convergence
     * @param vertices array containing all vertices of the simplex
     * @return the point/cost pairs giving the minimal cost
     * @exception CostException if the cost function throws one during
     * the search
     * @exception ConvergenceException if none of the starts did
     * converge (it is not thrown if at least one start did converge)
     */
    public PointCostPair minimize(CostFunction f, int maxEvaluations,
                                  ConvergenceChecker checker,
                                  double[][] vertices)
    throws CostException, ConvergenceException {

        // set up optimizer
        buildSimplex(vertices);
        setSingleStart();

        // compute minimum
        return minimize(f, maxEvaluations, checker);

    }

    /** Minimizes a cost function.
     * <p>The simplex is built from all its vertices.</p>
     * <p>The optimization is performed in multi-start mode.</p>
     * @param f cost function
     * @param maxEvaluations maximal number of function calls for each
     * start (note that the number will be checked <em>after</em>
     * complete simplices have been evaluated, this means that in some
     * cases this number will be exceeded by a few units, depending on
     * the dimension of the problem)
     * @param checker object to use to check for convergence
     * @param vertices array containing all vertices of the simplex
     * @param starts number of starts to perform (including the
     * first one), multi-start is disabled if value is less than or
     * equal to 1
     * @param seed seed for the random vector generator
     * @return the point/cost pairs giving the minimal cost
     * @exception NotPositiveDefiniteMatrixException if the vertices
     * array is degenerated
     * @exception CostException if the cost function throws one during
     * the search
     * @exception ConvergenceException if none of the starts did
     * converge (it is not thrown if at least one start did converge)
     */
    public PointCostPair minimize(CostFunction f, int maxEvaluations,
                                  ConvergenceChecker checker,
                                  double[][] vertices,
                                  int starts, long seed)
    throws NotPositiveDefiniteMatrixException,
    CostException, ConvergenceException {

        try {
            // store the points into the simplex
            buildSimplex(vertices);

            // compute the statistical properties of the simplex points
            VectorialMean meanStat = new VectorialMean(vertices[0].length);
            VectorialCovariance covStat = new VectorialCovariance(vertices[0].length, true);
            for (int i = 0; i < vertices.length; ++i) {
                meanStat.increment(vertices[i]);
                covStat.increment(vertices[i]);
            }
            double[] mean = meanStat.getResult();
            RealMatrix covariance = covStat.getResult();
            

            RandomGenerator rg = new JDKRandomGenerator();
            rg.setSeed(seed);
            RandomVectorGenerator rvg =
                new CorrelatedRandomVectorGenerator(mean,
                                                    covariance, 1.0e-12 * covariance.getNorm(),
                                                    new UniformRandomGenerator(rg));
            setMultiStart(starts, rvg);

            // compute minimum
            return minimize(f, maxEvaluations, checker);

        } catch (DimensionMismatchException dme) {
            // this should not happen
            throw new RuntimeException("internal error");
        }

    }

    /** Minimizes a cost function.
     * <p>The simplex is built randomly.</p>
     * <p>The optimization is performed in single-start mode.</p>
     * @param f cost function
     * @param maxEvaluations maximal number of function calls for each
     * start (note that the number will be checked <em>after</em>
     * complete simplices have been evaluated, this means that in some
     * cases this number will be exceeded by a few units, depending on
     * the dimension of the problem)
     * @param checker object to use to check for convergence
     * @param generator random vector generator
     * @return the point/cost pairs giving the minimal cost
     * @exception CostException if the cost function throws one during
     * the search
     * @exception ConvergenceException if none of the starts did
     * converge (it is not thrown if at least one start did converge)
     */
    public PointCostPair minimize(CostFunction f, int maxEvaluations,
                                  ConvergenceChecker checker,
                                  RandomVectorGenerator generator)
    throws CostException, ConvergenceException {

        // set up optimizer
        buildSimplex(generator);
        setSingleStart();

        // compute minimum