Actual source code: boundproj.c

      2: #include "boundproj.h"    /*I "tao_solver.h" I*/

      4: static int TaoSetOptions_LineSearch2(TAO_SOLVER,void*);
      5: static int TaoView_LineSearch2(TAO_SOLVER,void*);
      6: static int TaoApply_LineSearch2(TAO_SOLVER,TaoVec*,TaoVec*,TaoVec*,TaoVec*,
      7:                                 double*,double*,int*,void*);
      8: static int TaoDestroy_LineSearch2(TAO_SOLVER,void*);

     10: /* ---------------------------------------------------------- */
     13: static int TaoDestroy_LineSearch2(TAO_SOLVER tao,void *linectx)
     14: {
     15:   int  info;
     16:   TAO_LINESEARCH2 *ctx = (TAO_LINESEARCH2 *)linectx;

     18:   TaoFunctionBegin;
     19:   if (ctx->setupcalled==1){
     20:     info=TaoVecDestroy(ctx->W2);CHKERRQ(info);
     21:   }
     22:   info = TaoFree(ctx);CHKERRQ(info);
     23:   TaoFunctionReturn(0);
     24: }
     25: /*------------------------------------------------------------*/
     28: static int TaoSetOptions_LineSearch2(TAO_SOLVER tao,void *linectx)
     29: {
     30:   TAO_LINESEARCH2 *ctx = (TAO_LINESEARCH2 *)linectx;
     31:   int            info;

     33:   TaoFunctionBegin;
     34:   info = TaoOptionsHead("Projected line search options");CHKERRQ(info);
     35:   info = TaoOptionInt("-tao_ls_maxfev","max function evals in line search","",ctx->maxfev,&ctx->maxfev,0);CHKERRQ(info);
     36:   info = TaoOptionDouble("-tao_ls_ftol","tol for sufficient decrease","",ctx->ftol,&ctx->ftol,0);CHKERRQ(info);
     37:   info = TaoOptionDouble("-tao_ls_stepmin","lower bound for step","",ctx->stepmin,&ctx->stepmin,0);CHKERRQ(info);
     38:   info = TaoOptionDouble("-tao_ls_stepmax","upper bound for step","",ctx->stepmax,&ctx->stepmax,0);CHKERRQ(info);
     39:   info = TaoOptionsTail();CHKERRQ(info);

     41:   TaoFunctionReturn(0);
     42: }
     43: /*------------------------------------------------------------*/

     47: static int TaoView_LineSearch2(TAO_SOLVER tao,void *ctx)
     48: {
     49:   TAO_LINESEARCH2 *ls = (TAO_LINESEARCH2 *)ctx;
     50:   int            info;

     52:   TaoFunctionBegin;
     53:   info = TaoPrintInt(tao,"  Line search: maxf=%d,",ls->maxfev);CHKERRQ(info);
     54:   info = TaoPrintDouble(tao," ftol=%g\n",ls->ftol);CHKERRQ(info);
     55:   TaoFunctionReturn(0);
     56: }
     57: /* ---------------------------------------------------------- */
     60: /* @ TaoApply_LineSearch - This routine performs a line search. It
     61:    backtracks until the Armijo conditions are satisfied.

     63:    Input Parameters:
     64: +  tao - TAO_SOLVER context
     65: .  X - current iterate (on output X contains new iterate, X + step*S)
     66: .  S - search direction
     67: .  f - objective function evaluated at X
     68: .  G - gradient evaluated at X
     69: .  W - work vector
     70: .  gdx - inner product of gradient and the direction of the first linear manifold being searched
     71: -  step - initial estimate of step length

     73:    Output parameters:
     74: +  f - objective function evaluated at new iterate, X + step*S
     75: .  G - gradient evaluated at new iterate, X + step*S
     76: .  X - new iterate
     77: -  step - final step length

     79:    Info is set to one of:
     80: .   0 - the line search succeeds; the sufficient decrease
     81:    condition and the directional derivative condition hold

     83:    negative number if an input parameter is invalid
     84: .   -1 -  step < 0
     85: .   -2 -  ftol < 0
     86: -   -7 -  maxfev < 0

     88:    positive number > 1 if the line search otherwise terminates
     89: +    2 -  Relative width of the interval of uncertainty is
     90:          at most rtol.
     91: .    3 -  Maximum number of function evaluations (maxfev) has
     92:          been reached.
     93: .    4 -  Step is at the lower bound, stepmin.
     94: .    5 -  Step is at the upper bound, stepmax.
     95: .    6 -  Rounding errors may prevent further progress.
     96:          There may not be a step that satisfies the
     97:          sufficient decrease and curvature conditions.
     98:          Tolerances may be too small.
     99: +    7 -  Search direction is not a descent direction.

    101:    Notes:
    102:    This routine is used within the TAO_NLS method.
    103: @ */
    104: static int TaoApply_LineSearch2(TAO_SOLVER tao,TaoVec* X,TaoVec* G,TaoVec* S,TaoVec* Gold,double *f,
    105:                         double *step,int *info2,void*ctx)
    106: {
    107:   TAO_LINESEARCH2 *neP = (TAO_LINESEARCH2 *) ctx;
    108:   int       info, i;
    109:   double zero=0.0;
    110:   double finit,gdx,dginit,actred,prered,rho,d1=0,d2=0,d3=0;
    111:   TaoVec* XL, *XU, *Xold=neP->W2;
    112:   TaoTruth flag;

    114:   TaoFunctionBegin;
    115:   /* neP->stepmin - lower bound for step */
    116:   /* neP->stepmax - upper bound for step */
    117:   /* neP->rtol           - relative tolerance for an acceptable step */
    118:   /* neP->ftol           - tolerance for sufficient decrease condition */
    119:   /* neP->gtol           - tolerance for curvature condition */
    120:   /* neP->nfev           - number of function evaluations */
    121:   /* neP->maxfev  - maximum number of function evaluations */

    123:   /* Check input parameters for errors */

    125:   *info2 = 0;
    126:     /* After 2 failures, Check that search direction is a descent direction */
    127:     /* This test is not really sufficient.  */
    128:   if (neP->setupcalled){
    129:     info=X->Compatible(neP->W2,&flag); CHKERRQ(info);
    130:     if (flag==TAO_FALSE){
    131:       info=TaoVecDestroy(neP->W2); CHKERRQ(info);neP->W2=0;
    132:       neP->setupcalled=0;
    133:     }
    134:   }
    135:   if (neP->setupcalled==0){
    136:     info=X->Clone(&neP->W2); CHKERRQ(info);
    137:     Xold=neP->W2;
    138:     neP->setupcalled=1;
    139:   }
    140:   info = TaoGetVariableBounds(tao,&XL,&XU); CHKERRQ(info);
    141:   info = Gold->ScaleCopyFrom(-1.0,S); CHKERRQ(info);
    142:   info = Gold->BoundGradientProjection(S,XL,X,XU); CHKERRQ(info);
    143:   info = Gold->Dot(G,&dginit);CHKERRQ(info);  /* dginit = G^T S */
    144:   dginit*=-1;
    145:   gdx=dginit;
    146:   if (dginit >= zero) {
    147:     info = G->CopyFrom(Gold); CHKERRQ(info);
    148:     info = X->CopyFrom(Xold); CHKERRQ(info);
    149:     info = PetscLogInfo((tao,"TaoApply_LineSearch2:Search direction not a descent direction\n")); CHKERRQ(info);
    150:     *info2 = 7; TaoFunctionReturn(0);
    151:   }

    153:   info = Xold->CopyFrom(X); CHKERRQ(info);
    154:   info = Gold->CopyFrom(G); CHKERRQ(info);
    155:   info=TaoGetVariableBounds(tao,&XL,&XU); CHKERRQ(info);
    156:   if (*step < zero) {
    157:     info = PetscLogInfo((tao,"TaoApply_LineSearch2:Line search error: step (%g) < 0\n",*step)); CHKERRQ(info);
    158:     *info2 = -1; TaoFunctionReturn(0);
    159:   } else if (neP->ftol < zero) {
    160:     info = PetscLogInfo((tao,"TaoApply_LineSearch2:Line search error: ftol (%g) < 0\n",neP->ftol)); CHKERRQ(info);
    161:     *info2 = -2; TaoFunctionReturn(0);
    162:   } else if (neP->maxfev < zero) {
    163:     info = PetscLogInfo((tao,"TaoApply_LineSearch2:Line search error: maxfev (%d) < 0\n",neP->maxfev)); CHKERRQ(info);
    164:     *info2 = -7; TaoFunctionReturn(0);
    165:   }

    169:   /* Initialization */
    170:   neP->nfev = 0;
    171:   finit = *f;
    172:   for (i=0; i< neP->maxfev; i++) {
    173:
    174:     /* Force the step to be within the bounds */
    175:     *step = TaoMax(*step,neP->stepmin);
    176:     *step = TaoMin(*step,neP->stepmax);
    177:
    178:     if (i==0){
    179:       info = X->Axpy(*step,S);CHKERRQ(info);
    180:       info = X->Median(XL,X,XU);CHKERRQ(info);
    181:       info = TaoComputeMeritFunctionGradient(tao,X,f,G);
    182:       neP->nfev++;
    183:       if (info==0 && *f==*f){ /* Successful function evaluation */
    184:         actred = *f - finit;
    185:         rho = actred/( (*step) * (-TaoAbsScalar(gdx)) );
    186:         if (actred < 0 && rho > neP->ftol){
    187:           break;
    188:         } else{
    189:           info=Xold->StepBoundInfo(XL,XU,S,&d3,&d2,&d1);CHKERRQ(info);
    190:           info = PetscLogInfo((tao,"stepmax = %10f\n",d1)); CHKERRQ(info);

    192:           *step = TaoMin(*step,d1);
    193:           info = Gold->Dot(X,&d1); CHKERRQ(info);
    194:           info = Gold->Dot(Xold,&prered); CHKERRQ(info);
    195:           prered=d1-prered;
    196:           rho = actred/(-TaoAbsScalar(prered));
    197:           if (actred < 0 && rho > neP->ftol){
    198:             break;
    199:           } else{
    200:                 *step = (*step)/2;
    201:           }
    202:         }
    203:       } else { /* Function could not be evaluated at this point */
    204:         *step = (*step)*0.7;
    205:       }

    207:     } else {
    208:       info = X->Waxpby(*step,S,1.0,Xold);CHKERRQ(info);
    209:       info = X->Median(XL,X,XU);CHKERRQ(info);
    210:
    211:       info = G->Waxpby(-1,Xold,1.0,X);CHKERRQ(info);
    212:       info = G->Dot(Gold,&prered); CHKERRQ(info);
    213:       info = TaoComputeMeritFunctionGradient(tao,X,f,G); CHKERRQ(info);
    214:       neP->nfev++;
    215:       if (info==0 && *f==*f){ /* Successful function evaluation */
    216:         actred = *f - finit;
    217:         rho = actred/(-TaoAbsScalar(prered));
    218:         /*
    219:            If sufficient progress has been obtained, accept the
    220:            point.   Prered could be positive or negative.
    221:            Otherwise, backtrack.
    222:         */

    224:         if (actred < 0 && rho > neP->ftol){
    225:           info = PetscLogInfo((tao,"TaoApply_LineSearch2: steplength: %g, actual reduction: %g (hopefully positive), Predicted reduction: %g, rho: %g\n",*step,-actred,prered,rho)); CHKERRQ(info);
    226:           break;
    227:         } else {
    228:           info = PetscLogInfo((tao,"TaoApply_LineSearch2: steplength: %g, actual reduction: %g (hopefully positive), Predicted reduction: %g, rho: %g\n",*step,-actred,prered,rho)); CHKERRQ(info);
    229:           if (i<5){
    230:             *step = (*step)/2;
    231:           } else {
    232:             *step = (*step)/2;
    233:           }
    234:         }
    235:       } else { /* Function could not be evaluated at this point */
    236:           info = PetscLogInfo((tao,"TaoApply_LineSearch2: Problem in function evaluation\n")); CHKERRQ(info);
    237:         *step = (*step)*0.7;
    238:       }
    239:     }

    241:   }
    242:   /* Convergence testing */
    243:
    244:   if (*step <= neP->stepmin || *step >= neP->stepmax) {
    245:     *info2 = 6;
    246:     info = PetscLogInfo((tao,"TaoApply_LineSearch2:Rounding errors may prevent further progress.  May not be a step satisfying\n")); CHKERRQ(info);
    247:     info = PetscLogInfo((tao,"TaoApply_LineSearch2:sufficient decrease and curvature conditions. Tolerances may be too small.\n")); CHKERRQ(info);
    248:   }
    249:   if (*step == neP->stepmax) {
    250:     info = PetscLogInfo((tao,"TaoApply_LineSearch2:Step is at the upper bound, stepmax (%g)\n",neP->stepmax)); CHKERRQ(info);
    251:     *info2 = 5;
    252:   }
    253:   if (*step == neP->stepmin) {
    254:     info = PetscLogInfo((tao,"TaoApply_LineSearch2:Step is at the lower bound, stepmin (%g)\n",neP->stepmin)); CHKERRQ(info);
    255:     *info2 = 4;
    256:   }
    257:   if (neP->nfev >= neP->maxfev) {
    258:     info = PetscLogInfo((tao,"TaoApply_LineSearch2:Number of line search function evals (%d) > maximum (%d)\n",neP->nfev,neP->maxfev)); CHKERRQ(info);
    259:     *info2 = 3;
    260:   }
    261:   if ((neP->bracket) && (neP->stepmax - neP->stepmin <= neP->rtol*neP->stepmax)){
    262:     info = PetscLogInfo((tao,"TaoApply_LineSearch2:Relative width of interval of uncertainty is at most rtol (%g)\n",neP->rtol)); CHKERRQ(info);
    263:     *info2 = 2;
    264:   }
    265:   /*
    266:   if ((*f <= ftest1) && (TaoAbsDouble(dg) <= neP->gtol*(-dginit))) {
    267:     info = PetscLogInfo((tao,"TaoApply_LineSearch2:Line search success: Sufficient decrease and directional deriv conditions hold\n")); CHKERRQ(info);
    268:     *info2 = 1;
    269:   }
    270:   */
    271:
    272:   /* Finish computations */
    273:   /*
    274:   info = PetscLogInfo((tao,"TaoApply_LineSearch2:%d function evals in line search, step = %10.4f\n",neP->nfev,*step)); CHKERRQ(info);
    275:   */
    276:   TaoFunctionReturn(0);
    277: }

    280: /* ---------------------------------------------------------- */

    284: /*@C
    285:    TaoCreateProjectedLineSearch - Create a line search

    287:    Input Parameters:
    288: .  tao - TAO_SOLVER context

    290:    Note:
    291:    This routine is used within the following TAO bound constrained
    292:    minimization solvers: TRON (tao_tron) and limited memory variable metric
    293:    (tao_blmvm). This line search projects points onto the feasible, bound
    294:    constrained region.  It only enforces the Armijo descent condition.

    296:    Level: developer

    298: .keywords: TAO_SOLVER, linesearch
    299: @*/
    300: int TaoCreateProjectedLineSearch(TAO_SOLVER tao)
    301: {
    302:   int info;
    303:   TAO_LINESEARCH2 *neP;

    305:   TaoFunctionBegin;

    307:   info = TaoNew(TAO_LINESEARCH2,&neP); CHKERRQ(info);
    308:   info = PetscLogObjectMemory(tao,sizeof(TAO_LINESEARCH2)); CHKERRQ(info);

    310:   neP->ftol                  = 0.001;
    311:   neP->rtol                  = 0.0;
    312:   neP->gtol                  = 0.0;
    313:   neP->stepmin                  = 1.0e-30;
    314:   neP->stepmax                  = 1.0e+20;
    315:   neP->nfev                  = 0;
    316:   neP->bracket                  = 0;
    317:   neP->infoc              = 1;
    318:   neP->maxfev                  = 30;
    319:   neP->W2                 = TAO_NULL;
    320:   neP->setupcalled        = 0;

    322:   info = TaoSetLineSearch(tao,0,
    323:                           TaoSetOptions_LineSearch2,
    324:                           TaoApply_LineSearch2,
    325:                           TaoView_LineSearch2,
    326:                           TaoDestroy_LineSearch2,
    327:                           (void *) neP);CHKERRQ(info);

    329:   TaoFunctionReturn(0);
    330: }