ACLOCAL_AMFLAGS=-I ./m4
-SUBDIRS=lbfgs subplex direct stogo api
+SUBDIRS=lbfgs subplex direct stogo api . test
EXTRA_DIST=COPYRIGHT autogen.sh nlopt.pc.in m4
libnlopt_la_SOURCES =
#include "nlopt.h"
#include "config.h"
+static const char nlopt_algorithm_names[NLOPT_NUM_ALGORITHMS][128] = {
+ "DIRECT (global)",
+ "Subplex (local)",
+ "StoGO (global)",
+ "Low-storage BFGS (LBFGS) (local)"
+};
+
+const char *nlopt_algorithm_name(nlopt_algorithm a)
+{
+ if (a < 0 || a >= NLOPT_NUM_ALGORITHMS) return "UNKNOWN";
+ return nlopt_algorithm_names[a];
+}
+
static int my_isinf(double x) {
return x == HUGE_VAL
#ifdef HAVE_ISINF
#define MIN(a,b) ((a) < (b) ? (a) : (b))
nlopt_result nlopt_minimize(
- nlopt_method method,
+ nlopt_algorithm algorithm,
int n, nlopt_func f, void *f_data,
const double *lb, const double *ub, /* bounds */
double *x, /* in: initial guess, out: minimizer */
d.f = f;
d.f_data = f_data;
- switch (method) {
+ switch (algorithm) {
case NLOPT_GLOBAL_DIRECT:
switch (direct_optimize(f_direct, &d, n, lb, ub, x, fmin,
maxeval, 500, ftol_rel, ftol_abs,
void *func_data);
typedef enum {
- /* non-gradient methods */
- NLOPT_GLOBAL_DIRECT,
+ /* non-gradient algorithms */
+ NLOPT_GLOBAL_DIRECT = 0,
NLOPT_LOCAL_SUBPLEX,
- /* gradient-based methods */
+ /* gradient-based algorithms */
NLOPT_GLOBAL_STOGO,
- NLOPT_LOCAL_LBFGS
-} nlopt_method;
+ NLOPT_LOCAL_LBFGS,
+
+ NLOPT_NUM_ALGORITHMS /* not an algorithm, just the number of them */
+} nlopt_algorithm;
+
+extern const char *nlopt_algorithm_name(nlopt_algorithm a);
typedef enum {
NLOPT_FAILURE = -1, /* generic failure code */
} nlopt_result;
extern nlopt_result nlopt_minimize(
- nlopt_method method,
+ nlopt_algorithm algorithm,
int n, nlopt_func f, void *f_data,
const double *lb, const double *ub, /* bounds */
double *x, /* in: initial guess, out: minimizer */
dnl Checks for typedefs, structures, and compiler characteristics.
AC_HEADER_STDC
AC_HEADER_TIME
+AC_CHECK_HEADERS([unistd.h getopt.h])
AC_C_CONST
AC_C_INLINE
stogo/Makefile
lbfgs/Makefile
subplex/Makefile
+ test/Makefile
])
AC_OUTPUT
double MacEpsilon ;
int FC=0, GC=0 ;
+int stogo_verbose = 0; /* set to nonzero for verbose output */
+
Global::Global(RTBox D, Pobj o, Pgrad g, GlobalParams P): Domain(D) {
dim=Domain.GetDim();
box.AddTrial(tmpTrial) ;
if (tmpTrial.objval<=fbound+mu && tmpTrial.objval<=box.fmin+mu) {
- cout << "Found a candidate, x=" << tmpTrial.xvals;
- cout << " F=" <<tmpTrial.objval << " FC=" << FC << endl;
+ if (stogo_verbose) {
+ cout << "Found a candidate, x=" << tmpTrial.xvals;
+ cout << " F=" <<tmpTrial.objval << " FC=" << FC << endl;
+ }
SolSet.push_back(tmpTrial);
}
#ifdef GS_DEBUG
MacEpsilon=eps(); // Get machine precision
if (det_pnts>2*dim+1) {
det_pnts=2*dim+1;
- cout << "Warning: Reducing det_pnts to " << det_pnts << endl;
+ if (stogo_verbose)
+ cout << "Warning: Reducing det_pnts to " << det_pnts << endl;
}
// Initialize timer
if (!InTime()) {
done=TRUE;
- cout << "The program has run out of time or function evaluations\n";
+ if (stogo_verbose)
+ cout << "The program has run out of time or function evaluations\n";
break;
}
} // inner while-loop
- cout << endl << "*** Inner loop completed ***" << endl ;
+ if (stogo_verbose)
+ cout << endl << "*** Inner loop completed ***" << endl ;
// Reduce SolSet if necessary
SolSet.erase(remove_if(SolSet.begin(), SolSet.end(),
TrialGT(fbound+mu)),SolSet.end());
if (InTime()) {
- cout << "Current set of minimizers (" << SolSet.size() << ")" << endl ;
- DispMinimizers() ;
+ if (stogo_verbose) {
+ cout << "Current set of minimizers (" << SolSet.size() << ")" << endl ;
+ DispMinimizers() ;
+ }
while (!Garbage.empty()) {
box=Garbage.top() ;
}
} // Outer while-loop
- cout << "Number of outer iterations : " << outer_iter << endl;
- cout << "Number of unexplored boxes : " << CandSet.size() << endl;
- cout << "Number of boxes in garbage : " << Garbage.size() << endl;
- cout << "Number of elements in SolSet : " << SolSet.size() << endl;
- cout << "Number of function evaluations : " << FC << endl;
- cout << "Number of gradient evaluations : " << GC << endl;
+ if (stogo_verbose) {
+ cout << "Number of outer iterations : " << outer_iter << endl;
+ cout << "Number of unexplored boxes : " << CandSet.size() << endl;
+ cout << "Number of boxes in garbage : " << Garbage.size() << endl;
+ cout << "Number of elements in SolSet : " << SolSet.size() << endl;
+ cout << "Number of function evaluations : " << FC << endl;
+ cout << "Number of gradient evaluations : " << GC << endl;
+ }
if (axis != -1) {
// Return minimizer when doing the AV method
}
/************* Various utility functions ****************/
-long int Global::GetTime()
+double Global::GetTime()
{
time_t ctime; time(&ctime);
- return (long int)difftime(ctime,StartTime);
+ return difftime(ctime,StartTime);
}
bool Global::InTime()
{
- return (!maxtime || GetTime()<maxtime) && (!maxeval || numeval<maxeval);
+ return (maxtime <= 0.0 || GetTime()<maxtime) && (!maxeval || numeval<maxeval);
}
double Global::GetMinValue() {
#include "tools.h"
using namespace std;
+extern "C" int stogo_verbose;
+
typedef void dom(RTBox) ;
typedef dom* Pdom ;
class GlobalParams {
public:
- long int maxtime, maxeval;
+ double maxtime;
+ long int maxeval;
double eps_cl, mu, rshift;
int det_pnts, rnd_pnts;
};
void ClearSolSet();
void AddPoint(RCRVector, double);
- long int GetTime();
+ double GetTime();
bool InTime();
private:
objective_func fgrad, void *data,
double *x, double *fmin,
const double *l, const double *u,
- long int maxeval, long int maxtime)
+ long int maxeval, double maxtime)
{
GlobalParams params;
objective_func fgrad, void *data,
double *x, double *fmin,
const double *l, const double *u,
- long int maxeval, long int maxtime);
+ long int maxeval, double maxtime);
+
+extern int stogo_verbose; /* set to nonzero for verbose output */
#ifdef __cplusplus
}
--- /dev/null
+AM_CPPFLAGS = -I$(top_srcdir)/api
+
+noinst_PROGRAMS = testopt
+
+testopt_SOURCES = testfuncs.c testfuncs.h testopt.cpp
+testopt_LDADD = $(top_builddir)/libnlopt.la
--- /dev/null
+#include <stdio.h>
+#include <math.h>
+
+#include "testfuncs.h"
+#include "config.h"
+
+#define UNUSED(x) (void) x
+
+static double sqr(double x) { return x * x; }
+
+int testfuncs_verbose = 0;
+int testfuncs_counter = 0;
+
+static double testfuncs_status(int n, const double *x, double f)
+{
+ ++testfuncs_counter;
+ if (testfuncs_verbose) {
+ int i;
+ printf("f_%d (%g", testfuncs_counter, x[0]);
+ for (i = 1; i < n; ++i) printf(", %g", x[i]);
+ printf(") = %g\n", f);
+ }
+ return f;
+}
+
+#define RETURN(f) return testfuncs_status(n, x, f);
+
+/****************************************************************************/
+static double rosenbrock_f(int n, const double *x, double *grad, void *data)
+{
+ double a = x[1] - x[0] * x[0], b = 1 - x[0];
+ UNUSED(data);
+ if (grad) {
+ grad[0] = -400 * a * x[0] - 2*b;
+ grad[1] = -200 * a;
+ }
+ RETURN(100 * sqr(a) + sqr(b));
+}
+
+static const double rosenbrock_lb[2] = {-2, -2};
+static const double rosenbrock_ub[2] = {2, 2};
+static const double rosenbrock_xmin[2] = {1, 1};
+
+/****************************************************************************/
+static double mccormic_f(int n, const double *x, double *grad, void *data)
+{
+ double a = x[0] + x[1], b = x[0] - x[1];
+ UNUSED(data);
+ if (grad) {
+ grad[0] = cos(a) + 2*b - 1.5;
+ grad[1] = cos(a) - 2*b + 2.5;
+ }
+ RETURN(sin(a) + sqr(b) - 1.5*x[0] + 2.5*x[1] + 1);
+}
+
+static const double mccormic_lb[2] = {-1.5, -3};
+static const double mccormic_ub[2] = {4, 4};
+static const double mccormic_xmin[2] = {-0.54719, 1.54719};
+
+/****************************************************************************/
+static double boxbetts_f(int n, const double *x, double *grad, void *data)
+{
+ int i;
+ double f = 0;
+ UNUSED(data);
+ if (grad)
+ grad[0] = grad[1] = grad[2] = 0;
+ for (i = 1; i <= 10; ++i) {
+ double e0 = exp(-0.1*i*x[0]);
+ double e1 = exp(-0.1*i*x[1]);
+ double e2 = exp(-0.1*i) - exp((double) -i);
+ double g = e0 - e1 - e2 * x[2];
+ f += sqr(g);
+ if (grad) {
+ grad[0] += (2 * g) * (-0.1*i*e0);
+ grad[1] += (2 * g) * (0.1*i*e1);
+ grad[2] += -(2 * g) * e2;
+ }
+ }
+ RETURN(f);
+}
+
+static const double boxbetts_lb[3] = {0.9,9,0.9};
+static const double boxbetts_ub[3] = {1.2,11.2,1.2};
+static const double boxbetts_xmin[3] = {1,10,1};
+
+/****************************************************************************/
+static double paviani_f(int n, const double *x, double *grad, void *data)
+{
+ int i;
+ double f = 0, prod = 1;
+ UNUSED(data);
+ if (grad) for (i = 0; i < 10; ++i) grad[i] = 0;
+ for (i = 0; i < 10; ++i) {
+ double ln1 = log(x[i] - 2);
+ double ln2 = log(10 - x[i]);
+ f += sqr(ln1) + sqr(ln2);
+ if (grad)
+ grad[i] += 2 * ln1 / (x[i] - 2) - 2 * ln2 / (10 - x[i]);
+ prod *= x[i];
+ }
+ f -= (prod = pow(prod, 0.2));
+ if (grad)
+ for (i = 0; i < 10; ++i)
+ grad[i] -= 0.2 * prod / x[i];
+ RETURN(f);
+}
+
+static const double paviani_lb[10] = {2.001,2.001,2.001,2.001,2.001,2.001,2.001,2.001,2.001,2.001};
+static const double paviani_ub[10] = {9.999,9.999,9.999,9.999,9.999,9.999,9.999,9.999,9.999,9.999};
+static const double paviani_xmin[10] = {9.350266,9.350266,9.350266,9.350266,9.350266,9.350266,9.350266,9.350266,9.350266,9.350266};
+
+/****************************************************************************/
+static double grosenbrock_f(int n, const double *x, double *grad, void *data)
+{
+ int i;
+ double f = 0;
+ UNUSED(data);
+ if (grad) grad[0] = 0;
+ for (i = 0; i < 29; ++i) {
+ double a = x[i+1] - x[i] * x[i], b = 1 - x[i];
+ if (grad) {
+ grad[i] += -400 * a * x[i] - 2*b;
+ grad[i+1] = -200 * a;
+ }
+ f += 100 * sqr(a) + sqr(b);
+ }
+ RETURN(f);
+}
+
+static const double grosenbrock_lb[30] = {-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30};
+static const double grosenbrock_ub[30] = {30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30};
+static const double grosenbrock_xmin[30] = {1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1};
+
+/****************************************************************************/
+static double goldsteinprice_f(int n, const double *x, double *grad, void *data)
+{
+ double x0, x1, a1, a12, a2, b1, b12, b2;
+ UNUSED(data);
+ x0 = x[0]; x1 = x[1];
+ a1 = x0+x1+1; a12 = sqr(a1);
+ a2 = 19 - 14*x0 + 3*x0*x0 - 14*x1 + 6*x0*x1 + 3*x1*x1;
+ b1 = 2*x0-3*x1; b12 = sqr(b1);
+ b2 = 18 - 32*x0 + 12*x0*x0 + 48*x1 - 36*x0*x1 + 27*x1*x1;
+ if (grad) {
+ grad[0] = (1 + a12 * a2) * (2 * b1 * 2 * b2
+ + b12 * (-32 + 24*x0 - 36*x1))
+ + (2 * a1 * a2 + a12 * (-14 + 6*x0 + 6*x1)) * (30 + b12 * b2);
+ grad[1] = (1 + a12 * a2) * (2 * b1 * (-3) * b2
+ + b12 * (48 - 36*x0 + 54 * x1))
+ + (2 * a1 * a2 + a12 * (-14 + 6*x0 + 6*x1)) * (30 + b12 * b2);
+ }
+ RETURN((1 + a12 * a2) * (30 + b12 * b2));
+}
+
+static const double goldsteinprice_lb[2] = {-2, -2};
+static const double goldsteinprice_ub[2] = {2, 2};
+static const double goldsteinprice_xmin[2] = {0, -1};
+
+/****************************************************************************/
+/****************************************************************************/
+
+const testfunc testfuncs[NTESTFUNCS] = {
+ { rosenbrock_f, NULL, 1, 2,
+ rosenbrock_lb, rosenbrock_ub, rosenbrock_xmin,
+ 0.0, "Rosenbrock function" },
+ { mccormic_f, NULL, 1, 2,
+ mccormic_lb, mccormic_ub, mccormic_xmin,
+ -1.9133, "McCormic function" },
+ { boxbetts_f, NULL, 1, 3,
+ boxbetts_lb, boxbetts_ub, boxbetts_xmin,
+ 0.0, "Box and Betts exponential quadratic sum" },
+ { paviani_f, NULL, 1, 10,
+ paviani_lb, paviani_ub, paviani_xmin,
+ -45.778470, "Paviani function" },
+ { grosenbrock_f, NULL, 1, 30,
+ grosenbrock_lb, grosenbrock_ub, grosenbrock_xmin,
+ 0.0, "Generalized Rosenbrock function" },
+ { goldsteinprice_f, NULL, 1, 2,
+ goldsteinprice_lb, goldsteinprice_ub, goldsteinprice_xmin,
+ 3.0, "Goldstein and Price function" }
+};
--- /dev/null
+#ifndef TESTFUNCS_H
+#define TESTFUNCS_H
+
+#ifdef __cplusplus
+extern "C"
+{
+#endif /* __cplusplus */
+
+#include "nlopt.h"
+
+typedef struct {
+ nlopt_func f;
+ void *f_data;
+ int has_gradient;
+ int n;
+ const double *lb, *ub, *xmin;
+ double fmin;
+ const char *name;
+} testfunc;
+
+#define NTESTFUNCS 6
+extern const testfunc testfuncs[NTESTFUNCS];
+
+extern int testfuncs_verbose;
+extern int testfuncs_counter;
+
+#ifdef __cplusplus
+} /* extern "C" */
+#endif /* __cplusplus */
+
+#endif
--- /dev/null
+#include <stdlib.h>
+#include <stdio.h>
+#include <math.h>
+
+#include "config.h"
+
+#ifdef HAVE_UNISTD_H
+# include <unistd.h>
+#endif
+#ifdef HAVE_GETOPT_H
+# include <getopt.h>
+#endif
+#if TIME_WITH_SYS_TIME
+# include <sys/time.h>
+# include <time.h>
+#else
+# if HAVE_SYS_TIME_H
+# include <sys/time.h>
+# else
+# include <time.h>
+# endif
+#endif
+
+#include "nlopt.h"
+#include "testfuncs.h"
+
+static double urand(double a, double b)
+{
+ return a + (rand() * (b - a) / RAND_MAX);
+}
+
+static nlopt_algorithm algorithm = NLOPT_GLOBAL_DIRECT;
+static double ftol_rel = 0, ftol_abs = 0, xtol_rel = 0;
+static int maxeval = 1000;
+static double maxtime = 0.0;
+
+static int test_function(int ifunc)
+{
+ testfunc func;
+ int i;
+ double *x, fmin, f0;
+
+ if (ifunc < 0 || ifunc >= NTESTFUNCS) {
+ fprintf(stderr, "testopt: invalid function %d\n", ifunc);
+ return 0;
+ }
+ func = testfuncs[ifunc];
+ if (!func.has_gradient && algorithm >= NLOPT_GLOBAL_STOGO) {
+ fprintf(stderr,
+ "testopt: A function with gradients is required for %s\n",
+ nlopt_algorithm_name(algorithm));
+ return 0;
+ }
+ x = (double *) malloc(sizeof(double) * func.n * 2);
+ if (!x) { fprintf(stderr, "testopt: Out of memory!\n"); return 0; }
+
+
+ printf("-----------------------------------------------------------\n");
+ printf("Optimizing %s (%d dims) using %s algorithm\n",
+ func.name, func.n, nlopt_algorithm_name(algorithm));
+ printf("Starting guess x = [");
+ for (i = 0; i < func.n; ++i)
+ printf(" %g", x[i] = urand(func.lb[i], func.ub[i]));
+ printf("]\n");
+ f0 = func.f(func.n, x, x + func.n, func.f_data);
+ printf("Starting function value = %g\n", f0);
+
+ if (testfuncs_verbose && func.has_gradient) {
+ printf("checking gradient:\n");
+ for (i = 0; i < func.n; ++i) {
+ double f;
+ x[i] *= 1 + 1e-6;
+ f = func.f(func.n, x, NULL, func.f_data);
+ x[i] /= 1 + 1e-6;
+ printf(" grad[%d] = %g vs. numerical derivative %g\n",
+ i, x[i + func.n], (f - f0) / (x[i] * 1e-6));
+ }
+ }
+
+ testfuncs_counter = 0;
+ if (nlopt_minimize(algorithm,
+ func.n, func.f, func.f_data,
+ func.lb, func.ub,
+ x, &fmin,
+ HUGE_VAL, ftol_rel, ftol_abs, xtol_rel, NULL,
+ maxeval, maxtime) < 0) {
+ fprintf(stderr, "testopt: error in nlopt_minimize\n");
+ return 0;
+ }
+ printf("Found minimum f = %g after %d evaluations.\n",
+ fmin, testfuncs_counter);
+ printf("Minimum at x = [");
+ for (i = 0; i < func.n; ++i) printf(" %g", x[i]);
+ printf("]\n");
+ printf("vs. global minimum f = %g at x = [", func.fmin);
+ for (i = 0; i < func.n; ++i) printf(" %g", func.xmin[i]);
+ printf("]\n");
+ printf("|f - fmin| = %g, |f - fmin| / |fmin| = %e\n",
+ fabs(fmin - func.fmin), fabs(fmin - func.fmin) / fabs(func.fmin));
+
+ free(x);
+ return 1;
+}
+
+static void usage(FILE *f)
+{
+ fprintf(f, "Usage: testopt [OPTIONS]\n"
+ "Options:\n"
+ " -h : print this help\n"
+ " -v : verbose mode\n"
+ " -r <s> : use random seed <s> for starting guesses\n"
+ " -a <n> : use optimization algorithm <n>\n"
+ " -o <n> : use objective function <n>\n"
+ " -e <n> : use at most <n> evals (default: %d)\n",
+ maxeval);
+}
+
+int main(int argc, char **argv)
+{
+ int c;
+
+ srand((unsigned) time(NULL));
+ testfuncs_verbose = 0;
+
+ while ((c = getopt(argc, argv, "hvra:o:e:")) != -1)
+ switch (c) {
+ case 'h':
+ usage(stdout);
+ return EXIT_SUCCESS;
+ case 'v':
+ testfuncs_verbose = 1;
+ break;
+ case 'r':
+ srand((unsigned) atoi(optarg));
+ break;
+ case 'a':
+ c = atoi(optarg);
+ if (c < 0 || c >= NLOPT_NUM_ALGORITHMS) {
+ fprintf(stderr, "testopt: invalid algorithm %d\n", c);
+ return EXIT_FAILURE;
+ }
+ algorithm = (nlopt_algorithm) c;
+ break;
+ case 'o':
+ if (!test_function(atoi(optarg)))
+ return EXIT_FAILURE;
+ break;
+ case 'e':
+ maxeval = atoi(optarg);
+ break;
+ default:
+ fprintf(stderr, "harminv: invalid argument -%c\n", c);
+ usage(stderr);
+ return EXIT_FAILURE;
+ }
+
+ return EXIT_SUCCESS;
+}