Check bounds in multi-start vector optimizers.

JIRA: MATH-914

git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1454746 13f79535-47bb-0310-9956-ffa450edef68

src/changes/changes.xml

@@ -55,6 +55,9 @@ This is a minor release: It combines bug fixes and new features.
   Changes to existing features were made in a backwards-compatible
   way such as to allow drop-in replacement of the v3.1[.1] JAR file.
 ">
+      <action dev="luc" type="add" issue="MATH-914" >
+        Check bounds in multi-start vector optimizers.
+      </action>
       <action dev="luc" type="add" issue="MATH-941" due-to="Piotr Wydrych" >
         Added discrete distributions.
       </action>

src/main/java/org/apache/commons/math3/optim/BaseMultiStartMultivariateOptimizer.java

@@ -18,6 +18,7 @@ package org.apache.commons.math3.optim;
 
 import org.apache.commons.math3.exception.MathIllegalStateException;
 import org.apache.commons.math3.exception.NotStrictlyPositiveException;
+import org.apache.commons.math3.exception.TooManyEvaluationsException;
 import org.apache.commons.math3.random.RandomVectorGenerator;
 
 /**
@@ -59,7 +60,15 @@ public abstract class BaseMultiStartMultivariateOptimizer<PAIR>
 
     /**
      * Create a multi-start optimizer from a single-start optimizer.
-     *
+     * <p>
+     * Note that if there are bounds constraints (see {@link #getLowerBound()}
+     * and {@link #getUpperBound()}), then a simple rejection algorithm is used
+     * at each restart. This implies that the random vector generator should have
+     * a good probability to generate vectors in the bounded domain, otherwise the
+     * rejection algorithm will hit the {@link #getMaxEvaluations()} count without
+     * generating a proper restart point. Users must be take great care of the <a
+     * href="http://en.wikipedia.org/wiki/Curse_of_dimensionality">curse of dimensionality</a>.
+     * </p>
      * @param optimizer Single-start optimizer to wrap.
      * @param starts Number of starts to perform. If {@code starts == 1},
      * the {@link #optimize(OptimizationData[]) optimize} will return the
@@ -157,8 +166,8 @@ public abstract class BaseMultiStartMultivariateOptimizer<PAIR>
         clear();
 
         final int maxEval = getMaxEvaluations();
-        final double[] min = getLowerBound(); // XXX Should be used to enforce bounds (see below).
+        final double[] min = getLowerBound();
-        final double[] max = getUpperBound(); // XXX Should be used to enforce bounds (see below).
+        final double[] max = getUpperBound();
         final double[] startPoint = getStartPoint();
 
         // Multi-start loop.
@@ -168,9 +177,24 @@ public abstract class BaseMultiStartMultivariateOptimizer<PAIR>
                 // Decrease number of allowed evaluations.
                 optimData[maxEvalIndex] = new MaxEval(maxEval - totalEvaluations);
                 // New start value.
-                final double[] s = (i == 0) ?
+                double[] s = null;
-                    startPoint :
+                if (i == 0) {
-                    generator.nextVector(); // XXX This does not enforce bounds!
+                    s = startPoint;
+                } else {
+                    int attempts = 0;
+                    while (s == null) {
+                        if (attempts++ >= getMaxEvaluations()) {
+                            throw new TooManyEvaluationsException(getMaxEvaluations());
+                        }
+                        s = generator.nextVector();
+                        for (int k = 0; s != null && k < s.length; ++k) {
+                            if ((min != null && s[k] < min[k]) || (max != null && s[k] > max[k])) {
+                                // reject the vector
+                                s = null;
+                            }
+                        }
+                    }
+                }
                 optimData[initialGuessIndex] = new InitialGuess(s);
                 // Optimize.
                 final PAIR result = optimizer.optimize(optimData);

src/test/java/org/apache/commons/math3/optim/nonlinear/vector/MultiStartMultivariateVectorOptimizerTest.java

@@ -16,14 +16,16 @@
  */
 package org.apache.commons.math3.optim.nonlinear.vector;
 
-import org.apache.commons.math3.analysis.MultivariateVectorFunction;
 import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
