From f92619892167e8bf7d91bec6b617281fac49a28a Mon Sep 17 00:00:00 2001 From: stevenj Date: Sat, 30 Aug 2008 06:27:37 -0400 Subject: [PATCH] added NEWUOA, unified starting step-size for derivative-free algorithms darcs-hash:20080830102737-c8de0-2f48550a7c042db585110d7a05b566c02fc7b892.gz --- Makefile.am | 4 +- api/Makefile.am | 2 +- api/nlopt.c | 101 +++++++++++++++++-------------- api/nlopt.h | 2 + api/nlopt_minimize.3 | 10 ++- api/nlopt_minimize_constrained.3 | 10 ++- configure.ac | 1 + octave/Makefile.am | 2 +- octave/NLOPT_LN_NEWUOA.m | 5 ++ 9 files changed, 86 insertions(+), 51 deletions(-) create mode 100644 octave/NLOPT_LN_NEWUOA.m diff --git a/Makefile.am b/Makefile.am index d716f64..144845f 100644 --- 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 8fa22fd..1ec9af2 100644 --- 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 0d35d58..49fa534 100644 --- 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 a7871ff..0af459f 100644 --- 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 7dbb8a3..5aa75ac 100644 --- 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 b2b2158..064ec50 100644 --- 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 ea89584..2645ad2 100644 --- 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 55aeecf..2b2fb43 100644 --- 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 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; -- 2.30.2