chiark / gitweb /
added test program and some test objectives
authorstevenj <stevenj@alum.mit.edu>
Thu, 23 Aug 2007 20:17:09 +0000 (16:17 -0400)
committerstevenj <stevenj@alum.mit.edu>
Thu, 23 Aug 2007 20:17:09 +0000 (16:17 -0400)
darcs-hash:20070823201709-c8de0-e8894aceee594f8bd61b5d13d85c7ae5d0692dfa.gz

12 files changed:
Makefile.am
api/nlopt.c
api/nlopt.h
configure.ac
stogo/global.cc
stogo/global.h
stogo/stogo.cc
stogo/stogo.h
test/Makefile.am [new file with mode: 0644]
test/testfuncs.c [new file with mode: 0644]
test/testfuncs.h [new file with mode: 0644]
test/testopt.cpp [new file with mode: 0644]

index 817c623043b39d7cf950eaabb725f05e55b904f2..0fa9377287f11a0827e85c2e22cffd90c4b8db6b 100644 (file)
@@ -3,7 +3,7 @@ lib_LTLIBRARIES = libnlopt.la
 
 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 = 
index e63e3290a1f571ee8c5c05e4706de60601d6732a..ef3acf959e047088621ac7dede46a8598fcd2e72 100644 (file)
@@ -4,6 +4,19 @@
 #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
@@ -46,7 +59,7 @@ static double f_direct(int n, const double *x, int *undefined, void *data_)
 #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 */
@@ -59,7 +72,7 @@ nlopt_result nlopt_minimize(
      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,
index f8a8c30da98fc491428985c1de83c10689c7cc09..60c21a4cf3abfd388c24d818a950f4541c6e93ab 100644 (file)
@@ -11,14 +11,18 @@ typedef double (*nlopt_func)(int n, const double *x,
                             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 */
@@ -34,7 +38,7 @@ typedef enum {
 } 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 */
index 14580de0c2ba21cdfb1afe62b7a8d5a0e0beffa3..92fcebee11cd1293b8e497c1fd40bb36c6182016 100644 (file)
@@ -22,6 +22,7 @@ AC_PROG_LIBTOOL
 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
 
@@ -82,6 +83,7 @@ AC_CONFIG_FILES([
    stogo/Makefile
    lbfgs/Makefile
    subplex/Makefile
+   test/Makefile
 ])
 
 AC_OUTPUT
index 94963584b28ca4315dd8d957aa5ed198515204e3..80decbe876565b36368d8689f74720418c8a209d 100644 (file)
@@ -26,6 +26,8 @@ time_t   StartTime;
 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();
@@ -121,8 +123,10 @@ double Global::NewtonTest(RTBox box, int axis, RCRVector x_av, int *noutside) {
       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
@@ -187,7 +191,8 @@ void Global::Search(int axis, RCRVector x_av){
   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
@@ -221,19 +226,23 @@ void Global::Search(int axis, RCRVector x_av){
 
       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() ;
@@ -247,12 +256,14 @@ void Global::Search(int axis, RCRVector x_av){
     }
   } // 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
@@ -262,15 +273,15 @@ void Global::Search(int axis, RCRVector x_av){
 }
 
 /************* 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() {
index 26a6ab85a700e32c64a33848fbfd90d6f1cc727c..e9d8814c035e79fd901363b1bde03b5b16cd5555 100644 (file)
@@ -6,6 +6,8 @@
 #include "tools.h"
 using namespace std;
 
+extern "C" int stogo_verbose;
+
 typedef void dom(RTBox) ;
 typedef dom* Pdom ;
 
@@ -22,7 +24,8 @@ typedef objgrad* Pobjgrad ;
 
 class GlobalParams {
 public:
-  long int maxtime, maxeval;
+  double maxtime;
+  long int maxeval;
   double eps_cl, mu, rshift;
   int det_pnts, rnd_pnts;
 };
@@ -62,7 +65,7 @@ public:
   void ClearSolSet();
   void AddPoint(RCRVector, double);
 
-  long int GetTime();
+  double GetTime();
   bool InTime();
 
 private:
index eb9e73f3dee9825a5514f3a31193e8b40ec16082..97b0db352d38dea968607a8bfbe5f3939154d9b2 100644 (file)
@@ -30,7 +30,7 @@ int stogo_minimize(int n,
                   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;
 
index 2e395802fe1f1d7f1c8e03a10c6a28036b9e7f3c..4ae727e14c7afcda1c4271648ac2c1e072e86f0f 100644 (file)
@@ -51,7 +51,9 @@ int stogo_minimize(int n,
                    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
 }
diff --git a/test/Makefile.am b/test/Makefile.am
new file mode 100644 (file)
index 0000000..a7895c5
--- /dev/null
@@ -0,0 +1,6 @@
+AM_CPPFLAGS = -I$(top_srcdir)/api
+
+noinst_PROGRAMS = testopt
+
+testopt_SOURCES = testfuncs.c testfuncs.h testopt.cpp
+testopt_LDADD = $(top_builddir)/libnlopt.la
diff --git a/test/testfuncs.c b/test/testfuncs.c
new file mode 100644 (file)
index 0000000..39d6d79
--- /dev/null
@@ -0,0 +1,182 @@
+#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" }
+};
diff --git a/test/testfuncs.h b/test/testfuncs.h
new file mode 100644 (file)
index 0000000..df81446
--- /dev/null
@@ -0,0 +1,31 @@
+#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
diff --git a/test/testopt.cpp b/test/testopt.cpp
new file mode 100644 (file)
index 0000000..e2921ba
--- /dev/null
@@ -0,0 +1,158 @@
+#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;
+}