-import org.apache.commons.math3.exception.MathIllegalStateException;
+import org.apache.commons.math3.analysis.MultivariateVectorFunction;
+import org.apache.commons.math3.exception.TooManyEvaluationsException;
 import org.apache.commons.math3.linear.BlockRealMatrix;
 import org.apache.commons.math3.linear.RealMatrix;
-import org.apache.commons.math3.optim.MaxEval;
 import org.apache.commons.math3.optim.InitialGuess;
+import org.apache.commons.math3.optim.MaxEval;
+import org.apache.commons.math3.optim.OptimizationData;
 import org.apache.commons.math3.optim.PointVectorValuePair;
+import org.apache.commons.math3.optim.SimpleBounds;
 import org.apache.commons.math3.optim.SimpleVectorValueChecker;
 import org.apache.commons.math3.optim.nonlinear.vector.jacobian.GaussNewtonOptimizer;
 import org.apache.commons.math3.random.GaussianRandomGenerator;
@@ -96,10 +98,10 @@ import org.junit.Test;
  * @author Luc Maisonobe (non-minpack tests and minpack tests Java translation)
  */
 public class MultiStartMultivariateVectorOptimizerTest {
+
     @Test(expected=NullPointerException.class)
     public void testGetOptimaBeforeOptimize() {
-        LinearProblem problem
+
-            = new LinearProblem(new double[][] { { 2 } }, new double[] { 3 });
         JacobianMultivariateVectorOptimizer underlyingOptimizer
             = new GaussNewtonOptimizer(true, new SimpleVectorValueChecker(1e-6, 1e-6));
         JDKRandomGenerator g = new JDKRandomGenerator();
@@ -145,8 +147,45 @@ public class MultiStartMultivariateVectorOptimizerTest {
         Assert.assertEquals(100, optimizer.getMaxEvaluations());
     }
 
+    @Test
+    public void testIssue914() {
+        LinearProblem problem = new LinearProblem(new double[][] { { 2 } }, new double[] { 3 });
+        JacobianMultivariateVectorOptimizer underlyingOptimizer =
+                new GaussNewtonOptimizer(true, new SimpleVectorValueChecker(1e-6, 1e-6)) {
+            public PointVectorValuePair optimize(OptimizationData... optData) {
+                // filter out simple bounds, as they are not supported
+                // by the underlying optimizer, and we don't really care for this test
+                OptimizationData[] filtered = optData.clone();
+                for (int i = 0; i < filtered.length; ++i) {
+                    if (filtered[i] instanceof SimpleBounds) {
+                        filtered[i] = null;
+                    }
+                }
+                return super.optimize(filtered);
+            }
+        };
+        JDKRandomGenerator g = new JDKRandomGenerator();
+        g.setSeed(16069223052l);
+        RandomVectorGenerator generator =
+                new UncorrelatedRandomVectorGenerator(1, new GaussianRandomGenerator(g));
+        MultiStartMultivariateVectorOptimizer optimizer =
+                new MultiStartMultivariateVectorOptimizer(underlyingOptimizer, 10, generator);
+
+        optimizer.optimize(new MaxEval(100),
+                           problem.getModelFunction(),
+                           problem.getModelFunctionJacobian(),
+                           problem.getTarget(),
+                           new Weight(new double[] { 1 }),
+                           new InitialGuess(new double[] { 0 }),
+                           new SimpleBounds(new double[] { -1.0e-10 }, new double[] {  1.0e-10 }));
+        PointVectorValuePair[] optima = optimizer.getOptima();
+        // only the first start should have succeeded
+        Assert.assertEquals(1, optima.length);
+
+    }
+
     /**
-     * Test demonstrating that the user exception is fnally thrown if none
+     * Test demonstrating that the user exception is finally thrown if none
      * of the runs succeed.
      */
     @Test(expected=TestException.class)
@@ -170,7 +209,9 @@ public class MultiStartMultivariateVectorOptimizerTest {
                                }));
     }
 
-    private static class TestException extends RuntimeException {}
+    private static class TestException extends RuntimeException {
+
+    private static final long serialVersionUID = 1L;}
 
     private static class LinearProblem {
         private final RealMatrix factors;