Actual source code: linesearchbt.c

petsc-3.7.3 2016-08-01
Report Typos and Errors
  1: #include <petsc/private/linesearchimpl.h> /*I  "petscsnes.h"  I*/
  2: #include <petsc/private/snesimpl.h>

  4: typedef struct {
  5:   PetscReal alpha;        /* sufficient decrease parameter */
  6: } SNESLineSearch_BT;

 10: /*@
 11:    SNESLineSearchBTSetAlpha - Sets the descent parameter, alpha, in the BT linesearch variant.

 13:    Input Parameters:
 14: +  linesearch - linesearch context
 15: -  alpha - The descent parameter

 17:    Level: intermediate

 19: .seealso: SNESLineSearchSetLambda(), SNESLineSearchGetTolerances() SNESLINESEARCHBT
 20: @*/
 21: PetscErrorCode SNESLineSearchBTSetAlpha(SNESLineSearch linesearch, PetscReal alpha)
 22: {
 23:   SNESLineSearch_BT *bt;

 27:   bt        = (SNESLineSearch_BT*)linesearch->data;
 28:   bt->alpha = alpha;
 29:   return(0);
 30: }


 35: /*@
 36:    SNESLineSearchBTGetAlpha - Gets the descent parameter, alpha, in the BT linesearch variant.

 38:    Input Parameters:
 39: .  linesearch - linesearch context

 41:    Output Parameters:
 42: .  alpha - The descent parameter

 44:    Level: intermediate

 46: .seealso: SNESLineSearchGetLambda(), SNESLineSearchGetTolerances() SNESLINESEARCHBT
 47: @*/
 48: PetscErrorCode SNESLineSearchBTGetAlpha(SNESLineSearch linesearch, PetscReal *alpha)
 49: {
 50:   SNESLineSearch_BT *bt;

 54:   bt     = (SNESLineSearch_BT*)linesearch->data;
 55:   *alpha = bt->alpha;
 56:   return(0);
 57: }

 61: static PetscErrorCode  SNESLineSearchApply_BT(SNESLineSearch linesearch)
 62: {
 63:   PetscBool         changed_y,changed_w;
 64:   PetscErrorCode    ierr;
 65:   Vec               X,F,Y,W,G;
 66:   SNES              snes;
 67:   PetscReal         fnorm, xnorm, ynorm, gnorm;
 68:   PetscReal         lambda,lambdatemp,lambdaprev,minlambda,maxstep,initslope,alpha,stol;
 69:   PetscReal         t1,t2,a,b,d;
 70:   PetscReal         f;
 71:   PetscReal         g,gprev;
 72:   PetscViewer       monitor;
 73:   PetscInt          max_its,count;
 74:   SNESLineSearch_BT *bt;
 75:   Mat               jac;
 76:   PetscErrorCode    (*objective)(SNES,Vec,PetscReal*,void*);

 79:   SNESLineSearchGetVecs(linesearch, &X, &F, &Y, &W, &G);
 80:   SNESLineSearchGetNorms(linesearch, &xnorm, &fnorm, &ynorm);
 81:   SNESLineSearchGetLambda(linesearch, &lambda);
 82:   SNESLineSearchGetSNES(linesearch, &snes);
 83:   SNESLineSearchGetDefaultMonitor(linesearch, &monitor);
 84:   SNESLineSearchGetTolerances(linesearch,&minlambda,&maxstep,NULL,NULL,NULL,&max_its);
 85:   SNESGetTolerances(snes,NULL,NULL,&stol,NULL,NULL);
 86:   SNESGetObjective(snes,&objective,NULL);
 87:   bt   = (SNESLineSearch_BT*)linesearch->data;
 88:   alpha = bt->alpha;

 90:   SNESGetJacobian(snes, &jac, NULL, NULL, NULL);
 91:   if (!jac && !objective) SETERRQ(PetscObjectComm((PetscObject)linesearch), PETSC_ERR_USER, "SNESLineSearchBT requires a Jacobian matrix");

 93:   SNESLineSearchPreCheck(linesearch,X,Y,&changed_y);
 94:   SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_SUCCEEDED);

 96:   VecNormBegin(Y, NORM_2, &ynorm);
 97:   VecNormBegin(X, NORM_2, &xnorm);
 98:   VecNormEnd(Y, NORM_2, &ynorm);
 99:   VecNormEnd(X, NORM_2, &xnorm);

