chiark / gitweb /
Interpret NaN as a soft feasibility constraint
authorJoshua Nathaniel Pritikin <jpritikin@pobox.com>
Thu, 19 Mar 2015 12:07:46 +0000 (08:07 -0400)
committerJoshua Nathaniel Pritikin <jpritikin@pobox.com>
Fri, 20 Mar 2015 18:31:25 +0000 (14:31 -0400)
slsqp/slsqp.c

index ba27d7fe4fe22af4c6300b30a66882d69bbbb972..7a341d60785ff49f362407a31f9c8befe627fe99 100644 (file)
@@ -2044,12 +2044,16 @@ L190:
     *mode = line == 1 ? -2 : 1;
     goto L330;
 L200:
-    if (h1 <= h3 / ten || line > 10) {
-       goto L240;
+    if (nlopt_isfinite(h1)) {
+           if (h1 <= h3 / ten || line > 10) {
+                   goto L240;
+           }
+           /* Computing MAX */
+           d__1 = h3 / (two * (h3 - h1));
+           alpha = MAX2(d__1,alfmin);
+    } else {
+           alpha = MAX2(alpha*.5,alfmin);
     }
-/* Computing MAX */
-    d__1 = h3 / (two * (h3 - h1));
-    alpha = MAX2(d__1,alfmin);
     goto L190;
 /*   EXACT LINESEARCH */
 L210:
@@ -2508,84 +2512,86 @@ nlopt_result nlopt_slsqp(unsigned n, nlopt_func f, void *f_data,
              stop->nevals++;
              if (nlopt_stop_forced(stop)) {
                  fcur = HUGE_VAL; ret = NLOPT_FORCED_STOP; goto done; }
-             want_grad = 0;
-             ii = 0;
-             for (i = 0; i < p; ++i) {
-                 unsigned j, k;
-                 nlopt_eval_constraint(c+ii, newcgrad, h+i, n, xcur);
-                 if (nlopt_stop_forced(stop)) {
-                     ret = NLOPT_FORCED_STOP; goto done; }
-                 for (k = 0; k < h[i].m; ++k, ++ii) {
-                     infeasibility_cur =
-                         MAX2(infeasibility_cur, fabs(c[ii]));
-                     feasible_cur =
-                         feasible_cur && fabs(c[ii]) <= h[i].tol[k];
-                     if (newcgrad) {
-                         for (j = 0; j < n; ++ j)
-                             cgrad[j*U(mpi1) + ii] = cgradtmp[k*n + j];
+             if (nlopt_isfinite(fcur)) {
+                 want_grad = 0;
+                 ii = 0;
+                 for (i = 0; i < p; ++i) {
+                     unsigned j, k;
+                     nlopt_eval_constraint(c+ii, newcgrad, h+i, n, xcur);
+                     if (nlopt_stop_forced(stop)) {
+                         ret = NLOPT_FORCED_STOP; goto done; }
+                     for (k = 0; k < h[i].m; ++k, ++ii) {
+                         infeasibility_cur =
+                             MAX2(infeasibility_cur, fabs(c[ii]));
+                         feasible_cur =
+                             feasible_cur && fabs(c[ii]) <= h[i].tol[k];
+                         if (newcgrad) {
+                             for (j = 0; j < n; ++ j)
+                                 cgrad[j*U(mpi1) + ii] = cgradtmp[k*n + j];
+                         }
                      }
                  }
-             }
-             for (i = 0; i < m; ++i) {
-                 unsigned j, k;
-                 nlopt_eval_constraint(c+ii, newcgrad, fc+i, n, xcur);
-                 if (nlopt_stop_forced(stop)) {
-                     ret = NLOPT_FORCED_STOP; goto done; }
-                 for (k = 0; k < fc[i].m; ++k, ++ii) {
-                     infeasibility_cur =
-                         MAX2(infeasibility_cur, c[ii]);
-                     feasible_cur =
-                         feasible_cur && c[ii] <= fc[i].tol[k];
-                     if (newcgrad) {
-                         for (j = 0; j < n; ++ j)
-                             cgrad[j*U(mpi1) + ii] = -cgradtmp[k*n + j];
+                 for (i = 0; i < m; ++i) {
+                     unsigned j, k;
+                     nlopt_eval_constraint(c+ii, newcgrad, fc+i, n, xcur);
+                     if (nlopt_stop_forced(stop)) {
+                         ret = NLOPT_FORCED_STOP; goto done; }
+                     for (k = 0; k < fc[i].m; ++k, ++ii) {
+                         infeasibility_cur =
+                             MAX2(infeasibility_cur, c[ii]);
+                         feasible_cur =
+                             feasible_cur && c[ii] <= fc[i].tol[k];
+                         if (newcgrad) {
+                             for (j = 0; j < n; ++ j)
+                                 cgrad[j*U(mpi1) + ii] = -cgradtmp[k*n + j];
+                         }
+                         c[ii] = -c[ii]; /* slsqp sign convention */
                      }
-                     c[ii] = -c[ii]; /* slsqp sign convention */
                  }
              }
              break;}
-         case 0: /* required accuracy for solution obtained */
-             goto done;
-         case 8: /* positive directional derivative for linesearch */
-             /* relaxed convergence check for a feasible_cur point,
-                as in the SLSQP code (except xtol as well as ftol) */
-             ret = NLOPT_ROUNDOFF_LIMITED; /* usually why deriv>0 */
-             if (feasible_cur) {
-                 double save_ftol_rel = stop->ftol_rel;
-                 double save_xtol_rel = stop->xtol_rel;
-                 double save_ftol_abs = stop->ftol_abs;
-                 stop->ftol_rel *= 10;
-                 stop->ftol_abs *= 10;
-                 stop->xtol_rel *= 10;
-                 if (nlopt_stop_ftol(stop, fcur, state.f0))
-                     ret = NLOPT_FTOL_REACHED;
-                 else if (nlopt_stop_x(stop, xcur, state.x0))
-                     ret = NLOPT_XTOL_REACHED;
-                 stop->ftol_rel = save_ftol_rel;
-                 stop->ftol_abs = save_ftol_abs;
-                 stop->xtol_rel = save_xtol_rel;
-             }
-             goto done;
-         case 5: /* singular matrix E in LSQ subproblem */
-         case 6: /* singular matrix C in LSQ subproblem */
-         case 7: /* rank-deficient equality constraint subproblem HFTI */
-             ret = NLOPT_ROUNDOFF_LIMITED;
-             goto done;
-         case 4: /* inequality constraints incompatible */
-         case 3: /* more than 3*n iterations in LSQ subproblem */
-         case 9: /* more than iter iterations in SQP */
-             ret = NLOPT_FAILURE;
-             goto done;
-         case 2: /* number of equality constraints > n */
-         default: /* >= 10: working space w or jw too small */
-             ret = NLOPT_INVALID_ARGS;
-             goto done;
+             case 0: /* required accuracy for solution obtained */
+                 goto done;
+             case 8: /* positive directional derivative for linesearch */
+                 /* relaxed convergence check for a feasible_cur point,
+                    as in the SLSQP code (except xtol as well as ftol) */
+                 ret = NLOPT_ROUNDOFF_LIMITED; /* usually why deriv>0 */
+                 if (feasible_cur) {
+                     double save_ftol_rel = stop->ftol_rel;
+                     double save_xtol_rel = stop->xtol_rel;
+                     double save_ftol_abs = stop->ftol_abs;
+                     stop->ftol_rel *= 10;
+                     stop->ftol_abs *= 10;
+                     stop->xtol_rel *= 10;
+                     if (nlopt_stop_ftol(stop, fcur, state.f0))
+                         ret = NLOPT_FTOL_REACHED;
+                     else if (nlopt_stop_x(stop, xcur, state.x0))
+                         ret = NLOPT_XTOL_REACHED;
+                     stop->ftol_rel = save_ftol_rel;
+                     stop->ftol_abs = save_ftol_abs;
+                     stop->xtol_rel = save_xtol_rel;
+                 }
+                 goto done;
+             case 5: /* singular matrix E in LSQ subproblem */
+             case 6: /* singular matrix C in LSQ subproblem */
+             case 7: /* rank-deficient equality constraint subproblem HFTI */
+                 ret = NLOPT_ROUNDOFF_LIMITED;
+                 goto done;
+             case 4: /* inequality constraints incompatible */
+             case 3: /* more than 3*n iterations in LSQ subproblem */
+             case 9: /* more than iter iterations in SQP */
+                 ret = NLOPT_FAILURE;
+                 goto done;
+             case 2: /* number of equality constraints > n */
+             default: /* >= 10: working space w or jw too small */
+                 ret = NLOPT_INVALID_ARGS;
+                 goto done;
          }
          prev_mode = mode;
 
          /* update best point so far */
-         if ((fcur < *minf && (feasible_cur || !feasible))
-             || (!feasible && infeasibility_cur < infeasibility)) {
+         if (nlopt_isfinite(fcur) && ((fcur < *minf && (feasible_cur || !feasible))
+                                      || (!feasible && infeasibility_cur < infeasibility))) {
               *minf = fcur;
               feasible = feasible_cur;
               infeasibility = infeasibility_cur;