Actual source code: morethuente.c
1: /*$Id: morethuente.c 1.35 05/05/10 15:21:38-05:00 sarich@zorak.(none) $*/
3: #include "morethuente.h" /*I "tao_solver.h
4: " I*/
6: static int TaoApply_LineSearch(TAO_SOLVER,TaoVec*,TaoVec*,TaoVec*,TaoVec*,
7: double*,double*,int*,void*);
9: /*
10: The subroutine mcstep is taken from the work of Jorge Nocedal.
11: this is a variant of More' and Thuente's routine.
13: subroutine mcstep
15: the purpose of mcstep is to compute a safeguarded step for
16: a linesearch and to update an interval of uncertainty for
17: a minimizer of the function.
19: the parameter stx contains the step with the least function
20: value. the parameter stp contains the current step. it is
21: assumed that the derivative at stx is negative in the
22: direction of the step. if bracket is set true then a
23: minimizer has been bracketed in an interval of uncertainty
24: with endpoints stx and sty.
26: the subroutine statement is
28: subroutine mcstep(stx,fx,dx,sty,fy,dy,stp,fp,dp,bracket,
29: stpmin,stpmax,info)
31: where
33: stx, fx, and dx are variables which specify the step,
34: the function, and the derivative at the best step obtained
35: so far. the derivative must be negative in the direction
36: of the step, that is, dx and stp-stx must have opposite
37: signs. on output these parameters are updated appropriately.
39: sty, fy, and dy are variables which specify the step,
40: the function, and the derivative at the other endpoint of
41: the interval of uncertainty. on output these parameters are
42: updated appropriately.
44: stp, fp, and dp are variables which specify the step,
45: the function, and the derivative at the current step.
46: if bracket is set true then on input stp must be
47: between stx and sty. on output stp is set to the new step.
49: bracket is a logical variable which specifies if a minimizer
50: has been bracketed. if the minimizer has not been bracketed
51: then on input bracket must be set false. if the minimizer
52: is bracketed then on output bracket is set true.
54: stpmin and stpmax are input variables which specify lower
55: and upper bounds for the step.
57: info is an integer output variable set as follows:
58: if info = 1,2,3,4,5, then the step has been computed
59: according to one of the five cases below. otherwise
60: info = 0, and this indicates improper input parameters.
62: subprograms called
64: fortran-supplied ... abs,max,min,sqrt
66: argonne national laboratory. minpack project. june 1983
67: jorge j. more', david j. thuente
69: */
73: int TaoStep_LineSearch(TAO_SOLVER tao,double *stx,double *fx,double *dx,
74: double *sty,double *fy,double *dy,double *stp,
75: double *fp,double *dp)
76: {
77: TAO_LINESEARCH *neP = (TAO_LINESEARCH *) tao->linectx;
78: double gamma1, p, q, r, s, sgnd, stpc, stpf, stpq, theta;
79: double two = 2.0, zero = 0.0;
80: int bound;
82: TaoFunctionBegin;
83: /* Check the input parameters for errors */
84: neP->infoc = 0;
85: if (neP->bracket && (*stp <= TaoMin(*stx,*sty) || (*stp >= TaoMax(*stx,*sty))))
86: SETERRQ(1,"bad stp in bracket");
87: if (*dx * (*stp-*stx) >= zero) SETERRQ(1,"dx * (stp-stx) >= 0");
88: if (neP->stepmax < neP->stepmin) SETERRQ(1,"stepmax > stepmin");
90: /* Determine if the derivatives have opposite sign */
91: sgnd = *dp * (*dx/TaoAbsDouble(*dx));
93: /* Case 1: a higher function value.
94: the minimum is bracketed. if the cubic step is closer
95: to stx than the quadratic step, the cubic step is taken,
96: else the average of the cubic and quadratic steps is taken.
97: */
98: if (*fp > *fx) {
99: neP->infoc = 1;
100: bound = 1;
101: theta = 3 * (*fx - *fp) / (*stp - *stx) + *dx + *dp;
102: s = TaoMax(TaoAbsDouble(theta),TaoAbsDouble(*dx));
103: s = TaoMax(s,TaoAbsDouble(*dp));
104: gamma1 = s*sqrt(pow(theta/s,two) - (*dx/s)*(*dp/s));
105: if (*stp < *stx) gamma1 = -gamma1;
106: /* Can p be 0? Check */
107: p = (gamma1 - *dx) + theta;
108: q = ((gamma1 - *dx) + gamma1) + *dp;
109: r = p/q;
110: stpc = *stx + r*(*stp - *stx);
111: stpq = *stx + ((*dx/((*fx-*fp)/(*stp-*stx)+*dx))*0.5) * (*stp - *stx);
112: if (TaoAbsDouble(stpc-*stx) < TaoAbsDouble(stpq-*stx)) {
113: stpf = stpc;
114: } else {
115: stpf = stpc + 0.5*(stpq - stpc);
116: }
117: neP->bracket = 1;
118: }
119: /*
120: Case 2: A lower function value and derivatives of
121: opposite sign. the minimum is bracketed. if the cubic
122: step is closer to stx than the quadratic (secant) step,
123: the cubic step is taken, else the quadratic step is taken.
124: */
125: else if (sgnd < zero) {
126: neP->infoc = 2;
127: bound = 0;
128: theta = 3*(*fx - *fp)/(*stp - *stx) + *dx + *dp;
129: s = TaoMax(TaoAbsDouble(theta),TaoAbsDouble(*dx));
130: s = TaoMax(s,TaoAbsDouble(*dp));
131: gamma1 = s*sqrt(pow(theta/s,two) - (*dx/s)*(*dp/s));
132: if (*stp > *stx) gamma1 = -gamma1;
133: p = (gamma1 - *dp) + theta;
134: q = ((gamma1 - *dp) + gamma1) + *dx;
135: r = p/q;
136: stpc = *stp + r*(*stx - *stp);
137: stpq = *stp + (*dp/(*dp-*dx))*(*stx - *stp);
138: if (TaoAbsDouble(stpc-*stp) > TaoAbsDouble(stpq-*stp)) stpf = stpc;
139: else stpf = stpq;
140: neP->bracket = 1;
141: }
143: /* Case 3: A lower function value, derivatives of the
144: same sign, and the magnitude of the derivative decreases.
145: the cubic step is only used if the cubic tends to infinity
146: in the direction of the step or if the minimum of the cubic
147: is beyond stp. otherwise the cubic step is defined to be
148: either stepmin or stepmax. the quadratic (secant) step is also
149: computed and if the minimum is bracketed then the the step
150: closest to stx is taken, else the step farthest away is taken.
151: */
153: else if (TaoAbsDouble(*dp) < TaoAbsDouble(*dx)) {
154: neP->infoc = 3;
155: bound = 1;
156: theta = 3*(*fx - *fp)/(*stp - *stx) + *dx + *dp;
157: s = TaoMax(TaoAbsDouble(theta),TaoAbsDouble(*dx));
158: s = TaoMax(s,TaoAbsDouble(*dp));
160: /* The case gamma1 = 0 only arises if the cubic does not tend
161: to infinity in the direction of the step. */
162: gamma1 = s*sqrt(TaoMax(zero,pow(theta/s,two) - (*dx/s)*(*dp/s)));
163: if (*stp > *stx) gamma1 = -gamma1;
164: p = (gamma1 - *dp) + theta;
165: q = (gamma1 + (*dx - *dp)) + gamma1;
166: r = p/q;
167: if (r < zero && gamma1 != zero) stpc = *stp + r*(*stx - *stp);
168: else if (*stp > *stx) stpc = neP->stepmax;
169: else stpc = neP->stepmin;
170: stpq = *stp + (*dp/(*dp-*dx)) * (*stx - *stp);
171: if (neP->bracket) {
172: if (TaoAbsDouble(*stp-stpc) < TaoAbsDouble(*stp-stpq)) stpf = stpc;
173: else stpf = stpq;
174: }
175: else {
176: if (TaoAbsDouble(*stp-stpc) > TaoAbsDouble(*stp-stpq)) stpf = stpc;
177: else stpf = stpq;
178: }
179: }
181: /* Case 4: A lower function value, derivatives of the
182: same sign, and the magnitude of the derivative does
183: not decrease. if the minimum is not bracketed, the step
184: is either stpmin or stpmax, else the cubic step is taken.
185: */
186: else {
187: neP->infoc = 4;
188: bound = 0;
189: if (neP->bracket) {
190: theta = 3*(*fp - *fy)/(*sty - *stp) + *dy + *dp;
191: s = TaoMax(TaoAbsDouble(theta),TaoAbsDouble(*dy));
192: s = TaoMax(s,TaoAbsDouble(*dp));
193: gamma1 = s*sqrt(pow(theta/s,two) - (*dy/s)*(*dp/s));
194: if (*stp > *sty) gamma1 = -gamma1;
195: p = (gamma1 - *dp) + theta;
196: q = ((gamma1 - *dp) + gamma1) + *dy;
197: r = p/q;
198: stpc = *stp + r*(*sty - *stp);
199: stpf = stpc;
200: } else if (*stp > *stx) {
201: stpf = neP->stepmax;
202: } else {
203: stpf = neP->stepmin;
204: }
205: }
207: /* Update the interval of uncertainty. This update does not
208: depend on the new step or the case analysis above. */
210: if (*fp > *fx) {
211: *sty = *stp;
212: *fy = *fp;
213: *dy = *dp;
214: } else {
215: if (sgnd < zero) {
216: *sty = *stx;
217: *fy = *fx;
218: *dy = *dx;
219: }
220: *stx = *stp;
221: *fx = *fp;
222: *dx = *dp;
223: }
225: /* Compute the new step and safeguard it */
226: stpf = TaoMin(neP->stepmax,stpf);
227: stpf = TaoMax(neP->stepmin,stpf);
228: *stp = stpf;
229: if (neP->bracket && bound) {
230: if (*sty > *stx) *stp = TaoMin(*stx+0.66*(*sty-*stx),*stp);
231: else *stp = TaoMax(*stx+0.66*(*sty-*stx),*stp);
232: }
233: TaoFunctionReturn(0);
234: }
236: /* ---------------------------------------------------------- */
239: int TaoDestroy_LineSearch(TAO_SOLVER tao, void *ctx)
240: {
241: int info;
243: TaoFunctionBegin;
244: info = TaoFree(tao->linectx);CHKERRQ(info);
245: TaoFunctionReturn(0);
246: }
247: /*------------------------------------------------------------*/
250: int TaoSetOptions_LineSearch(TAO_SOLVER tao, void*linectx)
251: {
252: TAO_LINESEARCH *ctx = (TAO_LINESEARCH *)linectx;
253: int info;
255: TaoFunctionBegin;
256: info = TaoOptionsHead("More-Thuente line line search options for unconstrained minimization");CHKERRQ(info);
257: info = TaoOptionInt("-tao_ls_maxfev","max function evals in line search","",ctx->maxfev,&ctx->maxfev,0);CHKERRQ(info);
258: info = TaoOptionDouble("-tao_ls_ftol","tol for sufficient decrease","",ctx->ftol,&ctx->ftol,0);CHKERRQ(info);
259: info = TaoOptionDouble("-tao_ls_gtol","tol for curvature condition","",ctx->gtol,&ctx->gtol,0);CHKERRQ(info);
260: info = TaoOptionDouble("-tao_ls_rtol","relative tol for acceptable step","",ctx->rtol,&ctx->rtol,0);CHKERRQ(info);
261: info = TaoOptionDouble("-tao_ls_stepmin","lower bound for step","",ctx->stepmin,&ctx->stepmin,0);CHKERRQ(info);
262: info = TaoOptionDouble("-tao_ls_stepmax","upper bound for step","",ctx->stepmax,&ctx->stepmax,0);CHKERRQ(info);
263: info = TaoOptionsTail();CHKERRQ(info);
264: TaoFunctionReturn(0);
265: }
266: /*------------------------------------------------------------*/
270: int TaoView_LineSearch(TAO_SOLVER tao, void*ctx)
271: {
272: TAO_LINESEARCH *ls = (TAO_LINESEARCH *)ctx;
273: int info;
275: TaoFunctionBegin;
277: info=TaoPrintInt(tao," More'-Thuente Line search: maxf=%d,",ls->maxfev);CHKERRQ(info);
278: info=TaoPrintDouble(tao," ftol=%g,",ls->ftol);CHKERRQ(info);
279: info=TaoPrintDouble(tao," rtol=%g,",ls->rtol);CHKERRQ(info);
280: info=TaoPrintDouble(tao," gtol=%g\n",ls->gtol);CHKERRQ(info);
282: TaoFunctionReturn(0);
283: }
286: /* ---------------------------------------------------------- */
289: /* @
290: TaoApply_LineSearch - This routine performs a line search algorithm.
292: Input Parameters:
293: + tao - TAO_SOLVER context
294: . X - current iterate (on output X contains new iterate, X + step*S)
295: . S - search direction
296: . f - objective function evaluated at X
297: . G - gradient evaluated at X
298: . W - work vector
299: - step - initial estimate of step length
301: Output parameters:
302: + f - objective function evaluated at new iterate, X + step*S
303: . G - gradient evaluated at new iterate, X + step*S
304: . X - new iterate
305: . gnorm - 2-norm of G
306: - step - final step length
309: Info is set to one of:
310: + 0 - the line search succeeds; the sufficient decrease
311: condition and the directional derivative condition hold
313: negative number if an input parameter is invalid
314: . -1 - step < 0
315: . -2 - ftol < 0
316: . -3 - rtol < 0
317: . -4 - gtol < 0
318: . -5 - stepmin < 0
319: . -6 - stepmax < stepmin
320: - -7 - maxfev < 0
322: positive number > 1 if the line search otherwise terminates
323: + 2 - Relative width of the interval of uncertainty is
324: at most rtol.
325: . 3 - Maximum number of function evaluations (maxfev) has
326: been reached.
327: . 4 - Step is at the lower bound, stepmin.
328: . 5 - Step is at the upper bound, stepmax.
329: . 6 - Rounding errors may prevent further progress.
330: There may not be a step that satisfies the
331: sufficient decrease and curvature conditions.
332: Tolerances may be too small.
333: - 7 - Search direction is not a descent direction.
335: Notes:
336: This algorithm is taken from More' and Thuente, "Line search algorithms
337: with guaranteed sufficient decrease", Argonne National Laboratory,
338: Technical Report MCS-P330-1092.
340: Notes:
341: This routine is used within the following TAO unconstrained minimization
342: solvers: Newton linesearch (tao_nls), limited memory variable metric
343: (tao_lmvm), and conjugate gradient (tao_cg).
345: Level: advanced
347: .keywords: TAO_SOLVER, linesearch
348: @ */
349: static int TaoApply_LineSearch(TAO_SOLVER tao,TaoVec* X,TaoVec* G,TaoVec* S,TaoVec* W,double *f,
350: double *step,int *info2,void*ctx)
351: {
352: TAO_LINESEARCH *neP = (TAO_LINESEARCH *) tao->linectx;
353: double zero = 0.0, two = 2.0, p5 = 0.5, p66 = 0.66, xtrapf = 4.0;
354: double finit, width, width1, dginit,fm, fxm, fym, dgm, dgxm, dgym;
355: double dgx, dgy, dg, fx, fy, stx, sty, dgtest, ftest1=0.0;
356: int info, i, stage1;
358: #if defined(PETSC_USE_COMPLEX)
359: PetscScalar cdginit, cdg, cstep = 0.0;
360: #endif
362: TaoFunctionBegin;
363: /* neP->stepmin - lower bound for step */
364: /* neP->stepmax - upper bound for step */
365: /* neP->rtol - relative tolerance for an acceptable step */
366: /* neP->ftol - tolerance for sufficient decrease condition */
367: /* neP->gtol - tolerance for curvature condition */
368: /* neP->nfev - number of function evaluations */
369: /* neP->maxfev - maximum number of function evaluations */
371: /* Check input parameters for errors */
372: if (*step < zero) {
373: info = PetscLogInfo((tao,"TaoApply_LineSearch:Line search error: step (%g) < 0\n",*step)); CHKERRQ(info);
374: *info2 = -1; TaoFunctionReturn(0);
375: } else if (neP->ftol < zero) {
376: info = PetscLogInfo((tao,"TaoApply_LineSearch:Line search error: ftol (%g) < 0\n",neP->ftol)); CHKERRQ(info);
377: *info2 = -2; TaoFunctionReturn(0);
378: } else if (neP->rtol < zero) {
379: info = PetscLogInfo((tao,"TaoApply_LineSearch:Line search error: rtol (%g) < 0\n",neP->rtol)); CHKERRQ(info);
380: *info2 = -3; TaoFunctionReturn(0);
381: } else if (neP->gtol < zero) {
382: info = PetscLogInfo((tao,"TaoApply_LineSearch:Line search error: gtol (%g) < 0\n",neP->gtol)); CHKERRQ(info);
383: *info2 = -4; TaoFunctionReturn(0);
384: } else if (neP->stepmin < zero) {
385: info = PetscLogInfo((tao,"TaoApply_LineSearch:Line search error: stepmin (%g) < 0\n",neP->stepmin)); CHKERRQ(info);
386: *info2 = -5; TaoFunctionReturn(0);
387: } else if (neP->stepmax < neP->stepmin) {
388: info = PetscLogInfo((tao,"TaoApply_LineSearch:Line search error: stepmax (%g) < stepmin (%g)\n", neP->stepmax,neP->stepmin)); CHKERRQ(info);
389: *info2 = -6; TaoFunctionReturn(0);
390: } else if (neP->maxfev < zero) {
391: info = PetscLogInfo((tao,"TaoApply_LineSearch:Line search error: maxfev (%d) < 0\n",neP->maxfev)); CHKERRQ(info);
392: *info2 = -7; TaoFunctionReturn(0);
393: }
395: /* Check that search direction is a descent direction */
397: #if defined(PETSC_USE_COMPLEX)
398: info = G->Dot(S,&cdginit);CHKERRQ(info); dginit = TaoReal(cdginit);
399: #else
400: info = G->Dot(S,&dginit);CHKERRQ(info);
401: #endif
403: if (dginit >= zero) {
404: info = PetscLogInfo((tao,"TaoApply_LineSearch:Search direction not a descent direction\n")); CHKERRQ(info);
405: *info2 = 7; TaoFunctionReturn(0);
406: }
408: /* Initialization */
409: neP->bracket = 0;
410: *info2 = 0;
411: stage1 = 1;
412: finit = *f;
413: dgtest = neP->ftol * dginit;
414: width = neP->stepmax - neP->stepmin;
415: width1 = width * two;
416: info = W->CopyFrom(X);CHKERRQ(info);
417: /* Variable dictionary:
418: stx, fx, dgx - the step, function, and derivative at the best step
419: sty, fy, dgy - the step, function, and derivative at the other endpoint
420: of the interval of uncertainty
421: step, f, dg - the step, function, and derivative at the current step */
423: stx = zero;
424: fx = finit;
425: dgx = dginit;
426: sty = zero;
427: fy = finit;
428: dgy = dginit;
429:
430: neP->nfev = 0;
431: for (i=0; i< neP->maxfev; i++) {
432: /* Set min and max steps to correspond to the interval of uncertainty */
433: if (neP->bracket) {
434: neP->stepmin = TaoMin(stx,sty);
435: neP->stepmax = TaoMax(stx,sty);
436: } else {
437: neP->stepmin = stx;
438: neP->stepmax = *step + xtrapf * (*step - stx);
439: }
441: /* Force the step to be within the bounds */
442: *step = TaoMax(*step,neP->stepmin);
443: *step = TaoMin(*step,neP->stepmax);
445: /* If an unusual termination is to occur, then let step be the lowest
446: point obtained thus far */
447: if (((neP->bracket) && (*step <= neP->stepmin || *step >= neP->stepmax)) ||
448: ((neP->bracket) && (neP->stepmax - neP->stepmin <= neP->rtol * neP->stepmax)) ||
449: (neP->nfev >= neP->maxfev - 1) || (neP->infoc == 0))
450: *step = stx;
451:
452: #if defined(PETSC_USE_COMPLEX)
453: cstep = *step;
454: info = X->Waxpby(cstep,S,1.0,W);CHKERRQ(info);
455: #else
456: info = X->Waxpby(*step,S,1.0,W);CHKERRQ(info); /* X = W + step*S */
457: #endif
458: info = TaoComputeMeritFunctionGradient(tao,X,f,G);CHKERRQ(info);
459: neP->nfev++;
460: #if defined(PETSC_USE_COMPLEX)
461: info = G->Dot(S,&cdg);CHKERRQ(info); dg = TaoReal(cdg);
462: #else
463: info = G->Dot(S,&dg);CHKERRQ(info); /* dg = G^T S */
464: #endif
465: ftest1 = finit + *step * dgtest;
467: /* Convergence testing */
468: if (((neP->bracket) && (*step <= neP->stepmin||*step >= neP->stepmax)) || (!neP->infoc)) {
469: info = PetscLogInfo((tao,"TaoApply_LineSearch:Rounding errors may prevent further progress. May not be a step satisfying\n")); CHKERRQ(info);
470: info = PetscLogInfo((tao,"TaoApply_LineSearch:sufficient decrease and curvature conditions. Tolerances may be too small.\n")); CHKERRQ(info);
471: *info2 = 6; break;
472: }
473: if ((*step == neP->stepmax) && (*f <= ftest1) && (dg <= dgtest)) {
474: info = PetscLogInfo((tao,"TaoApply_LineSearch:Step is at the upper bound, stepmax (%g)\n",neP->stepmax)); CHKERRQ(info);
475: *info2 = 5; break;
476: }
477: if ((*step == neP->stepmin) && (*f >= ftest1) && (dg >= dgtest)) {
478: info = PetscLogInfo((tao,"TaoApply_LineSearch:Step is at the lower bound, stepmin (%g)\n",neP->stepmin)); CHKERRQ(info);
479: *info2 = 4; break;
480: }
481: if (neP->nfev >= neP->maxfev) {
482: info = PetscLogInfo((tao,"TaoApply_LineSearch:Number of line search function evals (%d) > maximum (%d)\n",neP->nfev,neP->maxfev)); CHKERRQ(info);
483: *info2 = 3; break;
484: }
485: if ((neP->bracket) && (neP->stepmax - neP->stepmin <= neP->rtol*neP->stepmax)){
486: info = PetscLogInfo((tao,"TaoApply_LineSearch:Relative width of interval of uncertainty is at most rtol (%g)\n",neP->rtol)); CHKERRQ(info);
487: *info2 = 2; break;
488: }
489: if ((*f <= ftest1) && (TaoAbsDouble(dg) <= neP->gtol*(-dginit))) {
490: info = PetscLogInfo((tao,"TaoApply_LineSearch:Line search success: Sufficient decrease and directional deriv conditions hold\n")); CHKERRQ(info);
491: *info2 = 0; break;
492: } else {
493: if (*f > ftest1){
494: info = PetscLogInfo((tao,"TaoApply_LineSearch: Sufficient decrease only if %g < %g\n",*f,ftest1)); CHKERRQ(info);
495: }
496: if (TaoAbsDouble(dg) > neP->gtol*(-dginit)){
497: info = PetscLogInfo((tao,"TaoApply_LineSearch: Sufficient directional deriv only if %g < %g * %g\n",TaoAbsDouble(dg),neP->gtol,(-dginit))); CHKERRQ(info);
498: }
499: }
501: /* In the first stage, we seek a step for which the modified function
502: has a nonpositive value and nonnegative derivative */
503: if ((stage1) && (*f <= ftest1) && (dg >= dginit * TaoMin(neP->ftol, neP->gtol)))
504: stage1 = 0;
506: /* A modified function is used to predict the step only if we
507: have not obtained a step for which the modified function has a
508: nonpositive function value and nonnegative derivative, and if a
509: lower function value has been obtained but the decrease is not
510: sufficient */
512: if ((stage1) && (*f <= fx) && (*f > ftest1)) {
513: fm = *f - *step * dgtest; /* Define modified function */
514: fxm = fx - stx * dgtest; /* and derivatives */
515: fym = fy - sty * dgtest;
516: dgm = dg - dgtest;
517: dgxm = dgx - dgtest;
518: dgym = dgy - dgtest;
520: /* Update the interval of uncertainty and compute the new step */
521: info = TaoStep_LineSearch(tao,&stx,&fxm,&dgxm,&sty,&fym,&dgym,step,&fm,&dgm);CHKERRQ(info);
523: fx = fxm + stx * dgtest; /* Reset the function and */
524: fy = fym + sty * dgtest; /* gradient values */
525: dgx = dgxm + dgtest;
526: dgy = dgym + dgtest;
527: } else {
528: /* Update the interval of uncertainty and compute the new step */
529: info = TaoStep_LineSearch(tao,&stx,&fx,&dgx,&sty,&fy,&dgy,step,f,&dg);CHKERRQ(info);
530: }
532: /* Force a sufficient decrease in the interval of uncertainty */
533: if (neP->bracket) {
534: if (TaoAbsDouble(sty - stx) >= p66 * width1) *step = stx + p5*(sty - stx);
535: width1 = width;
536: width = TaoAbsDouble(sty - stx);
537: }
538: }
540: /* Finish computations */
541: info = PetscLogInfo((tao,"TaoApply_LineSearch:%d function evals in line search, step = %10.4f\n",neP->nfev,*step)); CHKERRQ(info);
542: /*
543: info = G->Norm2(dginitt);CHKERRQ(info);
544: */
545: TaoFunctionReturn(0);
546: }
548: /* ---------------------------------------------------------- */
552: /*@C
553: TaoCreateMoreThuenteLineSearch - Create a line search
555: Input Parameters:
556: + tao - TAO_SOLVER context
557: . fftol - the sufficient descent parameter , greater than 0.
558: - ggtol - the curvature tolerance, greater than 0, less than 1.
561: Note:
562: If either fftol or ggtol is 0, default parameters will be used.
564: Note:
565: This algorithm is taken from More' and Thuente, "Line search algorithms
566: with guaranteed sufficient decrease", Argonne National Laboratory,
567: Technical Report MCS-P330-1092.
569: Note:
570: This line search enforces the strong Wolfe conditions for unconstrained
571: optimization. This routine is used within the following TAO unconstrained
572: minimization solvers: Newton linesearch (tao_nls), limited memory variable
573: metric (tao_lmvm), and nonlinear conjugate gradient methods.
575: Level: developer
577: .keywords: TAO_SOLVER, linesearch
578: @*/
579: int TaoCreateMoreThuenteLineSearch(TAO_SOLVER tao, double fftol, double ggtol)
580: {
581: int info;
582: TAO_LINESEARCH *neP;
584: TaoFunctionBegin;
586: info = TaoNew(TAO_LINESEARCH,&neP); CHKERRQ(info);
587: info = PetscLogObjectMemory(tao,sizeof(TAO_LINESEARCH)); CHKERRQ(info);
589: if (fftol<=0) fftol=0.001;
590: if (ggtol<=0) ggtol=0.9;
592: neP->ftol = fftol;
593: neP->rtol = 1.0e-10;
594: neP->gtol = ggtol;
595: neP->stepmin = 1.0e-20;
596: neP->stepmax = 1.0e+20;
597: neP->nfev = 0;
598: neP->bracket = 0;
599: neP->infoc = 1;
600: neP->maxfev = 30;
602: info = TaoSetLineSearch(tao,0,
603: TaoSetOptions_LineSearch,
604: TaoApply_LineSearch,
605: TaoView_LineSearch,
606: TaoDestroy_LineSearch,
607: (void *) neP);CHKERRQ(info);
609: TaoFunctionReturn(0);
610: }
615: /*
618: */