added NEWUOA, unified starting step-size for derivative-free algorithms

darcs-hash:20080830102737-c8de0-2f48550a7c042db585110d7a05b566c02fc7b892.gz

diff --git a/Makefile.am b/Makefile.am
index d716f6457d167be25eea2a5d497a57b81f2e0ff5..144845fe6c2e6cda0df8d3708abe89d3b2206ae7 100644 (file)
--- a/Makefile.am
+++ b/Makefile.am
@@ -8,7 +8,7 @@ CXX_DIRS = stogo
 CXX_LIBS = stogo/libstogo.la
 endif
 
-SUBDIRS= util subplex direct cdirect $(CXX_DIRS) praxis luksan crs mlsl mma cobyla lbfgs api . octave test
+SUBDIRS= util subplex direct cdirect $(CXX_DIRS) praxis luksan crs mlsl mma cobyla newuoa lbfgs api . octave test
 EXTRA_DIST=COPYRIGHT autogen.sh nlopt.pc.in m4
 
 if WITH_NOCEDAL
@@ -19,7 +19,7 @@ libnlopt@NLOPT_SUFFIX@_la_SOURCES =
 libnlopt@NLOPT_SUFFIX@_la_LIBADD = subplex/libsubplex.la               \
 direct/libdirect.la cdirect/libcdirect.la $(CXX_LIBS)                  \
 praxis/libpraxis.la $(NOCEDAL_LBFGS) luksan/libluksan.la crs/libcrs.la \
-mlsl/libmlsl.la mma/libmma.la cobyla/libcobyla.la api/libapi.la util/libutil.la
+mlsl/libmlsl.la mma/libmma.la cobyla/libcobyla.la newuoa/libnewuoa.la api/libapi.la util/libutil.la
 
 libnlopt@NLOPT_SUFFIX@_la_LDFLAGS = -no-undefined -version-info @SHARED_VERSION_INFO@
 

diff --git a/api/Makefile.am b/api/Makefile.am
index 8fa22fd8ec4ebaa01021496c814d93f244242c1a..1ec9af25e805f900872091ae87f4f53d7d310346 100644 (file)
--- a/api/Makefile.am
+++ b/api/Makefile.am
@@ -1,4 +1,4 @@
-AM_CPPFLAGS = -I$(top_srcdir)/cdirect -I$(top_srcdir)/direct -I$(top_srcdir)/stogo -I$(top_srcdir)/subplex -I$(top_srcdir)/praxis -I$(top_srcdir)/lbfgs -I$(top_srcdir)/luksan -I$(top_srcdir)/crs -I$(top_srcdir)/mlsl -I$(top_srcdir)/mma -I$(top_srcdir)/cobyla -I$(top_srcdir)/util
+AM_CPPFLAGS = -I$(top_srcdir)/cdirect -I$(top_srcdir)/direct -I$(top_srcdir)/stogo -I$(top_srcdir)/subplex -I$(top_srcdir)/praxis -I$(top_srcdir)/lbfgs -I$(top_srcdir)/luksan -I$(top_srcdir)/crs -I$(top_srcdir)/mlsl -I$(top_srcdir)/mma -I$(top_srcdir)/cobyla -I$(top_srcdir)/newuoa -I$(top_srcdir)/util
 
 include_HEADERS = nlopt.h nlopt.f
 noinst_LTLIBRARIES = libapi.la

diff --git a/api/nlopt.c b/api/nlopt.c
index 0d35d586e36ea060720387aa1ca6b28fd8d1a5d9..49fa53484045242063fa5dcbc0367cb2c71750e2 100644 (file)
--- a/api/nlopt.c
+++ b/api/nlopt.c
@@ -99,7 +99,8 @@ static const char nlopt_algorithm_names[NLOPT_NUM_ALGORITHMS][256] = {
      "Multi-level single-linkage (MLSL), quasi-random (global, no-derivative)",
      "Multi-level single-linkage (MLSL), quasi-random (global, derivative)",
      "Method of Moving Asymptotes (MMA) (local, derivative)",
-     "COBYLA (Constrained Optimization BY Linear Approximations) (local, no-derivative)"
+     "COBYLA (Constrained Optimization BY Linear Approximations) (local, no-derivative)",
+     "NEWUOA unconstrained optimization via quadratic models (local, no-derivative)"
 };
 
 const char *nlopt_algorithm_name(nlopt_algorithm a)
