Actual source code: gpcglinesearch.c

  2: #include "src/bound/impls/gpcg/gpcglinesearch.h"    /*I "tao_solver.h" I*/

  4: /* ---------------------------------------------------------- */
  7: static int TaoGPCGDestroyLineSearch(TAO_SOLVER tao, void*lsctx)
  8: {
  9:   int  info;
 10:   TAO_GPCGLINESEARCH *ctx = (TAO_GPCGLINESEARCH *)lsctx;

 12:   TaoFunctionBegin;
 13:   if (ctx->setupcalled==1){
 14:     info = TaoVecDestroy(ctx->W2);CHKERRQ(info);
 15:     info = TaoVecDestroy(ctx->Gold);CHKERRQ(info);
 16:   }
 17:   info = TaoFree(ctx);CHKERRQ(info);
 18:   TaoFunctionReturn(0);
 19: }
 20: /*------------------------------------------------------------*/
 23: static int TaoGPCGSetOptionsLineSearch(TAO_SOLVER tao, void*linectx)
 24: {
 25:   TAO_GPCGLINESEARCH *ctx = (TAO_GPCGLINESEARCH *)linectx;
 26:   double         tmp;
 27:   int            info;
 28:   TaoInt         itmp;
 29:   TaoTruth     flg;

 31:   TaoFunctionBegin;
 32:   info = TaoOptionsHead("GPCG line search options");CHKERRQ(info);

 34:   info = TaoOptionInt("-tao_nls_maxfev","max function evals in line search",0,ctx->maxfev,&itmp,&flg);CHKERRQ(info);
 35:   if (flg) {ctx->maxfev = itmp;}
 36:   info = TaoOptionDouble("-tao_nls_ftol","tol for sufficient decrease",0,ctx->ftol,&tmp,&flg);CHKERRQ(info);
 37:   if (flg) {ctx->ftol = tmp;}
 38:   info = TaoOptionDouble("-tao_nls_gtol","tol for curvature condition",0,ctx->gtol,&tmp,&flg);CHKERRQ(info);
 39:   if (flg) {ctx->gtol = tmp;}
 40:   info = TaoOptionDouble("-tao_nls_rtol","relative tol for acceptable step",0,ctx->rtol,&tmp,&flg);CHKERRQ(info);
 41:   if (flg) {ctx->rtol = tmp;}
 42:   info = TaoOptionDouble("-tao_nls_stepmin","lower bound for step",0,ctx->stepmin,&tmp,&flg);CHKERRQ(info);
 43:   if (flg) {ctx->stepmin = tmp;}
 44:   info = TaoOptionDouble("-tao_nls_stepmax","upper bound for step",0,ctx->stepmax,&tmp,&flg);CHKERRQ(info);
 45:   if (flg) {ctx->stepmax = tmp;}
 46:   info = TaoOptionsTail();CHKERRQ(info);

 48:   TaoFunctionReturn(0);
 49: }


 52: /*------------------------------------------------------------*/
 55: static int TaoGPCGViewLineSearch(TAO_SOLVER tao,void *ctx)
 56: {
 57:   TAO_GPCGLINESEARCH *ls = (TAO_GPCGLINESEARCH *)ctx;
 58:   int            info;

 60:   TaoFunctionBegin;
 61:   info = TaoPrintInt(tao,"  Line search: maxf=%d,",ls->maxfev);CHKERRQ(info);
 62:   info = TaoPrintDouble(tao," ftol=%g,",ls->ftol);CHKERRQ(info);
 63:   info = TaoPrintDouble(tao," rtol=%g,",ls->rtol);CHKERRQ(info);
 64:   info = TaoPrintDouble(tao," gtol=%g\n",ls->gtol);CHKERRQ(info);
 65:   TaoFunctionReturn(0);
 66: }

 68: /*------------------------------------------------------------*/
 71: static int TaoGPCGApplyLineSearch(TAO_SOLVER tao,TaoVec* X, 
 72:                            TaoVec* G,TaoVec* S,TaoVec* W,
 73:                            double *f, double *f_full, double *step, TaoInt *info2,
 74:                            void*ctx)
 75: {
 76:   TAO_GPCGLINESEARCH *neP = (TAO_GPCGLINESEARCH *) ctx;
 77:   int       info;
 78:   TaoInt i;
 79:   double zero=0.0;
 80:   double d1,finit,actred,prered,rho, gdx;
 81:   TaoVec* XL, *XU, *Xold=neP->W2,*Gold=neP->Gold;
 82:   TaoTruth flag;

 84:   TaoFunctionBegin;
 85:   /* neP->stepmin - lower bound for step */
 86:   /* neP->stepmax - upper bound for step */
 87:   /* neP->rtol           - relative tolerance for an acceptable step */
 88:   /* neP->ftol           - tolerance for sufficient decrease condition */
 89:   /* neP->gtol           - tolerance for curvature condition */
 90:   /* neP->nfev           - number of function evaluations */
 91:   /* neP->maxfev  - maximum number of function evaluations */

 93:   /* Check input parameters for errors */
 94:   *info2=0;
 95:   if (neP->setupcalled){
 96:     info=X->Compatible(neP->W2,&flag); CHKERRQ(info);
 97:     if (flag==TAO_FALSE){
 98:       info=TaoVecDestroy(neP->W2); CHKERRQ(info);neP->W2=0;
 99:       info=TaoVecDestroy(neP->Gold); CHKERRQ(info);neP->Gold=0;
100:       neP->setupcalled=0;
101:     }
102:   }

104:   if (neP->setupcalled==0){
105:     info = X->Clone(&neP->W2); CHKERRQ(info);
106:     Xold=neP->W2;
107:     info = X->Clone(&neP->Gold); CHKERRQ(info);
108:     Gold=neP->Gold;
109:     neP->setupcalled=1;
110:   }

112:   info = G->Dot(S,&gdx); CHKERRQ(info);
113:   info = Xold->CopyFrom(X); CHKERRQ(info);
114:   info = Gold->CopyFrom(G); CHKERRQ(info);
115:   info = TaoGetVariableBounds(tao,&XL,&XU); CHKERRQ(info);
116:   info = X->StepBoundInfo(XL,XU,S,&rho,&actred,&d1);CHKERRQ(info);
117:   rho=0; actred=0;
118:   *step = TaoMin(*step,d1);

120:   if (*step < zero) {
121:     info = PetscInfo1(tao,"TaoGPCGApplyLineSearch:Line search error: step (%g) < 0\n",*step); CHKERRQ(info);
122:     *info2 = -1; TaoFunctionReturn(0);
123:   } else if (neP->ftol < zero) {
124:     info = PetscInfo1(tao,"TaoGPCGApplyLineSearch:Line search error: ftol (%g) < 0\n",neP->ftol); CHKERRQ(info);
125:     *info2 = -2; TaoFunctionReturn(0);
126:   } else if (neP->rtol < zero) {
127:     info = PetscInfo1(tao,"TaoGPCGApplyLineSearch:Line search error: rtol (%g) < 0\n",neP->rtol); CHKERRQ(info);
128:     *info2 = -3; TaoFunctionReturn(0);
129:   } else if (neP->gtol < zero) {
130:     info = PetscInfo1(tao,"TaoGPCGApplyLineSearch:Line search error: gtol (%g) < 0\n",neP->gtol); CHKERRQ(info);
131:     *info2 = -4; TaoFunctionReturn(0);
132:   } else if (neP->stepmin < zero) {
133:     info = PetscInfo1(tao,"TaoGPCGApplyLineSearch:Line search error: stepmin (%g) < 0\n",neP->stepmin); CHKERRQ(info);
134:     *info2 = -5; TaoFunctionReturn(0);
135:   } else if (neP->stepmax < neP->stepmin) {
136:     info = PetscInfo2(tao,"TaoGPCGApplyLineSearch:Line search error: stepmax (%g) < stepmin (%g)\n",neP->stepmax,neP->stepmin); CHKERRQ(info);
137:     *info2 = -6; TaoFunctionReturn(0);
138:   } else if (neP->maxfev < zero) {
139:     info = PetscInfo1(tao,"TaoGPCGApplyLineSearch:Line search error: maxfev (%d) < 0\n",neP->maxfev); CHKERRQ(info);
140:     *info2 = -7; TaoFunctionReturn(0);
141:   }


144:   /* Check that search direction is a descent direction */
145:   /*
146:   info = VecDot(G,S,&dginit);CHKERRQ(info);  / * dginit = G^T S * /
147:   if (dginit >= zero) {
148:     info = PetscLogInfo((tao,"TaoGPCGApplyLineSearch:Search direction not a descent direction\n")); CHKERRQ(info);
149:     *info2 = 7; return(0);
150:   }
151:   */
152:   /* Initialization */
153:   neP->nfev = 0;
154:   finit = *f;
155:   for (i=0; i< neP->maxfev; i++) {
156:     
157:     /* Force the step to be within the bounds */
158:     *step = TaoMax(*step,neP->stepmin);
159:     *step = TaoMin(*step,neP->stepmax);
160:     
161:     info = X->Waxpby(*step,S,1.0,Xold); CHKERRQ(info);
162:     info = X->Median(XL,X,XU); CHKERRQ(info);

164:     info = TaoGPCGComputeFunctionGradient(tao, X, f, G); CHKERRQ(info);
165:     if (0 == i) {
166:       *f_full = *f;
167:     }

169:     actred = *f - finit;
170:     info = W->Waxpby(-1.0,Xold,1.0,X); CHKERRQ(info);
171:     info = W->Dot(Gold,&prered); CHKERRQ(info);
172:     if (fabs(prered)<1.0e-100) prered=1.0e-12;
173:     rho = actred/prered;
174:     /* 
175:        If sufficient progress has been obtained, accept the
176:        point.  Otherwise, backtrack. 
177:     */

179:     if (rho > neP->ftol){
180:       break;
181:     } else{
182:       *step = (*step)/2;
183:     }
184:   }

186:   /* Convergence testing */
187:   
188:   if (*step <= neP->stepmin || *step >= neP->stepmax) {
189:     *info2 = 6;
190:     info = PetscInfo(tao,"TaoGPCGApplyLineSearch:Rounding errors may prevent further progress.  May not be a step satisfying\n"); CHKERRQ(info);
191:     info = PetscInfo(tao,"TaoGPCGApplyLineSearch:sufficient decrease and curvature conditions. Tolerances may be too small.\n"); CHKERRQ(info);
192:   }
193:   if (*step == neP->stepmax) {
194:     info = PetscInfo1(tao,"TaoGPCGApplyLineSearch:Step is at the upper bound, stepmax (%g)\n",neP->stepmax); CHKERRQ(info);
195:     *info2 = 5;
196:   }
197:   if (*step == neP->stepmin) {
198:     info = PetscInfo1(tao,"TaoGPCGApplyLineSearch:Step is at the lower bound, stepmin (%g)\n",neP->stepmin); CHKERRQ(info);
199:     *info2 = 4;
200:   }
201:   if (neP->nfev >= neP->maxfev) {
202:     info = PetscInfo2(tao,"TaoGPCGApplyLineSearch:Number of line search function evals (%d) > maximum (%d)\n",neP->nfev,neP->maxfev); CHKERRQ(info);
203:     *info2 = 3;
204:   }
205:   if ((neP->bracket) && (neP->stepmax - neP->stepmin <= neP->rtol*neP->stepmax)){
206:     info = PetscInfo1(tao,"TaoGPCGApplyLineSearch:Relative width of interval of uncertainty is at most rtol (%g)\n",neP->rtol); CHKERRQ(info);
207:     *info2 = 2;
208:   }
209:   /*
210:   if ((*f <= ftest1) && (PetscAbsDouble(dg) <= neP->gtol*(-dginit))) {
211:     info = PetscLogInfo((tao,"TaoGPCGApplyLineSearch:Line search success: Sufficient decrease and directional deriv conditions hold\n")); CHKERRQ(info);
212:     *info2 = 1;
213:   }
214:   */
215:   
216:   /* Finish computations */
217:   info = PetscInfo2(tao,"TaoGPCGApplyLineSearch:%d function evals in line search, step = %10.4f\n",neP->nfev,*step); CHKERRQ(info);

219:   TaoFunctionReturn(0);
220: }

