chiark / gitweb /
use NLopt stop criteria in Luksan code: in particular, absolute tolerances were missi...
[nlopt.git] / luksan / plip.c
1 #include <limits.h>
2 #include <stdlib.h>
3 #include <math.h>
4 #include <string.h>
5 #include "luksan.h"
6
7 #define max(a,b) ((a) > (b) ? (a) : (b))
8 #define min(a,b) ((a) < (b) ? (a) : (b))
9
10 /* *********************************************************************** */
11 /* SUBROUTINE PLIP               ALL SYSTEMS                   01/09/22 */
12 /* PURPOSE : */
13 /* GENERAL SUBROUTINE FOR LARGE-SCALE BOX CONSTRAINED MINIMIZATION THAT */
14 /* USE THE SHIFTED LIMITED MEMORY VARIABLE METRIC METHOD BASED ON THE */
15 /* PRODUCT FORM UPDATES. */
16
17 /* PARAMETERS : */
18 /*  II  NF  NUMBER OF VARIABLES. */
19 /*  II  NB  CHOICE OF SIMPLE BOUNDS. NB=0-SIMPLE BOUNDS SUPPRESSED. */
20 /*         NB>0-SIMPLE BOUNDS ACCEPTED. */
21 /*  RI  X(NF)  VECTOR OF VARIABLES. */
22 /*  II  IX(NF)  VECTOR CONTAINING TYPES OF BOUNDS. IX(I)=0-VARIABLE */
23 /*         X(I) IS UNBOUNDED. IX(I)=1-LOVER BOUND XL(I).LE.X(I). */
24 /*         IX(I)=2-UPPER BOUND X(I).LE.XU(I). IX(I)=3-TWO SIDE BOUND */
25 /*         XL(I).LE.X(I).LE.XU(I). IX(I)=5-VARIABLE X(I) IS FIXED. */
26 /*  RI  XL(NF)  VECTOR CONTAINING LOWER BOUNDS FOR VARIABLES. */
27 /*  RI  XU(NF)  VECTOR CONTAINING UPPER BOUNDS FOR VARIABLES. */
28 /*  RA  GF(NF)  GRADIENT OF THE OBJECTIVE FUNCTION. */
29 /*  RO  S(NF)  DIRECTION VECTOR. */
30 /*  RU  XO(NF)  VECTORS OF VARIABLES DIFFERENCE. */
31 /*  RI  GO(NF)  GRADIENTS DIFFERENCE. */
32 /*  RA  SO(NF)  AUXILIARY VECTOR. */
33 /*  RA  XM(NF*MF)  AUXILIARY VECTOR. */
34 /*  RA  XR(MF)  AUXILIARY VECTOR. */
35 /*  RA  GR(MF)  AUXILIARY VECTOR. */
36 /*  RI  XMAX  MAXIMUM STEPSIZE. */
37 /*  RI  TOLX  TOLERANCE FOR CHANGE OF VARIABLES. */
38 /*  RI  TOLF  TOLERANCE FOR CHANGE OF FUNCTION VALUES. */
39 /*  RI  TOLB  TOLERANCE FOR THE FUNCTION VALUE. */
40 /*  RI  TOLG  TOLERANCE FOR THE GRADIENT NORM. */
41 /*  RI  MINF_EST  ESTIMATION OF THE MINIMUM FUNCTION VALUE. */
42 /*  RO  GMAX  MAXIMUM PARTIAL DERIVATIVE. */
43 /*  RO  F  VALUE OF THE OBJECTIVE FUNCTION. */
44 /*  II  MIT  MAXIMUM NUMBER OF ITERATIONS. */
45 /*  II  MFV  MAXIMUM NUMBER OF FUNCTION EVALUATIONS. */
46 /*  II  IEST  ESTIMATION INDICATOR. IEST=0-MINIMUM IS NOT ESTIMATED. */
47 /*         IEST=1-MINIMUM IS ESTIMATED BY THE VALUE MINF_EST. */
48 /*  II  MET  METHOD USED. MET=1-RANK-ONE METHOD. MET=2-RANK-TWO */
49 /*         METHOD. */
50 /*  II  MF  NUMBER OF LIMITED MEMORY STEPS. */
51 /*  IO  ITERM  VARIABLE THAT INDICATES THE CAUSE OF TERMINATION. */
52 /*         ITERM=1-IF ABS(X-XO) WAS LESS THAN OR EQUAL TO TOLX IN */
53 /*                   MTESX (USUALLY TWO) SUBSEQUEBT ITERATIONS. */
54 /*         ITERM=2-IF ABS(F-FO) WAS LESS THAN OR EQUAL TO TOLF IN */
55 /*                   MTESF (USUALLY TWO) SUBSEQUEBT ITERATIONS. */
56 /*         ITERM=3-IF F IS LESS THAN OR EQUAL TO TOLB. */
57 /*         ITERM=4-IF GMAX IS LESS THAN OR EQUAL TO TOLG. */
58 /*         ITERM=6-IF THE TERMINATION CRITERION WAS NOT SATISFIED, */
59 /*                   BUT THE SOLUTION OBTAINED IS PROBABLY ACCEPTABLE. */
60 /*         ITERM=11-IF NIT EXCEEDED MIT. ITERM=12-IF NFV EXCEEDED MFV. */
61 /*         ITERM=13-IF NFG EXCEEDED MFG. ITERM<0-IF THE METHOD FAILED. */
62
63 /* VARIABLES IN COMMON /STAT/ (STATISTICS) : */
64 /*  IO  NRES  NUMBER OF RESTARTS. */
65 /*  IO  NDEC  NUMBER OF MATRIX DECOMPOSITION. */
66 /*  IO  NIN  NUMBER OF INNER ITERATIONS. */
67 /*  IO  NIT  NUMBER OF ITERATIONS. */
68 /*  IO  NFV  NUMBER OF FUNCTION EVALUATIONS. */
69 /*  IO  NFG  NUMBER OF GRADIENT EVALUATIONS. */
70 /*  IO  NFH  NUMBER OF HESSIAN EVALUATIONS. */
71
72 /* SUBPROGRAMS USED : */
73 /*  S   PCBS04  ELIMINATION OF BOX CONSTRAINT VIOLATIONS. */
74 /*  S   PS1L01  STEPSIZE SELECTION USING LINE SEARCH. */
75 /*  S   PULSP3  SHIFTED VARIABLE METRIC UPDATE. */
76 /*  S   PULVP3  SHIFTED LIMITED-MEMORY VARIABLE METRIC UPDATE. */
77 /*  S   PYADC0  ADDITION OF A BOX CONSTRAINT. */
78 /*  S   PYFUT1  TEST ON TERMINATION. */
79 /*  S   PYRMC0  DELETION OF A BOX CONSTRAINT. */
80 /*  S   PYTRCD  COMPUTATION OF PROJECTED DIFFERENCES FOR THE VARIABLE METRIC */
81 /*         UPDATE. */
82 /*  S   PYTRCG  COMPUTATION OF THE PROJECTED GRADIENT. */
83 /*  S   PYTRCS  COMPUTATION OF THE PROJECTED DIRECTION VECTOR. */
84 /*  S   MXDRMM  MULTIPLICATION OF A ROWWISE STORED DENSE RECTANGULAR */
85 /*         MATRIX A BY A VECTOR X. */
86 /*  S   MXDCMD  MULTIPLICATION OF A COLUMNWISE STORED DENSE RECTANGULAR */
87 /*         MATRIX A BY A VECTOR X AND ADDITION OF THE SCALED VECTOR */
88 /*         ALF*Y. */
89 /*  S   MXUCOP  COPYING OF A VECTOR. */
90 /*  S   MXUDIR  VECTOR AUGMENTED BY THE SCALED VECTOR. */
91 /*  RF  MXUDOT  DOT PRODUCT OF TWO VECTORS. */
92 /*  S   MXUNEG  COPYING OF A VECTOR WITH CHANGE OF THE SIGN. */
93 /*  S   MXUZER  VECTOR ELEMENTS CORRESPONDING TO ACTIVE BOUNDS ARE SET */
94 /*         TO ZERO. */
95 /*  S   MXVCOP  COPYING OF A VECTOR. */
96
97 /* EXTERNAL SUBROUTINES : */
98 /*  SE  OBJ  COMPUTATION OF THE VALUE OF THE OBJECTIVE FUNCTION. */
99 /*         CALLING SEQUENCE: CALL OBJ(NF,X,FF) WHERE NF IS THE NUMBER */
100 /*         OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND FF IS */
101 /*         THE VALUE OF THE OBJECTIVE FUNCTION. */
102 /*  SE  DOBJ  COMPUTATION OF THE GRADIENT OF THE OBJECTIVE FUNCTION. */
103 /*         CALLING SEQUENCE: CALL DOBJ(NF,X,GF) WHERE NF IS THE NUMBER */
104 /*         OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND GF(NF) */
105 /*         IS THE GRADIENT OF THE OBJECTIVE FUNCTION. */
106 /* -- OBJ and DOBJ are replaced by a single function, objgrad, in NLopt */
107
108 /* METHOD : */
109 /* HYBRID METHOD WITH SPARSE MARWIL UPDATES FOR SPARSE LEAST SQUARES */
110 /* PROBLEMS. */
111
112 static void plip_(int *nf, int *nb, double *x, int *
113                   ix, double *xl, double *xu, double *gf, double *s, 
114                   double *xo, double *go, double *so, double *xm, 
115                   double *xr, double *gr, double *xmax, double *tolx, 
116                   double *tolf, double *tolb, double *tolg, 
117                   nlopt_stopping *stop, double *
118                   minf_est, double *gmax, double *f, int *mit, int *mfv, 
119                   int *iest, int *met, int *mf, 
120                   int *iterm, stat_common *stat_1,
121                   nlopt_func objgrad, void *objgrad_data)
122 {
123     /* System generated locals */
124     int i__1;
125     double d__1, d__2;
126
127     /* Builtin functions */
128
129     /* Local variables */
130     int i__, n;
131     double p, r__;
132     int kd, ld;
133     double fo, fp;
134     int nn;
135     double po, pp, ro, rp;
136     int kbf, mec, mfg;
137     double par;
138     int mes, kit;
139     double alf1, alf2, eta0, eta9, par1, par2;
140     int mes1, mes2, mes3, met3;
141     double eps8, eps9;
142     int meta, mred, nred, iold;
143     double maxf, dmax__;
144     int xstop = 0;
145     int inew;
146     double told;
147     int ites;
148     double rmin, rmax, umax, tolp, tols;
149     int isys;
150     int ires1, ires2;
151     int iterd, iterh, mtesf, ntesf;
152     double gnorm;
153     int iters, irest, inits, kters, maxst;
154     double snorm;
155     int mtesx, ntesx;
156
157 /*     INITIATION */
158
159     /* Parameter adjustments */
160     --gr;
161     --xr;
162     --xm;
163     --so;
164     --go;
165     --xo;
166     --s;
167     --gf;
168     --xu;
169     --xl;
170     --ix;
171     --x;
172
173     /* Function Body */
174     kbf = 0;
175     if (*nb > 0) {
176         kbf = 2;
177     }
178     stat_1->nres = 0;
179     stat_1->ndec = 0;
180     stat_1->nin = 0;
181     stat_1->nit = 0;
182     stat_1->nfg = 0;
183     stat_1->nfh = 0;
184     isys = 0;
185     ites = 1;
186     mtesx = 2;
187     mtesf = 2;
188     inits = 2;
189     *iterm = 0;
190     iterd = 0;
191     iters = 2;
192     iterh = 0;
193     kters = 3;
194     irest = 0;
195     ires1 = 999;
196     ires2 = 0;
197     mred = 10;
198     meta = 1;
199     met3 = 4;
200     mec = 4;
201     mes = 4;
202     mes1 = 2;
203     mes2 = 2;
204     mes3 = 2;
205     eta0 = 1e-15;
206     eta9 = 1e120;
207     eps8 = 1.;
208     eps9 = 1e-8;
209     alf1 = 1e-10;
210     alf2 = 1e10;
211     rmax = eta9;
212     dmax__ = eta9;
213     maxf = 1e20;
214     if (*iest <= 0) {
215          *minf_est = -HUGE_VAL; /* changed from -1e60 by SGJ */
216     }
217     if (*iest > 0) {
218         *iest = 1;
219     }
220     if (*xmax <= 0.) {
221         *xmax = 1e16;
222     }
223     if (*tolx <= 0.) {
224         *tolx = 1e-16;
225     }
226     if (*tolf <= 0.) {
227         *tolf = 1e-14;
228     }
229     if (*tolg <= 0.) {
230          *tolg = 1e-8; /* SGJ: was 1e-6, but this sometimes stops too soon */
231     }
232 #if 0
233     /* removed by SGJ: this check prevented us from using minf_max <= 0,
234        which doesn't make sense.  Instead, if you don't want to have a
235        lower limit, you should set minf_max = -HUGE_VAL */
236     if (*tolb <= 0.) {
237         *tolb = *minf_est + 1e-16;
238     }
239 #endif
240     told = 1e-4;
241     tols = 1e-4;
242     tolp = .9;
243     if (*met <= 0) {
244         *met = 2;
245     }
246     /* changed by SGJ: default is no limit (INT_MAX) on # iterations/fevals */
247     if (*mit <= 0) {
248         *mit = INT_MAX;
249     }
250     if (*mfv <= 0) {
251         *mfv = INT_MAX;
252     }
253     mfg = *mfv;
254     kd = 1;
255     ld = -1;
256     kit = -(ires1 * *nf + ires2);
257     fo = *minf_est;
258
259 /*     INITIAL OPERATIONS WITH SIMPLE BOUNDS */
260
261     if (kbf > 0) {
262         i__1 = *nf;
263         for (i__ = 1; i__ <= i__1; ++i__) {
264             if ((ix[i__] == 3 || ix[i__] == 4) && xu[i__] <= xl[i__]) {
265                 xu[i__] = xl[i__];
266                 ix[i__] = 5;
267             } else if (ix[i__] == 5 || ix[i__] == 6) {
268                 xl[i__] = x[i__];
269                 xu[i__] = x[i__];
270                 ix[i__] = 5;
271             }
272 /* L2: */
273         }
274         luksan_pcbs04__(nf, &x[1], &ix[1], &xl[1], &xu[1], &eps9, &kbf);
275         luksan_pyadc0__(nf, &n, &x[1], &ix[1], &xl[1], &xu[1], &inew);
276     }
277     if (*iterm != 0) {
278         goto L11190;
279     }
280     if (nlopt_stop_time(stop)) { *iterm = 100; goto L11190; }
281     *f = objgrad(*nf, &x[1], &gf[1], objgrad_data);
282     ++stop->nevals;
283     ++stat_1->nfg;
284 L11120:
285     luksan_pytrcg__(nf, nf, &ix[1], &gf[1], &umax, gmax, &kbf, &iold);
286     luksan_pyfut1__(nf, f, &fo, &umax, gmax, xstop, stop, tolg, 
287             &kd, &stat_1->nit, &kit, mit, &stat_1->nfg, &mfg, 
288             &ntesx, &mtesx, &ntesf, &mtesf, &ites, &ires1, &ires2, &irest, &
289             iters, iterm);
290     if (*iterm != 0) {
291         goto L11190;
292     }
293     if (nlopt_stop_time(stop)) { *iterm = 100; goto L11190; }
294     if (kbf > 0 && rmax > 0.) {
295         luksan_pyrmc0__(nf, &n, &ix[1], &gf[1], &eps8, &umax, gmax, &rmax, &
296                 iold, &irest);
297     }
298 L11130:
299     if (irest > 0) {
300         nn = 0;
301         par = 1.;
302         ld = min(ld,1);
303         if (kit < stat_1->nit) {
304             ++stat_1->nres;
305             kit = stat_1->nit;
306         } else {
307             *iterm = -10;
308             if (iters < 0) {
309                 *iterm = iters - 5;
310             }
311         }
312     }
313     if (*iterm != 0) {
314         goto L11190;
315     }
316     if (nlopt_stop_time(stop)) { *iterm = 100; goto L11190; }
317
318 /*     DIRECTION DETERMINATION */
319
320     gnorm = sqrt(luksan_mxudot__(nf, &gf[1], &gf[1], &ix[1], &kbf));
321
322 /*     NEWTON LIKE STEP */
323
324     luksan_mxuneg__(nf, &gf[1], &s[1], &ix[1], &kbf);
325     luksan_mxdrmm__(nf, &nn, &xm[1], &s[1], &gr[1]);
326     luksan_mxdcmd__(nf, &nn, &xm[1], &gr[1], &par, &s[1], &s[1]);
327     luksan_mxuzer__(nf, &s[1], &ix[1], &kbf);
328     iterd = 1;
329     snorm = sqrt(luksan_mxudot__(nf, &s[1], &s[1], &ix[1], &kbf));
330
331 /*     TEST ON DESCENT DIRECTION AND PREPARATION OF LINE SEARCH */
332
333     if (kd > 0) {
334         p = luksan_mxudot__(nf, &gf[1], &s[1], &ix[1], &kbf);
335     }
336     if (iterd < 0) {
337         *iterm = iterd;
338     } else {
339
340 /*     TEST ON DESCENT DIRECTION */
341
342         if (snorm <= 0.) {
343             irest = max(irest,1);
344         } else if (p + told * gnorm * snorm <= 0.) {
345             irest = 0;
346         } else {
347
348 /*     UNIFORM DESCENT CRITERION */
349
350             irest = max(irest,1);
351         }
352         if (irest == 0) {
353
354 /*     PREPARATION OF LINE SEARCH */
355
356             nred = 0;
357             rmin = alf1 * gnorm / snorm;
358 /* Computing MIN */
359             d__1 = alf2 * gnorm / snorm, d__2 = *xmax / snorm;
360             rmax = min(d__1,d__2);
361         }
362     }
363     if (*iterm != 0) {
364         goto L11190;
365     }
366     if (nlopt_stop_time(stop)) { *iterm = 100; goto L11190; }
367     if (irest != 0) {
368         goto L11130;
369     }
370     luksan_pytrcs__(nf, &x[1], &ix[1], &xo[1], &xl[1], &xu[1], &gf[1], &go[1],
371              &s[1], &ro, &fp, &fo, f, &po, &p, &rmax, &eta9, &kbf);
372     if (rmax == 0.) {
373         goto L11175;
374     }
375 L11170:
376     luksan_ps1l01__(&r__, &rp, f, &fo, &fp, &p, &po, &pp, minf_est, &maxf, &rmin, 
377             &rmax, &tols, &tolp, &par1, &par2, &kd, &ld, &stat_1->nit, &kit, &
378             nred, &mred, &maxst, iest, &inits, &iters, &kters, &mes, &isys);
379     if (isys == 0) {
380         goto L11174;
381     }
382     luksan_mxudir__(nf, &r__, &s[1], &xo[1], &x[1], &ix[1], &kbf);
383     luksan_pcbs04__(nf, &x[1], &ix[1], &xl[1], &xu[1], &eps9, &kbf);
384     *f = objgrad(*nf, &x[1], &gf[1], objgrad_data);
385     ++stop->nevals;
386     ++stat_1->nfg;
387     p = luksan_mxudot__(nf, &gf[1], &s[1], &ix[1], &kbf);
388     goto L11170;
389 L11174:
390     if (iters <= 0) {
391         r__ = 0.;
392         *f = fo;
393         p = po;
394         luksan_mxvcop__(nf, &xo[1], &x[1]);
395         luksan_mxvcop__(nf, &go[1], &gf[1]);
396         irest = max(irest,1);
397         ld = kd;
398         goto L11130;
399     }
400     luksan_mxuneg__(nf, &go[1], &s[1], &ix[1], &kbf);
401     luksan_pytrcd__(nf, &x[1], &ix[1], &xo[1], &gf[1], &go[1], &r__, f, &fo, &
402             p, &po, &dmax__, &kbf, &kd, &ld, &iters);
403     xstop = nlopt_stop_dx(stop, &x[1], &xo[1]);
404     luksan_mxucop__(nf, &gf[1], &so[1], &ix[1], &kbf);
405     if (nn < *mf) {
406         luksan_pulsp3__(nf, &nn, mf, &xm[1], &gr[1], &xo[1], &go[1], &r__, &
407                 po, &par, &iterh, &met3);
408     } else {
409         luksan_pulvp3__(nf, &nn, &xm[1], &xr[1], &gr[1], &s[1], &so[1], &xo[1]
410                 , &go[1], &r__, &po, &par, &iterh, &mec, &met3, met);
411     }
412 L11175:
413     if (iterh != 0) {
414         irest = max(irest,1);
415     }
416     if (kbf > 0) {
417         luksan_pyadc0__(nf, &n, &x[1], &ix[1], &xl[1], &xu[1], &inew);
418     }
419     goto L11120;
420 L11190:
421     return;
422 } /* plip_ */
423
424 /* NLopt wrapper around plip_, handling dynamic allocation etc. */
425 nlopt_result luksan_plip(int n, nlopt_func f, void *f_data,
426                          const double *lb, const double *ub, /* bounds */
427                          double *x, /* in: initial guess, out: minimizer */
428                          double *minf,
429                          nlopt_stopping *stop,
430                          int method) /* 1 or 2, see below */
431 {
432      int i, *ix, nb = 1;
433      double *work, *xl, *xu, *gf, *s, *xo, *go, *so, *xm, *xr, *gr;
434      double gmax, minf_est;
435      double xmax = 0; /* no maximum */
436      double tolg = 0; /* default gradient tolerance */
437      int iest = 0; /* we have no estimate of min function value */
438      int mit = 0; /* default no limit on #iterations */
439      int mfv = stop->maxeval;
440      stat_common stat;
441      int iterm;
442      int mf;
443
444      ix = (int*) malloc(sizeof(int) * n);
445      if (!ix) return NLOPT_OUT_OF_MEMORY;
446
447      /* FIXME: what should we set mf to?  The example program tlis.for
448         sets it to zero as far as I can tell, but it seems to greatly
449         improve convergence to make it > 0.  The computation time
450         per iteration, and of course the memory, seem to go as O(n * mf),
451         and we'll assume that the main limiting factor is the memory.
452         We'll assume that at least MEMAVAIL memory, or 4*n memory, whichever
453         is bigger, is available. */
454      mf = max(MEMAVAIL/n, 4);
455      if (stop->maxeval && stop->maxeval <= mf)
456           mf = max(stop->maxeval - 5, 1); /* mf > maxeval seems not good */
457
458  retry_alloc:
459      work = (double*) malloc(sizeof(double) * (n * 7 + max(n,n*mf) + 
460                                                max(n,mf)*2));
461      if (!work) { 
462           if (mf > 0) {
463                mf = 0; /* allocate minimal memory */
464                goto retry_alloc;
465           }
466           free(ix);
467           return NLOPT_OUT_OF_MEMORY;
468      }
469
470      xl = work; xu = xl + n;
471      gf = xu + n; s = gf + n; xo = s + n; go = xo + n; so = go + n;
472      xm = so + n;
473      xr = xm + max(n*mf,n); gr = xr + max(n,mf);
474
475      for (i = 0; i < n; ++i) {
476           int lbu = lb[i] <= -0.99 * HUGE_VAL; /* lb unbounded */
477           int ubu = ub[i] >= 0.99 * HUGE_VAL;  /* ub unbounded */
478           ix[i] = lbu ? (ubu ? 0 : 2) : (ubu ? 1 : (lb[i] == ub[i] ? 5 : 3));
479           xl[i] = lb[i];
480           xu[i] = ub[i];
481      }
482
483      /* ?  xo does not seem to be initialized in the
484         original Fortran code, but it is used upon
485         input to plip if mf > 0 ... perhaps ALLOCATE initializes
486         arrays to zero by default? */
487      memset(xo, 0, sizeof(double) * max(n,n*mf));
488
489      plip_(&n, &nb, x, ix, xl, xu, 
490            gf, s, xo, go, so, xm, xr, gr,
491            &xmax,
492
493            /* fixme: pass tol_rel and tol_abs and use NLopt check */
494            &stop->xtol_rel,
495            &stop->ftol_rel,
496            &stop->minf_max,
497            &tolg,
498            stop,
499
500            &minf_est, &gmax,
501            minf,
502            &mit, &mfv,
503            &iest,
504            &method, /* 1 == rank-one method VAR1, 2 == rank-two method VAR2 */
505            &mf,
506            &iterm, &stat,
507            f, f_data);
508
509      free(work);
510      free(ix);
511
512      switch (iterm) {
513          case 1: return NLOPT_XTOL_REACHED;
514          case 2: return NLOPT_FTOL_REACHED;
515          case 3: return NLOPT_MINF_MAX_REACHED;
516          case 4: return NLOPT_SUCCESS; /* gradient tolerance reached */
517          case 6: return NLOPT_SUCCESS;
518          case 12: case 13: return NLOPT_MAXEVAL_REACHED;
519          default: return NLOPT_FAILURE;
520      }
521 }