101:   if (ynorm == 0.0) {
102:     if (monitor) {
103:       PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);
104:       PetscViewerASCIIPrintf(monitor,"    Line search: Initial direction and size is 0\n");
105:       PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);
106:     }
107:     VecCopy(X,W);
108:     VecCopy(F,G);
109:     SNESLineSearchSetNorms(linesearch,xnorm,fnorm,ynorm);
110:     SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_REDUCT);
111:     return(0);
112:   }
113:   if (ynorm > maxstep) {        /* Step too big, so scale back */
114:     if (monitor) {
115:       PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);
116:       PetscViewerASCIIPrintf(monitor,"    Line search: Scaling step by %14.12e old ynorm %14.12e\n", (double)(maxstep/ynorm),(double)ynorm);
117:       PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);
118:     }
119:     VecScale(Y,maxstep/(ynorm));
120:     ynorm = maxstep;
121:   }

123:   /* if the SNES has an objective set, use that instead of the function value */
124:   if (objective) {
125:     SNESComputeObjective(snes,X,&f);
126:   } else {
127:     f = fnorm*fnorm;
128:   }

130:   /* compute the initial slope */
131:   if (objective) {
132:     /* slope comes from the function (assumed to be the gradient of the objective */
133:     VecDotRealPart(Y,F,&initslope);
134:   } else {
135:     /* slope comes from the normal equations */
136:     MatMult(jac,Y,W);
137:     VecDotRealPart(F,W,&initslope);
138:     if (initslope > 0.0)  initslope = -initslope;
139:     if (initslope == 0.0) initslope = -1.0;
140:   }

142:   VecWAXPY(W,-lambda,Y,X);
143:   if (linesearch->ops->viproject) {
144:     (*linesearch->ops->viproject)(snes, W);
145:   }
146:   if (snes->nfuncs >= snes->max_funcs) {
147:     PetscInfo(snes,"Exceeded maximum function evaluations, while checking full step length!\n");
148:     snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
149:     SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION);
150:     return(0);
151:   }

153:   if (objective) {
154:     SNESComputeObjective(snes,W,&g);
155:   } else {
156:     (*linesearch->ops->snesfunc)(snes,W,G);
157:     if (linesearch->ops->vinorm) {
158:       gnorm = fnorm;
159:       (*linesearch->ops->vinorm)(snes, G, W, &gnorm);
160:     } else {
161:       VecNorm(G,NORM_2,&gnorm);
162:     }
163:     g = PetscSqr(gnorm);
164:   }
165:   SNESLineSearchMonitor(linesearch);

167:   if (PetscIsInfOrNanReal(g)) {
168:     SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_NANORINF);
169:     PetscInfo(snes,"Aborted due to Nan or Inf in function evaluation\n");
170:     return(0);
171:   }
172:   if (!objective) {
173:     PetscInfo2(snes,"Initial fnorm %14.12e gnorm %14.12e\n", (double)fnorm, (double)gnorm);
174:   }
175:   if (.5*g <= .5*f + lambda*alpha*initslope) { /* Sufficient reduction or step tolerance convergence */
176:     if (monitor) {
177:       PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);
178:       if (!objective) {
179:         PetscViewerASCIIPrintf(monitor,"    Line search: Using full step: fnorm %14.12e gnorm %14.12e\n", (double)fnorm, (double)gnorm);
180:       } else {
181:         PetscViewerASCIIPrintf(monitor,"    Line search: Using full step: obj0 %14.12e obj %14.12e\n", (double)f, (double)g);
182:       }
183:       PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);
184:     }
185:   } else {
186:     /* Since the full step didn't work and the step is tiny, quit */
187:     if (stol*xnorm > ynorm) {
188:       SNESLineSearchSetNorms(linesearch,xnorm,fnorm,ynorm);
189:       SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_REDUCT);
190:       if (monitor) {
191:         PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);
192:         PetscViewerASCIIPrintf(monitor,"    Line search: Aborted due to ynorm < stol*xnorm (%14.12e < %14.12e) and inadequate full step.\n",(double)ynorm,(double)stol*xnorm);
193:         PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);
194:       }
195:       return(0);
196:     }
197:     /* Fit points with quadratic */
198:     lambdatemp = -initslope/(g - f - 2.0*lambda*initslope);
199:     lambdaprev = lambda;
200:     gprev      = g;
201:     if (lambdatemp > .5*lambda)  lambdatemp = .5*lambda;
202:     if (lambdatemp <= .1*lambda) lambda = .1*lambda;
203:     else                         lambda = lambdatemp;

