From 60837d8a4e304f302d9c890bfd16c892e15d8597 Mon Sep 17 00:00:00 2001 From: stevenj Date: Thu, 8 Jul 2010 17:23:53 -0400 Subject: [PATCH] in functions with constraints, make sure forced_stop stops immediately, before evaluating any more user callbacks darcs-hash:20100708212353-c8de0-5060b668f7e91e93ca182f3c2a8ee7618b4d5032.gz --- NEWS | 4 ++++ auglag/auglag.c | 19 +++++++++++++++++-- cobyla/cobyla.c | 8 ++++++-- isres/isres.c | 6 ++++++ mma/mma.c | 6 ++++++ 5 files changed, 39 insertions(+), 4 deletions(-) diff --git a/NEWS b/NEWS index 5014cc1..e8bd3b0 100644 --- a/NEWS +++ b/NEWS @@ -4,6 +4,10 @@ NLopt 2.1.2 (8 July 2010) now pass a 2-dimensional array for the gradient argument, rather than a flattened 1d array. +* Improved handling of exceptions and forced stops for constrained + optimization, making sure that no constraints are evaluated after + the stop. + * Return an NLOPT_INVALID_ARGS error if more than n equality constraints are added in an n-dimensional problem. diff --git a/auglag/auglag.c b/auglag/auglag.c index 1b59aa8..f43934e 100644 --- a/auglag/auglag.c +++ b/auglag/auglag.c @@ -34,9 +34,12 @@ static double auglag(unsigned n, const double *x, double *grad, void *data) unsigned j, k; L = d->f(n, x, grad, d->f_data); + d->stop->nevals++; + if (nlopt_stop_forced(d->stop)) return L; for (ii = i = 0; i < d->p; ++i) { nlopt_eval_constraint(restmp, gradtmp, d->h + i, n, x); + if (nlopt_stop_forced(d->stop)) return L; for (k = 0; k < d->h[i].m; ++k) { double h = restmp[k] + lambda[ii++] / rho; L += 0.5 * rho * h*h; @@ -47,6 +50,7 @@ static double auglag(unsigned n, const double *x, double *grad, void *data) for (ii = i = 0; i < d->m; ++i) { nlopt_eval_constraint(restmp, gradtmp, d->fc + i, n, x); + if (nlopt_stop_forced(d->stop)) return L; for (k = 0; k < d->fc[i].m; ++k) { double fc = restmp[k] + mu[ii++] / rho; if (fc > 0) { @@ -56,8 +60,6 @@ static double auglag(unsigned n, const double *x, double *grad, void *data) } } } - - d->stop->nevals++; return L; } @@ -142,10 +144,14 @@ nlopt_result auglag_minimize(int n, nlopt_func f, void *f_data, double con2 = 0; d.stop->nevals++; fcur = f(n, xcur, NULL, f_data); + if (nlopt_stop_forced(stop)) { + ret = NLOPT_FORCED_STOP; goto done; } penalty = 0; feasible = 1; for (i = 0; i < d.p; ++i) { nlopt_eval_constraint(d.restmp, NULL, d.h + i, n, xcur); + if (nlopt_stop_forced(stop)) { + ret = NLOPT_FORCED_STOP; goto done; } for (k = 0; k < d.h[i].m; ++k) { double hi = d.restmp[k]; penalty += fabs(hi); @@ -155,6 +161,8 @@ nlopt_result auglag_minimize(int n, nlopt_func f, void *f_data, } for (i = 0; i < d.m; ++i) { nlopt_eval_constraint(d.restmp, NULL, d.fc + i, n, xcur); + if (nlopt_stop_forced(stop)) { + ret = NLOPT_FORCED_STOP; goto done; } for (k = 0; k < d.fc[i].m; ++k) { double fci = d.restmp[k]; penalty += fci > 0 ? fci : 0; @@ -191,6 +199,8 @@ nlopt_result auglag_minimize(int n, nlopt_func f, void *f_data, d.stop->nevals++; fcur = f(n, xcur, NULL, f_data); + if (nlopt_stop_forced(stop)) { + ret = NLOPT_FORCED_STOP; goto done; } if (auglag_verbose) printf("auglag: fcur = %g\n", fcur); @@ -199,6 +209,8 @@ nlopt_result auglag_minimize(int n, nlopt_func f, void *f_data, feasible = 1; for (i = ii = 0; i < d.p; ++i) { nlopt_eval_constraint(d.restmp, NULL, d.h + i, n, xcur); + if (nlopt_stop_forced(stop)) { + ret = NLOPT_FORCED_STOP; goto done; } for (k = 0; k < d.h[i].m; ++k) { double hi = d.restmp[k]; double newlam = d.lambda[ii] + d.rho * hi; @@ -210,6 +222,8 @@ nlopt_result auglag_minimize(int n, nlopt_func f, void *f_data, } for (i = ii = 0; i < d.m; ++i) { nlopt_eval_constraint(d.restmp, NULL, d.fc + i, n, xcur); + if (nlopt_stop_forced(stop)) { + ret = NLOPT_FORCED_STOP; goto done; } for (k = 0; k < d.fc[i].m; ++k) { double fci = d.restmp[k]; double newmu = d.mu[ii] + d.rho * fci; @@ -269,6 +283,7 @@ nlopt_result auglag_minimize(int n, nlopt_func f, void *f_data, if (ICM == 0) return NLOPT_FTOL_REACHED; } while (1); +done: free(xcur); return ret; } diff --git a/cobyla/cobyla.c b/cobyla/cobyla.c index ea981cd..990f389 100644 --- a/cobyla/cobyla.c +++ b/cobyla/cobyla.c @@ -71,6 +71,7 @@ typedef struct { double *xtmp; const double *lb, *ub; double *con_tol; + nlopt_stopping *stop; } func_wrap_state; static int func_wrap(int n, int m, double *x, double *f, double *con, @@ -92,15 +93,18 @@ static int func_wrap(int n, int m, double *x, double *f, double *con, } *f = s->f(n, xtmp, NULL, s->f_data); + if (nlopt_stop_forced(s->stop)) return 1; i = 0; for (j = 0; j < s->m_orig; ++j) { nlopt_eval_constraint(con + i, NULL, s->fc+j, n, xtmp); + if (nlopt_stop_forced(s->stop)) return 1; for (k = 0; k < s->fc[j].m; ++k) con[i + k] = -con[i + k]; i += s->fc[j].m; } for (j = 0; j < s->p; ++j) { nlopt_eval_constraint(con + i, NULL, s->h+j, n, xtmp); + if (nlopt_stop_forced(s->stop)) return 1; for (k = 0; k < s->h[j].m; ++k) con[(i + s->h[j].m) + k] = -con[i + k]; i += 2 * s->h[j].m; @@ -111,7 +115,6 @@ static int func_wrap(int n, int m, double *x, double *f, double *con, if (!nlopt_isinf(ub[j])) con[i++] = ub[j] - x[j]; } - if (m != i) return 1; /* ... bug?? */ return 0; } @@ -185,6 +188,7 @@ nlopt_result cobyla_minimize(int n, nlopt_func f, void *f_data, s.p = p; s.h = h; s.lb = lb; s.ub = ub; + s.stop = stop; s.xtmp = (double *) malloc(sizeof(double) * n); if (!s.xtmp) return NLOPT_OUT_OF_MEMORY; @@ -548,7 +552,7 @@ L40: if (*iprint >= 1) { fprintf(stderr, "cobyla: user requested end of minimization.\n"); } - rc = NLOPT_FAILURE; + rc = NLOPT_FORCED_STOP; goto L600; } diff --git a/isres/isres.c b/isres/isres.c index 5b4fec3..02b6f3c 100644 --- a/isres/isres.c +++ b/isres/isres.c @@ -131,10 +131,14 @@ nlopt_result isres_minimize(int n, nlopt_func f, void *f_data, double gpenalty; stop->nevals++; fval[k] = f(n, xs + k*n, NULL, f_data); + if (nlopt_stop_forced(stop)) { + ret = NLOPT_FORCED_STOP; goto done; } penalty[k] = 0; for (c = 0; c < m; ++c) { /* inequality constraints */ nlopt_eval_constraint(results, NULL, fc + c, n, xs + k*n); + if (nlopt_stop_forced(stop)) { + ret = NLOPT_FORCED_STOP; goto done; } for (ires = 0; ires < fc[c].m; ++ires) { double gval = results[ires]; if (gval > fc[c].tol[ires]) feasible = 0; @@ -146,6 +150,8 @@ nlopt_result isres_minimize(int n, nlopt_func f, void *f_data, for (c = m; c < mp; ++c) { /* equality constraints */ nlopt_eval_constraint(results, NULL, h + (c-m), n, xs + k*n); + if (nlopt_stop_forced(stop)) { + ret = NLOPT_FORCED_STOP; goto done; } for (ires = 0; ires < h[c-m].m; ++ires) { double hval = results[ires]; if (fabs(hval) > h[c-m].tol[ires]) feasible = 0; diff --git a/mma/mma.c b/mma/mma.c index f9fa5f2..f3db384 100644 --- a/mma/mma.c +++ b/mma/mma.c @@ -210,12 +210,14 @@ nlopt_result mma_minimize(unsigned n, nlopt_func f, void *f_data, dd.fval = fcur = *minf = f(n, x, dfdx, f_data); stop->nevals++; memcpy(xcur, x, sizeof(double) * n); + if (nlopt_stop_forced(stop)) { ret = NLOPT_FORCED_STOP; goto done; } feasible = 1; infeasibility = 0; for (i = ifc = 0; ifc < mfc; ++ifc) { nlopt_eval_constraint(fcval + i, dfcdx + i*n, fc + ifc, n, x); i += fc[ifc].m; + if (nlopt_stop_forced(stop)) { ret = NLOPT_FORCED_STOP; goto done; } } for (i = 0; i < m; ++i) { feasible = feasible && (fcval[i] <= 0 || isnan(fcval[i])); @@ -285,6 +287,8 @@ nlopt_result mma_minimize(unsigned n, nlopt_func f, void *f_data, fcur = f(n, xcur, dfdx_cur, f_data); stop->nevals++; + if (nlopt_stop_forced(stop)) { + ret = NLOPT_FORCED_STOP; goto done; } feasible_cur = 1; infeasibility_cur = 0; new_infeasible_constraint = 0; inner_done = dd.gval >= fcur; @@ -292,6 +296,8 @@ nlopt_result mma_minimize(unsigned n, nlopt_func f, void *f_data, nlopt_eval_constraint(fcval_cur + i, dfcdx_cur + i*n, fc + ifc, n, xcur); i += fc[ifc].m; + if (nlopt_stop_forced(stop)) { + ret = NLOPT_FORCED_STOP; goto done; } } for (i = ifc = 0; ifc < mfc; ++ifc) { unsigned i0 = i, inext = i + fc[ifc].m; -- 2.30.2