Actual source code: snesut.c

 2:  #include src/snes/snesimpl.h

  6: /*@C
  7:    SNESVecViewMonitor - Monitors progress of the SNES solvers by calling 
  8:    VecView() for the approximate solution at each iteration.

 10:    Collective on SNES

 12:    Input Parameters:
 13: +  snes - the SNES context
 14: .  its - iteration number
 15: .  fgnorm - 2-norm of residual
 16: -  dummy - either a viewer or PETSC_NULL

 18:    Level: intermediate

 20: .keywords: SNES, nonlinear, vector, monitor, view

 22: .seealso: SNESSetMonitor(), SNESDefaultMonitor(), VecView()
 23: @*/
 24: PetscErrorCode SNESVecViewMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
 25: {
 27:   Vec            x;
 28:   PetscViewer    viewer = (PetscViewer) dummy;

 31:   SNESGetSolution(snes,&x);
 32:   if (!viewer) {
 33:     MPI_Comm comm;
 34:     PetscObjectGetComm((PetscObject)snes,&comm);
 35:     viewer = PETSC_VIEWER_DRAW_(comm);
 36:   }
 37:   VecView(x,viewer);

 39:   return(0);
 40: }

 44: /*@C
 45:    SNESVecViewResidualMonitor - Monitors progress of the SNES solvers by calling 
 46:    VecView() for the residual at each iteration.

 48:    Collective on SNES

 50:    Input Parameters:
 51: +  snes - the SNES context
 52: .  its - iteration number
 53: .  fgnorm - 2-norm of residual
 54: -  dummy - either a viewer or PETSC_NULL

 56:    Level: intermediate

 58: .keywords: SNES, nonlinear, vector, monitor, view

 60: .seealso: SNESSetMonitor(), SNESDefaultMonitor(), VecView()
 61: @*/
 62: PetscErrorCode SNESVecViewResidualMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
 63: {
 65:   Vec            x;
 66:   PetscViewer    viewer = (PetscViewer) dummy;

 69:   SNESGetFunction(snes,&x,0,0);
 70:   if (!viewer) {
 71:     MPI_Comm comm;
 72:     PetscObjectGetComm((PetscObject)snes,&comm);
 73:     viewer = PETSC_VIEWER_DRAW_(comm);
 74:   }
 75:   VecView(x,viewer);

 77:   return(0);
 78: }

 82: /*@C
 83:    SNESVecViewUpdateMonitor - Monitors progress of the SNES solvers by calling 
 84:    VecView() for the UPDATE to the solution at each iteration.

 86:    Collective on SNES

 88:    Input Parameters:
 89: +  snes - the SNES context
 90: .  its - iteration number
 91: .  fgnorm - 2-norm of residual
 92: -  dummy - either a viewer or PETSC_NULL

 94:    Level: intermediate

 96: .keywords: SNES, nonlinear, vector, monitor, view

 98: .seealso: SNESSetMonitor(), SNESDefaultMonitor(), VecView()
 99: @*/
100: PetscErrorCode SNESVecViewUpdateMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
101: {
103:   Vec            x;
104:   PetscViewer    viewer = (PetscViewer) dummy;

107:   SNESGetSolutionUpdate(snes,&x);
108:   if (!viewer) {
109:     MPI_Comm comm;
110:     PetscObjectGetComm((PetscObject)snes,&comm);
111:     viewer = PETSC_VIEWER_DRAW_(comm);
112:   }
113:   VecView(x,viewer);

115:   return(0);
116: }

120: /*@C
121:    SNESDefaultMonitor - Monitors progress of the SNES solvers (default).

123:    Collective on SNES

125:    Input Parameters:
126: +  snes - the SNES context
127: .  its - iteration number
128: .  fgnorm - 2-norm of residual
129: -  dummy - unused context

131:    Notes:
132:    This routine prints the residual norm at each iteration.

134:    Level: intermediate

136: .keywords: SNES, nonlinear, default, monitor, norm

138: .seealso: SNESSetMonitor(), SNESVecViewMonitor()
139: @*/
140: PetscErrorCode SNESDefaultMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
141: {
143:   PetscViewer    viewer = (PetscViewer) dummy;

146:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(snes->comm);
147:   PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %14.12e \n",its,fgnorm);
148:   return(0);
149: }

