Added:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/MultivariateFunctionPenaltyAdapter.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/MultivariateFunctionPenaltyAdapter.java?rev=1420684&view=auto
==============================================================================
---
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/MultivariateFunctionPenaltyAdapter.java
(added)
+++
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/MultivariateFunctionPenaltyAdapter.java
Wed Dec 12 14:10:38 2012
@@ -0,0 +1,187 @@
+/*
+ * 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;
+
+import org.apache.commons.math3.analysis.MultivariateFunction;
+import org.apache.commons.math3.exception.DimensionMismatchException;
+import org.apache.commons.math3.exception.NumberIsTooSmallException;
+import org.apache.commons.math3.util.FastMath;
+import org.apache.commons.math3.util.MathUtils;
+
+/**
+ * <p>Adapter extending bounded {@link MultivariateFunction} to an unbouded
+ * domain using a penalty function.</p>
+ *
+ * <p>
+ * This adapter can be used to wrap functions subject to simple bounds on
+ * parameters so they can be used by optimizers that do <em>not</em> directly
+ * support simple bounds.
+ * </p>
+ * <p>
+ * The principle is that the user function that will be wrapped will see its
+ * parameters bounded as required, i.e when its {@code value} method is called
+ * with argument array {@code point}, the elements array will fulfill
requirement
+ * {@code lower[i] <= point[i] <= upper[i]} for all i. Some of the components
+ * may be unbounded or bounded only on one side if the corresponding bound is
+ * set to an infinite value. The optimizer will not manage the user function by
+ * itself, but it will handle this adapter and it is this adapter that will
take
+ * care the bounds are fulfilled. The adapter {@link #value(double[])} method
will
+ * be called by the optimizer with unbound parameters, and the adapter will
check
+ * if the parameters is within range or not. If it is in range, then the
underlying
+ * user function will be called, and if it is not the value of a penalty
function
+ * will be returned instead.
+ * </p>
+ * <p>
+ * This adapter is only a poor-man's solution to simple bounds optimization
+ * constraints that can be used with simple optimizers like
+ * {@link
org.apache.commons.math3.optim.nonlinear.scalar.noderiv.SimplexOptimizer
+ * SimplexOptimizer}.
+ * A better solution is to use an optimizer that directly supports simple
bounds like
+ * {@link
org.apache.commons.math3.optim.nonlinear.scalar.noderiv.CMAESOptimizer
+ * CMAESOptimizer} or
+ * {@link
org.apache.commons.math3.optim.nonlinear.scalar.noderiv.BOBYQAOptimizer
+ * BOBYQAOptimizer}.
+ * One caveat of this poor-man's solution is that if start point or start
simplex
+ * is completely outside of the allowed range, only the penalty function is
used,
+ * and the optimizer may converge without ever entering the range.
+ * </p>
+ *
+ * @see MultivariateFunctionMappingAdapter
+ *
+ * @version $Id: MultivariateFunctionPenaltyAdapter.java 1416643 2012-12-03
19:37:14Z tn $
+ * @since 3.0
+ */
+public class MultivariateFunctionPenaltyAdapter
+ implements MultivariateFunction {
+ /** Underlying bounded function. */
+ private final MultivariateFunction bounded;
+ /** Lower bounds. */
+ private final double[] lower;
+ /** Upper bounds. */
+ private final double[] upper;
+ /** Penalty offset. */
+ private final double offset;
+ /** Penalty scales. */
+ private final double[] scale;
+
+ /**
+ * Simple constructor.
+ * <p>
+ * When the optimizer provided points are out of range, the value of the
+ * penalty function will be used instead of the value of the underlying
+ * function. In order for this penalty to be effective in rejecting this
+ * point during the optimization process, the penalty function value should
+ * be defined with care. This value is computed as:
+ * <pre>
+ * penalty(point) = offset + ∑<sub>i</sub>[scale[i] *
√|point[i]-boundary[i]|]
+ * </pre>
+ * where indices i correspond to all the components that violates their
boundaries.
+ * </p>
+ * <p>
+ * So when attempting a function minimization, offset should be larger than
+ * the maximum expected value of the underlying function and scale
components
+ * should all be positive. When attempting a function maximization, offset
+ * should be lesser than the minimum expected value of the underlying
function
+ * and scale components should all be negative.
+ * minimization, and lesser than the minimum expected value of the
underlying
+ * function when attempting maximization.
+ * </p>
+ * <p>
+ * These choices for the penalty function have two properties. First, all
out
+ * of range points will return a function value that is worse than the
value
+ * returned by any in range point. Second, the penalty is worse for large
+ * boundaries violation than for small violations, so the optimizer has an
hint
+ * about the direction in which it should search for acceptable points.
+ * </p>
+ * @param bounded bounded function
+ * @param lower lower bounds for each element of the input parameters array
+ * (some elements may be set to {@code Double.NEGATIVE_INFINITY} for
+ * unbounded values)
+ * @param upper upper bounds for each element of the input parameters array
+ * (some elements may be set to {@code Double.POSITIVE_INFINITY} for
+ * unbounded values)
+ * @param offset base offset of the penalty function
+ * @param scale scale of the penalty function
+ * @exception DimensionMismatchException if lower bounds, upper bounds and
+ * scales are not consistent, either according to dimension or to bounadary
+ * values
+ */
+ public MultivariateFunctionPenaltyAdapter(final MultivariateFunction
bounded,
+ final double[] lower, final
double[] upper,
+ final double offset, final
double[] scale) {
+
+ // safety checks
+ MathUtils.checkNotNull(lower);
+ MathUtils.checkNotNull(upper);
+ MathUtils.checkNotNull(scale);
+ if (lower.length != upper.length) {
+ throw new DimensionMismatchException(lower.length, upper.length);
+ }
+ if (lower.length != scale.length) {
+ throw new DimensionMismatchException(lower.length, scale.length);
+ }
+ for (int i = 0; i < lower.length; ++i) {
+ // note the following test is written in such a way it also fails
for NaN
+ if (!(upper[i] >= lower[i])) {
+ throw new NumberIsTooSmallException(upper[i], lower[i], true);
+ }
+ }
+
+ this.bounded = bounded;
+ this.lower = lower.clone();
+ this.upper = upper.clone();
+ this.offset = offset;
+ this.scale = scale.clone();
+ }
+
+ /**
+ * Computes the underlying function value from an unbounded point.
+ * <p>
+ * This method simply returns the value of the underlying function
+ * if the unbounded point already fulfills the bounds, and compute
+ * a replacement value using the offset and scale if bounds are
+ * violated, without calling the function at all.
+ * </p>
+ * @param point unbounded point
+ * @return either underlying function value or penalty function value
+ */
+ public double value(double[] point) {
+
+ for (int i = 0; i < scale.length; ++i) {
+ if ((point[i] < lower[i]) || (point[i] > upper[i])) {
+ // bound violation starting at this component
+ double sum = 0;
+ for (int j = i; j < scale.length; ++j) {
+ final double overshoot;
+ if (point[j] < lower[j]) {
+ overshoot = scale[j] * (lower[j] - point[j]);
+ } else if (point[j] > upper[j]) {
+ overshoot = scale[j] * (point[j] - upper[j]);
+ } else {
+ overshoot = 0;
+ }
+ sum += FastMath.sqrt(overshoot);
+ }
+ return offset + sum;
+ }
+ }
+
+ // all boundaries are fulfilled, we are in the expected
+ // domain of the underlying function
+ return bounded.value(point);
+ }
+}
Propchange:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/MultivariateFunctionPenaltyAdapter.java
------------------------------------------------------------------------------
svn:eol-style = native
Added:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/MultivariateOptimizer.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/MultivariateOptimizer.java?rev=1420684&view=auto
==============================================================================
---
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/MultivariateOptimizer.java
(added)
+++
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/MultivariateOptimizer.java
Wed Dec 12 14:10:38 2012
@@ -0,0 +1,119 @@
+/*
+ * 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;
+
+import org.apache.commons.math3.analysis.MultivariateFunction;
+import org.apache.commons.math3.optim.BaseMultivariateOptimizer;
+import org.apache.commons.math3.optim.OptimizationData;
+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.ObjectiveFunction;
+import org.apache.commons.math3.exception.TooManyEvaluationsException;
+
+/**
+ * Base class for a multivariate scalar function optimizer.
+ *
+ * @version $Id$
+ * @since 3.1
+ */
+public abstract class MultivariateOptimizer
+ extends BaseMultivariateOptimizer<PointValuePair> {
+ /** Objective function. */
+ private MultivariateFunction function;
+ /** Type of optimization. */
+ private GoalType goal;
+
+ /**
+ * @param checker Convergence checker.
+ */
+ protected MultivariateOptimizer(ConvergenceChecker<PointValuePair>
checker) {
+ super(checker);
+ }
+
+ /**
+ * {@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 ObjectiveFunction}</li>
+ * <li>{@link GoalType}</li>
+ * </ul>
+ * @return {@inheritDoc}
+ * @throws TooManyEvaluationsException if the maximal number of
+ * evaluations is exceeded.
+ */
+ @Override
+ public PointValuePair optimize(OptimizationData... optData)
+ throws TooManyEvaluationsException {
+ // 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 ObjectiveFunction}</li>
+ * <li>{@link GoalType}</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 GoalType) {
+ goal = (GoalType) data;
+ continue;
+ }
+ if (data instanceof ObjectiveFunction) {
+ function = ((ObjectiveFunction) data).getObjectiveFunction();
+ continue;
+ }
+ }
+ }
+
+ /**
+ * @return the optimization type.
+ */
+ public GoalType getGoalType() {
+ return goal;
+ }
+
+ /**
+ * 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 is exceeded.
+ */
+ protected double computeObjectiveValue(double[] params) {
+ super.incrementEvaluationCount();
+ return function.value(params);
+ }
+}
Propchange:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/MultivariateOptimizer.java
------------------------------------------------------------------------------
svn:eol-style = native
Added:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/ObjectiveFunctionGradient.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/ObjectiveFunctionGradient.java?rev=1420684&view=auto
==============================================================================
---
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/ObjectiveFunctionGradient.java
(added)
+++
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/ObjectiveFunctionGradient.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.scalar;
+
+import org.apache.commons.math3.analysis.MultivariateVectorFunction;
+import org.apache.commons.math3.optim.OptimizationData;
+
+/**
+ * Gradient of the scalar function to be optimized.
+ *
+ * @version $Id$
+ * @since 3.1
+ */
+public class ObjectiveFunctionGradient implements OptimizationData {
+ /** Function to be optimized. */
+ private final MultivariateVectorFunction gradient;
+
+ /**
+ * @param g Gradient of the function to be optimized.
+ */
+ public ObjectiveFunctionGradient(MultivariateVectorFunction g) {
+ gradient = g;
+ }
+
+ /**
+ * Gets the gradient of the function to be optimized.
+ *
+ * @return the objective function gradient.
+ */
+ public MultivariateVectorFunction getObjectiveFunctionGradient() {
+ return gradient;
+ }
+}
Propchange:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/ObjectiveFunctionGradient.java
------------------------------------------------------------------------------
svn:eol-style = native
Added:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/NonLinearConjugateGradientOptimizer.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/NonLinearConjugateGradientOptimizer.java?rev=1420684&view=auto
==============================================================================
---
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/NonLinearConjugateGradientOptimizer.java
(added)
+++
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/NonLinearConjugateGradientOptimizer.java
Wed Dec 12 14:10:38 2012
@@ -0,0 +1,393 @@
+/*
+ * 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.gradient;
+
+import org.apache.commons.math3.analysis.UnivariateFunction;
+import org.apache.commons.math3.analysis.solvers.BrentSolver;
+import org.apache.commons.math3.analysis.solvers.UnivariateSolver;
+import org.apache.commons.math3.exception.MathInternalError;
+import org.apache.commons.math3.exception.MathIllegalStateException;
+import org.apache.commons.math3.exception.TooManyEvaluationsException;
+import org.apache.commons.math3.exception.util.LocalizedFormats;
+import org.apache.commons.math3.optim.OptimizationData;
+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.GradientMultivariateOptimizer;
+import org.apache.commons.math3.util.FastMath;
+
+/**
+ * Non-linear conjugate gradient optimizer.
+ * <p>
+ * This class supports both the Fletcher-Reeves and the Polak-Ribière
+ * update formulas for the conjugate search directions.
+ * It also supports optional preconditioning.
+ * </p>
+ *
+ * @version $Id: NonLinearConjugateGradientOptimizer.java 1416643 2012-12-03
19:37:14Z tn $
+ * @since 2.0
+ */
+public class NonLinearConjugateGradientOptimizer
+ extends GradientMultivariateOptimizer {
+ /** Update formula for the beta parameter. */
+ private final Formula updateFormula;
+ /** Preconditioner (may be null). */
+ private final Preconditioner preconditioner;
+ /** solver to use in the line search (may be null). */
+ private final UnivariateSolver solver;
+ /** Initial step used to bracket the optimum in line search. */
+ private double initialStep = 1;
+
+ /**
+ * Constructor with default {@link BrentSolver line search solver} and
+ * {@link IdentityPreconditioner preconditioner}.
+ *
+ * @param updateFormula formula to use for updating the β parameter,
+ * must be one of {@link Formula#FLETCHER_REEVES} or
+ * {@link Formula#POLAK_RIBIERE}.
+ * @param checker Convergence checker.
+ */
+ public NonLinearConjugateGradientOptimizer(final Formula updateFormula,
+
ConvergenceChecker<PointValuePair> checker) {
+ this(updateFormula,
+ checker,
+ new BrentSolver(),
+ new IdentityPreconditioner());
+ }
+
+ /**
+ * Available choices of update formulas for the updating the parameter
+ * that is used to compute the successive conjugate search directions.
+ * For non-linear conjugate gradients, there are
+ * two formulas:
+ * <ul>
+ * <li>Fletcher-Reeves formula</li>
+ * <li>Polak-Ribière formula</li>
+ * </ul>
+ *
+ * On the one hand, the Fletcher-Reeves formula is guaranteed to converge
+ * if the start point is close enough of the optimum whether the
+ * Polak-Ribière formula may not converge in rare cases. On the
+ * other hand, the Polak-Ribière formula is often faster when it
+ * does converge. Polak-Ribière is often used.
+ *
+ * @since 2.0
+ */
+ public static enum Formula {
+ /** Fletcher-Reeves formula. */
+ FLETCHER_REEVES,
+ /** Polak-Ribière formula. */
+ POLAK_RIBIERE
+ }
+
+ /**
+ * The initial step is a factor with respect to the search direction
+ * (which itself is roughly related to the gradient of the function).
+ * <br/>
+ * It is used to find an interval that brackets the optimum in line
+ * search.
+ *
+ * @since 3.1
+ */
+ public static class BracketingStep implements OptimizationData {
+ /** Initial step. */
+ private final double initialStep;
+
+ /**
+ * @param step Initial step for the bracket search.
+ */
+ public BracketingStep(double step) {
+ initialStep = step;
+ }
+
+ /**
+ * Gets the initial step.
+ *
+ * @return the initial step.
+ */
+ public double getBracketingStep() {
+ return initialStep;
+ }
+ }
+
+ /**
+ * Constructor with default {@link IdentityPreconditioner preconditioner}.
+ *
+ * @param updateFormula formula to use for updating the β parameter,
+ * must be one of {@link Formula#FLETCHER_REEVES} or
+ * {@link Formula#POLAK_RIBIERE}.
+ * @param checker Convergence checker.
+ * @param lineSearchSolver Solver to use during line search.
+ */
+ public NonLinearConjugateGradientOptimizer(final Formula updateFormula,
+
ConvergenceChecker<PointValuePair> checker,
+ final UnivariateSolver
lineSearchSolver) {
+ this(updateFormula,
+ checker,
+ lineSearchSolver,
+ new IdentityPreconditioner());
+ }
+
+ /**
+ * @param updateFormula formula to use for updating the β parameter,
+ * must be one of {@link Formula#FLETCHER_REEVES} or
+ * {@link Formula#POLAK_RIBIERE}.
+ * @param checker Convergence checker.
+ * @param lineSearchSolver Solver to use during line search.
+ * @param preconditioner Preconditioner.
+ */
+ public NonLinearConjugateGradientOptimizer(final Formula updateFormula,
+
ConvergenceChecker<PointValuePair> checker,
+ final UnivariateSolver
lineSearchSolver,
+ final Preconditioner
preconditioner) {
+ super(checker);
+
+ this.updateFormula = updateFormula;
+ solver = lineSearchSolver;
+ this.preconditioner = preconditioner;
+ initialStep = 1;
+ }
+
+ /**
+ * {@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.GoalType}</li>
+ * <li>{@link org.apache.commons.math3.optim.ObjectiveFunction}</li>
+ * <li>{@link
org.apache.commons.math3.optim.nonlinear.scalar.ObjectiveFunctionGradient}</li>
+ * <li>{@link BracketingStep}</li>
+ * </ul>
+ * @return {@inheritDoc}
+ * @throws TooManyEvaluationsException if the maximal number of
+ * evaluations (of the objective function) is exceeded.
+ */
+ @Override
+ public PointValuePair optimize(OptimizationData... optData)
+ throws TooManyEvaluationsException {
+ // Retrieve settings.
+ parseOptimizationData(optData);
+ // Set up base class and perform computation.
+ return super.optimize(optData);
+ }
+
+ /** {@inheritDoc} */
+ @Override
+ protected PointValuePair doOptimize() {
+ final ConvergenceChecker<PointValuePair> checker =
getConvergenceChecker();
+ final double[] point = getStartPoint();
+ final GoalType goal = getGoalType();
+ final int n = point.length;
+ double[] r = computeObjectiveGradient(point);
+ if (goal == GoalType.MINIMIZE) {
+ for (int i = 0; i < n; i++) {
+ r[i] = -r[i];
+ }
+ }
+
+ // Initial search direction.
+ double[] steepestDescent = preconditioner.precondition(point, r);
+ double[] searchDirection = steepestDescent.clone();
+
+ double delta = 0;
+ for (int i = 0; i < n; ++i) {
+ delta += r[i] * searchDirection[i];
+ }
+
+ PointValuePair current = null;
+ int iter = 0;
+ int maxEval = getMaxEvaluations();
+ while (true) {
+ ++iter;
+
+ final double objective = computeObjectiveValue(point);
+ PointValuePair previous = current;
+ current = new PointValuePair(point, objective);
+ if (previous != null) {
+ if (checker.converged(iter, previous, current)) {
+ // We have found an optimum.
+ return current;
+ }
+ }
+
+ // Find the optimal step in the search direction.
+ final UnivariateFunction lsf = new LineSearchFunction(point,
searchDirection);
+ final double uB = findUpperBound(lsf, 0, initialStep);
+ // XXX Last parameters is set to a value close to zero in order to
+ // work around the divergence problem in the "testCircleFitting"
+ // unit test (see MATH-439).
+ final double step = solver.solve(maxEval, lsf, 0, uB, 1e-15);
+ maxEval -= solver.getEvaluations(); // Subtract used up
evaluations.
+
+ // Validate new point.
+ for (int i = 0; i < point.length; ++i) {
+ point[i] += step * searchDirection[i];
+ }
+
+ r = computeObjectiveGradient(point);
+ if (goal == GoalType.MINIMIZE) {
+ for (int i = 0; i < n; ++i) {
+ r[i] = -r[i];
+ }
+ }
+
+ // Compute beta.
+ final double deltaOld = delta;
+ final double[] newSteepestDescent =
preconditioner.precondition(point, r);
+ delta = 0;
+ for (int i = 0; i < n; ++i) {
+ delta += r[i] * newSteepestDescent[i];
+ }
+
+ final double beta;
+ switch (updateFormula) {
+ case FLETCHER_REEVES:
+ beta = delta / deltaOld;
+ break;
+ case POLAK_RIBIERE:
+ double deltaMid = 0;
+ for (int i = 0; i < r.length; ++i) {
+ deltaMid += r[i] * steepestDescent[i];
+ }
+ beta = (delta - deltaMid) / deltaOld;
+ break;
+ default:
+ // Should never happen.
+ throw new MathInternalError();
+ }
+ steepestDescent = newSteepestDescent;
+
+ // Compute conjugate search direction.
+ if (iter % n == 0 ||
+ beta < 0) {
+ // Break conjugation: reset search direction.
+ searchDirection = steepestDescent.clone();
+ } else {
+ // Compute new conjugate search direction.
+ for (int i = 0; i < n; ++i) {
+ searchDirection[i] = steepestDescent[i] + beta *
searchDirection[i];
+ }
+ }
+ }
+ }
+
+ /**
+ * 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 InitialStep}</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 BracketingStep) {
+ initialStep = ((BracketingStep) data).getBracketingStep();
+ // If more data must be parsed, this statement _must_ be
+ // changed to "continue".
+ break;
+ }
+ }
+ }
+
+ /**
+ * Finds the upper bound b ensuring bracketing of a root between a and b.
+ *
+ * @param f function whose root must be bracketed.
+ * @param a lower bound of the interval.
+ * @param h initial step to try.
+ * @return b such that f(a) and f(b) have opposite signs.
+ * @throws MathIllegalStateException if no bracket can be found.
+ */
+ private double findUpperBound(final UnivariateFunction f,
+ final double a, final double h) {
+ final double yA = f.value(a);
+ double yB = yA;
+ for (double step = h; step < Double.MAX_VALUE; step *= FastMath.max(2,
yA / yB)) {
+ final double b = a + step;
+ yB = f.value(b);
+ if (yA * yB <= 0) {
+ return b;
+ }
+ }
+ throw new
MathIllegalStateException(LocalizedFormats.UNABLE_TO_BRACKET_OPTIMUM_IN_LINE_SEARCH);
+ }
+
+ /** Default identity preconditioner. */
+ public static class IdentityPreconditioner implements Preconditioner {
+ /** {@inheritDoc} */
+ public double[] precondition(double[] variables, double[] r) {
+ return r.clone();
+ }
+ }
+
+ /**
+ * Internal class for line search.
+ * <p>
+ * The function represented by this class is the dot product of
+ * the objective function gradient and the search direction. Its
+ * value is zero when the gradient is orthogonal to the search
+ * direction, i.e. when the objective function value is a local
+ * extremum along the search direction.
+ * </p>
+ */
+ private class LineSearchFunction implements UnivariateFunction {
+ /** Current point. */
+ private final double[] currentPoint;
+ /** Search direction. */
+ private final double[] searchDirection;
+
+ /**
+ * @param point Current point.
+ * @param direction Search direction.
+ */
+ public LineSearchFunction(double[] point,
+ double[] direction) {
+ currentPoint = point.clone();
+ searchDirection = direction.clone();
+ }
+
+ /** {@inheritDoc} */
+ public double value(double x) {
+ // current point in the search direction
+ final double[] shiftedPoint = currentPoint.clone();
+ for (int i = 0; i < shiftedPoint.length; ++i) {
+ shiftedPoint[i] += x * searchDirection[i];
+ }
+
+ // gradient of the objective function
+ final double[] gradient = computeObjectiveGradient(shiftedPoint);
+
+ // dot product with the search direction
+ double dotProduct = 0;
+ for (int i = 0; i < gradient.length; ++i) {
+ dotProduct += gradient[i] * searchDirection[i];
+ }
+
+ return dotProduct;
+ }
+ }
+}
Propchange:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/NonLinearConjugateGradientOptimizer.java
------------------------------------------------------------------------------
svn:eol-style = native
Added:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/Preconditioner.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/Preconditioner.java?rev=1420684&view=auto
==============================================================================
---
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/Preconditioner.java
(added)
+++
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/Preconditioner.java
Wed Dec 12 14:10:38 2012
@@ -0,0 +1,45 @@
+/*
+ * 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.gradient;
+
+/**
+ * This interface represents a preconditioner for differentiable scalar
+ * objective function optimizers.
+ * @version $Id: Preconditioner.java 1416643 2012-12-03 19:37:14Z tn $
+ * @since 2.0
+ */
+public interface Preconditioner {
+ /**
+ * Precondition a search direction.
+ * <p>
+ * The returned preconditioned search direction must be computed fast or
+ * the algorithm performances will drop drastically. A classical approach
+ * is to compute only the diagonal elements of the hessian and to divide
+ * the raw search direction by these elements if they are all positive.
+ * If at least one of them is negative, it is safer to return a clone of
+ * the raw search direction as if the hessian was the identity matrix. The
+ * rationale for this simplified choice is that a negative diagonal element
+ * means the current point is far from the optimum and preconditioning will
+ * not be efficient anyway in this case.
+ * </p>
+ * @param point current point at which the search direction was computed
+ * @param r raw search direction (i.e. opposite of the gradient)
+ * @return approximation of H<sup>-1</sup>r where H is the objective
function hessian
+ */
+ double[] precondition(double[] point, double[] r);
+}
Propchange:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/Preconditioner.java
------------------------------------------------------------------------------
svn:eol-style = native
Added:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/package-info.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/package-info.java?rev=1420684&view=auto
==============================================================================
---
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/package-info.java
(added)
+++
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/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.gradient;
+
+/**
+ * This package provides optimization algorithms that require derivatives.
+ */
Propchange:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/gradient/package-info.java
------------------------------------------------------------------------------
svn:eol-style = native
Added:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/AbstractSimplex.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/AbstractSimplex.java?rev=1420684&view=auto
==============================================================================
---
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/AbstractSimplex.java
(added)
+++
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/AbstractSimplex.java
Wed Dec 12 14:10:38 2012
@@ -0,0 +1,346 @@
+/*
+ * 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.Arrays;
+import java.util.Comparator;
+
+import org.apache.commons.math3.analysis.MultivariateFunction;
+import org.apache.commons.math3.exception.NotStrictlyPositiveException;
+import org.apache.commons.math3.exception.DimensionMismatchException;
+import org.apache.commons.math3.exception.ZeroException;
+import org.apache.commons.math3.exception.OutOfRangeException;
+import org.apache.commons.math3.exception.NullArgumentException;
+import org.apache.commons.math3.exception.MathIllegalArgumentException;
+import org.apache.commons.math3.exception.util.LocalizedFormats;
+import org.apache.commons.math3.optim.PointValuePair;
+import org.apache.commons.math3.optim.OptimizationData;
+
+/**
+ * This class implements the simplex concept.
+ * It is intended to be used in conjunction with {@link SimplexOptimizer}.
+ * <br/>
+ * The initial configuration of the simplex is set by the constructors
+ * {@link #AbstractSimplex(double[])} or {@link #AbstractSimplex(double[][])}.
+ * The other {@link #AbstractSimplex(int) constructor} will set all steps
+ * to 1, thus building a default configuration from a unit hypercube.
+ * <br/>
+ * Users <em>must</em> call the {@link #build(double[]) build} method in order
+ * to create the data structure that will be acted on by the other methods of
+ * this class.
+ *
+ * @see SimplexOptimizer
+ * @version $Id: AbstractSimplex.java 1397759 2012-10-13 01:12:58Z erans $
+ * @since 3.0
+ */
+public abstract class AbstractSimplex implements OptimizationData {
+ /** Simplex. */
+ private PointValuePair[] simplex;
+ /** Start simplex configuration. */
+ private double[][] startConfiguration;
+ /** Simplex dimension (must be equal to {@code simplex.length - 1}). */
+ private final int dimension;
+
+ /**
+ * Build a unit hypercube simplex.
+ *
+ * @param n Dimension of the simplex.
+ */
+ protected AbstractSimplex(int n) {
+ this(n, 1d);
+ }
+
+ /**
+ * Build a hypercube simplex with the given side length.
+ *
+ * @param n Dimension of the simplex.
+ * @param sideLength Length of the sides of the hypercube.
+ */
+ protected AbstractSimplex(int n,
+ double sideLength) {
+ this(createHypercubeSteps(n, sideLength));
+ }
+
+ /**
+ * The start configuration for simplex is built from a box parallel to
+ * the canonical axes of the space. The simplex is the subset of vertices
+ * of a box parallel to the canonical axes. It is built as the path
followed
+ * while traveling from one vertex of the box to the diagonally opposite
+ * vertex moving only along the box edges. The first vertex of the box will
+ * be located at the start point of the optimization.
+ * As an example, in dimension 3 a simplex has 4 vertices. Setting the
+ * steps to (1, 10, 2) and the start point to (1, 1, 1) would imply the
+ * start simplex would be: { (1, 1, 1), (2, 1, 1), (2, 11, 1), (2, 11, 3)
}.
+ * The first vertex would be set to the start point at (1, 1, 1) and the
+ * last vertex would be set to the diagonally opposite vertex at (2, 11,
3).
+ *
+ * @param steps Steps along the canonical axes representing box edges. They
+ * may be negative but not zero.
+ * @throws NullArgumentException if {@code steps} is {@code null}.
+ * @throws ZeroException if one of the steps is zero.
+ */
+ protected AbstractSimplex(final double[] steps) {
+ if (steps == null) {
+ throw new NullArgumentException();
+ }
+ if (steps.length == 0) {
+ throw new ZeroException();
+ }
+ dimension = steps.length;
+
+ // Only the relative position of the n final vertices with respect
+ // to the first one are stored.
+ startConfiguration = new double[dimension][dimension];
+ for (int i = 0; i < dimension; i++) {
+ final double[] vertexI = startConfiguration[i];
+ for (int j = 0; j < i + 1; j++) {
+ if (steps[j] == 0) {
+ throw new
ZeroException(LocalizedFormats.EQUAL_VERTICES_IN_SIMPLEX);
+ }
+ System.arraycopy(steps, 0, vertexI, 0, j + 1);
+ }
+ }
+ }
+
+ /**
+ * The real initial simplex will be set up by moving the reference
+ * simplex such that its first point is located at the start point of the
+ * optimization.
+ *
+ * @param referenceSimplex Reference simplex.
+ * @throws NotStrictlyPositiveException if the reference simplex does not
+ * contain at least one point.
+ * @throws DimensionMismatchException if there is a dimension mismatch
+ * in the reference simplex.
+ * @throws IllegalArgumentException if one of its vertices is duplicated.
+ */
+ protected AbstractSimplex(final double[][] referenceSimplex) {
+ if (referenceSimplex.length <= 0) {
+ throw new
NotStrictlyPositiveException(LocalizedFormats.SIMPLEX_NEED_ONE_POINT,
+ referenceSimplex.length);
+ }
+ dimension = referenceSimplex.length - 1;
+
+ // Only the relative position of the n final vertices with respect
+ // to the first one are stored.
+ startConfiguration = new double[dimension][dimension];
+ final double[] ref0 = referenceSimplex[0];
+
+ // Loop over vertices.
+ for (int i = 0; i < referenceSimplex.length; i++) {
+ final double[] refI = referenceSimplex[i];
+
+ // Safety checks.
+ if (refI.length != dimension) {
+ throw new DimensionMismatchException(refI.length, dimension);
+ }
+ for (int j = 0; j < i; j++) {
+ final double[] refJ = referenceSimplex[j];
+ boolean allEquals = true;
+ for (int k = 0; k < dimension; k++) {
+ if (refI[k] != refJ[k]) {
+ allEquals = false;
+ break;
+ }
+ }
+ if (allEquals) {
+ throw new
MathIllegalArgumentException(LocalizedFormats.EQUAL_VERTICES_IN_SIMPLEX,
+ i, j);
+ }
+ }
+
+ // Store vertex i position relative to vertex 0 position.
+ if (i > 0) {
+ final double[] confI = startConfiguration[i - 1];
+ for (int k = 0; k < dimension; k++) {
+ confI[k] = refI[k] - ref0[k];
+ }
+ }
+ }
+ }
+
+ /**
+ * Get simplex dimension.
+ *
+ * @return the dimension of the simplex.
+ */
+ public int getDimension() {
+ return dimension;
+ }
+
+ /**
+ * Get simplex size.
+ * After calling the {@link #build(double[]) build} method, this method
will
+ * will be equivalent to {@code getDimension() + 1}.
+ *
+ * @return the size of the simplex.
+ */
+ public int getSize() {
+ return simplex.length;
+ }
+
+ /**
+ * Compute the next simplex of the algorithm.
+ *
+ * @param evaluationFunction Evaluation function.
+ * @param comparator Comparator to use to sort simplex vertices from best
+ * to worst.
+ * @throws org.apache.commons.math3.exception.TooManyEvaluationsException
+ * if the algorithm fails to converge.
+ */
+ public abstract void iterate(final MultivariateFunction evaluationFunction,
+ final Comparator<PointValuePair> comparator);
+
+ /**
+ * Build an initial simplex.
+ *
+ * @param startPoint First point of the simplex.
+ * @throws DimensionMismatchException if the start point does not match
+ * simplex dimension.
+ */
+ public void build(final double[] startPoint) {
+ if (dimension != startPoint.length) {
+ throw new DimensionMismatchException(dimension, startPoint.length);
+ }
+
+ // Set first vertex.
+ simplex = new PointValuePair[dimension + 1];
+ simplex[0] = new PointValuePair(startPoint, Double.NaN);
+
+ // Set remaining vertices.
+ for (int i = 0; i < dimension; i++) {
+ final double[] confI = startConfiguration[i];
+ final double[] vertexI = new double[dimension];
+ for (int k = 0; k < dimension; k++) {
+ vertexI[k] = startPoint[k] + confI[k];
+ }
+ simplex[i + 1] = new PointValuePair(vertexI, Double.NaN);
+ }
+ }
+
+ /**
+ * Evaluate all the non-evaluated points of the simplex.
+ *
+ * @param evaluationFunction Evaluation function.
+ * @param comparator Comparator to use to sort simplex vertices from best
to worst.
+ * @throws org.apache.commons.math3.exception.TooManyEvaluationsException
+ * if the maximal number of evaluations is exceeded.
+ */
+ public void evaluate(final MultivariateFunction evaluationFunction,
+ final Comparator<PointValuePair> comparator) {
+ // Evaluate the objective function at all non-evaluated simplex points.
+ for (int i = 0; i < simplex.length; i++) {
+ final PointValuePair vertex = simplex[i];
+ final double[] point = vertex.getPointRef();
+ if (Double.isNaN(vertex.getValue())) {
+ simplex[i] = new PointValuePair(point,
evaluationFunction.value(point), false);
+ }
+ }
+
+ // Sort the simplex from best to worst.
+ Arrays.sort(simplex, comparator);
+ }
+
+ /**
+ * Replace the worst point of the simplex by a new point.
+ *
+ * @param pointValuePair Point to insert.
+ * @param comparator Comparator to use for sorting the simplex vertices
+ * from best to worst.
+ */
+ protected void replaceWorstPoint(PointValuePair pointValuePair,
+ final Comparator<PointValuePair>
comparator) {
+ for (int i = 0; i < dimension; i++) {
+ if (comparator.compare(simplex[i], pointValuePair) > 0) {
+ PointValuePair tmp = simplex[i];
+ simplex[i] = pointValuePair;
+ pointValuePair = tmp;
+ }
+ }
+ simplex[dimension] = pointValuePair;
+ }
+
+ /**
+ * Get the points of the simplex.
+ *
+ * @return all the simplex points.
+ */
+ public PointValuePair[] getPoints() {
+ final PointValuePair[] copy = new PointValuePair[simplex.length];
+ System.arraycopy(simplex, 0, copy, 0, simplex.length);
+ return copy;
+ }
+
+ /**
+ * Get the simplex point stored at the requested {@code index}.
+ *
+ * @param index Location.
+ * @return the point at location {@code index}.
+ */
+ public PointValuePair getPoint(int index) {
+ if (index < 0 ||
+ index >= simplex.length) {
+ throw new OutOfRangeException(index, 0, simplex.length - 1);
+ }
+ return simplex[index];
+ }
+
+ /**
+ * Store a new point at location {@code index}.
+ * Note that no deep-copy of {@code point} is performed.
+ *
+ * @param index Location.
+ * @param point New value.
+ */
+ protected void setPoint(int index, PointValuePair point) {
+ if (index < 0 ||
+ index >= simplex.length) {
+ throw new OutOfRangeException(index, 0, simplex.length - 1);
+ }
+ simplex[index] = point;
+ }
+
+ /**
+ * Replace all points.
+ * Note that no deep-copy of {@code points} is performed.
+ *
+ * @param points New Points.
+ */
+ protected void setPoints(PointValuePair[] points) {
+ if (points.length != simplex.length) {
+ throw new DimensionMismatchException(points.length,
simplex.length);
+ }
+ simplex = points;
+ }
+
+ /**
+ * Create steps for a unit hypercube.
+ *
+ * @param n Dimension of the hypercube.
+ * @param sideLength Length of the sides of the hypercube.
+ * @return the steps.
+ */
+ private static double[] createHypercubeSteps(int n,
+ double sideLength) {
+ final double[] steps = new double[n];
+ for (int i = 0; i < n; i++) {
+ steps[i] = sideLength;
+ }
+ return steps;
+ }
+}
Propchange:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optim/nonlinear/scalar/noderiv/AbstractSimplex.java
------------------------------------------------------------------------------
svn:eol-style = native
