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 2012/12/12 15:11:04 UTC
svn commit: r1420684 [7/15] - in /commons/proper/math/trunk/src:
main/java/org/apache/commons/math3/exception/
main/java/org/apache/commons/math3/exception/util/
main/java/org/apache/commons/math3/fitting/
main/java/org/apache/commons/math3/optim/ main...
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/PowellOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/PowellOptimizer.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/PowellOptimizer.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/PowellOptimizer.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,356 @@
+/*
+ * 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.analysis.UnivariateFunction;
+import org.apache.commons.math3.exception.NumberIsTooSmallException;
+import org.apache.commons.math3.exception.NotStrictlyPositiveException;
+import org.apache.commons.math3.optim.MaxEval;
+import org.apache.commons.math3.optim.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.univariate.BracketFinder;
+import org.apache.commons.math3.optim.univariate.BrentOptimizer;
+import org.apache.commons.math3.optim.univariate.UnivariatePointValuePair;
+import org.apache.commons.math3.optim.univariate.SimpleUnivariateValueChecker;
+import org.apache.commons.math3.optim.univariate.SearchInterval;
+import org.apache.commons.math3.optim.univariate.UnivariateObjectiveFunction;
+
+/**
+ * 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 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/>
+ * The internal line search optimizer is a {@link BrentOptimizer} with a
+ * convergence checker set to {@link SimpleUnivariateValueChecker}.
+ *
+ * @version $Id: PowellOptimizer.java 1413594 2012-11-26 13:16:39Z erans $
+ * @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(lineRel,
+ lineAbs);
+ }
+
+ /**
+ * 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() {
+ 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();
+ int iter = 0;
+ while (true) {
+ ++iter;
+
+ 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) { // User-defined stopping criteria.
+ if (checker != null) {
+ stop = checker.converged(iter, 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;
+ }
+
+ /**
+ * Class for finding the minimum of the objective function along a given
+ * direction.
+ */
+ private class LineSearch extends BrentOptimizer {
+ /**
+ * Value that will pass the precondition check for {@link BrentOptimizer}
+ * but will not pass the convergence check, so that the custom checker
+ * will always decide when to stop the line search.
+ */
+ private static final double REL_TOL_UNUSED = 1e-15;
+ /**
+ * Value that will pass the precondition check for {@link BrentOptimizer}
+ * but will not pass the convergence check, so that the custom checker
+ * will always decide when to stop the line search.
+ */
+ private static final double ABS_TOL_UNUSED = Double.MIN_VALUE;
+ /**
+ * Automatic bracketing.
+ */
+ private final BracketFinder bracket = new BracketFinder();
+
+ /**
+ * The "BrentOptimizer" default stopping criterion uses the tolerances
+ * to check the domain (point) values, not the function values.
+ * We thus create a custom checker to use function values.
+ *
+ * @param rel Relative threshold.
+ * @param abs Absolute threshold.
+ */
+ LineSearch(double rel,
+ double abs) {
+ super(REL_TOL_UNUSED,
+ ABS_TOL_UNUSED,
+ new SimpleUnivariateValueChecker(rel, abs));
+ }
+
+ /**
+ * Find the minimum of the function {@code f(p + alpha * d)}.
+ *
+ * @param p Starting point.
+ * @param d Search direction.
+ * @return the optimum.
+ * @throws org.apache.commons.math3.exception.TooManyEvaluationsException
+ * if the number of evaluations is exceeded.
+ */
+ public UnivariatePointValuePair search(final double[] p, final double[] d) {
+ final int n = p.length;
+ final UnivariateFunction f = new UnivariateFunction() {
+ public double value(double alpha) {
+ 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);
+ // Passing "MAX_VALUE" as a dummy value because it is the enclosing
+ // class that counts the number of evaluations (and will eventually
+ // generate the exception).
+ return optimize(new MaxEval(Integer.MAX_VALUE),
+ new UnivariateObjectiveFunction(f),
+ goal,
+ new SearchInterval(bracket.getLo(),
+ bracket.getHi(),
+ bracket.getMid()));
+ }
+ }
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/PowellOptimizer.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/SimplexOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/SimplexOptimizer.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/SimplexOptimizer.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/SimplexOptimizer.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,201 @@
+/*
+ * 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 java.util.Comparator;
+import org.apache.commons.math3.analysis.MultivariateFunction;
+import org.apache.commons.math3.exception.NullArgumentException;
+import org.apache.commons.math3.optim.GoalType;
+import org.apache.commons.math3.optim.ConvergenceChecker;
+import org.apache.commons.math3.optim.PointValuePair;
+import org.apache.commons.math3.optim.SimpleValueChecker;
+import org.apache.commons.math3.optim.OptimizationData;
+import org.apache.commons.math3.optim.nonlinear.scalar.MultivariateOptimizer;
+
+/**
+ * This class implements simplex-based direct search optimization.
+ *
+ * <p>
+ * Direct search methods only use objective function values, they do
+ * not 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 discontinuities) 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 objective 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>
+ * The simplex update procedure ({@link NelderMeadSimplex} or
+ * {@link MultiDirectionalSimplex}) must be passed to the
+ * {@code optimize} method.
+ * </p>
+ * <p>
+ * Each call to {@code optimize} will re-use the start configuration of
+ * the current simplex and move it such that its first vertex is at the
+ * provided start point of the optimization.
+ * If the {@code optimize} method is called to solve a different problem
+ * and the number of parameters change, the simplex must be re-initialized
+ * to one with the appropriate dimensions.
+ * </p>
+ * <p>
+ * Convergence is checked by providing the <em>worst</em> points of
+ * previous and current simplex to the convergence checker, not the best
+ * ones.
+ * </p>
+ * <p>
+ * This simplex optimizer implementation does not directly support constrained
+ * optimization with simple bounds; so, for such optimizations, either a more
+ * dedicated algorithm must be used like
+ * {@link CMAESOptimizer} or {@link BOBYQAOptimizer}, or 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}.
+ * </p>
+ *
+ * @version $Id: SimplexOptimizer.java 1397759 2012-10-13 01:12:58Z erans $
+ * @since 3.0
+ */
+public class SimplexOptimizer extends MultivariateOptimizer {
+ /** Simplex update rule. */
+ private AbstractSimplex simplex;
+
+ /**
+ * @param checker Convergence checker.
+ */
+ public SimplexOptimizer(ConvergenceChecker<PointValuePair> checker) {
+ super(checker);
+ }
+
+ /**
+ * @param rel Relative threshold.
+ * @param abs Absolute threshold.
+ */
+ public SimplexOptimizer(double rel, double abs) {
+ this(new SimpleValueChecker(rel, abs));
+ }
+
+ /**
+ * {@inheritDoc}
+ *
+ * @param optData Optimization data.
+ * The following data will be looked for:
+ * <ul>
+ * <li>{@link org.apache.commons.math3.optim.MaxEval}</li>
+ * <li>{@link org.apache.commons.math3.optim.InitialGuess}</li>
+ * <li>{@link org.apache.commons.math3.optim.SimpleBounds}</li>
+ * <li>{@link AbstractSimplex}</li>
+ * </ul>
+ * @return {@inheritDoc}
+ */
+ @Override
+ public PointValuePair optimize(OptimizationData... optData) {
+ // Retrieve settings
+ parseOptimizationData(optData);
+ // Set up base class and perform computation.
+ return super.optimize(optData);
+ }
+
+ /** {@inheritDoc} */
+ @Override
+ protected PointValuePair doOptimize() {
+ if (simplex == null) {
+ throw new NullArgumentException();
+ }
+
+ // Indirect call to "computeObjectiveValue" in order to update the
+ // evaluations counter.
+ final MultivariateFunction evalFunc
+ = new MultivariateFunction() {
+ public double value(double[] point) {
+ return computeObjectiveValue(point);
+ }
+ };
+
+ final boolean isMinim = getGoalType() == GoalType.MINIMIZE;
+ final Comparator<PointValuePair> comparator
+ = new Comparator<PointValuePair>() {
+ public int compare(final PointValuePair o1,
+ final PointValuePair o2) {
+ final double v1 = o1.getValue();
+ final double v2 = o2.getValue();
+ return isMinim ? Double.compare(v1, v2) : Double.compare(v2, v1);
+ }
+ };
+
+ // Initialize search.
+ simplex.build(getStartPoint());
+ simplex.evaluate(evalFunc, comparator);
+
+ PointValuePair[] previous = null;
+ int iteration = 0;
+ final ConvergenceChecker<PointValuePair> checker = getConvergenceChecker();
+ while (true) {
+ if (iteration > 0) {
+ boolean converged = true;
+ for (int i = 0; i < simplex.getSize(); i++) {
+ PointValuePair prev = previous[i];
+ converged = converged &&
+ checker.converged(iteration, prev, simplex.getPoint(i));
+ }
+ if (converged) {
+ // We have found an optimum.
+ return simplex.getPoint(0);
+ }
+ }
+
+ // We still need to search.
+ previous = simplex.getPoints();
+ simplex.iterate(evalFunc, comparator);
+ ++iteration;
+ }
+ }
+
+ /**
+ * Scans the list of (required and optional) optimization data that
+ * characterize the problem.
+ *
+ * @param optData Optimization data.
+ * The following data will be looked for:
+ * <ul>
+ * <li>{@link AbstractSimplex}</li>
+ * </ul>
+ */
+ private void parseOptimizationData(OptimizationData... optData) {
+ // The existing values (as set by the previous call) are reused if
+ // not provided in the argument list.
+ for (OptimizationData data : optData) {
+ if (data instanceof AbstractSimplex) {
+ simplex = (AbstractSimplex) data;
+ // If more data must be parsed, this statement _must_ be
+ // changed to "continue".
+ break;
+ }
+ }
+ }
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/SimplexOptimizer.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/package-info.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/package-info.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/package-info.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/package-info.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,21 @@
+/*
+ * 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;
+
+/**
+ * This package provides optimization algorithms that do not require derivatives.
+ */
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/package-info.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/package-info.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/package-info.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/package-info.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/package-info.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,21 @@
+/*
+ * 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;
+
+/**
+ * Algorithms for optimizing a scalar function.
+ */
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/package-info.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/JacobianMultivariateVectorOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/JacobianMultivariateVectorOptimizer.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/JacobianMultivariateVectorOptimizer.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/JacobianMultivariateVectorOptimizer.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,114 @@
+/*
+ * 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.vector;
+
+import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
+import org.apache.commons.math3.optim.ConvergenceChecker;
+import org.apache.commons.math3.optim.OptimizationData;
+import org.apache.commons.math3.optim.PointVectorValuePair;
+import org.apache.commons.math3.exception.TooManyEvaluationsException;
+import org.apache.commons.math3.exception.DimensionMismatchException;
+
+/**
+ * Base class for implementing optimizers for multivariate vector
+ * differentiable functions.
+ * It contains boiler-plate code for dealing with Jacobian evaluation.
+ * It assumes that the rows of the Jacobian matrix iterate on the model
+ * functions while the columns iterate on the parameters; thus, the numbers
+ * of rows is equal to the dimension of the {@link Target} while the
+ * number of columns is equal to the dimension of the
+ * {@link org.apache.commons.math3.optim.InitialGuess InitialGuess}.
+ *
+ * @version $Id$
+ * @since 3.1
+ */
+public abstract class JacobianMultivariateVectorOptimizer
+ extends MultivariateVectorOptimizer {
+ /**
+ * Jacobian of the model function.
+ */
+ private MultivariateMatrixFunction jacobian;
+
+ /**
+ * @param checker Convergence checker.
+ */
+ protected JacobianMultivariateVectorOptimizer(ConvergenceChecker<PointVectorValuePair> checker) {
+ super(checker);
+ }
+
+ /**
+ * Computes the Jacobian matrix.
+ *
+ * @param params Point at which the Jacobian must be evaluated.
+ * @return the Jacobian at the specified point.
+ */
+ protected double[][] computeJacobian(final double[] params) {
+ return jacobian.value(params);
+ }
+
+ /**
+ * {@inheritDoc}
+ *
+ * @param optData Optimization data. The following data will be looked for:
+ * <ul>
+ * <li>{@link org.apache.commons.math3.optim.MaxEval}</li>
+ * <li>{@link org.apache.commons.math3.optim.InitialGuess}</li>
+ * <li>{@link org.apache.commons.math3.optim.SimpleBounds}</li>
+ * <li>{@link Target}</li>
+ * <li>{@link Weight}</li>
+ * <li>{@link ModelFunction}</li>
+ * <li>{@link ModelFunctionJacobian}</li>
+ * </ul>
+ * @return {@inheritDoc}
+ * @throws TooManyEvaluationsException if the maximal number of
+ * evaluations is exceeded.
+ * @throws DimensionMismatchException if the initial guess, target, and weight
+ * arguments have inconsistent dimensions.
+ */
+ @Override
+ public PointVectorValuePair optimize(OptimizationData... optData)
+ throws TooManyEvaluationsException,
+ DimensionMismatchException {
+ // Retrieve settings.
+ parseOptimizationData(optData);
+ // Set up base class and perform computation.
+ return super.optimize(optData);
+ }
+
+ /**
+ * Scans the list of (required and optional) optimization data that
+ * characterize the problem.
+ *
+ * @param optData Optimization data.
+ * The following data will be looked for:
+ * <ul>
+ * <li>{@link ModelFunctionJacobian}</li>
+ * </ul>
+ */
+ private void parseOptimizationData(OptimizationData... optData) {
+ // The existing values (as set by the previous call) are reused if
+ // not provided in the argument list.
+ for (OptimizationData data : optData) {
+ if (data instanceof ModelFunctionJacobian) {
+ jacobian = ((ModelFunctionJacobian) data).getModelFunctionJacobian();
+ // If more data must be parsed, this statement _must_ be
+ // changed to "continue".
+ break;
+ }
+ }
+ }
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/JacobianMultivariateVectorOptimizer.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/ModelFunction.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/ModelFunction.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/ModelFunction.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/ModelFunction.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,47 @@
+/*
+ * 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.vector;
+
+import org.apache.commons.math3.analysis.MultivariateVectorFunction;
+import org.apache.commons.math3.optim.OptimizationData;
+
+/**
+ * Model (vector) function to be optimized.
+ *
+ * @version $Id$
+ * @since 3.1
+ */
+public class ModelFunction implements OptimizationData {
+ /** Function to be optimized. */
+ private final MultivariateVectorFunction model;
+
+ /**
+ * @param m Model function to be optimized.
+ */
+ public ModelFunction(MultivariateVectorFunction m) {
+ model = m;
+ }
+
+ /**
+ * Gets the model function to be optimized.
+ *
+ * @return the model function.
+ */
+ public MultivariateVectorFunction getModelFunction() {
+ return model;
+ }
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/ModelFunction.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/ModelFunctionJacobian.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/ModelFunctionJacobian.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/ModelFunctionJacobian.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/ModelFunctionJacobian.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,47 @@
+/*
+ * 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.vector;
+
+import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
+import org.apache.commons.math3.optim.OptimizationData;
+
+/**
+ * Jacobian of the model (vector) function to be optimized.
+ *
+ * @version $Id$
+ * @since 3.1
+ */
+public class ModelFunctionJacobian implements OptimizationData {
+ /** Function to be optimized. */
+ private final MultivariateMatrixFunction jacobian;
+
+ /**
+ * @param j Jacobian of the model function to be optimized.
+ */
+ public ModelFunctionJacobian(MultivariateMatrixFunction j) {
+ jacobian = j;
+ }
+
+ /**
+ * Gets the Jacobian of the model function to be optimized.
+ *
+ * @return the model function Jacobian.
+ */
+ public MultivariateMatrixFunction getModelFunctionJacobian() {
+ return jacobian;
+ }
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/ModelFunctionJacobian.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/MultiStartMultivariateVectorOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/MultiStartMultivariateVectorOptimizer.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/MultiStartMultivariateVectorOptimizer.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/MultiStartMultivariateVectorOptimizer.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,121 @@
+/*
+ * 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.vector;
+
+import java.util.Collections;
+import java.util.List;
+import java.util.ArrayList;
+import java.util.Comparator;
+import org.apache.commons.math3.exception.NotStrictlyPositiveException;
+import org.apache.commons.math3.exception.NullArgumentException;
+import org.apache.commons.math3.linear.RealMatrix;
+import org.apache.commons.math3.linear.RealVector;
+import org.apache.commons.math3.linear.ArrayRealVector;
+import org.apache.commons.math3.random.RandomVectorGenerator;
+import org.apache.commons.math3.optim.BaseMultiStartMultivariateOptimizer;
+import org.apache.commons.math3.optim.PointVectorValuePair;
+
+/**
+ * Multi-start optimizer for a (vector) model function.
+ *
+ * This class wraps an optimizer in order to use it several times in
+ * turn with different starting points (trying to avoid being trapped
+ * in a local extremum when looking for a global one).
+ *
+ * @version $Id$
+ * @since 3.0
+ */
+public class MultiStartMultivariateVectorOptimizer
+ extends BaseMultiStartMultivariateOptimizer<PointVectorValuePair> {
+ /** Underlying optimizer. */
+ private final MultivariateVectorOptimizer optimizer;
+ /** Found optima. */
+ private final List<PointVectorValuePair> optima = new ArrayList<PointVectorValuePair>();
+
+ /**
+ * Create a multi-start optimizer from a single-start optimizer.
+ *
+ * @param optimizer Single-start optimizer to wrap.
+ * @param starts Number of starts to perform.
+ * If {@code starts == 1}, the result will be same as if {@code optimizer}
+ * is called directly.
+ * @param generator Random vector generator to use for restarts.
+ * @throws NullArgumentException if {@code optimizer} or {@code generator}
+ * is {@code null}.
+ * @throws NotStrictlyPositiveException if {@code starts < 1}.
+ */
+ public MultiStartMultivariateVectorOptimizer(final MultivariateVectorOptimizer optimizer,
+ final int starts,
+ final RandomVectorGenerator generator)
+ throws NullArgumentException,
+ NotStrictlyPositiveException {
+ super(optimizer, starts, generator);
+ this.optimizer = optimizer;
+ }
+
+ /**
+ * {@inheritDoc}
+ */
+ @Override
+ public PointVectorValuePair[] getOptima() {
+ Collections.sort(optima, getPairComparator());
+ return optima.toArray(new PointVectorValuePair[0]);
+ }
+
+ /**
+ * {@inheritDoc}
+ */
+ @Override
+ protected void store(PointVectorValuePair optimum) {
+ optima.add(optimum);
+ }
+
+ /**
+ * {@inheritDoc}
+ */
+ @Override
+ protected void clear() {
+ optima.clear();
+ }
+
+ /**
+ * @return a comparator for sorting the optima.
+ */
+ private Comparator<PointVectorValuePair> getPairComparator() {
+ return new Comparator<PointVectorValuePair>() {
+ private final RealVector target = new ArrayRealVector(optimizer.getTarget(), false);
+ private final RealMatrix weight = optimizer.getWeight();
+
+ public int compare(final PointVectorValuePair o1,
+ final PointVectorValuePair o2) {
+ if (o1 == null) {
+ return (o2 == null) ? 0 : 1;
+ } else if (o2 == null) {
+ return -1;
+ }
+ return Double.compare(weightedResidual(o1),
+ weightedResidual(o2));
+ }
+
+ private double weightedResidual(final PointVectorValuePair pv) {
+ final RealVector v = new ArrayRealVector(pv.getValueRef(), false);
+ final RealVector r = target.subtract(v);
+ return r.dotProduct(weight.operate(r));
+ }
+ };
+ }
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/MultiStartMultivariateVectorOptimizer.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/MultivariateVectorOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/MultivariateVectorOptimizer.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/MultivariateVectorOptimizer.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/MultivariateVectorOptimizer.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,164 @@
+/*
+ * 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.vector;
+
+import org.apache.commons.math3.exception.TooManyEvaluationsException;
+import org.apache.commons.math3.exception.DimensionMismatchException;
+import org.apache.commons.math3.analysis.MultivariateVectorFunction;
+import org.apache.commons.math3.optim.OptimizationData;
+import org.apache.commons.math3.optim.BaseMultivariateOptimizer;
+import org.apache.commons.math3.optim.ConvergenceChecker;
+import org.apache.commons.math3.optim.PointVectorValuePair;
+import org.apache.commons.math3.linear.RealMatrix;
+
+/**
+ * Base class for a multivariate vector function optimizer.
+ *
+ * @version $Id$
+ * @since 3.1
+ */
+public abstract class MultivariateVectorOptimizer
+ extends BaseMultivariateOptimizer<PointVectorValuePair> {
+ /** Target values for the model function at optimum. */
+ private double[] target;
+ /** Weight matrix. */
+ private RealMatrix weightMatrix;
+ /** Model function. */
+ private MultivariateVectorFunction model;
+
+ /**
+ * @param checker Convergence checker.
+ */
+ protected MultivariateVectorOptimizer(ConvergenceChecker<PointVectorValuePair> checker) {
+ super(checker);
+ }
+
+ /**
+ * Computes the objective function value.
+ * This method <em>must</em> be called by subclasses to enforce the
+ * evaluation counter limit.
+ *
+ * @param params Point at which the objective function must be evaluated.
+ * @return the objective function value at the specified point.
+ * @throws TooManyEvaluationsException if the maximal number of evaluations
+ * (of the model vector function) is exceeded.
+ */
+ protected double[] computeObjectiveValue(double[] params) {
+ super.incrementEvaluationCount();
+ return model.value(params);
+ }
+
+ /**
+ * {@inheritDoc}
+ *
+ * @param optData Optimization data. The following data will be looked for:
+ * <ul>
+ * <li>{@link org.apache.commons.math3.optim.MaxEval}</li>
+ * <li>{@link org.apache.commons.math3.optim.InitialGuess}</li>
+ * <li>{@link org.apache.commons.math3.optim.SimpleBounds}</li>
+ * <li>{@link Target}</li>
+ * <li>{@link Weight}</li>
+ * <li>{@link ModelFunction}</li>
+ * </ul>
+ * @return {@inheritDoc}
+ * @throws TooManyEvaluationsException if the maximal number of
+ * evaluations is exceeded.
+ * @throws DimensionMismatchException if the initial guess, target, and weight
+ * arguments have inconsistent dimensions.
+ */
+ public PointVectorValuePair optimize(OptimizationData... optData)
+ throws TooManyEvaluationsException,
+ DimensionMismatchException {
+ // Retrieve settings.
+ parseOptimizationData(optData);
+ // Check input consistency.
+ checkParameters();
+ // Set up base class and perform computation.
+ return super.optimize(optData);
+ }
+
+ /**
+ * Gets the weight matrix of the observations.
+ *
+ * @return the weight matrix.
+ */
+ public RealMatrix getWeight() {
+ return weightMatrix.copy();
+ }
+ /**
+ * Gets the observed values to be matched by the objective vector
+ * function.
+ *
+ * @return the target values.
+ */
+ public double[] getTarget() {
+ return target.clone();
+ }
+
+ /**
+ * Gets the number of observed values.
+ *
+ * @return the length of the target vector.
+ */
+ public int getTargetSize() {
+ return target.length;
+ }
+
+ /**
+ * Scans the list of (required and optional) optimization data that
+ * characterize the problem.
+ *
+ * @param optData Optimization data. The following data will be looked for:
+ * <ul>
+ * <li>{@link Target}</li>
+ * <li>{@link Weight}</li>
+ * <li>{@link ModelFunction}</li>
+ * </ul>
+ */
+ private void parseOptimizationData(OptimizationData... optData) {
+ // The existing values (as set by the previous call) are reused if
+ // not provided in the argument list.
+ for (OptimizationData data : optData) {
+ if (data instanceof ModelFunction) {
+ model = ((ModelFunction) data).getModelFunction();
+ continue;
+ }
+ if (data instanceof Target) {
+ target = ((Target) data).getTarget();
+ continue;
+ }
+ if (data instanceof Weight) {
+ weightMatrix = ((Weight) data).getWeight();
+ continue;
+ }
+ }
+ }
+
+ /**
+ * Check parameters consistency.
+ *
+ * @throws DimensionMismatchException if {@link #target} and
+ * {@link #weightMatrix} have inconsistent dimensions.
+ */
+ private void checkParameters() {
+ if (target.length != weightMatrix.getColumnDimension()) {
+ throw new DimensionMismatchException(target.length,
+ weightMatrix.getColumnDimension());
+ }
+ }
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/MultivariateVectorOptimizer.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/Target.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/Target.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/Target.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/Target.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,50 @@
+/*
+ * 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.vector;
+
+import org.apache.commons.math3.optim.OptimizationData;
+
+/**
+ * Target of the optimization procedure.
+ * They are the values which the objective vector function must reproduce
+ * When the parameters of the model have been optimized.
+ * <br/>
+ * Immutable class.
+ *
+ * @version $Id: Target.java 1416643 2012-12-03 19:37:14Z tn $
+ * @since 3.1
+ */
+public class Target implements OptimizationData {
+ /** Target values (of the objective vector function). */
+ private final double[] target;
+
+ /**
+ * @param observations Target values.
+ */
+ public Target(double[] observations) {
+ target = observations.clone();
+ }
+
+ /**
+ * Gets the initial guess.
+ *
+ * @return the initial guess.
+ */
+ public double[] getTarget() {
+ return target.clone();
+ }
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/Target.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/Weight.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/Weight.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/Weight.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/Weight.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,71 @@
+/*
+ * 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.vector;
+
+import org.apache.commons.math3.optim.OptimizationData;
+import org.apache.commons.math3.linear.RealMatrix;
+import org.apache.commons.math3.linear.MatrixUtils;
+import org.apache.commons.math3.linear.NonSquareMatrixException;
+
+/**
+ * Weight matrix of the residuals between model and observations.
+ * <br/>
+ * Immutable class.
+ *
+ * @version $Id: Weight.java 1416643 2012-12-03 19:37:14Z tn $
+ * @since 3.1
+ */
+public class Weight implements OptimizationData {
+ /** Weight matrix. */
+ private final RealMatrix weightMatrix;
+
+ /**
+ * Creates a diagonal weight matrix.
+ *
+ * @param weight List of the values of the diagonal.
+ */
+ public Weight(double[] weight) {
+ final int dim = weight.length;
+ weightMatrix = MatrixUtils.createRealMatrix(dim, dim);
+ for (int i = 0; i < dim; i++) {
+ weightMatrix.setEntry(i, i, weight[i]);
+ }
+ }
+
+ /**
+ * @param weight Weight matrix.
+ * @throws NonSquareMatrixException if the argument is not
+ * a square matrix.
+ */
+ public Weight(RealMatrix weight) {
+ if (weight.getColumnDimension() != weight.getRowDimension()) {
+ throw new NonSquareMatrixException(weight.getColumnDimension(),
+ weight.getRowDimension());
+ }
+
+ weightMatrix = weight.copy();
+ }
+
+ /**
+ * Gets the initial guess.
+ *
+ * @return the initial guess.
+ */
+ public RealMatrix getWeight() {
+ return weightMatrix.copy();
+ }
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/Weight.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/jacobian/AbstractLeastSquaresOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/jacobian/AbstractLeastSquaresOptimizer.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/jacobian/AbstractLeastSquaresOptimizer.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/jacobian/AbstractLeastSquaresOptimizer.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,269 @@
+/*
+ * 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.vector.jacobian;
+
+import org.apache.commons.math3.exception.DimensionMismatchException;
+import org.apache.commons.math3.exception.TooManyEvaluationsException;
+import org.apache.commons.math3.linear.ArrayRealVector;
+import org.apache.commons.math3.linear.RealMatrix;
+import org.apache.commons.math3.linear.DecompositionSolver;
+import org.apache.commons.math3.linear.MatrixUtils;
+import org.apache.commons.math3.linear.QRDecomposition;
+import org.apache.commons.math3.linear.EigenDecomposition;
+import org.apache.commons.math3.optim.OptimizationData;
+import org.apache.commons.math3.optim.ConvergenceChecker;
+import org.apache.commons.math3.optim.PointVectorValuePair;
+import org.apache.commons.math3.optim.nonlinear.vector.Weight;
+import org.apache.commons.math3.optim.nonlinear.vector.JacobianMultivariateVectorOptimizer;
+import org.apache.commons.math3.util.FastMath;
+
+/**
+ * Base class for implementing least-squares optimizers.
+ * It provides methods for error estimation.
+ *
+ * @version $Id$
+ * @since 3.1
+ */
+public abstract class AbstractLeastSquaresOptimizer
+ extends JacobianMultivariateVectorOptimizer {
+ /** Square-root of the weight matrix. */
+ private RealMatrix weightMatrixSqrt;
+ /** Cost value (square root of the sum of the residuals). */
+ private double cost;
+
+ /**
+ * @param checker Convergence checker.
+ */
+ protected AbstractLeastSquaresOptimizer(ConvergenceChecker<PointVectorValuePair> checker) {
+ super(checker);
+ }
+
+ /**
+ * Computes the weighted Jacobian matrix.
+ *
+ * @param params Model parameters at which to compute the Jacobian.
+ * @return the weighted Jacobian: W<sup>1/2</sup> J.
+ * @throws DimensionMismatchException if the Jacobian dimension does not
+ * match problem dimension.
+ */
+ protected RealMatrix computeWeightedJacobian(double[] params) {
+ return weightMatrixSqrt.multiply(MatrixUtils.createRealMatrix(computeJacobian(params)));
+ }
+
+ /**
+ * Computes the cost.
+ *
+ * @param residuals Residuals.
+ * @return the cost.
+ * @see #computeResiduals(double[])
+ */
+ protected double computeCost(double[] residuals) {
+ final ArrayRealVector r = new ArrayRealVector(residuals);
+ return FastMath.sqrt(r.dotProduct(getWeight().operate(r)));
+ }
+
+ /**
+ * Gets the root-mean-square (RMS) value.
+ *
+ * The RMS the root of the arithmetic mean of the square of all weighted
+ * residuals.
+ * This is related to the criterion that is minimized by the optimizer
+ * as follows: If <em>c</em> if the criterion, and <em>n</em> is the
+ * number of measurements, then the RMS is <em>sqrt (c/n)</em>.
+ *
+ * @return the RMS value.
+ */
+ public double getRMS() {
+ return FastMath.sqrt(getChiSquare() / getTargetSize());
+ }
+
+ /**
+ * Get a Chi-Square-like value assuming the N residuals follow N
+ * distinct normal distributions centered on 0 and whose variances are
+ * the reciprocal of the weights.
+ * @return chi-square value
+ */
+ public double getChiSquare() {
+ return cost * cost;
+ }
+
+ /**
+ * Gets the square-root of the weight matrix.
+ *
+ * @return the square-root of the weight matrix.
+ */
+ public RealMatrix getWeightSquareRoot() {
+ return weightMatrixSqrt.copy();
+ }
+
+ /**
+ * Sets the cost.
+ *
+ * @param cost Cost value.
+ */
+ protected void setCost(double cost) {
+ this.cost = cost;
+ }
+
+ /**
+ * Get the covariance matrix of the optimized parameters.
+ * <br/>
+ * Note that this operation involves the inversion of the
+ * <code>J<sup>T</sup>J</code> matrix, where {@code J} is the
+ * Jacobian matrix.
+ * The {@code threshold} parameter is a way for the caller to specify
+ * that the result of this computation should be considered meaningless,
+ * and thus trigger an exception.
+ *
+ * @param params Model parameters.
+ * @param threshold Singularity threshold.
+ * @return the covariance matrix.
+ * @throws org.apache.commons.math3.linear.SingularMatrixException
+ * if the covariance matrix cannot be computed (singular problem).
+ */
+ public double[][] computeCovariances(double[] params,
+ double threshold) {
+ // Set up the Jacobian.
+ final RealMatrix j = computeWeightedJacobian(params);
+
+ // Compute transpose(J)J.
+ final RealMatrix jTj = j.transpose().multiply(j);
+
+ // Compute the covariances matrix.
+ final DecompositionSolver solver
+ = new QRDecomposition(jTj, threshold).getSolver();
+ return solver.getInverse().getData();
+ }
+
+ /**
+ * Computes an estimate of the standard deviation of the parameters. The
+ * returned values are the square root of the diagonal coefficients of the
+ * covariance matrix, {@code sd(a[i]) ~= sqrt(C[i][i])}, where {@code a[i]}
+ * is the optimized value of the {@code i}-th parameter, and {@code C} is
+ * the covariance matrix.
+ *
+ * @param params Model parameters.
+ * @param covarianceSingularityThreshold Singularity threshold (see
+ * {@link #computeCovariances(double[],double) computeCovariances}).
+ * @return an estimate of the standard deviation of the optimized parameters
+ * @throws org.apache.commons.math3.linear.SingularMatrixException
+ * if the covariance matrix cannot be computed.
+ */
+ public double[] computeSigma(double[] params,
+ double covarianceSingularityThreshold) {
+ final int nC = params.length;
+ final double[] sig = new double[nC];
+ final double[][] cov = computeCovariances(params, covarianceSingularityThreshold);
+ for (int i = 0; i < nC; ++i) {
+ sig[i] = FastMath.sqrt(cov[i][i]);
+ }
+ return sig;
+ }
+
+ /**
+ * {@inheritDoc}
+ *
+ * @param optData Optimization data. The following data will be looked for:
+ * <ul>
+ * <li>{@link org.apache.commons.math3.optim.MaxEval}</li>
+ * <li>{@link org.apache.commons.math3.optim.InitialGuess}</li>
+ * <li>{@link org.apache.commons.math3.optim.SimpleBounds}</li>
+ * <li>{@link org.apache.commons.math3.optim.nonlinear.vector.Target}</li>
+ * <li>{@link org.apache.commons.math3.optim.nonlinear.vector.Weight}</li>
+ * <li>{@link org.apache.commons.math3.optim.nonlinear.vector.ModelFunction}</li>
+ * <li>{@link org.apache.commons.math3.optim.nonlinear.vector.ModelFunctionJacobian}</li>
+ * </ul>
+ * @return {@inheritDoc}
+ * @throws TooManyEvaluationsException if the maximal number of
+ * evaluations is exceeded.
+ * @throws DimensionMismatchException if the initial guess, target, and weight
+ * arguments have inconsistent dimensions.
+ */
+ @Override
+ public PointVectorValuePair optimize(OptimizationData... optData)
+ throws TooManyEvaluationsException {
+ // Retrieve settings.
+ parseOptimizationData(optData);
+ // Set up base class and perform computation.
+ return super.optimize(optData);
+ }
+
+ /**
+ * Computes the residuals.
+ * The residual is the difference between the observed (target)
+ * values and the model (objective function) value.
+ * There is one residual for each element of the vector-valued
+ * function.
+ *
+ * @param objectiveValue Value of the the objective function. This is
+ * the value returned from a call to
+ * {@link #computeObjectiveValue(double[]) computeObjectiveValue}
+ * (whose array argument contains the model parameters).
+ * @return the residuals.
+ * @throws DimensionMismatchException if {@code params} has a wrong
+ * length.
+ */
+ protected double[] computeResiduals(double[] objectiveValue) {
+ final double[] target = getTarget();
+ if (objectiveValue.length != target.length) {
+ throw new DimensionMismatchException(target.length,
+ objectiveValue.length);
+ }
+
+ final double[] residuals = new double[target.length];
+ for (int i = 0; i < target.length; i++) {
+ residuals[i] = target[i] - objectiveValue[i];
+ }
+
+ return residuals;
+ }
+
+ /**
+ * Scans the list of (required and optional) optimization data that
+ * characterize the problem.
+ * If the weight matrix is specified, the {@link #weightMatrixSqrt}
+ * field is recomputed.
+ *
+ * @param optData Optimization data. The following data will be looked for:
+ * <ul>
+ * <li>{@link Weight}</li>
+ * </ul>
+ */
+ private void parseOptimizationData(OptimizationData... optData) {
+ // The existing values (as set by the previous call) are reused if
+ // not provided in the argument list.
+ for (OptimizationData data : optData) {
+ if (data instanceof Weight) {
+ weightMatrixSqrt = squareRoot(((Weight) data).getWeight());
+ // If more data must be parsed, this statement _must_ be
+ // changed to "continue".
+ break;
+ }
+ }
+ }
+
+ /**
+ * Computes the square-root of the weight matrix.
+ *
+ * @param m Symmetric, positive-definite (weight) matrix.
+ * @return the square-root of the weight matrix.
+ */
+ private RealMatrix squareRoot(RealMatrix m) {
+ final EigenDecomposition dec = new EigenDecomposition(m);
+ return dec.getSquareRoot();
+ }
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/jacobian/AbstractLeastSquaresOptimizer.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/jacobian/GaussNewtonOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/jacobian/GaussNewtonOptimizer.java?rev=1420684&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/jacobian/GaussNewtonOptimizer.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/jacobian/GaussNewtonOptimizer.java Wed Dec 12 14:10:38 2012
@@ -0,0 +1,162 @@
+/*
+ * 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.vector.jacobian;
+
+import org.apache.commons.math3.exception.ConvergenceException;
+import org.apache.commons.math3.exception.NullArgumentException;
+import org.apache.commons.math3.exception.MathInternalError;
+import org.apache.commons.math3.exception.util.LocalizedFormats;
+import org.apache.commons.math3.linear.ArrayRealVector;
+import org.apache.commons.math3.linear.BlockRealMatrix;
+import org.apache.commons.math3.linear.DecompositionSolver;
+import org.apache.commons.math3.linear.LUDecomposition;
+import org.apache.commons.math3.linear.QRDecomposition;
+import org.apache.commons.math3.linear.RealMatrix;
+import org.apache.commons.math3.linear.SingularMatrixException;
+import org.apache.commons.math3.optim.ConvergenceChecker;
+import org.apache.commons.math3.optim.PointVectorValuePair;
+
+/**
+ * Gauss-Newton least-squares solver.
+ * <p>
+ * This class solve a least-square problem by solving the normal equations
+ * of the linearized problem at each iteration. Either LU decomposition or
+ * QR decomposition can be used to solve the normal equations. LU decomposition
+ * is faster but QR decomposition is more robust for difficult problems.
+ * </p>
+ *
+ * @version $Id: GaussNewtonOptimizer.java 1416643 2012-12-03 19:37:14Z tn $
+ * @since 2.0
+ *
+ */
+public class GaussNewtonOptimizer extends AbstractLeastSquaresOptimizer {
+ /** Indicator for using LU decomposition. */
+ private final boolean useLU;
+
+ /**
+ * Simple constructor with default settings.
+ * The normal equations will be solved using LU decomposition.
+ *
+ * @param checker Convergence checker.
+ */
+ public GaussNewtonOptimizer(ConvergenceChecker<PointVectorValuePair> checker) {
+ this(true, checker);
+ }
+
+ /**
+ * @param useLU If {@code true}, the normal equations will be solved
+ * using LU decomposition, otherwise they will be solved using QR
+ * decomposition.
+ * @param checker Convergence checker.
+ */
+ public GaussNewtonOptimizer(final boolean useLU,
+ ConvergenceChecker<PointVectorValuePair> checker) {
+ super(checker);
+ this.useLU = useLU;
+ }
+
+ /** {@inheritDoc} */
+ @Override
+ public PointVectorValuePair doOptimize() {
+ final ConvergenceChecker<PointVectorValuePair> checker
+ = getConvergenceChecker();
+
+ // Computation will be useless without a checker (see "for-loop").
+ if (checker == null) {
+ throw new NullArgumentException();
+ }
+
+ final double[] targetValues = getTarget();
+ final int nR = targetValues.length; // Number of observed data.
+
+ final RealMatrix weightMatrix = getWeight();
+ // Diagonal of the weight matrix.
+ final double[] residualsWeights = new double[nR];
+ for (int i = 0; i < nR; i++) {
+ residualsWeights[i] = weightMatrix.getEntry(i, i);
+ }
+
+ final double[] currentPoint = getStartPoint();
+ final int nC = currentPoint.length;
+
+ // iterate until convergence is reached
+ PointVectorValuePair current = null;
+ int iter = 0;
+ for (boolean converged = false; !converged;) {
+ ++iter;
+
+ // evaluate the objective function and its jacobian
+ PointVectorValuePair previous = current;
+ // Value of the objective function at "currentPoint".
+ final double[] currentObjective = computeObjectiveValue(currentPoint);
+ final double[] currentResiduals = computeResiduals(currentObjective);
+ final RealMatrix weightedJacobian = computeWeightedJacobian(currentPoint);
+ current = new PointVectorValuePair(currentPoint, currentObjective);
+
+ // build the linear problem
+ final double[] b = new double[nC];
+ final double[][] a = new double[nC][nC];
+ for (int i = 0; i < nR; ++i) {
+
+ final double[] grad = weightedJacobian.getRow(i);
+ final double weight = residualsWeights[i];
+ final double residual = currentResiduals[i];
+
+ // compute the normal equation
+ final double wr = weight * residual;
+ for (int j = 0; j < nC; ++j) {
+ b[j] += wr * grad[j];
+ }
+
+ // build the contribution matrix for measurement i
+ for (int k = 0; k < nC; ++k) {
+ double[] ak = a[k];
+ double wgk = weight * grad[k];
+ for (int l = 0; l < nC; ++l) {
+ ak[l] += wgk * grad[l];
+ }
+ }
+ }
+
+ try {
+ // solve the linearized least squares problem
+ RealMatrix mA = new BlockRealMatrix(a);
+ DecompositionSolver solver = useLU ?
+ new LUDecomposition(mA).getSolver() :
+ new QRDecomposition(mA).getSolver();
+ final double[] dX = solver.solve(new ArrayRealVector(b, false)).toArray();
+ // update the estimated parameters
+ for (int i = 0; i < nC; ++i) {
+ currentPoint[i] += dX[i];
+ }
+ } catch (SingularMatrixException e) {
+ throw new ConvergenceException(LocalizedFormats.UNABLE_TO_SOLVE_SINGULAR_PROBLEM);
+ }
+
+ // Check convergence.
+ if (previous != null) {
+ converged = checker.converged(iter, previous, current);
+ if (converged) {
+ setCost(computeCost(currentResiduals));
+ return current;
+ }
+ }
+ }
+ // Must never happen.
+ throw new MathInternalError();
+ }
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/vector/jacobian/GaussNewtonOptimizer.java
------------------------------------------------------------------------------
svn:eol-style = native