handle case of negative rescalings (from negative dx) in COBYLA and BOBYQA; thanks...

Ignore-this: 2eb3c38f99b5d23bc5cf82aba25f1305

darcs-hash:20120730155158-c8de0-97cceb8446c78f4ed375222d1026ebf29f646506.gz

diff --git a/bobyqa/bobyqa.c b/bobyqa/bobyqa.c
index c4ddaf01fdabb2344a09d0bd3af3f146448e01de..8fc454b67442ef885242e8efda5b230f4ab93eb0 100644 (file)
--- a/bobyqa/bobyqa.c
+++ b/bobyqa/bobyqa.c
@@ -3108,8 +3108,9 @@ nlopt_result bobyqa(int n, int npt, double *x,
     sxu = nlopt_new_rescaled(U(n), s, xu);
     if (!sxu) { ret = NLOPT_OUT_OF_MEMORY; goto done; }
     xu = sxu;
+    nlopt_reorder_bounds(n, sxl, sxu);
 
-    rhobeg = dx[0] / s[0]; /* equals all other dx[i] after rescaling */
+    rhobeg = fabs(dx[0] / s[0]); /* equals all other dx[i] after rescaling */
 
     calfun_data.s = s;
     calfun_data.xs = xs;
@@ -3119,8 +3120,8 @@ nlopt_result bobyqa(int n, int npt, double *x,
     /* SGJ, 2009: compute rhoend from NLopt stop info */
     rhoend = stop->xtol_rel * (rhobeg);
     for (j = 0; j < n; ++j)
-        if (rhoend < stop->xtol_abs[j] / s[j])
-             rhoend = stop->xtol_abs[j] / s[j];
+        if (rhoend < stop->xtol_abs[j] / fabs(s[j]))
+             rhoend = stop->xtol_abs[j] / fabs(s[j]);
 
 
 /*     This subroutine seeks the least value of a function of many variables, */

diff --git a/cobyla/cobyla.c b/cobyla/cobyla.c
index 6a9a80efdf33698c1334a91a313dcaf01fbfc18e..763e01cbba52e1a89805cbd12b1b21fc8a8b5404 100644 (file)
--- a/cobyla/cobyla.c
+++ b/cobyla/cobyla.c
@@ -203,16 +203,17 @@ nlopt_result cobyla_minimize(unsigned n, nlopt_func f, void *f_data,
      if (!s.lb) { ret = NLOPT_OUT_OF_MEMORY; goto done; }
      s.ub = nlopt_new_rescaled(n, s.scale, ub);
      if (!s.ub) { ret = NLOPT_OUT_OF_MEMORY; goto done; }
+     nlopt_reorder_bounds(n, s.lb, s.ub);
 
      s.xtmp = (double *) malloc(sizeof(double) * n);
      if (!s.xtmp) { ret = NLOPT_OUT_OF_MEMORY; goto done; }
 
      /* SGJ, 2008: compute rhoend from NLopt stop info */
-     rhobeg = dx[0] / s.scale[0];
+     rhobeg = fabs(dx[0] / s.scale[0]);
      rhoend = stop->xtol_rel * (rhobeg);
      for (j = 0; j < n; ++j)
-         if (rhoend < stop->xtol_abs[j] / s.scale[j])
-              rhoend = stop->xtol_abs[j] / s.scale[j];
+         if (rhoend < stop->xtol_abs[j] / fabs(s.scale[j]))
+              rhoend = stop->xtol_abs[j] / fabs(s.scale[j]);
 
      /* each equality constraint gives two inequality constraints */
      m = nlopt_count_constraints(m, fc) + 2 * nlopt_count_constraints(p, h);

diff --git a/util/nlopt-util.h b/util/nlopt-util.h
index fa60d2bb97c8d90ab2c74d1527b3016588d02896..5e2745ed09972aae608d08b8fa521e10f0c3a913 100644 (file)
--- a/util/nlopt-util.h
+++ b/util/nlopt-util.h
@@ -127,6 +127,7 @@ double *nlopt_compute_rescaling(unsigned n, const double *dx);
 double *nlopt_new_rescaled(unsigned n, const double *s, const double *x);
 void nlopt_rescale(unsigned n, const double *s, const double *x, double *xs);
 void nlopt_unscale(unsigned n, const double *s, const double *x, double *xs);
+void nlopt_reorder_bounds(unsigned n, double *lb, double *ub);
 
 #ifdef __cplusplus
 }  /* extern "C" */

diff --git a/util/rescale.c b/util/rescale.c
index 11f2608e8f9e15020de2f048ff7e78ffd717755c..58004e808791462ecf93e5d14c17f7cc365c96d7 100644 (file)
--- a/util/rescale.c
+++ b/util/rescale.c
@@ -66,3 +66,17 @@ double *nlopt_new_rescaled(unsigned n, const double *s, const double *x)
      nlopt_rescale(n, s, x, xs);
      return xs;
 }
+
+/* since rescaling can flip the signs of the x components and the bounds,
+   we may have to re-order the bounds in order to ensure that they
+   remain in the correct order */
+void nlopt_reorder_bounds(unsigned n, double *lb, double *ub)
+{
+     unsigned i;
+     for (i = 0; i < n; ++i)
+         if (lb[i] > ub[i]) {
+              double t = lb[i];
+              lb[i] = ub[i];
+              ub[i] = t;
+         }
+}