153: /*@C
154:    SNESRatioMonitor - Monitors progress of the SNES solvers by printing the ratio
155:    of residual norm at each iteration to the previous.

157:    Collective on SNES

159:    Input Parameters:
160: +  snes - the SNES context
161: .  its - iteration number
162: .  fgnorm - 2-norm of residual (or gradient)
163: -  dummy - unused context

165:    Level: intermediate

167: .keywords: SNES, nonlinear, monitor, norm

169: .seealso: SNESSetMonitor(), SNESVecViewMonitor()
170: @*/
171: PetscErrorCode SNESRatioMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
172: {
174:   PetscInt       len;
175:   PetscReal      *history;
176:   PetscViewer    viewer;

179:   viewer = PETSC_VIEWER_STDOUT_(snes->comm);

181:   SNESGetConvergenceHistory(snes,&history,PETSC_NULL,&len);
182:   if (!its || !history || its > len) {
183:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %14.12e \n",its,fgnorm);
184:   } else {
185:     PetscReal ratio = fgnorm/history[its-1];
186:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %14.12e %g \n",its,fgnorm,ratio);
187:   }
188:   return(0);
189: }

191: /*
192:    If the we set the history monitor space then we must destroy it
193: */
196: PetscErrorCode SNESRatioMonitorDestroy(void *history)
197: {

201:   PetscFree(history);
202:   return(0);
203: }

207: /*@C
208:    SNESSetRatioMonitor - Sets SNES to use a monitor that prints the 
209:    ratio of the function norm at each iteration.

211:    Collective on SNES

213:    Input Parameters:
214: .   snes - the SNES context

216:    Level: intermediate

218: .keywords: SNES, nonlinear, monitor, norm

220: .seealso: SNESSetMonitor(), SNESVecViewMonitor(), SNESDefaultMonitor()
221: @*/
222: PetscErrorCode SNESSetRatioMonitor(SNES snes)
223: {
225:   PetscReal      *history;


229:   SNESGetConvergenceHistory(snes,&history,PETSC_NULL,PETSC_NULL);
230:   if (!history) {
231:     PetscMalloc(100*sizeof(double),&history);
232:     SNESSetConvergenceHistory(snes,history,0,100,PETSC_TRUE);
233:     SNESSetMonitor(snes,SNESRatioMonitor,history,SNESRatioMonitorDestroy);
234:   } else {
235:     SNESSetMonitor(snes,SNESRatioMonitor,0,0);
236:   }
237:   return(0);
238: }