205:     VecWAXPY(W,-lambda,Y,X);
206:     if (linesearch->ops->viproject) {
207:       (*linesearch->ops->viproject)(snes, W);
208:     }
209:     if (snes->nfuncs >= snes->max_funcs) {
210:       PetscInfo1(snes,"Exceeded maximum function evaluations, while attempting quadratic backtracking! %D \n",snes->nfuncs);
211:       snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
212:       SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION);
213:       return(0);
214:     }
215:     if (objective) {
216:       SNESComputeObjective(snes,W,&g);
217:     } else {
218:       (*linesearch->ops->snesfunc)(snes,W,G);
219:       if (linesearch->ops->vinorm) {
220:         gnorm = fnorm;
221:         (*linesearch->ops->vinorm)(snes, G, W, &gnorm);
222:       } else {
223:         VecNorm(G,NORM_2,&gnorm);
224:       }
225:       g = PetscSqr(gnorm);
226:     }
227:     if (PetscIsInfOrNanReal(g)) {
228:       SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_NANORINF);
229:       PetscInfo(snes,"Aborted due to Nan or Inf in function evaluation\n");
230:       return(0);
231:     }
232:     if (monitor) {
233:       PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);
234:       if (!objective) {
235:         PetscViewerASCIIPrintf(monitor,"    Line search: gnorm after quadratic fit %14.12e\n",(double)gnorm);
236:       } else {
237:         PetscViewerASCIIPrintf(monitor,"    Line search: obj after quadratic fit %14.12e\n",(double)g);
238:       }
239:       PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);
240:     }
241:     if (.5*g < .5*f + lambda*alpha*initslope) { /* sufficient reduction */
242:       if (monitor) {
243:         PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);
244:         PetscViewerASCIIPrintf(monitor,"    Line search: Quadratically determined step, lambda=%18.16e\n",(double)lambda);
245:         PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);
246:       }
247:     } else {
248:       /* Fit points with cubic */
249:       for (count = 0; count < max_its; count++) {
250:         if (lambda <= minlambda) {
251:           if (monitor) {
252:             PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);
253:             PetscViewerASCIIPrintf(monitor,"    Line search: unable to find good step length! After %D tries \n",count);
254:             if (!objective) {
255:               PetscViewerASCIIPrintf(monitor,
256:                                             "    Line search: fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, minlambda=%18.16e, lambda=%18.16e, initial slope=%18.16e\n",
257:                                             (double)fnorm, (double)gnorm, (double)ynorm, (double)minlambda, (double)lambda, (double)initslope);
258:             } else {
259:               PetscViewerASCIIPrintf(monitor,
260:                                             "    Line search: obj(0)=%18.16e, obj=%18.16e, ynorm=%18.16e, minlambda=%18.16e, lambda=%18.16e, initial slope=%18.16e\n",
261:                                             (double)f, (double)g, (double)ynorm, (double)minlambda, (double)lambda, (double)initslope);
262:             }
263:             PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);
264:           }
265:           SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_REDUCT);
266:           return(0);
267:         }
268:         if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
269:           t1 = .5*(g - f) - lambda*initslope;
270:           t2 = .5*(gprev  - f) - lambdaprev*initslope;
271:           a  = (t1/(lambda*lambda) - t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev);
272:           b  = (-lambdaprev*t1/(lambda*lambda) + lambda*t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev);
273:           d  = b*b - 3*a*initslope;
274:           if (d < 0.0) d = 0.0;
275:           if (a == 0.0) lambdatemp = -initslope/(2.0*b);
276:           else lambdatemp = (-b + PetscSqrtReal(d))/(3.0*a);