@@ -122,6 +123,45 @@ void nlopt_srand_time() {
 
 /*************************************************************************/
 
+/* crude heuristics for initial step size of nonderivative algorithms */
+static double initial_step(int n, 
+                          const double *lb, const double *ub, const double *x)
+{
+     int i;
+     double step = HUGE_VAL;
+     for (i = 0; i < n; ++i) {
+         if (!nlopt_isinf(ub[i]) && !nlopt_isinf(lb[i])
+             && (ub[i] - lb[i]) * 0.25 < step && ub[i] > lb[i])
+              step = (ub[i] - lb[i]) * 0.25;
+         if (!nlopt_isinf(ub[i]) 
+             && ub[i] - x[i] < step && ub[i] > x[i])
+              step = ub[i] - x[i];
+         if (!nlopt_isinf(lb[i]) 
+             && x[i] - lb[i] < step && x[i] > lb[i])
+              step = x[i] - lb[i];
+     }
+     if (nlopt_isinf(step))
+         for (i = 0; i < n; ++i) {
+              if (!nlopt_isinf(ub[i]) 
+                  && ub[i] - x[i] < step && ub[i] > x[i] + 1e-4)
+                   step = ub[i] - x[i];
+              if (!nlopt_isinf(lb[i]) 
+                  && x[i] - lb[i] < step && x[i] > lb[i] + 1e-4)
+                   step = x[i] - lb[i];
+         }
+     if (nlopt_isinf(step)) {
+         step = 0;
+         for (i = 0; i < n; ++i)
+              if (fabs(x[i]) * 0.25 > step)
+                   step = fabs(x[i]) * 0.25;
+         if (step == 0)
+              step = 1;
+     }
+     return step;
+}
+             
+/*************************************************************************/
+
 typedef struct {
      nlopt_func f;
      void *f_data;
@@ -144,7 +184,7 @@ static double f_subplex(int n, const double *x, void *data_)
               return MY_INF;
 
      f = data->f(n, x, NULL, data->f_data);
-     return (isnan(f) ? MY_INF : f);
+     return (isnan(f) || nlopt_isinf(f) ? MY_INF : f);
 }
 
 #include "direct.h"
@@ -175,6 +215,7 @@ static double f_direct(int n, const double *x, int *undefined, void *data_)
 #include "mlsl.h"
 #include "mma.h"
 #include "cobyla.h"
+#include "newuoa.h"
 
 /*************************************************************************/
 
@@ -332,16 +373,8 @@ static nlopt_result nlopt_minimize_(
              int iret;
              double *scale = (double *) malloc(sizeof(double) * n);
              if (!scale) return NLOPT_OUT_OF_MEMORY;
-             for (i = 0; i < n; ++i) {
-                  if (!nlopt_isinf(ub[i]) && !nlopt_isinf(lb[i]))
-                       scale[i] = (ub[i] - lb[i]) * 0.01;
-                  else if (!nlopt_isinf(lb[i]) && x[i] > lb[i])
-                       scale[i] = (x[i] - lb[i]) * 0.01;
-                  else if (!nlopt_isinf(ub[i]) && x[i] < ub[i])
-                       scale[i] = (ub[i] - x[i]) * 0.01;
-                  else
-                       scale[i] = 0.01 * x[i] + 0.0001;
-             }
+             for (i = 0; i < n; ++i)
+                  scale[i] = initial_step(1, lb+i, ub+i, x+i);
              iret = nlopt_subplex(f_subplex, minf, x, n, &d, &stop, scale);
              free(scale);
              switch (iret) {
@@ -358,21 +391,10 @@ static nlopt_result nlopt_minimize_(
              break;
         }
 
-        case NLOPT_LN_PRAXIS: {
-             double h0 = HUGE_VAL;
-             for (i = 0; i < n; ++i) {
-                  if (!nlopt_isinf(ub[i]) && !nlopt_isinf(lb[i]))
-                       h0 = MIN(h0, (ub[i] - lb[i]) * 0.01);
-                  else if (!nlopt_isinf(lb[i]) && x[i] > lb[i])
-                       h0 = MIN(h0, (x[i] - lb[i]) * 0.01);
-                  else if (!nlopt_isinf(ub[i]) && x[i] < ub[i])
-                       h0 = MIN(h0, (ub[i] - x[i]) * 0.01);
-                  else
-                       h0 = MIN(h0, 0.01 * x[i] + 0.0001);
-             }
-             return praxis_(0.0, DBL_EPSILON, h0, n, x, f_subplex, &d,
+        case NLOPT_LN_PRAXIS:
+             return praxis_(0.0, DBL_EPSILON, 
+                            initial_step(n, lb, ub, x), n, x, f_subplex, &d,
                             &stop, minf);
-        }
 
 #ifdef WITH_NOCEDAL
         case NLOPT_LD_LBFGS_NOCEDAL: {
@@ -443,27 +465,16 @@ static nlopt_result nlopt_minimize_(
                                  lb, ub, x, minf, &stop,
                                  local_search_alg_deriv, 1e-12, 100000);
 
-        case NLOPT_LN_COBYLA: {
-             /* crude heuristics for initial step */
-             double rhobegin = HUGE_VAL;
-             for (i = 0; i < n; ++i) {
-                  if (!nlopt_isinf(ub[i]) && !nlopt_isinf(lb[i])
-                       && (ub[i] - lb[i]) * 0.25 < rhobegin && ub[i] != lb[i])
-                       rhobegin = (ub[i] - lb[i]) * 0.25;
-             }
-             if (nlopt_isinf(rhobegin)) {
-                  rhobegin = 0;
-                  for (i = 0; i < n; ++i)
-                       if (fabs(x[i]) * 0.25 > rhobegin)
-                            rhobegin = fabs(x[i]) * 0.25;
-                  if (rhobegin == 0)
-                       rhobegin = 1;
-             }
-             
+        case NLOPT_LN_COBYLA:
              return cobyla_minimize(n, f, f_data, 
                                     m, fc, fc_data, fc_datum_size,
-                                    lb, ub, x, minf, &stop, rhobegin);
-        }
+                                    lb, ub, x, minf, &stop,
+                                    initial_step(n, lb, ub, x));
+                                    
+        case NLOPT_LN_NEWUOA:
+             return newuoa(n, 2*n+1, x, initial_step(n, lb, ub, x),
+                           &stop, minf, f_subplex, &d);
+                                    
 
         default:
              return NLOPT_INVALID_ARGS;

diff --git a/api/nlopt.h b/api/nlopt.h
index a7871ff25932a53d033d5a57ea95ca7b7ccad808..0af459fe1f3db784bd1c6d3d83c6a64b91f47d07 100644 (file)
--- a/api/nlopt.h
+++ b/api/nlopt.h
@@ -87,6 +87,8 @@ typedef enum {
 
      NLOPT_LN_COBYLA,
 
+     NLOPT_LN_NEWUOA,
+
      NLOPT_NUM_ALGORITHMS /* not an algorithm, just the number of them */
 } nlopt_algorithm;
 

diff --git a/api/nlopt_minimize.3 b/api/nlopt_minimize.3
index 7dbb8a3f3adde6c4ed73756fa2e4ef03b33a1974..5aa75acd56d2972d2d1638fa35885645d51f54f0 100644 (file)
--- a/api/nlopt_minimize.3
+++ b/api/nlopt_minimize.3
@@ -328,12 +328,20 @@ function.
 Local (L) derivative-free (N) optimization using the COBYLA algorithm
 of Powell (Constrained Optimization BY Linear Approximations).
 The
-.B NLOPT_LD_COBYLA
+.B NLOPT_LN_COBYLA
 algorithm supports both bound-constrained and unconstrained optimization,
 and also supports an arbitrary number (\fIm\fR) of nonlinear constraints
 via the
 .BR nlopt_minimize_constrained ()
 function.
+.TP
+.B NLOPT_LN_NEWUOA
+Local (L) derivative-free (N) optimization using the NEWUOA algorithm
+of Powell, based on successive quadratic approximations of the objective
+function. The
+.B NLOPT_LN_NEWUOA
+algorithm was originally designed only for unconstrained optimization,
+and we only support bound constraints by an inefficient algorithm.
 .SH STOPPING CRITERIA
 Multiple stopping criteria for the optimization are supported, as
 specified by the following arguments to

diff --git a/api/nlopt_minimize_constrained.3 b/api/nlopt_minimize_constrained.3
index b2b21589a00d707af6b431e5b02cefe3c647dbbc..064ec50f9e7736aeb1f0def2a37605720a768d99 100644 (file)
--- a/api/nlopt_minimize_constrained.3
+++ b/api/nlopt_minimize_constrained.3
@@ -428,10 +428,18 @@ as described above.
 Local (L) derivative-free (N) optimization using the COBYLA algorithm
 of Powell (Constrained Optimization BY Linear Approximations).
 The
-.B NLOPT_LD_COBYLA
+.B NLOPT_LN_COBYLA
 algorithm supports both bound-constrained and unconstrained optimization,
 and also supports an arbitrary number (\fIm\fR) of nonlinear constraints
 as described above.
+.TP
+.B NLOPT_LN_NEWUOA
+Local (L) derivative-free (N) optimization using the NEWUOA algorithm
+of Powell, based on successive quadratic approximations of the objective
+function. The
+.B NLOPT_LN_NEWUOA
+algorithm was originally designed only for unconstrained optimization,
+and we only support bound constraints by an inefficient algorithm.
 .SH STOPPING CRITERIA
 Multiple stopping criteria for the optimization are supported, as
 specified by the following arguments to

diff --git a/configure.ac b/configure.ac
index ea895842d40e268d6f32140675ca4af3d7ca9304..2645ad227f27d1d11cf798730d87b82b8ac5f36f 100644 (file)
--- a/configure.ac
+++ b/configure.ac
@@ -212,6 +212,7 @@ AC_CONFIG_FILES([
    mlsl/Makefile
    mma/Makefile
    cobyla/Makefile
+   newuoa/Makefile
    test/Makefile
 ])
 

diff --git a/octave/Makefile.am b/octave/Makefile.am
index 55aeecf016bf2f8f7df4d83db7ce1c7d22b30c91..2b2fb43e9066e7a8c4d4519bff1d6d49f60d2ac5 100644 (file)
--- a/octave/Makefile.am
+++ b/octave/Makefile.am
@@ -1,6 +1,6 @@
 AM_CPPFLAGS = -I$(top_srcdir)/api 
 
-MFILES = NLOPT_GN_DIRECT.m NLOPT_GN_DIRECT_L.m NLOPT_GN_DIRECT_L_RAND.m NLOPT_GN_DIRECT_NOSCAL.m NLOPT_GN_DIRECT_L_NOSCAL.m NLOPT_GN_DIRECT_L_RAND_NOSCAL.m NLOPT_GN_ORIG_DIRECT.m NLOPT_GN_ORIG_DIRECT_L.m NLOPT_LN_SUBPLEX.m NLOPT_GD_STOGO.m NLOPT_GD_STOGO_RAND.m NLOPT_LD_LBFGS_NOCEDAL.m NLOPT_LD_LBFGS.m NLOPT_LN_PRAXIS.m NLOPT_LD_VAR1.m NLOPT_LD_VAR2.m NLOPT_LD_TNEWTON.m NLOPT_LD_TNEWTON_RESTART.m NLOPT_LD_TNEWTON_PRECOND.m NLOPT_LD_TNEWTON_PRECOND_RESTART.m NLOPT_GN_CRS2_LM.m NLOPT_GN_MLSL.m NLOPT_GD_MLSL.m NLOPT_GN_MLSL_LDS.m NLOPT_GD_MLSL_LDS.m NLOPT_LD_MMA.m NLOPT_LN_COBYLA.m 
+MFILES = NLOPT_GN_DIRECT.m NLOPT_GN_DIRECT_L.m NLOPT_GN_DIRECT_L_RAND.m NLOPT_GN_DIRECT_NOSCAL.m NLOPT_GN_DIRECT_L_NOSCAL.m NLOPT_GN_DIRECT_L_RAND_NOSCAL.m NLOPT_GN_ORIG_DIRECT.m NLOPT_GN_ORIG_DIRECT_L.m NLOPT_LN_SUBPLEX.m NLOPT_GD_STOGO.m NLOPT_GD_STOGO_RAND.m NLOPT_LD_LBFGS_NOCEDAL.m NLOPT_LD_LBFGS.m NLOPT_LN_PRAXIS.m NLOPT_LD_VAR1.m NLOPT_LD_VAR2.m NLOPT_LD_TNEWTON.m NLOPT_LD_TNEWTON_RESTART.m NLOPT_LD_TNEWTON_PRECOND.m NLOPT_LD_TNEWTON_PRECOND_RESTART.m NLOPT_GN_CRS2_LM.m NLOPT_GN_MLSL.m NLOPT_GD_MLSL.m NLOPT_GN_MLSL_LDS.m NLOPT_GD_MLSL_LDS.m NLOPT_LD_MMA.m NLOPT_LN_COBYLA.m NLOPT_LN_NEWUOA.m 
 
 #######################################################################
 octdir = $(OCT_INSTALL_DIR)

diff --git a/octave/NLOPT_LN_NEWUOA.m b/octave/NLOPT_LN_NEWUOA.m
new file mode 100644 (file)
index 0000000..8f5dbcd
--- /dev/null
+++ b/octave/NLOPT_LN_NEWUOA.m
@@ -0,0 +1,5 @@
+% NLOPT_LN_NEWUOA: NEWUOA unconstrained optimization via quadratic models (local, no-derivative)
+%
+% See nlopt_minimize for more information.
+function val = NLOPT_LN_NEWUOA
+  val = 27;