240: /* ---------------------------------------------------------------- */
243: /*
244:      Default (short) SNES Monitor, same as SNESDefaultMonitor() except
245:   it prints fewer digits of the residual as the residual gets smaller.
246:   This is because the later digits are meaningless and are often 
247:   different on different machines; by using this routine different 
248:   machines will usually generate the same output.
249: */
250: PetscErrorCode SNESDefaultSMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
251: {

255:   if (fgnorm > 1.e-9) {
256:     PetscPrintf(snes->comm,"%3D SNES Function norm %g \n",its,fgnorm);
257:   } else if (fgnorm > 1.e-11){
258:     PetscPrintf(snes->comm,"%3D SNES Function norm %5.3e \n",its,fgnorm);
259:   } else {
260:     PetscPrintf(snes->comm,"%3D SNES Function norm < 1.e-11\n",its);
261:   }
262:   return(0);
263: }
264: /* ---------------------------------------------------------------- */
267: /*@C 
268:    SNESConverged_LS - Monitors the convergence of the solvers for
269:    systems of nonlinear equations (default).

271:    Collective on SNES

273:    Input Parameters:
274: +  snes - the SNES context
275: .  xnorm - 2-norm of current iterate
276: .  pnorm - 2-norm of current step 
277: .  fnorm - 2-norm of function
278: -  dummy - unused context

280:    Output Parameter:
281: .   reason  - one of
282: $  SNES_CONVERGED_FNORM_ABS       - (fnorm < abstol),
283: $  SNES_CONVERGED_PNORM_RELATIVE  - (pnorm < xtol*xnorm),
284: $  SNES_CONVERGED_FNORM_RELATIVE  - (fnorm < rtol*fnorm0),
285: $  SNES_DIVERGED_FUNCTION_COUNT   - (nfct > maxf),
286: $  SNES_DIVERGED_FNORM_NAN        - (fnorm == NaN),
287: $  SNES_CONVERGED_ITERATING       - (otherwise),

289:    where
290: +    maxf - maximum number of function evaluations,
291:             set with SNESSetTolerances()
292: .    nfct - number of function evaluations,
293: .    abstol - absolute function norm tolerance,
294:             set with SNESSetTolerances()
295: -    rtol - relative function norm tolerance, set with SNESSetTolerances()

297:    Level: intermediate

299: .keywords: SNES, nonlinear, default, converged, convergence

301: .seealso: SNESSetConvergenceTest(), SNESEisenstatWalkerConverged()
302: @*/
303: PetscErrorCode SNESConverged_LS(SNES snes,PetscReal xnorm,PetscReal pnorm,PetscReal fnorm,SNESConvergedReason *reason,void *dummy)
304: {
306:   if (fnorm != fnorm) {
307:     PetscLogInfo(snes,"SNESConverged_LS:Failed to converged, function norm is NaN\n");
308:     *reason = SNES_DIVERGED_FNORM_NAN;
309:   } else if (fnorm <= snes->ttol) {
310:     PetscLogInfo(snes,"SNESConverged_LS:Converged due to function norm %g < %g (relative tolerance)\n",fnorm,snes->ttol);
311:     *reason = SNES_CONVERGED_FNORM_RELATIVE;
312:   } else if (fnorm < snes->abstol) {
313:     PetscLogInfo(snes,"SNESConverged_LS:Converged due to function norm %g < %g\n",fnorm,snes->abstol);
314:     *reason = SNES_CONVERGED_FNORM_ABS;
315:   } else if (pnorm < snes->xtol*xnorm) {
316:     PetscLogInfo(snes,"SNESConverged_LS:Converged due to small update length: %g < %g * %g\n",pnorm,snes->xtol,xnorm);
317:     *reason = SNES_CONVERGED_PNORM_RELATIVE;
318:   } else if (snes->nfuncs > snes->max_funcs) {
319:     PetscLogInfo(snes,"SNESConverged_LS:Exceeded maximum number of function evaluations: %D > %D\n",snes->nfuncs,snes->max_funcs);
320:     *reason = SNES_DIVERGED_FUNCTION_COUNT ;
321:   } else {
322:     *reason = SNES_CONVERGED_ITERATING;
323:   }
324:   return(0);
325: }
326: /* ------------------------------------------------------------ */
329: /*@
330:    SNES_KSP_SetConvergenceTestEW - Sets alternative convergence test
331:    for the linear solvers within an inexact Newton method.  

333:    Collective on SNES

335:    Input Parameter:
336: .  snes - SNES context

338:    Notes:
339:    Currently, the default is to use a constant relative tolerance for 
340:    the inner linear solvers.  Alternatively, one can use the 
341:    Eisenstat-Walker method, where the relative convergence tolerance 
342:    is reset at each Newton iteration according progress of the nonlinear 
343:    solver. 

345:    Level: advanced

347:    Reference:
348:    S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an 
349:    inexact Newton method", SISC 17 (1), pp.16-32, 1996.

351: .keywords: SNES, KSP, Eisenstat, Walker, convergence, test, inexact, Newton
352: @*/
353: PetscErrorCode SNES_KSP_SetConvergenceTestEW(SNES snes)
354: {
356:   snes->ksp_ewconv = PETSC_TRUE;
357:   return(0);
358: }

362: /*@
363:    SNES_KSP_SetParametersEW - Sets parameters for Eisenstat-Walker
364:    convergence criteria for the linear solvers within an inexact
365:    Newton method.

367:    Collective on SNES
368:  
369:    Input Parameters:
370: +    snes - SNES context
371: .    version - version 1 or 2 (default is 2)
372: .    rtol_0 - initial relative tolerance (0 <= rtol_0 < 1)
373: .    rtol_max - maximum relative tolerance (0 <= rtol_max < 1)
374: .    alpha - power for version 2 rtol computation (1 < alpha <= 2)
375: .    alpha2 - power for safeguard
376: .    gamma2 - multiplicative factor for version 2 rtol computation
377:               (0 <= gamma2 <= 1)
378: -    threshold - threshold for imposing safeguard (0 < threshold < 1)

380:    Note:
381:    Use PETSC_DEFAULT to retain the default for any of the parameters.

383:    Level: advanced

385:    Reference:
386:    S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an 
387:    inexact Newton method", Utah State University Math. Stat. Dept. Res. 
388:    Report 6/94/75, June, 1994, to appear in SIAM J. Sci. Comput. 

390: .keywords: SNES, KSP, Eisenstat, Walker, set, parameters

392: .seealso: SNES_KSP_SetConvergenceTestEW()
393: @*/
394: PetscErrorCode SNES_KSP_SetParametersEW(SNES snes,PetscInt version,PetscReal rtol_0,PetscReal rtol_max,PetscReal gamma2,PetscReal alpha,
395:                                         PetscReal alpha2,PetscReal threshold)
396: {
397:   SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx*)snes->kspconvctx;

