#include <stdlib.h>
#include <math.h>
+#include <string.h>
#include "nlopt.h"
#include "nlopt-util.h"
typedef struct {
nlopt_func f;
void *f_data;
- const double *lb, *ub;
+ const double *lb, *ub, *x0;
+ double *xtmp;
} nlopt_data;
+#define RECENTER 1 /* 0 to disable recentering */
+
+/* for global-search algorithms that ignore the starting guess,
+ but always check the center of the search box, we perform a
+ coordinate transformation to put the initial guess x0 at the
+ center of the box, and store the transformed x in xtmp. */
+static void recenter_x(int n, const double *x,
+ const double *lb, const double *ub,
+ const double *x0, double *xtmp)
+{
+ int i;
+ for (i = 0; i < n; ++i) {
+#if RECENTER
+ /* Lagrange interpolating polynomial */
+ double xm = 0.5 * (lb[i] + ub[i]);
+ double dlu = 1. / (lb[i] - ub[i]);
+ double dlm = 1. / (lb[i] - xm);
+ double dum = 1. / (ub[i] - xm);
+ double dxu = x[i] - ub[i];
+ double dxl = x[i] - lb[i];
+ double dxm = x[i] - xm;
+ xtmp[i] = (lb[i] * (dxu * dlu) * (dxm * dlm)
+ - ub[i] * (dxl * dlu) * (dxm * dum)
+ + x0[i] * (dxl * dlm) * (dxu * dum));
+#else
+ xtmp[i] = x[i];
+#endif
+ }
+}
+
+/* transform grad from df/dxtmp to df/dx */
+static void recenter_grad(int n, const double *x,
+ const double *lb, const double *ub,
+ const double *x0,
+ double *grad)
+{
+#if RECENTER
+ if (grad) {
+ int i;
+ for (i = 0; i < n; ++i) {
+ double xm = 0.5 * (lb[i] + ub[i]);
+ double dlu = 1. / (lb[i] - ub[i]);
+ double dlm = 1. / (lb[i] - xm);
+ double dum = 1. / (ub[i] - xm);
+ double dxu = x[i] - ub[i];
+ double dxl = x[i] - lb[i];
+ double dxm = x[i] - xm;
+ grad[i] *= (lb[i] * dlu*dlm * (dxm + dxu)
+ - ub[i] * dum*dlu * (dxm + dxl)
+ + x0[i] * dlm*dum * (dxu + dxl));
+ }
+ }
+#endif
+}
+
+static double f_recenter(int n, const double *x, double *grad, void *data_)
+{
+ nlopt_data *data = (nlopt_data *) data_;
+ double f;
+ recenter_x(n, x, data->lb, data->ub, data->x0, data->xtmp);
+ f = data->f(n, data->xtmp, grad, data->f_data);
+ recenter_grad(n, x, data->lb, data->ub, data->x0, grad);
+ return f;
+}
+
#include "subplex.h"
static double f_subplex(int n, const double *x, void *data_)
static double f_direct(int n, const double *x, int *undefined, void *data_)
{
nlopt_data *data = (nlopt_data *) data_;
- double f = data->f(n, x, NULL, data->f_data);
+ double f;
+ recenter_x(n, x, data->lb, data->ub, data->x0, data->xtmp);
+ f = data->f(n, data->xtmp, NULL, data->f_data);
*undefined = isnan(f) || my_isinf(f);
return f;
}
#include "stogo.h"
+
#include "l-bfgs-b.h"
/*************************************************************************/
d.f_data = f_data;
d.lb = lb;
d.ub = ub;
+ d.x0 = d.xtmp = NULL;
/* make sure rand generator is inited */
if (!nlopt_srand_called)
switch (algorithm) {
case NLOPT_GLOBAL_DIRECT:
- case NLOPT_GLOBAL_DIRECT_L:
- switch (direct_optimize(f_direct, &d, n, lb, ub, x, fmin,
- maxeval, 500, ftol_rel, ftol_abs,
- xtol_rel, xtol_rel,
- DIRECT_UNKNOWN_FGLOBAL, -1.0,
- NULL,
- algorithm == NLOPT_GLOBAL_DIRECT
- ? DIRECT_ORIGINAL
- : DIRECT_GABLONSKY)) {
+ case NLOPT_GLOBAL_DIRECT_L: {
+ int iret;
+ d.xtmp = (double *) malloc(sizeof(double) * n*2);
+ if (!d.xtmp) return NLOPT_OUT_OF_MEMORY;
+ memcpy(d.xtmp + n, x, sizeof(double) * n); d.x0 = d.xtmp + n;
+ iret = direct_optimize(f_direct, &d, n, lb, ub, x, fmin,
+ maxeval, 500, ftol_rel, ftol_abs,
+ xtol_rel, xtol_rel,
+ DIRECT_UNKNOWN_FGLOBAL, -1.0,
+ NULL,
+ algorithm == NLOPT_GLOBAL_DIRECT
+ ? DIRECT_ORIGINAL
+ : DIRECT_GABLONSKY);
+ recenter_x(n, x, lb, ub, d.x0, x);
+ free(d.xtmp);
+ switch (iret) {
case DIRECT_INVALID_BOUNDS:
case DIRECT_MAXFEVAL_TOOBIG:
case DIRECT_INVALID_ARGS:
return NLOPT_OUT_OF_MEMORY;
}
break;
+ }
case NLOPT_GLOBAL_STOGO:
- case NLOPT_GLOBAL_STOGO_RANDOMIZED:
- if (!stogo_minimize(n, f, f_data, x, fmin, lb, ub,
- maxeval, maxtime,
- algorithm == NLOPT_GLOBAL_STOGO
- ? 0 : 2*n))
- return NLOPT_FAILURE;
+ case NLOPT_GLOBAL_STOGO_RANDOMIZED: {
+ int iret;
+ d.xtmp = (double *) malloc(sizeof(double) * n*2);
+ if (!d.xtmp) return NLOPT_OUT_OF_MEMORY;
+ memcpy(d.xtmp + n, x, sizeof(double) * n); d.x0 = d.xtmp + n;
+ iret = stogo_minimize(n, f_recenter, &d, x, fmin, lb, ub,
+ maxeval, maxtime,
+ algorithm == NLOPT_GLOBAL_STOGO
+ ? 0 : 2*n);
+ recenter_x(n, x, lb, ub, d.x0, x);
+ free(d.xtmp);
+ if (!iret) return NLOPT_FAILURE;
break;
+ }
case NLOPT_LOCAL_SUBPLEX: {
int iret;
.IR algorithm .
The minimum function value found is returned in
.IR fmin ,
-with the corresponding design variable values stored in the array
+with the corresponding design variable values returned in the array
.I x
of length
.IR n .
+The input values in
+.I x
+should be a starting guess for the optimum.
The inputs
.I lb
and
.IR f ,
and other algorithms do not require derivatives. Some of the
algorithms attempt to find a global minimum within the given bounds,
-and others find only a local minimum (with the initial value of
-.I x
-as a starting guess).
+and others find only a local minimum.
.PP
The
.B nlopt_minimize
RVector m(dim), x(dim);
if (det_pnts>0) {
- // Add midpoint
box.Midpoint(m) ;
- tmpTrial.xvals=m ; tmpTrial.objval=DBL_MAX ;
- SampleBox.AddTrial(tmpTrial) ;
+ tmpTrial.objval=DBL_MAX ;
// Add the rest
i=1 ; flag=1 ; dir=0 ;
x=m ;
}
i++ ;
}
+ // Add midpoint
+ tmpTrial.xvals=m ;
+ SampleBox.AddTrial(tmpTrial) ;
}
}
double maxgrad=0 ;
// Create sampling points
- FillRegular(SampleBox, box);
FillRandom(SampleBox, box);
+ FillRegular(SampleBox, box);
// Perform the actual sampling
while ( !SampleBox.EmptyBox() ) {
static nlopt_algorithm algorithm = NLOPT_GLOBAL_DIRECT;
static double ftol_rel = 0, ftol_abs = 0, xtol_rel = 0, xtol_abs = 0, fmin_max = -HUGE_VAL;
-static int maxeval = 1000;
+static int maxeval = 1000, iterations = 1;
static double maxtime = 0.0;
+static double xinit_tol = -1;
static void listalgs(FILE *f)
{
static int test_function(int ifunc)
{
testfunc func;
- int i;
+ int i, iter;
double *x, fmin, f0, *xtabs;
nlopt_result ret;
double start = nlopt_seconds();
printf("-----------------------------------------------------------\n");
printf("Optimizing %s (%d dims) using %s algorithm\n",
func.name, func.n, nlopt_algorithm_name(algorithm));
+ printf("lower bounds at lb = [");
+ for (i = 0; i < func.n; ++i) printf(" %g", func.lb[i]);
+ printf("]\n");
+ printf("upper bounds at ub = [");
+ for (i = 0; i < func.n; ++i) printf(" %g", func.ub[i]);
+ printf("]\n");
+
+
printf("Starting guess x = [");
- for (i = 0; i < func.n; ++i)
- printf(" %g", x[i] = nlopt_urand(func.lb[i], func.ub[i]));
+ for (i = 0; i < func.n; ++i) {
+ if (xinit_tol < 0) { /* random starting point near center of box */
+ double dx = (func.ub[i] - func.lb[i]) * 0.25;
+ double xm = 0.5 * (func.ub[i] + func.lb[i]);
+ x[i] = nlopt_urand(xm - dx, xm + dx);
+ }
+ else
+ x[i] = nlopt_urand(-xinit_tol, xinit_tol)
+ + (1 + nlopt_urand(-xinit_tol, xinit_tol)) * func.xmin[i];
+ printf(" %g", x[i]);
+ }
printf("]\n");
f0 = func.f(func.n, x, x + func.n, func.f_data);
printf("Starting function value = %g\n", f0);
}
}
- testfuncs_counter = 0;
- ret = nlopt_minimize(algorithm,
- func.n, func.f, func.f_data,
- func.lb, func.ub,
- x, &fmin,
- fmin_max, ftol_rel, ftol_abs, xtol_rel, xtabs,
- maxeval, maxtime);
- printf("finished after %g seconds.\n", nlopt_seconds() - start);
- printf("return code %d from nlopt_minimize\n", ret);
- if (ret < 0) {
- fprintf(stderr, "testopt: error in nlopt_minimize\n");
- return 0;
+ for (iter = 0; iter < iterations; ++iter) {
+ testfuncs_counter = 0;
+ ret = nlopt_minimize(algorithm,
+ func.n, func.f, func.f_data,
+ func.lb, func.ub,
+ x, &fmin,
+ fmin_max, ftol_rel, ftol_abs, xtol_rel, xtabs,
+ maxeval, maxtime);
+ printf("finished after %g seconds.\n", nlopt_seconds() - start);
+ printf("return code %d from nlopt_minimize\n", ret);
+ if (ret < 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("|f - fmin| = %g, |f - fmin| / |fmin| = %e\n",
+ fabs(fmin - func.fmin), fabs(fmin - func.fmin) / fabs(func.fmin));
}
- 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;
" -h : print this help\n"
" -L : list available algorithms and objective functions\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"
+ " -0 <x> : starting guess within <x> + (1+<x>) * optimum\n"
" -e <n> : use at most <n> evals (default: %d, 0 to disable)\n"
" -t <t> : use at most <t> seconds (default: disabled)\n"
" -x <t> : relative tolerance <t> on x (default: disabled)\n"
" -f <t> : relative tolerance <t> on f (default: disabled)\n"
" -F <t> : absolute tolerance <t> on f (default: disabled)\n"
" -m <m> : minimize f until <m> is reached (default: disabled)\n"
+ " -i <n> : iterate optimization <n> times (default: 1)\n"
+ " -r <s> : use random seed <s> for starting guesses\n"
, maxeval);
}
if (argc <= 1)
usage(stdout);
- while ((c = getopt(argc, argv, "hLvr:a:o:e:t:x:X:f:F:m:")) != -1)
+ while ((c = getopt(argc, argv, "hLv0:r:a:o:i:e:t:x:X:f:F:m:")) != -1)
switch (c) {
case 'h':
usage(stdout);
case 'e':
maxeval = atoi(optarg);
break;
+ case 'i':
+ iterations = atoi(optarg);
+ break;
case 't':
maxtime = atof(optarg);
break;
case 'm':
fmin_max = atof(optarg);
break;
+ case '0':
+ xinit_tol = atof(optarg);
+ break;
default:
fprintf(stderr, "harminv: invalid argument -%c\n", c);
usage(stderr);
~ianmdlvl / nlopt.git / commitdiff