/*************************************************************************/
-nlopt_result nlopt_minimize(
+/* same as nlopt_minimize, but xtol_abs is required to be non-NULL */
+static nlopt_result nlopt_minimize_(
nlopt_algorithm algorithm,
int n, nlopt_func f, void *f_data,
const double *lb, const double *ub, /* bounds */
{
int i;
nlopt_data d;
+ nlopt_stopping stop;
+
d.f = f;
d.f_data = f_data;
d.lb = lb;
if (lb[i] > ub[i] || x[i] < lb[i] || x[i] > ub[i])
return NLOPT_INVALID_ARGS;
+ stop.n = n;
+ stop.fmin_max = (isnan(fmin_max) || (my_isinf(fmin_max) && fmin_max < 0))
+ ? -MY_INF : fmin_max;
+ stop.ftol_rel = ftol_rel;
+ stop.ftol_abs = ftol_abs;
+ stop.xtol_rel = xtol_rel;
+ stop.xtol_abs = xtol_abs;
+ stop.nevals = 0;
+ stop.maxeval = maxeval;
+ stop.maxtime = maxtime;
+ stop.start = nlopt_seconds();
+
switch (algorithm) {
case NLOPT_GLOBAL_DIRECT:
case NLOPT_GLOBAL_DIRECT_L:
else
scale[i] = 0.01 * x[i] + 0.0001;
}
- iret = subplex(f_subplex, fmin, x, n, &d, xtol_rel, maxeval,
- fmin_max, !my_isinf(fmin_max), scale);
+ iret = subplex(f_subplex, fmin, x, n, &d, &stop, scale);
free(scale);
switch (iret) {
case -2: return NLOPT_INVALID_ARGS;
+ case -10: return NLOPT_MAXTIME_REACHED;
case -1: return NLOPT_MAXEVAL_REACHED;
case 0: return NLOPT_XTOL_REACHED;
case 1: return NLOPT_SUCCESS;
case 2: return NLOPT_FMIN_MAX_REACHED;
+ case 20: return NLOPT_FTOL_REACHED;
+ case -200: return NLOPT_OUT_OF_MEMORY;
+ default: return NLOPT_FAILURE; /* unknown return code */
}
break;
}
return NLOPT_SUCCESS;
}
+
+nlopt_result nlopt_minimize(
+ 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 */
+ double *fmin, /* out: minimum */
+ double fmin_max, double ftol_rel, double ftol_abs,
+ double xtol_rel, const double *xtol_abs,
+ int maxeval, double maxtime)
+{
+ nlopt_result ret;
+ if (xtol_abs)
+ ret = nlopt_minimize_(algorithm, n, f, f_data, lb, ub,
+ x, fmin, fmin_max, ftol_rel, ftol_abs,
+ xtol_rel, xtol_abs, maxeval, maxtime);
+ else {
+ int i;
+ double *xtol = (double *) malloc(sizeof(double) * n);
+ if (!xtol) return NLOPT_OUT_OF_MEMORY;
+ for (i = 0; i < n; ++i) xtol[i] = -1;
+ ret = nlopt_minimize_(algorithm, n, f, f_data, lb, ub,
+ x, fmin, fmin_max, ftol_rel, ftol_abs,
+ xtol_rel, xtol, maxeval, maxtime);
+ free(xtol);
+ }
+ return ret;
+}
.B fmin_max
Stop when a function value less than or equal to
.I fmin_max
-is found. Set to Inf (see constraints section above) to disable.
+is found. Set to -Inf or NaN (see constraints section above) to disable.
.TP
.B ftol_rel
Relative tolerance on function value: stop when an optimization step
+AM_CPPFLAGS = -I$(top_srcdir)/util
+
noinst_LTLIBRARIES = libsubplex.la
libsubplex_la_SOURCES = subplex.c subplex.h
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
+#include <string.h>
#include "subplex.h"
*/
static int simplx_(D_fp f, void *fdata, integer *n, doublereal *step, integer *
- ns, integer *ips, integer *maxnfe, logical *cmode, doublereal *x,
+ ns, integer *ips, nlopt_stopping *stop, logical *cmode, doublereal *x,
doublereal *fx, integer *nfe, doublereal *s, doublereal *fs, integer *
iflag)
{
if (usubc_1.irepl > 0) {
isubc_1.new__ = FALSE_;
evalf_((D_fp)f,fdata, ns, &ips[1], &s[s_dim1 + 1], n, &x[1], &fs[1], nfe);
+ stop->nevals++;
} else {
fs[1] = *fx;
}
for (j = 2; j <= i__1; ++j) {
evalf_((D_fp)f, fdata,ns, &ips[1], &s[j * s_dim1 + 1], n, &x[1], &fs[j],
nfe);
+ stop->nevals++;
/* L10: */
}
il = 1;
goto L40;
}
evalf_((D_fp)f,fdata, ns, &ips[1], &s[itemp * s_dim1 + 1], n, &x[1], &fr, nfe);
+ stop->nevals++;
if (fr < fs[il]) {
/* expand */
goto L40;
}
evalf_((D_fp)f,fdata, ns, &ips[1], &s[ih * s_dim1 + 1], n, &x[1], &fe, nfe);
+ stop->nevals++;
if (fe < fr) {
fs[ih] = fe;
} else {
}
evalf_((D_fp)f,fdata, ns, &ips[1], &s[itemp * s_dim1 + 1], n, &x[1], &fc,
nfe);
+ stop->nevals++;
/* Computing MIN */
d__1 = fr, d__2 = fs[ih];
if (fc < min(d__1,d__2)) {
}
evalf_((D_fp)f,fdata, ns, &ips[1], &s[j * s_dim1 + 1], n, &x[1],
&fs[j], nfe);
+ stop->nevals++;
}
/* L30: */
}
*fx = isubc_1.sfbest;
}
L50:
- if (usubc_1.nfstop > 0 && *fx <= isubc_1.sfstop && usubc_1.nfxe >=
- usubc_1.nfstop) {
- *iflag = 2;
- } else if (*nfe >= *maxnfe) {
- *iflag = -1;
- } else if (dist_(ns, &s[ih * s_dim1 + 1], &s[il * s_dim1 + 1]) <= tol ||
- small) {
- *iflag = 0;
- } else {
- goto L20;
- }
+ if (*fx < stop->fmin_max)
+ *iflag = 2;
+ else if (nlopt_stop_evals(stop))
+ *iflag = -1;
+ else if (nlopt_stop_time(stop))
+ *iflag = -10;
+ else if (dist_(ns, &s[ih * s_dim1 + 1], &s[il * s_dim1 + 1]) <= tol
+ || small)
+ *iflag = 0;
+ else
+ goto L20;
/* end main loop, return best point */
-lf2c -lm (in that order)
*/
-static int subplx_(D_fp f, void *fdata, integer *n, doublereal *tol, integer *
- maxnfe, integer *mode, const doublereal *scale, doublereal *x, doublereal *
- fx, integer *nfe, doublereal *work, integer *iwork, integer *iflag)
+static int subplx_(D_fp f, void *fdata, integer *n,
+ nlopt_stopping *stop, integer *mode,
+ const doublereal *scale, doublereal *x, doublereal *fx,
+ integer *nfe, doublereal *work, integer *iwork,
+ integer *iflag)
{
/* Initialized data */
static integer istptr;
static doublereal scl, dum;
static integer ins;
- static doublereal sfx;
+ static doublereal sfx, sfx_old, *x_old;
--scale;
--work;
--iwork;
+ x_old = work;
+ work += *n;
/* Function Body */
/* ----------------------------------------------------------- */
if (*n < 1) {
goto L120;
}
- if (*tol < 0.) {
- goto L120;
- }
- if (*maxnfe < 1) {
- goto L120;
- }
if (scale[1] >= 0.) {
i__1 = *n;
for (i__ = 1; i__ <= i__1; ++i__) {
isubc_1.new__ = TRUE_;
usubc_1.initx = TRUE_;
evalf_((D_fp)f, fdata, &c__0, &iwork[1], &dum, n, &x[1], &sfx, nfe);
+ stop->nevals++;
usubc_1.initx = FALSE_;
} else {
ipptr = 1;
/* simplex loop */
-
+ sfx_old = sfx;
+ memcpy(&x_old[1], &x[1], sizeof(doublereal) * *n);
L60:
ns = iwork[ins];
L70:
- simplx_((D_fp)f, fdata, n, &work[istptr], &ns, &iwork[ipptr], maxnfe, &cmode, &x[
- 1], &sfx, nfe, &work[isptr], &work[ifsptr], iflag);
+ simplx_((D_fp)f, fdata, n, &work[istptr], &ns, &iwork[ipptr], stop, &cmode, &x[1], &sfx, nfe, &work[isptr], &work[ifsptr], iflag);
cmode = FALSE_;
if (*iflag != 0) {
goto L110;
/* check termination */
+ if (nlopt_stop_ftol(stop, sfx, sfx_old)) {
+ *iflag = 20;
+ goto L110;
+ }
+ if (nlopt_stop_x(stop, &x[1], &x_old[1])) {
+ *iflag = 0;
+ goto L110;
+ }
+
L90:
istep = istptr;
i__1 = *n;
d__1)) * usubc_1.psi;
/* Computing MAX */
d__6 = (d__3 = x[i__], abs(d__3));
- if (max(d__4,d__5) / max(d__6,1.) > *tol) {
+ if (max(d__4,d__5) / max(d__6,1.) > stop->xtol_rel) {
setstp_(&nsubs, n, &work[1], &work[istptr]);
goto L40;
}
fmin: (output) value of f at minimum
x[n]: (input) starting guess position, (output) computed minimum
fdata: data pointer passed to f
+
+ old args:
tol: relative error tolerance for x
maxnfe: maximum number of function evaluations
fmin_max, use_fmin_max: if use_fmin_max, stop when f <= fmin_max
+
+ new args: nlopt_stopping *stop (stopping criteria)
+
scale[n]: (input) scale & initial stepsizes for components of x
(if *scale < 0, |*scale| is used for all components)
Return value:
= -2 : invalid input
- = -1 : maxnfe exceeded
+ = -10 : maxtime exceeded
+ = -1 : maxevals exceeded
= 0 : tol satisfied
= 1 : limit of machine precision
= 2 : fstop reached (fstop usage is determined by values of
options minf, nfstop, and irepl. default: f(x) not
tested against fstop)
+ = 20 : ftol reached
+ = -200 : out of memory
*/
int subplex(subplex_func f, double *fmin, double *x, int n, void *fdata,
- double tol, int maxnfe,
- double fmin_max, int use_fmin_max,
- const double *scale)
+ nlopt_stopping *stop,
+ const double *scale)
{
int mode = 0, *iwork, nsmax, nsmin, errflag, nfe;
double *work;
nsmax = min(5,n);
nsmin = min(2,n);
- work = (double*) malloc(sizeof(double) * (2*n + nsmax*(nsmax+4) + 1));
+ work = (double*) malloc(sizeof(double) * (3*n + nsmax*(nsmax+4) + 1));
+ if (!work)
+ return -200;
iwork = (int*) malloc(sizeof(int) * (n + n/nsmin + 1));
- if (!work || !iwork) {
- fprintf(stderr, "subplex: error, out of memory!\n");
- exit(EXIT_FAILURE);
- }
-
- if (use_fmin_max) { /* stop when fmin_max is reached */
- subopt_(&n);
- usubc_1.nfstop = 1;
- usubc_1.fstop = fmin_max;
- mode = 2;
+ if (!iwork) {
+ free(work);
+ return -200;
}
subplx_(f,fdata, &n,
- &tol, &maxnfe, &mode,
+ stop, &mode,
scale, x,
fmin, &nfe,
work, iwork, &errflag);
#ifndef SUBPLEX_H
#define SUBPLEX_H
+#include "nlopt-util.h"
+
#ifdef __cplusplus
extern "C"
{
typedef double (*subplex_func)(int n, const double *x, void *func_data);
extern int subplex(subplex_func f, double *fmin, double *x, int n, void *fdata,
- double tol, int maxnfe,
- double fmin_max, int use_fmin_max,
+ nlopt_stopping *stop,
const double *scale);
#ifdef __cplusplus
#include "testfuncs.h"
static nlopt_algorithm algorithm = NLOPT_GLOBAL_DIRECT;
-static double ftol_rel = 0, ftol_abs = 0, xtol_rel = 0;
+static double ftol_rel = 0, ftol_abs = 0, xtol_rel = 0, xtol_abs = 0, fmin_max = -HUGE_VAL;
static int maxeval = 1000;
static double maxtime = 0.0;
{
testfunc func;
int i;
- double *x, fmin, f0;
+ double *x, fmin, f0, *xtabs;
nlopt_result ret;
double start = nlopt_seconds();
nlopt_algorithm_name(algorithm));
return 0;
}
- x = (double *) malloc(sizeof(double) * func.n * 2);
+ x = (double *) malloc(sizeof(double) * func.n * 3);
if (!x) { fprintf(stderr, "testopt: Out of memory!\n"); return 0; }
+ xtabs = x + func.n * 2;
+ for (i = 0; i < func.n; ++i) xtabs[i] = xtol_abs;
+
printf("-----------------------------------------------------------\n");
printf("Optimizing %s (%d dims) using %s algorithm\n",
func.n, func.f, func.f_data,
func.lb, func.ub,
x, &fmin,
- HUGE_VAL, ftol_rel, ftol_abs, xtol_rel, NULL,
+ 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);
" -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);
+ " -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"
+ " -X <t> : absolute 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"
+ , maxeval);
}
int main(int argc, char **argv)
if (argc <= 1)
usage(stdout);
- while ((c = getopt(argc, argv, "hLvra:o:e:")) != -1)
+ while ((c = getopt(argc, argv, "hLvr:a:o:e:t:x:X:f:F:m:")) != -1)
switch (c) {
case 'h':
usage(stdout);
testfuncs_verbose = 1;
break;
case 'r':
- srand((unsigned) atoi(optarg));
+ nlopt_srand((unsigned long) atoi(optarg));
break;
case 'a':
c = atoi(optarg);
case 'e':
maxeval = atoi(optarg);
break;
+ case 't':
+ maxtime = atof(optarg);
+ break;
+ case 'x':
+ xtol_rel = atof(optarg);
+ break;
+ case 'X':
+ xtol_abs = atof(optarg);
+ break;
+ case 'f':
+ ftol_rel = atof(optarg);
+ break;
+ case 'F':
+ ftol_abs = atof(optarg);
+ break;
+ case 'm':
+ fmin_max = atof(optarg);
+ break;
default:
fprintf(stderr, "harminv: invalid argument -%c\n", c);
usage(stderr);
noinst_LTLIBRARIES = libutil.la
-libutil_la_SOURCES = mt19937ar.c timer.c nlopt-util.h
+libutil_la_SOURCES = mt19937ar.c timer.c stop.c nlopt-util.h
extern void nlopt_init_genrand(unsigned long s);
extern double nlopt_urand(double a, double b);
+/* stopping criteria */
+typedef struct {
+ int n;
+ double fmin_max;
+ double ftol_rel;
+ double ftol_abs;
+ double xtol_rel;
+ const double *xtol_abs;
+ int nevals, maxeval;
+ double maxtime, start;
+} nlopt_stopping;
+extern int nlopt_stop_f(const nlopt_stopping *stop, double f, double oldf);
+extern int nlopt_stop_ftol(const nlopt_stopping *stop, double f, double oldf);
+extern int nlopt_stop_x(const nlopt_stopping *stop,
+ const double *x, const double *oldx);
+extern int nlopt_stop_xs(const nlopt_stopping *stop,
+ const double *xs, const double *oldxs,
+ const double *scale_min, const double *scale_max);
+extern int nlopt_stop_evals(const nlopt_stopping *stop);
+extern int nlopt_stop_time(const nlopt_stopping *stop);
+
#ifdef __cplusplus
} /* extern "C" */
#endif /* __cplusplus */
--- /dev/null
+#include <math.h>
+#include "nlopt-util.h"
+
+/* utility routines to implement the various stopping criteria */
+
+static int relstop(double old, double new, double reltol, double abstol)
+{
+ return(fabs(new - old) < abstol
+ || fabs(new - old) < reltol * (fabs(new) + fabs(old)) * 0.5);
+}
+
+int nlopt_stop_ftol(const nlopt_stopping *s, const double f, double oldf)
+{
+ return (relstop(oldf, f, s->ftol_rel, s->ftol_abs));
+}
+
+int nlopt_stop_f(const nlopt_stopping *s, const double f, double oldf)
+{
+ return (f <= s->fmin_max || nlopt_stop_ftol(s, f, oldf));
+}
+
+int nlopt_stop_x(const nlopt_stopping *s, const double *x, const double *oldx)
+{
+ int i;
+ for (i = 0; i < s->n; ++i)
+ if (relstop(oldx[i], x[i], s->xtol_rel, s->xtol_abs[i]))
+ return 1;
+ return 0;
+}
+
+static double sc(double x, double smin, double smax)
+{
+ return smin + x * (smax - smin);
+}
+
+/* some of the algorithms rescale x to a unit hypercube, so we need to
+ scale back before we can compare to the tolerances */
+int nlopt_stop_xs(const nlopt_stopping *s,
+ const double *xs, const double *oldxs,
+ const double *scale_min, const double *scale_max)
+{
+ int i;
+ for (i = 0; i < s->n; ++i)
+ if (relstop(sc(oldxs[i], scale_min[i], scale_max[i]),
+ sc(xs[i], scale_min[i], scale_max[i]),
+ s->xtol_rel, s->xtol_abs[i]))
+ return 1;
+ return 0;
+}
+
+int nlopt_stop_evals(const nlopt_stopping *s)
+{
+ return (s->maxeval > 0 && s->nevals >= s->maxeval);
+}
+
+int nlopt_stop_time(const nlopt_stopping *s)
+{
+ return (s->maxtime > 0 && nlopt_seconds() - s->start >= s->maxtime);
+}