1 /* Copyright (c) 2007-2008 Massachusetts Institute of Technology
3 * Permission is hereby granted, free of charge, to any person obtaining
4 * a copy of this software and associated documentation files (the
5 * "Software"), to deal in the Software without restriction, including
6 * without limitation the rights to use, copy, modify, merge, publish,
7 * distribute, sublicense, and/or sell copies of the Software, and to
8 * permit persons to whom the Software is furnished to do so, subject to
9 * the following conditions:
11 * The above copyright notice and this permission notice shall be
12 * included in all copies or substantial portions of the Software.
14 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
15 * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
16 * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
17 * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE
18 * LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION
19 * OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
20 * WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
29 #include "nlopt-util.h"
31 #define MIN(a,b) ((a) < (b) ? (a) : (b))
33 /*************************************************************************/
35 int nlopt_isinf(double x) {
36 return fabs(x) >= HUGE_VAL * 0.99
44 static int my_isnan(double x) { return x != x; }
45 # define isnan my_isnan
48 /*************************************************************************/
50 void nlopt_version(int *major, int *minor, int *bugfix)
52 *major = MAJOR_VERSION;
53 *minor = MINOR_VERSION;
54 *bugfix = BUGFIX_VERSION;
57 /*************************************************************************/
59 static const char nlopt_algorithm_names[NLOPT_NUM_ALGORITHMS][256] = {
60 "DIRECT (global, no-derivative)",
61 "DIRECT-L (global, no-derivative)",
62 "Randomized DIRECT-L (global, no-derivative)",
63 "Unscaled DIRECT (global, no-derivative)",
64 "Unscaled DIRECT-L (global, no-derivative)",
65 "Unscaled Randomized DIRECT-L (global, no-derivative)",
66 "Original DIRECT version (global, no-derivative)",
67 "Original DIRECT-L version (global, no-derivative)",
69 "StoGO (global, derivative-based)",
70 "StoGO with randomized search (global, derivative-based)",
72 "StoGO (NOT COMPILED)",
73 "StoGO randomized (NOT COMPILED)",
75 #ifdef WITH_NOCEDAL_LBFGS
76 "original NON-FREE L-BFGS code by Nocedal et al. (local, deriv.-based)",
78 "original NON-FREE L-BFGS code by Nocedal et al. (NOT COMPILED)",
80 "Limited-memory BFGS (L-BFGS) (local, derivative-based)",
81 "Principal-axis, praxis (local, no-derivative)",
82 "Limited-memory variable-metric, rank 1 (local, derivative-based)",
83 "Limited-memory variable-metric, rank 2 (local, derivative-based)",
84 "Truncated Newton (local, derivative-based)",
85 "Truncated Newton with restarting (local, derivative-based)",
86 "Preconditioned truncated Newton (local, derivative-based)",
87 "Preconditioned truncated Newton with restarting (local, derivative-based)",
88 "Controlled random search (CRS2) with local mutation (global, no-derivative)",
89 "Multi-level single-linkage (MLSL), random (global, no-derivative)",
90 "Multi-level single-linkage (MLSL), random (global, derivative)",
91 "Multi-level single-linkage (MLSL), quasi-random (global, no-derivative)",
92 "Multi-level single-linkage (MLSL), quasi-random (global, derivative)",
93 "Method of Moving Asymptotes (MMA) (local, derivative)",
94 "COBYLA (Constrained Optimization BY Linear Approximations) (local, no-derivative)",
95 "NEWUOA unconstrained optimization via quadratic models (local, no-derivative)",
96 "Bound-constrained optimization via NEWUOA-based quadratic models (local, no-derivative)",
97 "Nelder-Mead simplex algorithm (local, no-derivative)",
98 "Sbplx variant of Nelder-Mead (re-implementation of Rowan's Subplex) (local, no-derivative)",
99 "Augmented Lagrangian method (local, no-derivative)",
100 "Augmented Lagrangian method (local, derivative)",
101 "Augmented Lagrangian method for equality constraints (local, no-derivative)",
102 "Augmented Lagrangian method for equality constraints (local, derivative)",
103 "BOBYQA bound-constrained optimization via quadratic models (local, no-derivative)",
104 "ISRES evolutionary constrained optimization (global, no-derivative)",
107 const char *nlopt_algorithm_name(nlopt_algorithm a)
109 if (a < 0 || a >= NLOPT_NUM_ALGORITHMS) return "UNKNOWN";
110 return nlopt_algorithm_names[a];
113 /*************************************************************************/
115 static int nlopt_srand_called = 0;
116 void nlopt_srand(unsigned long seed) {
117 nlopt_srand_called = 1;
118 nlopt_init_genrand(seed);
121 void nlopt_srand_time() {
122 nlopt_srand(nlopt_time_seed());
125 /*************************************************************************/
127 /* crude heuristics for initial step size of nonderivative algorithms */
128 static double initial_step(int n,
129 const double *lb, const double *ub, const double *x)
132 double step = HUGE_VAL;
133 for (i = 0; i < n; ++i) {
134 if (!nlopt_isinf(ub[i]) && !nlopt_isinf(lb[i])
135 && (ub[i] - lb[i]) * 0.25 < step && ub[i] > lb[i])
136 step = (ub[i] - lb[i]) * 0.25;
137 if (!nlopt_isinf(ub[i])
138 && ub[i] - x[i] < step && ub[i] > x[i])
140 if (!nlopt_isinf(lb[i])
141 && x[i] - lb[i] < step && x[i] > lb[i])
144 if (nlopt_isinf(step))
145 for (i = 0; i < n; ++i) {
146 if (!nlopt_isinf(ub[i])
147 && ub[i] - x[i] < step && ub[i] > x[i] + 1e-4)
149 if (!nlopt_isinf(lb[i])
150 && x[i] - lb[i] < step && x[i] > lb[i] + 1e-4)
153 if (nlopt_isinf(step)) {
155 for (i = 0; i < n; ++i)
156 if (fabs(x[i]) * 0.25 > step)
157 step = fabs(x[i]) * 0.25;
164 /*************************************************************************/
169 const double *lb, *ub;
174 static double f_bound(int n, const double *x, void *data_)
177 nlopt_data *data = (nlopt_data *) data_;
180 /* some methods do not support bound constraints, but support
181 discontinuous objectives so we can just return Inf for invalid x */
182 for (i = 0; i < n; ++i)
183 if (x[i] < data->lb[i] || x[i] > data->ub[i])
186 f = data->f(n, x, NULL, data->f_data);
187 return (isnan(f) || nlopt_isinf(f) ? HUGE_VAL : f);
190 static double f_noderiv(int n, const double *x, void *data_)
192 nlopt_data *data = (nlopt_data *) data_;
193 return data->f(n, x, NULL, data->f_data);
198 static double f_direct(int n, const double *x, int *undefined, void *data_)
200 nlopt_data *data = (nlopt_data *) data_;
202 f = data->f(n, x, NULL, data->f_data);
203 *undefined = isnan(f) || nlopt_isinf(f);
214 # include "l-bfgs-b.h"
225 #include "neldermead.h"
230 #define AUGLAG_ALG(a) ((a) == NLOPT_LN_AUGLAG || \
231 (a) == NLOPT_LN_AUGLAG_EQ || \
232 (a) == NLOPT_LD_AUGLAG || \
233 (a) == NLOPT_LD_AUGLAG_EQ)
235 /*************************************************************************/
237 /* for "hybrid" algorithms that combine global search with some
238 local search algorithm, most of the time we anticipate that people
239 will stick with the default local search algorithm, so we
240 don't add this as a parameter to nlopt_minimize. Instead, the user
241 can call nlopt_{set/get}_hybrid_local_algorithm to get/set the defaults. */
243 /* default local-search algorithms */
244 static nlopt_algorithm local_search_alg_deriv = NLOPT_LD_MMA;
245 static nlopt_algorithm local_search_alg_nonderiv = NLOPT_LN_COBYLA;
247 static int local_search_maxeval = -1; /* no maximum by default */
249 void nlopt_get_local_search_algorithm(nlopt_algorithm *deriv,
250 nlopt_algorithm *nonderiv,
253 *deriv = local_search_alg_deriv;
254 *nonderiv = local_search_alg_nonderiv;
255 *maxeval = local_search_maxeval;
258 void nlopt_set_local_search_algorithm(nlopt_algorithm deriv,
259 nlopt_algorithm nonderiv,
262 local_search_alg_deriv = deriv;
263 local_search_alg_nonderiv = nonderiv;
264 local_search_maxeval = maxeval;
267 /*************************************************************************/
269 /* For stochastic algorithms, there is often an initial "population"
270 size to seed the search. Like above, we let the user
271 call nlopt_{set/get}_stochastic_population in order to get/set the
272 defaults. The special stochastic population size of "0" means
273 that the optimization algorithm should pick its default population. */
275 static int stochastic_population = 0;
277 int nlopt_get_stochastic_population(void) { return stochastic_population; }
278 void nlopt_set_stochastic_population(int pop) { stochastic_population = pop; }
280 /*************************************************************************/
282 /* same as nlopt_minimize_econstrained,
283 but xtol_abs is required to be non-NULL */
284 static nlopt_result nlopt_minimize_(
285 nlopt_algorithm algorithm,
286 int n, nlopt_func f, void *f_data,
287 int m, nlopt_func fc, void *fc_data, ptrdiff_t fc_datum_size,
288 int p, nlopt_func h, void *h_data, ptrdiff_t h_datum_size,
289 const double *lb, const double *ub, /* bounds */
290 double *x, /* in: initial guess, out: minimizer */
291 double *minf, /* out: minimum */
292 double minf_max, double ftol_rel, double ftol_abs,
293 double xtol_rel, const double *xtol_abs,
294 double htol_rel, double htol_abs,
295 int maxeval, double maxtime)
301 /* some basic argument checks */
302 if (!minf || !f || n < 0 || m < 0 || p < 0 || (p > 0 && !h)
303 || (m > 0 && !fc)) return NLOPT_INVALID_ARGS;
304 if (n == 0) { /* trivial case: no degrees of freedom */
305 *minf = f(n, x, NULL, f_data);
306 return NLOPT_SUCCESS;
308 else if (n < 0 || !lb || !ub || !x)
309 return NLOPT_INVALID_ARGS;
311 /* nonlinear constraints are only supported with some algorithms */
312 if (m != 0 && algorithm != NLOPT_LD_MMA && algorithm != NLOPT_LN_COBYLA
313 && !AUGLAG_ALG(algorithm) && algorithm != NLOPT_GN_ISRES)
314 return NLOPT_INVALID_ARGS;
316 /* nonlinear equality constraints (h(x) = 0) only via some algorithms */
317 if (p != 0 && !AUGLAG_ALG(algorithm) && algorithm != NLOPT_GN_ISRES)
318 return NLOPT_INVALID_ARGS;
327 /* make sure rand generator is inited */
328 if (!nlopt_srand_called)
329 nlopt_srand_time(); /* default is non-deterministic */
331 /* check bound constraints */
332 for (i = 0; i < n; ++i)
333 if (lb[i] > ub[i] || x[i] < lb[i] || x[i] > ub[i])
334 return NLOPT_INVALID_ARGS;
337 stop.minf_max = (isnan(minf_max)
338 || (nlopt_isinf(minf_max) && minf_max < 0))
339 ? -HUGE_VAL : minf_max;
340 stop.ftol_rel = ftol_rel;
341 stop.ftol_abs = ftol_abs;
342 stop.xtol_rel = xtol_rel;
343 stop.xtol_abs = xtol_abs;
344 stop.htol_rel = htol_rel;
345 stop.htol_abs = htol_abs;
347 stop.maxeval = maxeval;
348 stop.maxtime = maxtime;
349 stop.start = nlopt_seconds();
352 case NLOPT_GN_DIRECT:
353 case NLOPT_GN_DIRECT_L:
354 case NLOPT_GN_DIRECT_L_RAND:
355 return cdirect(n, f, f_data, lb, ub, x, minf, &stop, 0.0,
356 (algorithm != NLOPT_GN_DIRECT)
357 + 3 * (algorithm == NLOPT_GN_DIRECT_L_RAND ? 2 : (algorithm != NLOPT_GN_DIRECT))
358 + 9 * (algorithm == NLOPT_GN_DIRECT_L_RAND ? 1 : (algorithm != NLOPT_GN_DIRECT)));
360 case NLOPT_GN_DIRECT_NOSCAL:
361 case NLOPT_GN_DIRECT_L_NOSCAL:
362 case NLOPT_GN_DIRECT_L_RAND_NOSCAL:
363 return cdirect_unscaled(n, f, f_data, lb, ub, x, minf,
365 (algorithm != NLOPT_GN_DIRECT)
366 + 3 * (algorithm == NLOPT_GN_DIRECT_L_RAND ? 2 : (algorithm != NLOPT_GN_DIRECT))
367 + 9 * (algorithm == NLOPT_GN_DIRECT_L_RAND ? 1 : (algorithm != NLOPT_GN_DIRECT)));
369 case NLOPT_GN_ORIG_DIRECT:
370 case NLOPT_GN_ORIG_DIRECT_L:
371 switch (direct_optimize(f_direct, &d, n, lb, ub, x, minf,
372 maxeval, -1, 0.0, 0.0,
373 pow(xtol_rel, (double) n), -1.0,
376 algorithm == NLOPT_GN_ORIG_DIRECT
378 : DIRECT_GABLONSKY)) {
379 case DIRECT_INVALID_BOUNDS:
380 case DIRECT_MAXFEVAL_TOOBIG:
381 case DIRECT_INVALID_ARGS:
382 return NLOPT_INVALID_ARGS;
383 case DIRECT_INIT_FAILED:
384 case DIRECT_SAMPLEPOINTS_FAILED:
385 case DIRECT_SAMPLE_FAILED:
386 return NLOPT_FAILURE;
387 case DIRECT_MAXFEVAL_EXCEEDED:
388 case DIRECT_MAXITER_EXCEEDED:
389 return NLOPT_MAXEVAL_REACHED;
390 case DIRECT_GLOBAL_FOUND:
391 return NLOPT_MINF_MAX_REACHED;
393 case DIRECT_SIGMATOL:
394 return NLOPT_XTOL_REACHED;
395 case DIRECT_OUT_OF_MEMORY:
396 return NLOPT_OUT_OF_MEMORY;
401 case NLOPT_GD_STOGO_RAND:
403 if (!stogo_minimize(n, f, f_data, x, minf, lb, ub, &stop,
404 algorithm == NLOPT_GD_STOGO
406 return NLOPT_FAILURE;
409 return NLOPT_FAILURE;
413 /* lacking a free/open-source license, we no longer use
414 Rowan's code, and instead use by "sbplx" re-implementation */
415 case NLOPT_LN_SUBPLEX: {
417 double *scale = (double *) malloc(sizeof(double) * n);
418 if (!scale) return NLOPT_OUT_OF_MEMORY;
419 for (i = 0; i < n; ++i)
420 scale[i] = initial_step(1, lb+i, ub+i, x+i);
421 iret = nlopt_subplex(f_bound, minf, x, n, &d, &stop, scale);
424 case -2: return NLOPT_INVALID_ARGS;
425 case -10: return NLOPT_MAXTIME_REACHED;
426 case -1: return NLOPT_MAXEVAL_REACHED;
427 case 0: return NLOPT_XTOL_REACHED;
428 case 1: return NLOPT_SUCCESS;
429 case 2: return NLOPT_MINF_MAX_REACHED;
430 case 20: return NLOPT_FTOL_REACHED;
431 case -200: return NLOPT_OUT_OF_MEMORY;
432 default: return NLOPT_FAILURE; /* unknown return code */
438 case NLOPT_LN_PRAXIS:
439 return praxis_(0.0, DBL_EPSILON,
440 initial_step(n, lb, ub, x), n, x, f_bound, &d,
444 case NLOPT_LD_LBFGS_NOCEDAL: {
445 int iret, *nbd = (int *) malloc(sizeof(int) * n);
446 if (!nbd) return NLOPT_OUT_OF_MEMORY;
447 for (i = 0; i < n; ++i) {
448 int linf = nlopt_isinf(lb[i]) && lb[i] < 0;
449 int uinf = nlopt_isinf(ub[i]) && ub[i] > 0;
450 nbd[i] = linf && uinf ? 0 : (uinf ? 1 : (linf ? 3 : 2));
452 iret = lbfgsb_minimize(n, f, f_data, x, nbd, lb, ub,
453 MIN(n, 5), 0.0, ftol_rel,
454 xtol_abs ? *xtol_abs : xtol_rel,
459 case -1: return NLOPT_INVALID_ARGS;
460 case -2: default: return NLOPT_FAILURE;
464 *minf = f(n, x, NULL, f_data);
466 case 5: return NLOPT_MAXEVAL_REACHED;
467 case 2: return NLOPT_XTOL_REACHED;
468 case 1: return NLOPT_FTOL_REACHED;
469 default: return NLOPT_SUCCESS;
477 return luksan_plis(n, f, f_data, lb, ub, x, minf, &stop);
481 return luksan_plip(n, f, f_data, lb, ub, x, minf, &stop,
482 algorithm == NLOPT_LD_VAR1 ? 1 : 2);
484 case NLOPT_LD_TNEWTON:
485 case NLOPT_LD_TNEWTON_RESTART:
486 case NLOPT_LD_TNEWTON_PRECOND:
487 case NLOPT_LD_TNEWTON_PRECOND_RESTART:
488 return luksan_pnet(n, f, f_data, lb, ub, x, minf, &stop,
489 1 + (algorithm - NLOPT_LD_TNEWTON) % 2,
490 1 + (algorithm - NLOPT_LD_TNEWTON) / 2);
492 case NLOPT_GN_CRS2_LM:
493 return crs_minimize(n, f, f_data, lb, ub, x, minf, &stop,
494 stochastic_population, 0);
498 case NLOPT_GN_MLSL_LDS:
499 case NLOPT_GD_MLSL_LDS:
500 for (i = 0; i < n && stop.xtol_abs[i] <= 0; ++i) ;
501 if (stop.ftol_rel <= 0 && stop.ftol_abs <= 0 &&
502 stop.xtol_rel <= 0 && i == n) {
503 /* it is not sensible to call MLSL without *some*
504 nonzero tolerance for the local search */
505 stop.ftol_rel = 1e-15;
506 stop.xtol_rel = 1e-7;
508 return mlsl_minimize(n, f, f_data, lb, ub, x, minf, &stop,
509 (algorithm == NLOPT_GN_MLSL ||
510 algorithm == NLOPT_GN_MLSL_LDS)
511 ? local_search_alg_nonderiv
512 : local_search_alg_deriv,
513 local_search_maxeval,
514 stochastic_population,
515 algorithm >= NLOPT_GN_MLSL_LDS);
518 return mma_minimize(n, f, f_data, m, fc, fc_data, fc_datum_size,
519 lb, ub, x, minf, &stop,
520 local_search_alg_deriv, 1e-12, 100000);
522 case NLOPT_LN_COBYLA:
523 return cobyla_minimize(n, f, f_data,
524 m, fc, fc_data, fc_datum_size,
525 lb, ub, x, minf, &stop,
526 initial_step(n, lb, ub, x));
528 case NLOPT_LN_NEWUOA:
529 return newuoa(n, 2*n+1, x, 0, 0, initial_step(n, lb, ub, x),
530 &stop, minf, f_noderiv, &d);
532 case NLOPT_LN_NEWUOA_BOUND:
533 return newuoa(n, 2*n+1, x, lb, ub, initial_step(n, lb, ub, x),
534 &stop, minf, f_noderiv, &d);
536 case NLOPT_LN_BOBYQA:
537 return bobyqa(n, 2*n+1, x, lb, ub, initial_step(n, lb, ub, x),
538 &stop, minf, f_noderiv, &d);
540 case NLOPT_LN_NELDERMEAD:
544 double *scale = (double *) malloc(sizeof(double) * n);
545 if (!scale) return NLOPT_OUT_OF_MEMORY;
546 for (i = 0; i < n; ++i)
547 scale[i] = initial_step(1, lb+i, ub+i, x+i);
548 if (algorithm == NLOPT_LN_NELDERMEAD)
549 ret = nldrmd_minimize(n,f,f_data,lb,ub,x,minf,scale,&stop);
551 ret = sbplx_minimize(n,f,f_data,lb,ub,x,minf,scale,&stop);
556 case NLOPT_LN_AUGLAG:
557 case NLOPT_LN_AUGLAG_EQ:
558 return auglag_minimize(n, f, f_data,
559 m, fc, fc_data, fc_datum_size,
560 p, h, h_data, h_datum_size,
561 lb, ub, x, minf, &stop,
562 local_search_alg_nonderiv,
563 algorithm == NLOPT_LN_AUGLAG_EQ);
565 case NLOPT_LD_AUGLAG:
566 case NLOPT_LD_AUGLAG_EQ:
567 return auglag_minimize(n, f, f_data,
568 m, fc, fc_data, fc_datum_size,
569 p, h, h_data, h_datum_size,
570 lb, ub, x, minf, &stop,
571 local_search_alg_deriv,
572 algorithm == NLOPT_LD_AUGLAG_EQ);
575 return isres_minimize(n, f, f_data,
576 m, fc, fc_data, fc_datum_size,
577 p, h, h_data, h_datum_size,
578 lb, ub, x, minf, &stop,
579 stochastic_population);
582 return NLOPT_INVALID_ARGS;
585 return NLOPT_SUCCESS;
588 nlopt_result nlopt_minimize_econstrained(
589 nlopt_algorithm algorithm,
590 int n, nlopt_func f, void *f_data,
591 int m, nlopt_func fc, void *fc_data, ptrdiff_t fc_datum_size,
592 int p, nlopt_func h, void *h_data, ptrdiff_t h_datum_size,
593 const double *lb, const double *ub, /* bounds */
594 double *x, /* in: initial guess, out: minimizer */
595 double *minf, /* out: minimum */
596 double minf_max, double ftol_rel, double ftol_abs,
597 double xtol_rel, const double *xtol_abs,
598 double htol_rel, double htol_abs,
599 int maxeval, double maxtime)
603 ret = nlopt_minimize_(algorithm, n, f, f_data,
604 m, fc, fc_data, fc_datum_size,
605 p, h, h_data, h_datum_size,
607 x, minf, minf_max, ftol_rel, ftol_abs,
608 xtol_rel, xtol_abs, htol_rel, htol_abs,
612 double *xtol = (double *) malloc(sizeof(double) * n);
613 if (!xtol) return NLOPT_OUT_OF_MEMORY;
614 for (i = 0; i < n; ++i) xtol[i] = -1;
615 ret = nlopt_minimize_(algorithm, n, f, f_data,
616 m, fc, fc_data, fc_datum_size,
617 p, h, h_data, h_datum_size,
619 x, minf, minf_max, ftol_rel, ftol_abs,
620 xtol_rel, xtol, htol_rel, htol_abs,
627 nlopt_result nlopt_minimize_constrained(
628 nlopt_algorithm algorithm,
629 int n, nlopt_func f, void *f_data,
630 int m, nlopt_func fc, void *fc_data, ptrdiff_t fc_datum_size,
631 const double *lb, const double *ub, /* bounds */
632 double *x, /* in: initial guess, out: minimizer */
633 double *minf, /* out: minimum */
634 double minf_max, double ftol_rel, double ftol_abs,
635 double xtol_rel, const double *xtol_abs,
636 int maxeval, double maxtime)
638 return nlopt_minimize_econstrained(
639 algorithm, n, f, f_data,
640 m, fc, fc_data, fc_datum_size, 0, NULL, NULL, 0,
641 lb, ub, x, minf, minf_max, ftol_rel, ftol_abs,
642 xtol_rel, xtol_abs, ftol_rel, ftol_abs, maxeval, maxtime);
645 nlopt_result nlopt_minimize(
646 nlopt_algorithm algorithm,
647 int n, nlopt_func f, void *f_data,
648 const double *lb, const double *ub, /* bounds */
649 double *x, /* in: initial guess, out: minimizer */
650 double *minf, /* out: minimum */
651 double minf_max, double ftol_rel, double ftol_abs,
652 double xtol_rel, const double *xtol_abs,
653 int maxeval, double maxtime)
655 return nlopt_minimize_constrained(
656 algorithm, n, f, f_data, 0, NULL, NULL, 0,
657 lb, ub, x, minf, minf_max, ftol_rel, ftol_abs,
658 xtol_rel, xtol_abs, maxeval, maxtime);