278:         } else if (linesearch->order == SNES_LINESEARCH_ORDER_QUADRATIC) {
279:           lambdatemp = -initslope/(g - f - 2.0*initslope);
280:         } else SETERRQ(PetscObjectComm((PetscObject)linesearch), PETSC_ERR_SUP, "unsupported line search order for type bt");
281:         lambdaprev = lambda;
282:         gprev      = g;
283:         if (lambdatemp > .5*lambda)  lambdatemp = .5*lambda;
284:         if (lambdatemp <= .1*lambda) lambda     = .1*lambda;
285:         else                         lambda     = lambdatemp;
286:         VecWAXPY(W,-lambda,Y,X);
287:         if (linesearch->ops->viproject) {
288:           (*linesearch->ops->viproject)(snes,W);
289:         }
290:         if (snes->nfuncs >= snes->max_funcs) {
291:           PetscInfo1(snes,"Exceeded maximum function evaluations, while looking for good step length! %D \n",count);
292:           if (!objective) {
293:             PetscInfo5(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lambda=%18.16e, initial slope=%18.16e\n",
294:                               (double)fnorm,(double)gnorm,(double)ynorm,(double)lambda,(double)initslope);
295:           }
296:           SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION);
297:           snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
298:           return(0);
299:         }
300:         if (objective) {
301:           SNESComputeObjective(snes,W,&g);
302:         } else {
303:           (*linesearch->ops->snesfunc)(snes,W,G);
304:           if (linesearch->ops->vinorm) {
305:             gnorm = fnorm;
306:             (*linesearch->ops->vinorm)(snes, G, W, &gnorm);
307:           } else {
308:             VecNorm(G,NORM_2,&gnorm);
309:           }
310:           g = PetscSqr(gnorm);
311:         }
312:         if (PetscIsInfOrNanReal(gnorm)) {
313:           SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_NANORINF);
314:           PetscInfo(snes,"Aborted due to Nan or Inf in function evaluation\n");
315:           return(0);
316:         }
317:         if (.5*g < .5*f + lambda*alpha*initslope) { /* is reduction enough? */
318:           if (monitor) {
319:             PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);
320:             if (!objective) {
321:               if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
322:                 PetscViewerASCIIPrintf(monitor,"    Line search: Cubically determined step, current gnorm %14.12e lambda=%18.16e\n",(double)gnorm,(double)lambda);
323:               } else {
324:                 PetscViewerASCIIPrintf(monitor,"    Line search: Quadratically determined step, current gnorm %14.12e lambda=%18.16e\n",(double)gnorm,(double)lambda);
325:               }
326:               PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);
327:             } else {
328:               if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
329:                 PetscViewerASCIIPrintf(monitor,"    Line search: Cubically determined step, obj %14.12e lambda=%18.16e\n",(double)g,(double)lambda);
330:               } else {
331:                 PetscViewerASCIIPrintf(monitor,"    Line search: Quadratically determined step, obj %14.12e lambda=%18.16e\n",(double)g,(double)lambda);
332:               }
333:               PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);
334:             }
335:           }
336:           break;
337:         } else if (monitor) {
338:           PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);
339:           if (!objective) {
340:             if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
341:               PetscViewerASCIIPrintf(monitor,"    Line search: Cubic step no good, shrinking lambda, current gnorm %12.12e lambda=%18.16e\n",(double)gnorm,(double)lambda);
342:             } else {
343:               PetscViewerASCIIPrintf(monitor,"    Line search: Quadratic step no good, shrinking lambda, current gnorm %12.12e lambda=%18.16e\n",(double)gnorm,(double)lambda);
344:             }
345:             PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);
346:           } else {
347:             if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
348:               PetscViewerASCIIPrintf(monitor,"    Line search: Cubic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n",(double)g,(double)lambda);
349:             } else {
350:               PetscViewerASCIIPrintf(monitor,"    Line search: Quadratic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n",(double)g,(double)lambda);
351:             }
352:             PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);
353:           }
354:         }
355:       }
356:     }
357:   }

359:   /* postcheck */
360:   SNESLineSearchPostCheck(linesearch,X,Y,W,&changed_y,&changed_w);
361:   if (changed_y) {
362:     VecWAXPY(W,-lambda,Y,X);
363:     if (linesearch->ops->viproject) {
364:       (*linesearch->ops->viproject)(snes, W);
365:     }
366:   }
367:   if (changed_y || changed_w || objective) { /* recompute the function norm if the step has changed or the objective isn't the norm */
368:     (*linesearch->ops->snesfunc)(snes,W,G);
369:     if (linesearch->ops->vinorm) {
370:       gnorm = fnorm;
371:       (*linesearch->ops->vinorm)(snes, G, W, &gnorm);
372:     } else {
373:       VecNorm(G,NORM_2,&gnorm);
374:     }
375:     VecNorm(Y,NORM_2,&ynorm);
376:     if (PetscIsInfOrNanReal(gnorm)) {
377:       SNESLineSearchSetReason(linesearch,SNES_LINESEARCH_FAILED_NANORINF);
378:       PetscInfo(snes,"Aborted due to Nan or Inf in function evaluation\n");
379:       return(0);
380:     }
381:   }

