From 6cd45b41ceb297e5affb8aba89f838e46473ff8c Mon Sep 17 00:00:00 2001 From: stevenj Date: Mon, 5 Apr 2010 01:41:36 -0400 Subject: [PATCH] change new API to use unsigned where it makes sense darcs-hash:20100405054136-c8de0-a106e736e68aa365986a1cc91b9faf8492dd6ea1.gz --- api/deprecated.c | 25 ++++++++------- api/general.c | 2 +- api/nlopt-in.hpp | 26 ++++++++-------- api/nlopt-internal.h | 14 ++++----- api/nlopt.h | 26 +++++++++------- api/optimize.c | 70 +++++++++++++++++++++--------------------- api/options.c | 72 +++++++++++++++++++++----------------------- test/testfuncs.c | 60 ++++++++++++++++++------------------ util/nlopt-util.h | 2 +- util/stop.c | 6 ++-- 10 files changed, 155 insertions(+), 148 deletions(-) diff --git a/api/deprecated.c b/api/deprecated.c index 11539bd..6bf0768 100644 --- a/api/deprecated.c +++ b/api/deprecated.c @@ -53,15 +53,15 @@ int nlopt_stochastic_population = 0; int nlopt_get_stochastic_population(void) { return nlopt_stochastic_population; } void nlopt_set_stochastic_population(int pop) { - nlopt_stochastic_population = pop; } + nlopt_stochastic_population = pop <= 0 ? 0 : (unsigned) pop; } /*************************************************************************/ nlopt_result nlopt_minimize_econstrained( nlopt_algorithm algorithm, - int n, nlopt_func f, void *f_data, - int m, nlopt_func fc, void *fc_data_, ptrdiff_t fc_datum_size, - int p, nlopt_func h, void *h_data_, ptrdiff_t h_datum_size, + int n, nlopt_func_old f, void *f_data, + int m, nlopt_func_old fc, void *fc_data_, ptrdiff_t fc_datum_size, + int p, nlopt_func_old h, void *h_data_, ptrdiff_t h_datum_size, const double *lb, const double *ub, /* bounds */ double *x, /* in: initial guess, out: minimizer */ double *minf, /* out: minimum */ @@ -72,17 +72,20 @@ nlopt_result nlopt_minimize_econstrained( { char *fc_data = (char *) fc_data_; char *h_data = (char *) h_data_; - nlopt_opt opt = nlopt_create(algorithm, n); + nlopt_opt opt; nlopt_result ret; int i; + if (n < 0 || m < 0 || p < 0) return NLOPT_INVALID_ARGS; + + opt = nlopt_create(algorithm, (unsigned) n); if (!opt) return NLOPT_INVALID_ARGS; - ret = nlopt_set_min_objective(opt, f, f_data); + ret = nlopt_set_min_objective(opt, (nlopt_func) f, f_data); if (ret != NLOPT_SUCCESS) { nlopt_destroy(opt); return ret; } for (i = 0; i < m; ++i) { - ret = nlopt_add_inequality_constraint(opt, fc, + ret = nlopt_add_inequality_constraint(opt, (nlopt_func) fc, fc_data + i*fc_datum_size, 0.0); if (ret != NLOPT_SUCCESS) { nlopt_destroy(opt); return ret; } @@ -90,7 +93,7 @@ nlopt_result nlopt_minimize_econstrained( (void) htol_rel; /* unused */ for (i = 0; i < p; ++i) { - ret = nlopt_add_equality_constraint(opt, h, + ret = nlopt_add_equality_constraint(opt, (nlopt_func) h, h_data + i*h_datum_size, htol_abs); if (ret != NLOPT_SUCCESS) { nlopt_destroy(opt); return ret; } @@ -128,8 +131,8 @@ nlopt_result nlopt_minimize_econstrained( nlopt_result nlopt_minimize_constrained( nlopt_algorithm algorithm, - int n, nlopt_func f, void *f_data, - int m, nlopt_func fc, void *fc_data, ptrdiff_t fc_datum_size, + int n, nlopt_func_old f, void *f_data, + int m, nlopt_func_old fc, void *fc_data, ptrdiff_t fc_datum_size, const double *lb, const double *ub, /* bounds */ double *x, /* in: initial guess, out: minimizer */ double *minf, /* out: minimum */ @@ -146,7 +149,7 @@ nlopt_result nlopt_minimize_constrained( nlopt_result nlopt_minimize( nlopt_algorithm algorithm, - int n, nlopt_func f, void *f_data, + int n, nlopt_func_old f, void *f_data, const double *lb, const double *ub, /* bounds */ double *x, /* in: initial guess, out: minimizer */ double *minf, /* out: minimum */ diff --git a/api/general.c b/api/general.c index 3cb66c0..f52932f 100644 --- a/api/general.c +++ b/api/general.c @@ -96,7 +96,7 @@ static const char nlopt_algorithm_names[NLOPT_NUM_ALGORITHMS][256] = { const char *nlopt_algorithm_name(nlopt_algorithm a) { - if (a < 0 || a >= NLOPT_NUM_ALGORITHMS) return "UNKNOWN"; + if (((int) a) < 0 || a >= NLOPT_NUM_ALGORITHMS) return "UNKNOWN"; return nlopt_algorithm_names[a]; } diff --git a/api/nlopt-in.hpp b/api/nlopt-in.hpp index 1e76826..1aaa27d 100644 --- a/api/nlopt-in.hpp +++ b/api/nlopt-in.hpp @@ -47,10 +47,10 @@ namespace nlopt { public: // should return function value, and set grad to gradient // (x and grad are length n) - virtual double operator()(int n, const double *x, double *grad) = 0; + virtual double operator()(unsigned n, const double *x, double *grad) = 0; // should return function value (x is length n) - virtual double operator()(int n, const double *x) = 0; + virtual double operator()(unsigned n, const double *x) = 0; }; // (Note: it is inefficient to use std::vector for the arguments, @@ -73,7 +73,7 @@ namespace nlopt { } // nlopt_func wrapper around C++ "functional" - static double myfunc(int n, const double *x, double *grad, void *f_) { + static double myfunc(unsigned n, const double *x, double *grad, void *f_) { func *f = reinterpret_cast(f_); return grad ? (*f)(n, x, grad) : (*f)(n, x); } @@ -81,10 +81,10 @@ namespace nlopt { public: // Constructors etc. opt() : o(NULL) {} - opt(nlopt_algorithm a, int n) : o(nlopt_create(a, n)) { + opt(nlopt_algorithm a, unsigned n) : o(nlopt_create(a, n)) { if (!o) throw std::bad_alloc(); } - opt(algorithm a, int n) : o(nlopt_create(nlopt_algorithm(a), n)) { + opt(algorithm a, unsigned n) : o(nlopt_create(nlopt_algorithm(a), n)) { if (!o) throw std::bad_alloc(); } opt(const nlopt_opt o0) : o(nlopt_copy(o0)) { @@ -117,7 +117,7 @@ namespace nlopt { if (!o) throw std::invalid_argument("uninitialized nlopt::opt"); return algorithm(nlopt_get_algorithm(o)); } - int get_dimension() const { + unsigned get_dimension() const { if (!o) throw std::invalid_argument("uninitialized nlopt::opt"); return nlopt_get_dimension(o); } @@ -171,18 +171,18 @@ namespace nlopt { mythrow(ret); \ } \ void get_##name(std::vector &v) const { \ - if (o && unsigned(nlopt_get_dimension(o)) != v.size()) \ + if (o && nlopt_get_dimension(o) != v.size()) \ throw std::invalid_argument("dimension mismatch"); \ get_##name(v.empty() ? NULL : &v[0]); \ } \ std::vector get_##name(void) const { \ if (!o) throw std::invalid_argument("uninitialized nlopt::opt"); \ - std::vector v(unsigned(nlopt_get_dimension(o))); \ + std::vector v(nlopt_get_dimension(o)); \ get_##name(v); \ return v; \ } \ void set_##name(const std::vector &v) { \ - if (o && unsigned(nlopt_get_dimension(o)) != v.size()) \ + if (o && nlopt_get_dimension(o) != v.size()) \ throw std::invalid_argument("dimension mismatch"); \ set_##name(v.empty() ? NULL : &v[0]); \ } @@ -219,7 +219,7 @@ namespace nlopt { set_local_optimizer(lo.o); } - NLOPT_GETSET(int, population) + NLOPT_GETSET(unsigned, population) NLOPT_GETSET_VEC(initial_step) void set_default_initial_step(const double *x) { @@ -234,15 +234,15 @@ namespace nlopt { mythrow(ret); } void get_initial_step(const std::vector &x, std::vector &dx) const { - if (o && (unsigned(nlopt_get_dimension(o)) != x.size() - || unsigned(nlopt_get_dimension(o)) != dx.size())) + if (o && (nlopt_get_dimension(o) != x.size() + || nlopt_get_dimension(o) != dx.size())) throw std::invalid_argument("dimension mismatch"); get_initial_step(x.empty() ? NULL : &x[0], dx.empty() ? NULL : &dx[0]); } std::vector get_initial_step(const std::vector &x) const { if (!o) throw std::invalid_argument("uninitialized nlopt::opt"); - std::vector v(unsigned(nlopt_get_dimension(o))); + std::vector v(nlopt_get_dimension(o)); get_initial_step(x, v); return v; } diff --git a/api/nlopt-internal.h b/api/nlopt-internal.h index e74ea54..b2aab30 100644 --- a/api/nlopt-internal.h +++ b/api/nlopt-internal.h @@ -35,18 +35,18 @@ extern "C" struct nlopt_opt_s { nlopt_algorithm algorithm; /* the optimization algorithm (immutable) */ - int n; /* the dimension of the problem (immutable) */ + unsigned n; /* the dimension of the problem (immutable) */ nlopt_func f; void *f_data; /* objective function to minimize */ double *lb, *ub; /* lower and upper bounds (length n) */ - int m; /* number of inequality constraints */ - int m_alloc; /* number of inequality constraints allocated */ + unsigned m; /* number of inequality constraints */ + unsigned m_alloc; /* number of inequality constraints allocated */ nlopt_constraint *fc; /* inequality constraints, length m_alloc */ - int p; /* number of equality constraints */ - int p_alloc; /* number of inequality constraints allocated */ + unsigned p; /* number of equality constraints */ + unsigned p_alloc; /* number of inequality constraints allocated */ nlopt_constraint *h; /* equality constraints, length p_alloc */ /* stopping criteria */ @@ -58,7 +58,7 @@ struct nlopt_opt_s { /* algorithm-specific parameters */ nlopt_opt local_opt; /* local optimizer */ - int stochastic_population; /* population size for stochastic algs */ + unsigned stochastic_population; /* population size for stochastic algs */ double *dx; /* initial step sizes (length n) for nonderivative algs */ }; @@ -71,7 +71,7 @@ extern int nlopt_srand_called; /* whether the random seed is initialized */ extern nlopt_algorithm nlopt_local_search_alg_deriv; extern nlopt_algorithm nlopt_local_search_alg_nonderiv; extern int nlopt_local_search_maxeval; -extern int nlopt_stochastic_population; +extern unsigned nlopt_stochastic_population; /*********************************************************************/ diff --git a/api/nlopt.h b/api/nlopt.h index 39c56d1..c7dcdc9 100644 --- a/api/nlopt.h +++ b/api/nlopt.h @@ -55,7 +55,7 @@ extern "C" { #endif /* __cplusplus */ -typedef double (*nlopt_func)(int n, const double *x, +typedef double (*nlopt_func)(unsigned n, const double *x, double *gradient, /* NULL if not needed */ void *func_data); @@ -164,7 +164,7 @@ typedef struct nlopt_opt_s *nlopt_opt; /* the only immutable parameters of an optimization are the algorithm and the dimension n of the problem, since changing either of these could have side-effects on lots of other parameters */ -NLOPT_EXTERN nlopt_opt nlopt_create(nlopt_algorithm algorithm, int n); +NLOPT_EXTERN nlopt_opt nlopt_create(nlopt_algorithm algorithm, unsigned n); NLOPT_EXTERN void nlopt_destroy(nlopt_opt opt); NLOPT_EXTERN nlopt_opt nlopt_copy(const nlopt_opt opt); @@ -175,7 +175,7 @@ NLOPT_EXTERN nlopt_result nlopt_set_min_objective(nlopt_opt opt, nlopt_func f, void *f_data); NLOPT_EXTERN nlopt_algorithm nlopt_get_algorithm(const nlopt_opt opt); -NLOPT_EXTERN int nlopt_get_dimension(const nlopt_opt opt); +NLOPT_EXTERN unsigned nlopt_get_dimension(const nlopt_opt opt); /* constraints: */ @@ -230,8 +230,8 @@ NLOPT_EXTERN double nlopt_get_maxtime(nlopt_opt opt); NLOPT_EXTERN nlopt_result nlopt_set_local_optimizer(nlopt_opt opt, const nlopt_opt local_opt); -NLOPT_EXTERN nlopt_result nlopt_set_population(nlopt_opt opt, int pop); -NLOPT_EXTERN int nlopt_get_population(const nlopt_opt opt); +NLOPT_EXTERN nlopt_result nlopt_set_population(nlopt_opt opt, unsigned pop); +NLOPT_EXTERN unsigned nlopt_get_population(const nlopt_opt opt); NLOPT_EXTERN nlopt_result nlopt_set_default_initial_step(nlopt_opt opt, const double *x); @@ -255,9 +255,13 @@ NLOPT_EXTERN nlopt_result nlopt_get_initial_step(const nlopt_opt opt, # define NLOPT_DEPRECATED #endif +typedef double (*nlopt_func_old)(int n, const double *x, + double *gradient, /* NULL if not needed */ + void *func_data); + NLOPT_EXTERN nlopt_result nlopt_minimize( nlopt_algorithm algorithm, - int n, nlopt_func f, void *f_data, + int n, nlopt_func_old f, void *f_data, const double *lb, const double *ub, /* bounds */ double *x, /* in: initial guess, out: minimizer */ double *minf, /* out: minimum */ @@ -267,8 +271,8 @@ NLOPT_EXTERN nlopt_result nlopt_minimize( NLOPT_EXTERN nlopt_result nlopt_minimize_constrained( nlopt_algorithm algorithm, - int n, nlopt_func f, void *f_data, - int m, nlopt_func fc, void *fc_data, ptrdiff_t fc_datum_size, + int n, nlopt_func_old f, void *f_data, + int m, nlopt_func_old fc, void *fc_data, ptrdiff_t fc_datum_size, const double *lb, const double *ub, /* bounds */ double *x, /* in: initial guess, out: minimizer */ double *minf, /* out: minimum */ @@ -278,9 +282,9 @@ NLOPT_EXTERN nlopt_result nlopt_minimize_constrained( NLOPT_EXTERN nlopt_result nlopt_minimize_econstrained( nlopt_algorithm algorithm, - int n, nlopt_func f, void *f_data, - int m, nlopt_func fc, void *fc_data, ptrdiff_t fc_datum_size, - int p, nlopt_func h, void *h_data, ptrdiff_t h_datum_size, + int n, nlopt_func_old f, void *f_data, + int m, nlopt_func_old fc, void *fc_data, ptrdiff_t fc_datum_size, + int p, nlopt_func_old h, void *h_data, ptrdiff_t h_datum_size, const double *lb, const double *ub, /* bounds */ double *x, /* in: initial guess, out: minimizer */ double *minf, /* out: minimum */ diff --git a/api/optimize.c b/api/optimize.c index 595f44e..c2e596c 100644 --- a/api/optimize.c +++ b/api/optimize.c @@ -83,21 +83,21 @@ static double f_bound(int n, const double *x, void *data_) if (x[i] < data->lb[i] || x[i] > data->ub[i]) return HUGE_VAL; - f = data->f(n, x, NULL, data->f_data); + f = data->f((unsigned) n, x, NULL, data->f_data); return (isnan(f) || nlopt_isinf(f) ? HUGE_VAL : f); } static double f_noderiv(int n, const double *x, void *data_) { nlopt_data *data = (nlopt_data *) data_; - return data->f(n, x, NULL, data->f_data); + return data->f((unsigned) n, x, NULL, data->f_data); } static double f_direct(int n, const double *x, int *undefined, void *data_) { nlopt_data *data = (nlopt_data *) data_; double f; - f = data->f(n, x, NULL, data->f_data); + f = data->f((unsigned) n, x, NULL, data->f_data); *undefined = isnan(f) || nlopt_isinf(f); return f; } @@ -107,7 +107,7 @@ static double f_direct(int n, const double *x, int *undefined, void *data_) /* get min(dx) for algorithms requiring a scalar initial step size */ static nlopt_result initial_step(nlopt_opt opt, const double *x, double *step) { - int freedx = 0, i; + unsigned freedx = 0, i; if (!opt->dx) { freedx = 1; @@ -127,9 +127,9 @@ static nlopt_result initial_step(nlopt_opt opt, const double *x, double *step) /*********************************************************************/ /* return true if [lb,ub] is finite in every dimension (n dimensions) */ -static int finite_domain(int n, const double *lb, const double *ub) +static int finite_domain(unsigned n, const double *lb, const double *ub) { - int i; + unsigned i; for (i = 0; i < n; ++i) if (nlopt_isinf(ub[i] - lb[i])) return 0; return 1; @@ -147,7 +147,8 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) const double *lb, *ub; nlopt_algorithm algorithm; nlopt_func f; void *f_data; - int n, i; + unsigned n, i; + int ni; nlopt_data d; nlopt_stopping stop; @@ -155,6 +156,7 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) /* copy a few params to local vars for convenience */ n = opt->n; + ni = (int) n; /* most of the subroutines take "int" arg */ lb = opt->lb; ub = opt->ub; algorithm = opt->algorithm; f = opt->f; f_data = opt->f_data; @@ -196,7 +198,7 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) case NLOPT_GN_DIRECT_L: case NLOPT_GN_DIRECT_L_RAND: if (!finite_domain(n, lb, ub)) return NLOPT_INVALID_ARGS; - return cdirect(n, f, f_data, + return cdirect(ni, f, f_data, lb, ub, x, minf, &stop, 0.0, (algorithm != NLOPT_GN_DIRECT) + 3 * (algorithm == NLOPT_GN_DIRECT_L_RAND @@ -208,7 +210,7 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) case NLOPT_GN_DIRECT_L_NOSCAL: case NLOPT_GN_DIRECT_L_RAND_NOSCAL: if (!finite_domain(n, lb, ub)) return NLOPT_INVALID_ARGS; - return cdirect_unscaled(n, f, f_data, lb, ub, x, minf, + return cdirect_unscaled(ni, f, f_data, lb, ub, x, minf, &stop, 0.0, (algorithm != NLOPT_GN_DIRECT) + 3 * (algorithm == NLOPT_GN_DIRECT_L_RAND ? 2 : (algorithm != NLOPT_GN_DIRECT)) @@ -217,7 +219,7 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) case NLOPT_GN_ORIG_DIRECT: case NLOPT_GN_ORIG_DIRECT_L: if (!finite_domain(n, lb, ub)) return NLOPT_INVALID_ARGS; - switch (direct_optimize(f_direct, &d, n, lb, ub, x, minf, + switch (direct_optimize(f_direct, &d, ni, lb, ub, x, minf, stop.maxeval, -1, 0.0, 0.0, pow(stop.xtol_rel, (double) n), -1.0, stop.minf_max, 0.0, @@ -250,9 +252,9 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) case NLOPT_GD_STOGO_RAND: #ifdef WITH_CXX if (!finite_domain(n, lb, ub)) return NLOPT_INVALID_ARGS; - if (!stogo_minimize(n, f, f_data, x, minf, lb, ub, &stop, + if (!stogo_minimize(ni, f, f_data, x, minf, lb, ub, &stop, algorithm == NLOPT_GD_STOGO - ? 0 : POP(2*n))) + ? 0 : (int) POP(2*n))) return NLOPT_FAILURE; break; #else @@ -291,7 +293,7 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) if (initial_step(opt, x, &step) != NLOPT_SUCCESS) return NLOPT_OUT_OF_MEMORY; return praxis_(0.0, DBL_EPSILON, - step, n, x, f_bound, &d, &stop, minf); + step, ni, x, f_bound, &d, &stop, minf); } #ifdef WITH_NOCEDAL @@ -303,8 +305,8 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) int uinf = nlopt_isinf(ub[i]) && ub[i] > 0; nbd[i] = linf && uinf ? 0 : (uinf ? 1 : (linf ? 3 : 2)); } - iret = lbfgsb_minimize(n, f, f_data, x, nbd, lb, ub, - MIN(n, 5), 0.0, stop.ftol_rel, + iret = lbfgsb_minimize(ni, f, f_data, x, nbd, lb, ub, + MIN(ni, 5), 0.0, stop.ftol_rel, stop.xtol_abs[0] > 0 ? stop.xtol_abs[0] : stop.xtol_rel, stop.maxeval); @@ -329,25 +331,25 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) #endif case NLOPT_LD_LBFGS: - return luksan_plis(n, f, f_data, lb, ub, x, minf, &stop); + return luksan_plis(ni, f, f_data, lb, ub, x, minf, &stop); case NLOPT_LD_VAR1: case NLOPT_LD_VAR2: - return luksan_plip(n, f, f_data, lb, ub, x, minf, &stop, + return luksan_plip(ni, f, f_data, lb, ub, x, minf, &stop, algorithm == NLOPT_LD_VAR1 ? 1 : 2); case NLOPT_LD_TNEWTON: case NLOPT_LD_TNEWTON_RESTART: case NLOPT_LD_TNEWTON_PRECOND: case NLOPT_LD_TNEWTON_PRECOND_RESTART: - return luksan_pnet(n, f, f_data, lb, ub, x, minf, &stop, + return luksan_pnet(ni, f, f_data, lb, ub, x, minf, &stop, 1 + (algorithm - NLOPT_LD_TNEWTON) % 2, 1 + (algorithm - NLOPT_LD_TNEWTON) / 2); case NLOPT_GN_CRS2_LM: if (!finite_domain(n, lb, ub)) return NLOPT_INVALID_ARGS; - return crs_minimize(n, f, f_data, lb, ub, x, minf, &stop, - POP(0), 0); + return crs_minimize(ni, f, f_data, lb, ub, x, minf, &stop, + (int) POP(0), 0); case NLOPT_GN_MLSL: case NLOPT_GD_MLSL: @@ -384,8 +386,8 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) nlopt_set_ftol_rel(local_opt, 1e-15); nlopt_set_xtol_rel(local_opt, 1e-7); } - ret = mlsl_minimize(n, f, f_data, lb, ub, x, minf, &stop, - local_opt, POP(0), + ret = mlsl_minimize(ni, f, f_data, lb, ub, x, minf, &stop, + local_opt, (int) POP(0), algorithm >= NLOPT_GN_MLSL_LDS); if (!opt->local_opt) nlopt_destroy(local_opt); return ret; @@ -404,7 +406,7 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) nlopt_set_maxeval(dual_opt, LO(maxeval, 100000)); #undef LO - ret = mma_minimize(n, f, f_data, opt->m, opt->fc, + ret = mma_minimize(ni, f, f_data, (int) (opt->m), opt->fc, lb, ub, x, minf, &stop, dual_opt); nlopt_destroy(dual_opt); return ret; @@ -414,7 +416,7 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) double step; if (initial_step(opt, x, &step) != NLOPT_SUCCESS) return NLOPT_OUT_OF_MEMORY; - return cobyla_minimize(n, f, f_data, + return cobyla_minimize(ni, f, f_data, opt->m, opt->fc, lb, ub, x, minf, &stop, step); @@ -424,7 +426,7 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) double step; if (initial_step(opt, x, &step) != NLOPT_SUCCESS) return NLOPT_OUT_OF_MEMORY; - return newuoa(n, 2*n+1, x, 0, 0, step, + return newuoa(ni, 2*n+1, x, 0, 0, step, &stop, minf, f_noderiv, &d); } @@ -432,7 +434,7 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) double step; if (initial_step(opt, x, &step) != NLOPT_SUCCESS) return NLOPT_OUT_OF_MEMORY; - return newuoa(n, 2*n+1, x, lb, ub, step, + return newuoa(ni, 2*n+1, x, lb, ub, step, &stop, minf, f_noderiv, &d); } @@ -440,7 +442,7 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) double step; if (initial_step(opt, x, &step) != NLOPT_SUCCESS) return NLOPT_OUT_OF_MEMORY; - return bobyqa(n, 2*n+1, x, lb, ub, step, + return bobyqa(ni, 2*n+1, x, lb, ub, step, &stop, minf, f_noderiv, &d); } @@ -455,9 +457,9 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) return NLOPT_OUT_OF_MEMORY; } if (algorithm == NLOPT_LN_NELDERMEAD) - ret= nldrmd_minimize(n,f,f_data,lb,ub,x,minf,opt->dx,&stop); + ret= nldrmd_minimize(ni,f,f_data,lb,ub,x,minf,opt->dx,&stop); else - ret= sbplx_minimize(n,f,f_data,lb,ub,x,minf,opt->dx,&stop); + ret= sbplx_minimize(ni,f,f_data,lb,ub,x,minf,opt->dx,&stop); if (freedx) { free(opt->dx); opt->dx = NULL; } return ret; } @@ -482,7 +484,7 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) nlopt_set_maxeval(local_opt, nlopt_local_search_maxeval); nlopt_set_initial_step(local_opt, opt->dx); } - ret = auglag_minimize(n, f, f_data, + ret = auglag_minimize(ni, f, f_data, opt->m, opt->fc, opt->p, opt->h, lb, ub, x, minf, &stop, @@ -495,11 +497,11 @@ nlopt_result nlopt_optimize(nlopt_opt opt, double *x, double *minf) case NLOPT_GN_ISRES: if (!finite_domain(n, lb, ub)) return NLOPT_INVALID_ARGS; - return isres_minimize(n, f, f_data, - opt->m, opt->fc, - opt->p, opt->h, + return isres_minimize(ni, f, f_data, + (int) (opt->m), opt->fc, + (int) (opt->p), opt->h, lb, ub, x, minf, &stop, - POP(0)); + (int) POP(0)); default: return NLOPT_INVALID_ARGS; diff --git a/api/options.c b/api/options.c index 6e7e1a9..465405f 100644 --- a/api/options.c +++ b/api/options.c @@ -29,8 +29,6 @@ #include "nlopt-internal.h" #include "nlopt-util.h" -#define U(n) ((unsigned int) (n)) - /*************************************************************************/ void nlopt_destroy(nlopt_opt opt) @@ -46,11 +44,11 @@ void nlopt_destroy(nlopt_opt opt) } } -nlopt_opt nlopt_create(nlopt_algorithm algorithm, int n) +nlopt_opt nlopt_create(nlopt_algorithm algorithm, unsigned n) { nlopt_opt opt; - if (n < 0 || algorithm < 0 || algorithm >= NLOPT_NUM_ALGORITHMS) + if (((int) algorithm) < 0 || algorithm >= NLOPT_NUM_ALGORITHMS) return NULL; opt = (nlopt_opt) malloc(sizeof(struct nlopt_opt_s)); @@ -76,11 +74,11 @@ nlopt_opt nlopt_create(nlopt_algorithm algorithm, int n) opt->dx = NULL; if (n > 0) { - opt->lb = (double *) malloc(sizeof(double) * U(n)); + opt->lb = (double *) malloc(sizeof(double) * (n)); if (!opt->lb) goto oom; - opt->ub = (double *) malloc(sizeof(double) * U(n)); + opt->ub = (double *) malloc(sizeof(double) * (n)); if (!opt->ub) goto oom; - opt->xtol_abs = (double *) malloc(sizeof(double) * U(n)); + opt->xtol_abs = (double *) malloc(sizeof(double) * (n)); if (!opt->xtol_abs) goto oom; nlopt_set_lower_bounds1(opt, -HUGE_VAL); nlopt_set_upper_bounds1(opt, +HUGE_VAL); @@ -108,32 +106,32 @@ nlopt_opt nlopt_copy(const nlopt_opt opt) nopt->dx = NULL; if (opt->n > 0) { - nopt->lb = (double *) malloc(sizeof(double) * U(opt->n)); + nopt->lb = (double *) malloc(sizeof(double) * (opt->n)); if (!opt->lb) goto oom; - nopt->ub = (double *) malloc(sizeof(double) * U(opt->n)); + nopt->ub = (double *) malloc(sizeof(double) * (opt->n)); if (!opt->ub) goto oom; - nopt->xtol_abs = (double *) malloc(sizeof(double) * U(opt->n)); + nopt->xtol_abs = (double *) malloc(sizeof(double) * (opt->n)); if (!opt->xtol_abs) goto oom; - memcpy(nopt->lb, opt->lb, sizeof(double) * U(opt->n)); - memcpy(nopt->ub, opt->ub, sizeof(double) * U(opt->n)); - memcpy(nopt->xtol_abs, opt->xtol_abs, sizeof(double) * U(opt->n)); + memcpy(nopt->lb, opt->lb, sizeof(double) * (opt->n)); + memcpy(nopt->ub, opt->ub, sizeof(double) * (opt->n)); + memcpy(nopt->xtol_abs, opt->xtol_abs, sizeof(double) * (opt->n)); } if (opt->m) { nopt->m_alloc = opt->m; nopt->fc = (nlopt_constraint *) malloc(sizeof(nlopt_constraint) - * U(opt->m)); + * (opt->m)); if (!nopt->fc) goto oom; - memcpy(nopt->fc, opt->fc, sizeof(nlopt_constraint) * U(opt->m)); + memcpy(nopt->fc, opt->fc, sizeof(nlopt_constraint) * (opt->m)); } if (opt->p) { nopt->p_alloc = opt->p; nopt->h = (nlopt_constraint *) malloc(sizeof(nlopt_constraint) - * U(opt->p)); + * (opt->p)); if (!nopt->h) goto oom; - memcpy(nopt->h, opt->h, sizeof(nlopt_constraint) * U(opt->p)); + memcpy(nopt->h, opt->h, sizeof(nlopt_constraint) * (opt->p)); } if (opt->local_opt) { @@ -142,9 +140,9 @@ nlopt_opt nlopt_copy(const nlopt_opt opt) } if (opt->dx) { - nopt->dx = (double *) malloc(sizeof(double) * U(opt->n)); + nopt->dx = (double *) malloc(sizeof(double) * (opt->n)); if (!nopt->dx) goto oom; - memcpy(nopt->dx, opt->dx, sizeof(double) * U(opt->n)); + memcpy(nopt->dx, opt->dx, sizeof(double) * (opt->n)); } } return nopt; @@ -170,7 +168,7 @@ nlopt_result nlopt_set_min_objective(nlopt_opt opt, nlopt_func f, void *f_data) nlopt_result nlopt_set_lower_bounds(nlopt_opt opt, const double *lb) { if (opt && (opt->n == 0 || lb)) { - memcpy(opt->lb, lb, sizeof(double) * U(opt->n)); + memcpy(opt->lb, lb, sizeof(double) * (opt->n)); return NLOPT_SUCCESS; } return NLOPT_INVALID_ARGS; @@ -179,7 +177,7 @@ nlopt_result nlopt_set_lower_bounds(nlopt_opt opt, const double *lb) nlopt_result nlopt_set_lower_bounds1(nlopt_opt opt, double lb) { if (opt) { - int i; + unsigned i; for (i = 0; i < opt->n; ++i) opt->lb[i] = lb; return NLOPT_SUCCESS; @@ -190,7 +188,7 @@ nlopt_result nlopt_set_lower_bounds1(nlopt_opt opt, double lb) nlopt_result nlopt_get_lower_bounds(nlopt_opt opt, double *lb) { if (opt && (opt->n == 0 || lb)) { - memcpy(lb, opt->lb, sizeof(double) * U(opt->n)); + memcpy(lb, opt->lb, sizeof(double) * (opt->n)); return NLOPT_SUCCESS; } return NLOPT_INVALID_ARGS; @@ -199,7 +197,7 @@ nlopt_result nlopt_get_lower_bounds(nlopt_opt opt, double *lb) nlopt_result nlopt_set_upper_bounds(nlopt_opt opt, const double *ub) { if (opt && (opt->n == 0 || ub)) { - memcpy(opt->ub, ub, sizeof(double) * U(opt->n)); + memcpy(opt->ub, ub, sizeof(double) * (opt->n)); return NLOPT_SUCCESS; } return NLOPT_INVALID_ARGS; @@ -208,7 +206,7 @@ nlopt_result nlopt_set_upper_bounds(nlopt_opt opt, const double *ub) nlopt_result nlopt_set_upper_bounds1(nlopt_opt opt, double ub) { if (opt) { - int i; + unsigned i; for (i = 0; i < opt->n; ++i) opt->ub[i] = ub; return NLOPT_SUCCESS; @@ -219,7 +217,7 @@ nlopt_result nlopt_set_upper_bounds1(nlopt_opt opt, double ub) nlopt_result nlopt_get_upper_bounds(nlopt_opt opt, double *ub) { if (opt && (opt->n == 0 || ub)) { - memcpy(ub, opt->ub, sizeof(double) * U(opt->n)); + memcpy(ub, opt->ub, sizeof(double) * (opt->n)); return NLOPT_SUCCESS; } return NLOPT_INVALID_ARGS; @@ -240,7 +238,7 @@ nlopt_result nlopt_remove_inequality_constraints(nlopt_opt opt) return NLOPT_SUCCESS; } -static nlopt_result add_constraint(int *m, int *m_alloc, +static nlopt_result add_constraint(unsigned *m, unsigned *m_alloc, nlopt_constraint **c, nlopt_func fc, void *fc_data, double tol) @@ -252,7 +250,7 @@ static nlopt_result add_constraint(int *m, int *m_alloc, *m_alloc = 2 * (*m); *c = (nlopt_constraint *) realloc(*c, sizeof(nlopt_constraint) - * U(*m_alloc)); + * (*m_alloc)); if (!*c) { *m_alloc = *m = 0; return NLOPT_OUT_OF_MEMORY; @@ -343,7 +341,7 @@ nlopt_result nlopt_set_xtol_abs(nlopt_opt opt, const double *xtol_abs) nlopt_result nlopt_set_xtol_abs1(nlopt_opt opt, const double xtol_abs) { if (opt) { - int i; + unsigned i; for (i = 0; i < opt->n; ++i) opt->xtol_abs[i] = xtol_abs; return NLOPT_SUCCESS; @@ -364,7 +362,7 @@ GETSET(maxtime, double, maxtime) /*************************************************************************/ GET(algorithm, nlopt_algorithm, algorithm) -GET(dimension, int, n) +GET(dimension, unsigned, n) /*************************************************************************/ @@ -389,16 +387,16 @@ nlopt_result nlopt_set_local_optimizer(nlopt_opt opt, /*************************************************************************/ -GETSET(population, int, stochastic_population) +GETSET(population, unsigned, stochastic_population) /*************************************************************************/ nlopt_result nlopt_set_initial_step1(nlopt_opt opt, double dx) { - int i; + unsigned i; if (!opt || dx == 0) return NLOPT_INVALID_ARGS; if (!opt->dx && opt->n > 0) { - opt->dx = (double *) malloc(sizeof(double) * U(opt->n)); + opt->dx = (double *) malloc(sizeof(double) * (opt->n)); if (!opt->dx) return NLOPT_OUT_OF_MEMORY; } for (i = 0; i < opt->n; ++i) opt->dx[i] = dx; @@ -407,12 +405,12 @@ nlopt_result nlopt_set_initial_step1(nlopt_opt opt, double dx) nlopt_result nlopt_set_initial_step(nlopt_opt opt, const double *dx) { - int i; + unsigned i; if (!opt || !dx) return NLOPT_INVALID_ARGS; for (i = 0; i < opt->n; ++i) if (dx[i] == 0) return NLOPT_INVALID_ARGS; if (!opt->dx && nlopt_set_initial_step1(opt, 1) == NLOPT_OUT_OF_MEMORY) return NLOPT_OUT_OF_MEMORY; - memcpy(opt->dx, dx, sizeof(double) * U(opt->n)); + memcpy(opt->dx, dx, sizeof(double) * (opt->n)); return NLOPT_SUCCESS; } @@ -425,18 +423,18 @@ nlopt_result nlopt_get_initial_step(const nlopt_opt opt, const double *x, nlopt_opt o = (nlopt_opt) opt; /* discard const temporarily */ nlopt_result ret = nlopt_set_default_initial_step(o, x); if (ret != NLOPT_SUCCESS) return ret; - memcpy(dx, o->dx, sizeof(double) * U(opt->n)); + memcpy(dx, o->dx, sizeof(double) * (opt->n)); free(o->dx); o->dx = NULL; /* don't save, since x-dependent */ } else - memcpy(dx, opt->dx, sizeof(double) * U(opt->n)); + memcpy(dx, opt->dx, sizeof(double) * (opt->n)); return NLOPT_SUCCESS; } nlopt_result nlopt_set_default_initial_step(nlopt_opt opt, const double *x) { const double *lb, *ub; - int i; + unsigned i; if (!opt || !x) return NLOPT_INVALID_ARGS; lb = opt->lb; ub = opt->ub; diff --git a/test/testfuncs.c b/test/testfuncs.c index e589d82..6c89ee1 100644 --- a/test/testfuncs.c +++ b/test/testfuncs.c @@ -11,11 +11,11 @@ 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) +static double testfuncs_status(unsigned n, const double *x, double f) { ++testfuncs_counter; if (testfuncs_verbose) { - int i; + unsigned i; printf("f_%d (%g", testfuncs_counter, x[0]); for (i = 1; i < n; ++i) printf(", %g", x[i]); printf(") = %g\n", f); @@ -30,7 +30,7 @@ static double testfuncs_status(int n, const double *x, double f) #define PI4 12.5663706143592 /* 4*pi */ /****************************************************************************/ -static double rosenbrock_f(int n, const double *x, double *grad, void *data) +static double rosenbrock_f(unsigned n, const double *x, double *grad, void *data) { double a = x[1] - x[0] * x[0], b = 1 - x[0]; UNUSED(data); @@ -46,7 +46,7 @@ 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) +static double mccormic_f(unsigned n, const double *x, double *grad, void *data) { double a = x[0] + x[1], b = x[0] - x[1]; UNUSED(data); @@ -62,9 +62,9 @@ static const double mccormic_ub[2] = {4, 4}; static const double mccormic_xmin[2] = {-0.547197553, -1.54719756}; /****************************************************************************/ -static double boxbetts_f(int n, const double *x, double *grad, void *data) +static double boxbetts_f(unsigned n, const double *x, double *grad, void *data) { - int i; + unsigned i; double f = 0; UNUSED(data); if (grad) @@ -89,9 +89,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) +static double paviani_f(unsigned n, const double *x, double *grad, void *data) { - int i; + unsigned i; double f = 0, prod = 1; UNUSED(data); if (grad) for (i = 0; i < 10; ++i) grad[i] = 0; @@ -115,9 +115,9 @@ static const double paviani_ub[10] = {9.999,9.999,9.999,9.999,9.999,9.999,9.999, static const double paviani_xmin[10] = {9.35026583,9.35026583,9.35026583,9.35026583,9.35026583,9.35026583,9.35026583,9.35026583,9.35026583,9.35026583}; /****************************************************************************/ -static double grosenbrock_f(int n, const double *x, double *grad, void *data) +static double grosenbrock_f(unsigned n, const double *x, double *grad, void *data) { - int i; + unsigned i; double f = 0; UNUSED(data); if (grad) grad[0] = 0; @@ -137,7 +137,7 @@ static const double grosenbrock_ub[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) +static double goldsteinprice_f(unsigned n, const double *x, double *grad, void *data) { double x0, x1, a1, a12, a2, b1, b12, b2; UNUSED(data); @@ -162,7 +162,7 @@ static const double goldsteinprice_ub[2] = {2, 2}; static const double goldsteinprice_xmin[2] = {0, -1}; /****************************************************************************/ -static double shekel_f(int n, const double *x, double *grad, void *data) +static double shekel_f(unsigned n, const double *x, double *grad, void *data) { static const double A[10][4] = { {4,4,4,4}, {1,1,1,1}, @@ -175,10 +175,10 @@ static double shekel_f(int n, const double *x, double *grad, void *data) {6,2,6,2}, {7,3.6,7,3.6} }; static const double c[10] = {.1,.2,.2,.4,.4,.6,.3,.7,.5,.5}; - int i; + unsigned i; double f = 0; if (grad) for (i = 0; i < n; ++i) grad[i] = 0; - int m = *((int *) data); + unsigned m = *((unsigned *) data); for (i = 0; i < m; ++i) { double fi = 1.0 / (c[i] + sqr(x[0]-A[i][0]) @@ -196,7 +196,7 @@ static double shekel_f(int n, const double *x, double *grad, void *data) RETURN(f); } -static int shekel_m[3] = {5,7,10}; +static unsigned shekel_m[3] = {5,7,10}; static const double shekel_lb[4] = {0,0,0,0}; static const double shekel_ub[4] = {10,10,10,10}; static const double shekel0_xmin[4] = {4.000037154,4.000133276,4.000037154,4.000133276}; @@ -204,10 +204,10 @@ static const double shekel1_xmin[4] = {4.000572917,4.000689366,3.999489709,3.999 static const double shekel2_xmin[4] = {4.000746531,4.000592935,3.999663399,3.999509801}; /****************************************************************************/ -static double levy_f(int n, const double *x, double *grad, void *data) +static double levy_f(unsigned n, const double *x, double *grad, void *data) { UNUSED(data); - int i; + unsigned i; double a = x[n-1] - 1, b = 1 + sqr(sin(PI2*x[n-1])); double f = sqr(sin(PI3*x[0])) + a * b; if (grad) { @@ -235,9 +235,9 @@ static const double levy4_ub[4] = {10,10,10,10}; static const double levy4_xmin[4] = {1,1,1,-9.75235596}; /****************************************************************************/ -static double griewank_f(int n, const double *x, double *grad, void *data) +static double griewank_f(unsigned n, const double *x, double *grad, void *data) { - int i; + unsigned i; double f = 1, p = 1; UNUSED(data); for (i = 0; i < n; ++i) { @@ -257,7 +257,7 @@ static const double griewank_ub[10] = {600,600,600,600,600,600,600,600,600,600}; static const double griewank_xmin[10] = {0,0,0,0,0,0,0,0,0,0}; /****************************************************************************/ -static double sixhumpcamel_f(int n, const double *x, double *grad, void *data) +static double sixhumpcamel_f(unsigned n, const double *x, double *grad, void *data) { UNUSED(data); if (grad) { @@ -273,9 +273,9 @@ static const double sixhumpcamel_ub[2] = {5,5}; static const double sixhumpcamel_xmin[2] = {0.08984201317, -0.7126564032}; /****************************************************************************/ -static double convexcosh_f(int n, const double *x, double *grad, void *data) +static double convexcosh_f(unsigned n, const double *x, double *grad, void *data) { - int i; + unsigned i; double f = 1; UNUSED(data); for (i = 0; i < n; ++i) @@ -291,7 +291,7 @@ static const double convexcosh_ub[10] = {2,3,6,7,8,10,11,13,14,16}; static const double convexcosh_xmin[10] = {0,1,2,3,4,5,6,7,8,9}; /****************************************************************************/ -static double branin_f(int n, const double *x, double *grad, void *data) +static double branin_f(unsigned n, const double *x, double *grad, void *data) { double a = 1 - 2*x[1] + 0.05 * sin(PI4 * x[1]) - x[0]; double b = x[1] - 0.5 * sin(PI2 * x[0]); @@ -308,10 +308,10 @@ static const double branin_ub[2] = {10,10}; static const double branin_xmin[2] = {1,0}; /****************************************************************************/ -static double shubert_f(int n, const double *x, double *grad, void *data) +static double shubert_f(unsigned n, const double *x, double *grad, void *data) { UNUSED(data); - int i, j; + unsigned i, j; double f = 0; for (j = 1; j <= 5; ++j) for (i = 0; i < n; ++i) @@ -331,9 +331,9 @@ static const double shubert_ub[2] = {10,10}; static const double shubert_xmin[2] = {-6.774576143, -6.774576143}; /****************************************************************************/ -static double hansen_f(int n, const double *x, double *grad, void *data) +static double hansen_f(unsigned n, const double *x, double *grad, void *data) { - int i; + unsigned i; double a = 0, b = 0; UNUSED(data); for (i = 1; i <= 5; ++i) @@ -358,7 +358,7 @@ static const double hansen_ub[2] = {10,10}; static const double hansen_xmin[2] = {-1.306707704,-1.425128429}; /****************************************************************************/ -static double osc1d_f(int n, const double *x, double *grad, void *data) +static double osc1d_f(unsigned n, const double *x, double *grad, void *data) { double y = *x - 1.23456789; UNUSED(data); @@ -371,7 +371,7 @@ static const double osc1d_ub[1] = {5}; static const double osc1d_xmin[1] = {1.23456789}; /****************************************************************************/ -static double corner4d_f(int n, const double *x, double *grad, void *data) +static double corner4d_f(unsigned n, const double *x, double *grad, void *data) { UNUSED(data); UNUSED(n); @@ -391,7 +391,7 @@ static const double corner4d_ub[4] = {1,1,1,1}; static const double corner4d_xmin[4] = {0,0,0,0}; /****************************************************************************/ -static double side4d_f(int n, const double *x, double *grad, void *data) +static double side4d_f(unsigned n, const double *x, double *grad, void *data) { UNUSED(data); UNUSED(n); diff --git a/util/nlopt-util.h b/util/nlopt-util.h index e52710b..93d4cef 100644 --- a/util/nlopt-util.h +++ b/util/nlopt-util.h @@ -67,7 +67,7 @@ extern void nlopt_sobol_skip(nlopt_sobol s, unsigned n, double *x); /* stopping criteria */ typedef struct { - int n; + unsigned n; double minf_max; double ftol_rel; double ftol_abs; diff --git a/util/stop.c b/util/stop.c index 17ddbf5..9800d90 100644 --- a/util/stop.c +++ b/util/stop.c @@ -45,7 +45,7 @@ int nlopt_stop_f(const nlopt_stopping *s, const double f, double oldf) int nlopt_stop_x(const nlopt_stopping *s, const double *x, const double *oldx) { - int i; + unsigned i; for (i = 0; i < s->n; ++i) if (!relstop(oldx[i], x[i], s->xtol_rel, s->xtol_abs[i])) return 0; @@ -54,7 +54,7 @@ int nlopt_stop_x(const nlopt_stopping *s, const double *x, const double *oldx) int nlopt_stop_dx(const nlopt_stopping *s, const double *x, const double *dx) { - int i; + unsigned i; for (i = 0; i < s->n; ++i) if (!relstop(x[i] - dx[i], x[i], s->xtol_rel, s->xtol_abs[i])) return 0; @@ -72,7 +72,7 @@ int nlopt_stop_xs(const nlopt_stopping *s, const double *xs, const double *oldxs, const double *scale_min, const double *scale_max) { - int i; + unsigned 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]), -- 2.30.2