400:   if (!kctx) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context existing");
401:   if (version != PETSC_DEFAULT)   kctx->version   = version;
402:   if (rtol_0 != PETSC_DEFAULT)    kctx->rtol_0    = rtol_0;
403:   if (rtol_max != PETSC_DEFAULT)  kctx->rtol_max  = rtol_max;
404:   if (gamma2 != PETSC_DEFAULT)    kctx->gamma     = gamma2;
405:   if (alpha != PETSC_DEFAULT)     kctx->alpha     = alpha;
406:   if (alpha2 != PETSC_DEFAULT)    kctx->alpha2    = alpha2;
407:   if (threshold != PETSC_DEFAULT) kctx->threshold = threshold;
408:   if (kctx->rtol_0 < 0.0 || kctx->rtol_0 >= 1.0) {
409:     SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"0.0 <= rtol_0 < 1.0: %g",kctx->rtol_0);
410:   }
411:   if (kctx->rtol_max < 0.0 || kctx->rtol_max >= 1.0) {
412:     SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"0.0 <= rtol_max (%g) < 1.0\n",kctx->rtol_max);
413:   }
414:   if (kctx->threshold <= 0.0 || kctx->threshold >= 1.0) {
415:     SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"0.0 < threshold (%g) < 1.0\n",kctx->threshold);
416:   }
417:   if (kctx->gamma < 0.0 || kctx->gamma > 1.0) {
418:     SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"0.0 <= gamma (%g) <= 1.0\n",kctx->gamma);
419:   }
420:   if (kctx->alpha <= 1.0 || kctx->alpha > 2.0) {
421:     SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"1.0 < alpha (%g) <= 2.0\n",kctx->alpha);
422:   }
423:   if (kctx->version != 1 && kctx->version !=2) {
424:     SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Only versions 1 and 2 are supported: %D",kctx->version);
425:   }
426:   return(0);
427: }

431: PetscErrorCode SNES_KSP_EW_ComputeRelativeTolerance_Private(SNES snes,KSP ksp)
432: {
433:   SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx*)snes->kspconvctx;
434:   PetscReal           rtol = 0.0,stol;
435:   PetscErrorCode      ierr;

438:   if (!kctx) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context exists");
439:   if (!snes->iter) { /* first time in, so use the original user rtol */
440:     rtol = kctx->rtol_0;
441:   } else {
442:     if (kctx->version == 1) {
443:       rtol = (snes->norm - kctx->lresid_last)/kctx->norm_last;
444:       if (rtol < 0.0) rtol = -rtol;
445:       stol = pow(kctx->rtol_last,kctx->alpha2);
446:       if (stol > kctx->threshold) rtol = PetscMax(rtol,stol);
447:     } else if (kctx->version == 2) {
448:       rtol = kctx->gamma * pow(snes->norm/kctx->norm_last,kctx->alpha);
449:       stol = kctx->gamma * pow(kctx->rtol_last,kctx->alpha);
450:       if (stol > kctx->threshold) rtol = PetscMax(rtol,stol);
451:     } else SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Only versions 1 or 2 are supported: %D",kctx->version);
452:   }
453:   rtol = PetscMin(rtol,kctx->rtol_max);
454:   kctx->rtol_last = rtol;
455:   PetscLogInfo(snes,"SNES_KSP_EW_ComputeRelativeTolerance_Private: iter %D, Eisenstat-Walker (version %D) KSP rtol = %g\n",snes->iter,kctx->version,rtol);
456:   KSPSetTolerances(ksp,rtol,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
457:   kctx->norm_last = snes->norm;
458:   return(0);
459: }

463: PetscErrorCode SNES_KSP_EW_Converged_Private(KSP ksp,PetscInt n,PetscReal rnorm,KSPConvergedReason *reason,void *ctx)
464: {
465:   SNES                snes = (SNES)ctx;
466:   SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx*)snes->kspconvctx;
467:   PetscErrorCode      ierr;

470:   if (!kctx) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context set");
471:   if (!n) {SNES_KSP_EW_ComputeRelativeTolerance_Private(snes,ksp);}
472:   KSPDefaultConverged(ksp,n,rnorm,reason,ctx);
473:   kctx->lresid_last = rnorm;
474:   if (*reason) {
475:     PetscLogInfo(snes,"SNES_KSP_EW_Converged_Private: KSP iterations=%D, rnorm=%g\n",n,rnorm);
476:   }
477:   return(0);
478: }