383:   /* copy the solution over */
384:   VecCopy(W, X);
385:   VecCopy(G, F);
386:   VecNorm(X, NORM_2, &xnorm);
387:   SNESLineSearchSetLambda(linesearch, lambda);
388:   SNESLineSearchSetNorms(linesearch, xnorm, gnorm, ynorm);
389:   return(0);
390: }

394: PetscErrorCode SNESLineSearchView_BT(SNESLineSearch linesearch, PetscViewer viewer)
395: {
396:   PetscErrorCode    ierr;
397:   PetscBool         iascii;
398:   SNESLineSearch_BT *bt;

401:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
402:   bt   = (SNESLineSearch_BT*)linesearch->data;
403:   if (iascii) {
404:     if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
405:       PetscViewerASCIIPrintf(viewer, "  interpolation: cubic\n");
406:     } else if (linesearch->order == SNES_LINESEARCH_ORDER_QUADRATIC) {
407:       PetscViewerASCIIPrintf(viewer, "  interpolation: quadratic\n");
408:     }
409:     PetscViewerASCIIPrintf(viewer, "  alpha=%e\n", (double)bt->alpha);
410:   }
411:   return(0);
412: }


417: static PetscErrorCode SNESLineSearchDestroy_BT(SNESLineSearch linesearch)
418: {

422:   PetscFree(linesearch->data);
423:   return(0);
424: }


429: static PetscErrorCode SNESLineSearchSetFromOptions_BT(PetscOptionItems *PetscOptionsObject,SNESLineSearch linesearch)
430: {

432:   PetscErrorCode    ierr;
433:   SNESLineSearch_BT *bt = (SNESLineSearch_BT*)linesearch->data;

436:   PetscOptionsHead(PetscOptionsObject,"SNESLineSearch BT options");
437:   PetscOptionsReal("-snes_linesearch_alpha",   "Descent tolerance",        "SNESLineSearchBT", bt->alpha, &bt->alpha, NULL);
438:   PetscOptionsTail();
439:   return(0);
440: }


445: /*MC
446:    SNESLINESEARCHBT - Backtracking line search.

448:    This line search finds the minimum of a polynomial fitting of the L2 norm of the
449:    function or the objective function if it is provided with SNESSetObjective(). If this fit does not satisfy the conditions for progress, the interval shrinks
450:    and the fit is reattempted at most max_it times or until lambda is below minlambda.

452:    Options Database Keys:
453: +  -snes_linesearch_alpha<1e-4> - slope descent parameter
454: .  -snes_linesearch_damping<1.0> - initial step length
455: .  -snes_linesearch_maxstep <length> - if the length the initial step is larger than this then the
456:                                        step is scaled back to be of this length at the beginning of the line search
457: .  -snes_linesearch_max_it<40> - maximum number of shrinking step
458: .  -snes_linesearch_minlambda<1e-12> - minimum step length allowed
459: -  -snes_linesearch_order<cubic,quadratic> - order of the approximation

461:    Level: advanced

463:    Notes:
464:    This line search is taken from "Numerical Methods for Unconstrained
465:    Optimization and Nonlinear Equations" by Dennis and Schnabel, page 325.

467: .keywords: SNES, SNESLineSearch, damping

469: .seealso: SNESLineSearchCreate(), SNESLineSearchSetType()
470: M*/
471: PETSC_EXTERN PetscErrorCode SNESLineSearchCreate_BT(SNESLineSearch linesearch)
472: {

474:   SNESLineSearch_BT *bt;
475:   PetscErrorCode    ierr;

478:   linesearch->ops->apply          = SNESLineSearchApply_BT;
479:   linesearch->ops->destroy        = SNESLineSearchDestroy_BT;
480:   linesearch->ops->setfromoptions = SNESLineSearchSetFromOptions_BT;
481:   linesearch->ops->reset          = NULL;
482:   linesearch->ops->view           = SNESLineSearchView_BT;
483:   linesearch->ops->setup          = NULL;

485:   PetscNewLog(linesearch,&bt);

487:   linesearch->data    = (void*)bt;
488:   linesearch->max_its = 40;
489:   linesearch->order   = SNES_LINESEARCH_ORDER_CUBIC;
490:   bt->alpha           = 1e-4;
491:   return(0);
492: }