222: /* ---------------------------------------------------------- */
225: int TaoGPCGCreateLineSearch(TAO_SOLVER tao)
226: {
227:   int info;
228:   TAO_GPCGLINESEARCH *neP;

230:   TaoFunctionBegin;

232:   info = TaoNew(TAO_GPCGLINESEARCH,&neP);CHKERRQ(info);
233:   info = PetscLogObjectMemory(tao,sizeof(TAO_GPCGLINESEARCH)); CHKERRQ(info);
234:   neP->ftol                  = 0.05;
235:   neP->rtol                  = 0.0;
236:   neP->gtol                  = 0.0;
237:   neP->stepmin                  = 1.0e-20;
238:   neP->stepmax                  = 1.0e+20;
239:   neP->nfev                  = 0; 
240:   neP->bracket                  = 0; 
241:   neP->infoc              = 1;
242:   neP->maxfev                  = 30;
243:   neP->setupcalled        = 0;

245:   info = TaoSetLineSearch(tao,0,
246:                           TaoGPCGSetOptionsLineSearch,
247:                           TaoGPCGApplyLineSearch,
248:                           TaoGPCGViewLineSearch,
249:                           TaoGPCGDestroyLineSearch,
250:                           (void *) neP);CHKERRQ(info);

252:   TaoFunctionReturn(0);
253: }