-/* DIRect.f -- translated by f2c (version 20050501).
+/* DIRect-transp.f -- translated by f2c (version 20050501).
f2c output hand-cleaned by SGJ (August 2007).
*/
/* Table of constant values */
-static integer c__90000 = 90000;
-static integer c__600 = 600;
-static integer c__64 = 64;
-static integer c__3000 = 3000;
-
/* +-----------------------------------------------------------------------+ */
/* | Program : Direct.f | */
/* | Last modified : 07-16-2001 | */
/* | Lipschitz continues. However, DIRECT has proven to be effective on | */
/* | more complex problems than these. | */
/* +-----------------------------------------------------------------------+ */
-/* Subroutine */ void direct_direct_(fp fcn, doublereal *x, integer *n, doublereal *eps, doublereal epsabs, integer *maxf, integer *maxt, doublereal *fmin, doublereal *l,
+/* Subroutine */ void direct_direct_(fp fcn, doublereal *x, integer *n, doublereal *eps, doublereal epsabs, integer *maxf, integer *maxt, double starttime, double maxtime, int *force_stop, doublereal *minf, doublereal *l,
doublereal *u, integer *algmethod, integer *ierror, FILE *logfile,
doublereal *fglobal, doublereal *fglper, doublereal *volper,
doublereal *sigmaper, void *fcn_data)
integer i__1, i__2;
doublereal d__1;
+ /* changed by SGJ to be dynamically allocated ... would be
+ even better to use realloc, below, to grow these as needed */
+ integer MAXFUNC = *maxf <= 0 ? 101000 : (*maxf + 1000 + *maxf / 2);
+ integer MAXDEEP = *maxt <= 0 ? MAXFUNC/5: *maxt + 1000;
+ const integer MAXDIV = 5000;
+
/* Local variables */
integer increase;
doublereal *c__ = 0 /* was [90000][64] */, *f = 0 /*
/* FIXME: change sizes dynamically? */
#define MY_ALLOC(p, t, n) p = (t *) malloc(sizeof(t) * (n)); \
if (!(p)) { *ierror = -100; goto cleanup; }
- MY_ALLOC(c__, doublereal, 5760000);
- MY_ALLOC(f, doublereal, 180000);
- MY_ALLOC(s, integer, 6000);
- MY_ALLOC(w, doublereal, 64);
- MY_ALLOC(oldl, doublereal, 64);
- MY_ALLOC(oldu, doublereal, 64);
- MY_ALLOC(list2, integer, 128);
- MY_ALLOC(point, integer, 90000);
- MY_ALLOC(anchor, integer, 602);
- MY_ALLOC(length, integer, 5760000);
- MY_ALLOC(arrayi, integer, 64);
- MY_ALLOC(levels, doublereal, 601);
- MY_ALLOC(thirds, doublereal, 601);
+
+ /* Note that I've transposed c__, length, and f relative to the
+ original Fortran code. e.g. length was length(maxfunc,n)
+ in Fortran [ or actually length(maxfunc, maxdims), but by
+ using malloc I can just allocate n ], corresponding to
+ length[n][maxfunc] in C, but I've changed the code to access
+ it as length[maxfunc][n]. That is, the maxfunc direction
+ is the discontiguous one. This makes it easier to resize
+ dynamically (by adding contiguous rows) using realloc, without
+ having to move data around manually. */
+ MY_ALLOC(c__, doublereal, MAXFUNC * (*n));
+ MY_ALLOC(length, integer, MAXFUNC * (*n));
+ MY_ALLOC(f, doublereal, MAXFUNC * 2);
+ MY_ALLOC(point, integer, MAXFUNC);
+ if (*maxf <= 0) *maxf = MAXFUNC - 1000;
+
+ MY_ALLOC(s, integer, MAXDIV * 2);
+
+ MY_ALLOC(anchor, integer, MAXDEEP + 2);
+ MY_ALLOC(levels, doublereal, MAXDEEP + 1);
+ MY_ALLOC(thirds, doublereal, MAXDEEP + 1);
+ if (*maxt <= 0) *maxt = MAXDEEP;
+
+ MY_ALLOC(w, doublereal, (*n));
+ MY_ALLOC(oldl, doublereal, (*n));
+ MY_ALLOC(oldu, doublereal, (*n));
+ MY_ALLOC(list2, integer, (*n) * 2);
+ MY_ALLOC(arrayi, integer, (*n));
/* +-----------------------------------------------------------------------+ */
/* | SUBROUTINE Direct | */
/* | eps -- Exceeding value. If eps > 0, we use the same epsilon for | */
/* | all iterations. If eps < 0, we use the update formula from | */
/* | Jones: | */
-/* | eps = max(1.D-4*abs(fmin),epsfix), | */
+/* | eps = max(1.D-4*abs(minf),epsfix), | */
/* | where epsfix = abs(eps), the absolute value of eps which is| */
/* | passed to the function. | */
/* | maxf -- The maximum number of function evaluations. | */
/* | fglper -- Terminate the optimization when the percent error | */
/* | 100(f_min - fglobal)/max(1,abs(fglobal)) < fglper. | */
/* | volper -- Terminate the optimization when the volume of the | */
-/* | hyperrectangle S with f(c(S)) = fmin is less then volper | */
+/* | hyperrectangle S with f(c(S)) = minf is less then volper | */
/* | percent of the volume of the original hyperrectangle. | */
/* |sigmaper -- Terminate the optimization when the measure of the | */
-/* | hyperrectangle S with f(c(S)) = fmin is less then sigmaper.| */
+/* | hyperrectangle S with f(c(S)) = minf is less then sigmaper.| */
/* | | */
/* | User data that is passed through without being changed: | */
/* | fcn_data - opaque pointer to any user data | */
/* | X should be a good approximation to the global minimum | */
/* | for the function within the hyper-box. | */
/* | | */
-/* | fmin -- The value of the function at x. | */
+/* | minf -- The value of the function at x. | */
/* | Ierror -- Error flag. If Ierror is lower 0, an error has occured. The| */
/* | values of Ierror mean | */
/* | Fatal errors : | */
/* | 2 Number of iterations is equal to maxT. | */
/* | 3 The best function value found is within fglper of | */
/* | the (known) global optimum, that is | */
-/* | 100(fmin - fglobal/max(1,|fglobal|)) < fglper. | */
+/* | 100(minf - fglobal/max(1,|fglobal|)) < fglper. | */
/* | Note that this termination signal only occurs when | */
/* | the global optimal value is known, that is, a test | */
/* | function is optimized. | */
-/* | 4 The volume of the hyperrectangle with fmin at its | */
+/* | 4 The volume of the hyperrectangle with minf at its | */
/* | center is less than volper percent of the volume of | */
/* | the original hyperrectangle. | */
-/* | 5 The measure of the hyperrectangle with fmin at its | */
+/* | 5 The measure of the hyperrectangle with minf at its | */
/* | center is less than sigmaper. | */
/* | | */
/* | SUBROUTINEs used : | */
/* +-----------------------------------------------------------------------+ */
cheat = 0;
kmax = 1e10;
- mdeep = 600;
+ mdeep = MAXDEEP;
/* +-----------------------------------------------------------------------+ */
/* | Write the header of the logfile. | */
/* +-----------------------------------------------------------------------+ */
direct_dirheader_(logfile, &version, &x[1], n, eps, maxf, maxt, &l[1], &u[1],
- algmethod, &c__90000, &c__600, fglobal, fglper, ierror, &epsfix, &
+ algmethod, &MAXFUNC, &MAXDEEP, fglobal, fglper, ierror, &epsfix, &
iepschange, volper, sigmaper);
/* +-----------------------------------------------------------------------+ */
/* | If an error has occured while writing the header (we do some checking | */
/* +-----------------------------------------------------------------------+ */
/* | Initialiase the lists. | */
/* +-----------------------------------------------------------------------+ */
- direct_dirinitlist_(anchor, &ifree, point, f, &c__90000, &c__600);
+ direct_dirinitlist_(anchor, &ifree, point, f, &MAXFUNC, &MAXDEEP);
/* +-----------------------------------------------------------------------+ */
/* | Call the routine to initialise the mapping of x from the n-dimensional| */
/* | unit cube to the hypercube given by u and l. If an error occured, | */
/* +-----------------------------------------------------------------------+ */
direct_dirinit_(f, fcn, c__, length, &actdeep, point, anchor, &ifree,
logfile, arrayi, &maxi, list2, w, &x[1], &l[1], &u[1],
- fmin, &minpos, thirds, levels, &c__90000, &c__600, n, &c__64, &
- fmax, &ifeasiblef, &iinfesiblef, ierror, fcn_data, jones);
+ minf, &minpos, thirds, levels, &MAXFUNC, &MAXDEEP, n, n, &
+ fmax, &ifeasiblef, &iinfesiblef, ierror, fcn_data, jones,
+ starttime, maxtime, force_stop);
/* +-----------------------------------------------------------------------+ */
/* | Added error checking. | */
/* +-----------------------------------------------------------------------+ */
fprintf(logfile, "WARNING: Error occured in routine DIRsamplef..\n");
goto cleanup;
}
+ if (*ierror == -102) goto L100;
}
+ else if (*ierror == DIRECT_MAXTIME_EXCEEDED) goto L100;
numfunc = maxi + 1 + maxi;
actmaxdeep = 1;
oldpos = 0;
/* +-----------------------------------------------------------------------+ */
/* | If no feasible point has been found, give out the iteration, the | */
/* | number of function evaluations and a warning. Otherwise, give out | */
-/* | the iteration, the number of function evaluations done and fmin. | */
+/* | the iteration, the number of function evaluations done and minf. | */
/* +-----------------------------------------------------------------------+ */
if (ifeasiblef > 0) {
if (logfile)
} else {
if (logfile)
fprintf(logfile, "%d, %d, %g, %g\n",
- tstart-1, numfunc, *fmin, fmax);
+ tstart-1, numfunc, *minf, fmax);
}
/* +-----------------------------------------------------------------------+ */
/* +-----------------------------------------------------------------------+ */
/* | in the list S. | */
/* +-----------------------------------------------------------------------+ */
actdeep = actmaxdeep;
- direct_dirchoose_(anchor, s, &c__600, f, fmin, *eps, epsabs, levels, &maxpos, length,
- &c__90000, &c__600, &c__3000, n, logfile, &cheat, &
+ direct_dirchoose_(anchor, s, &MAXDEEP, f, minf, *eps, epsabs, levels, &maxpos, length,
+ &MAXFUNC, &MAXDEEP, &MAXDIV, n, logfile, &cheat, &
kmax, &ifeasiblef, jones);
/* +-----------------------------------------------------------------------+ */
/* | Add other hyperrectangles to S, which have the same level and the same| */
/* | JG 07/16/01 Added Errorflag. | */
/* +-----------------------------------------------------------------------+ */
if (*algmethod == 0) {
- direct_dirdoubleinsert_(anchor, s, &maxpos, point, f, &c__600, &c__90000,
- &c__3000, ierror);
+ direct_dirdoubleinsert_(anchor, s, &maxpos, point, f, &MAXDEEP, &MAXFUNC,
+ &MAXDIV, ierror);
if (*ierror == -6) {
if (logfile)
fprintf(logfile,
newtosample = 0;
i__2 = maxpos;
for (j = 1; j <= i__2; ++j) {
- actdeep = s[j + 2999];
+ actdeep = s[j + MAXDIV-1];
/* +-----------------------------------------------------------------------+ */
/* | If the actual index is a point to sample, do it. | */
/* +-----------------------------------------------------------------------+ */
/* +-----------------------------------------------------------------------+ */
/* | JG 09/24/00 Calculate the value delta used for sampling points. | */
/* +-----------------------------------------------------------------------+ */
- actdeep_div__ = direct_dirgetmaxdeep_(&s[j - 1], length, &c__90000,
+ actdeep_div__ = direct_dirgetmaxdeep_(&s[j - 1], length, &MAXFUNC,
n);
delta = thirds[actdeep_div__ + 1];
- actdeep = s[j + 2999];
+ actdeep = s[j + MAXDIV-1];
/* +-----------------------------------------------------------------------+ */
/* | If the current dept of division is only one under the maximal allowed | */
/* | dept, stop the computation. | */
anchor[actdeep + 1] = point[help - 1];
}
if (actdeep < 0) {
- actdeep = (integer) f[help - 1];
+ actdeep = (integer) f[(help << 1) - 2];
}
/* +-----------------------------------------------------------------------+ */
/* | Get the Directions in which to decrease the intervall-length. | */
/* +-----------------------------------------------------------------------+ */
- direct_dirget_i__(length, &help, arrayi, &maxi, n, &c__90000);
+ direct_dirget_i__(length, &help, arrayi, &maxi, n, &MAXFUNC);
/* +-----------------------------------------------------------------------+ */
/* | Sample the function. To do this, we first calculate the points where | */
/* | we need to sample the function. After checking for errors, we then do | */
/* +-----------------------------------------------------------------------+ */
direct_dirsamplepoints_(c__, arrayi, &delta, &help, &start, length,
logfile, f, &ifree, &maxi, point, &x[
- 1], &l[1], fmin, &minpos, &u[1], n, &c__90000, &
- c__600, &oops);
+ 1], &l[1], minf, &minpos, &u[1], n, &MAXFUNC, &
+ MAXDEEP, &oops);
if (oops > 0) {
if (logfile)
fprintf(logfile, "WARNING: Error occured in routine DIRsamplepoints.\n");
/* +-----------------------------------------------------------------------+ */
direct_dirsamplef_(c__, arrayi, &delta, &help, &start, length,
logfile, f, &ifree, &maxi, point, fcn, &x[
- 1], &l[1], fmin, &minpos, &u[1], n, &c__90000, &
- c__600, &oops, &fmax, &ifeasiblef, &iinfesiblef,
- fcn_data);
+ 1], &l[1], minf, &minpos, &u[1], n, &MAXFUNC, &
+ MAXDEEP, &oops, &fmax, &ifeasiblef, &iinfesiblef,
+ fcn_data, force_stop);
+ if (force_stop && *force_stop) {
+ *ierror = -102;
+ goto L100;
+ }
+ if (nlopt_stop_time_(starttime, maxtime)) {
+ *ierror = DIRECT_MAXTIME_EXCEEDED;
+ goto L100;
+ }
if (oops > 0) {
if (logfile)
fprintf(logfile, "WARNING: Error occured in routine DIRsamplef.\n");
/* | Divide the intervalls. | */
/* +-----------------------------------------------------------------------+ */
direct_dirdivide_(&start, &actdeep_div__, length, point, arrayi, &
- help, list2, w, &maxi, f, &c__90000, &c__600, n);
+ help, list2, w, &maxi, f, &MAXFUNC, &MAXDEEP, n);
/* +-----------------------------------------------------------------------+ */
/* | Insert the new intervalls into the list (sorted). | */
/* +-----------------------------------------------------------------------+ */
direct_dirinsertlist_(&start, anchor, point, f, &maxi, length, &
- c__90000, &c__600, n, &help, jones);
+ MAXFUNC, &MAXDEEP, n, &help, jones);
/* +-----------------------------------------------------------------------+ */
/* | Increase the number of function evaluations. | */
/* +-----------------------------------------------------------------------+ */
if (oldpos < minpos) {
if (logfile)
fprintf(logfile, "%d, %d, %g, %g\n",
- t, numfunc, *fmin, fmax);
+ t, numfunc, *minf, fmax);
}
/* +-----------------------------------------------------------------------+ */
/* | If no feasible point has been found, give out the iteration, the | */
/* +-----------------------------------------------------------------------+ */
/* +-----------------------------------------------------------------------+ */
/* | JG 01/22/01 Calculate the index for the hyperrectangle at which | */
-/* | fmin is assumed. We then calculate the volume of this | */
+/* | minf is assumed. We then calculate the volume of this | */
/* | hyperrectangle and store it in delta. This delta can be | */
/* | used to stop DIRECT once the volume is below a certain | */
/* | percentage of the original volume. Since the original | */
/* +-----------------------------------------------------------------------+ */
*ierror = jones;
jones = 0;
- actdeep_div__ = direct_dirgetlevel_(&minpos, length, &c__90000, n, jones);
+ actdeep_div__ = direct_dirgetlevel_(&minpos, length, &MAXFUNC, n, jones);
jones = *ierror;
/* +-----------------------------------------------------------------------+ */
/* | JG 07/16/01 Use precalculated values to calculate volume. | */
}
/* +-----------------------------------------------------------------------+ */
/* | JG 01/23/01 Calculate the measure for the hyperrectangle at which | */
-/* | fmin is assumed. If this measure is smaller then sigmaper,| */
+/* | minf is assumed. If this measure is smaller then sigmaper,| */
/* | we stop DIRECT. | */
/* +-----------------------------------------------------------------------+ */
- actdeep_div__ = direct_dirgetlevel_(&minpos, length, &c__90000, n, jones);
+ actdeep_div__ = direct_dirgetlevel_(&minpos, length, &MAXFUNC, n, jones);
delta = levels[actdeep_div__];
if (delta <= *sigmaper) {
*ierror = 5;
/* | global minimum value, terminate. This only makes sense if this optimal| */
/* | value is known, that is, in test problems. | */
/* +-----------------------------------------------------------------------+ */
- if ((*fmin - *fglobal) * 100 / divfactor <= *fglper) {
+ if ((*minf - *fglobal) * 100 / divfactor <= *fglper) {
*ierror = 3;
if (logfile)
- fprintf(logfile, "DIRECT stopped: fmin within fglper of global minimum.\n");
+ fprintf(logfile, "DIRECT stopped: minf within fglper of global minimum.\n");
goto L100;
}
/* +-----------------------------------------------------------------------+ */
/* +-----------------------------------------------------------------------+ */
if (iinfesiblef > 0) {
direct_dirreplaceinf_(&ifree, &ifreeold, f, c__, thirds, length, anchor,
- point, &u[1], &l[1], &c__90000, &c__600, &c__64, n,
+ point, &u[1], &l[1], &MAXFUNC, &MAXDEEP, n, n,
logfile, &fmax, jones);
}
ifreeold = ifree;
/* +-----------------------------------------------------------------------+ */
if (iepschange == 1) {
/* Computing MAX */
- d__1 = fabs(*fmin) * 1e-4;
+ d__1 = fabs(*minf) * 1e-4;
*eps = MAX(d__1,epsfix);
}
/* +-----------------------------------------------------------------------+ */
/* +-----------------------------------------------------------------------+ */
i__1 = *n;
for (i__ = 1; i__ <= i__1; ++i__) {
- x[i__] = c__[minpos + i__ * 90000 - 90001] * l[i__] + l[i__] * u[i__];
+ x[i__] = c__[i__ + minpos * i__1 - i__1-1] * l[i__] + l[i__] * u[i__];
u[i__] = oldu[i__ - 1];
l[i__] = oldl[i__ - 1];
/* L50: */
/* +-----------------------------------------------------------------------+ */
/* | Give out a summary of the run. | */
/* +-----------------------------------------------------------------------+ */
- direct_dirsummary_(logfile, &x[1], &l[1], &u[1], n, fmin, fglobal, &numfunc,
+ direct_dirsummary_(logfile, &x[1], &l[1], &u[1], n, minf, fglobal, &numfunc,
ierror);
/* +-----------------------------------------------------------------------+ */
/* | Format statements. | */