Actual source code: snesut.c
1: #define PETSCSNES_DLL
3: #include src/snes/snesimpl.h
7: /*@C
8: SNESVecViewMonitor - Monitors progress of the SNES solvers by calling
9: VecView() for the approximate solution at each iteration.
11: Collective on SNES
13: Input Parameters:
14: + snes - the SNES context
15: . its - iteration number
16: . fgnorm - 2-norm of residual
17: - dummy - either a viewer or PETSC_NULL
19: Level: intermediate
21: .keywords: SNES, nonlinear, vector, monitor, view
23: .seealso: SNESSetMonitor(), SNESDefaultMonitor(), VecView()
24: @*/
25: PetscErrorCode SNESVecViewMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
26: {
28: Vec x;
29: PetscViewer viewer = (PetscViewer) dummy;
32: SNESGetSolution(snes,&x);
33: if (!viewer) {
34: MPI_Comm comm;
35: PetscObjectGetComm((PetscObject)snes,&comm);
36: viewer = PETSC_VIEWER_DRAW_(comm);
37: }
38: VecView(x,viewer);
40: return(0);
41: }
45: /*@C
46: SNESVecViewResidualMonitor - Monitors progress of the SNES solvers by calling
47: VecView() for the residual at each iteration.
49: Collective on SNES
51: Input Parameters:
52: + snes - the SNES context
53: . its - iteration number
54: . fgnorm - 2-norm of residual
55: - dummy - either a viewer or PETSC_NULL
57: Level: intermediate
59: .keywords: SNES, nonlinear, vector, monitor, view
61: .seealso: SNESSetMonitor(), SNESDefaultMonitor(), VecView()
62: @*/
63: PetscErrorCode SNESVecViewResidualMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
64: {
66: Vec x;
67: PetscViewer viewer = (PetscViewer) dummy;
70: SNESGetFunction(snes,&x,0,0);
71: if (!viewer) {
72: MPI_Comm comm;
73: PetscObjectGetComm((PetscObject)snes,&comm);
74: viewer = PETSC_VIEWER_DRAW_(comm);
75: }
76: VecView(x,viewer);
78: return(0);
79: }
83: /*@C
84: SNESVecViewUpdateMonitor - Monitors progress of the SNES solvers by calling
85: VecView() for the UPDATE to the solution at each iteration.
87: Collective on SNES
89: Input Parameters:
90: + snes - the SNES context
91: . its - iteration number
92: . fgnorm - 2-norm of residual
93: - dummy - either a viewer or PETSC_NULL
95: Level: intermediate
97: .keywords: SNES, nonlinear, vector, monitor, view
99: .seealso: SNESSetMonitor(), SNESDefaultMonitor(), VecView()
100: @*/
101: PetscErrorCode SNESVecViewUpdateMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
102: {
104: Vec x;
105: PetscViewer viewer = (PetscViewer) dummy;
108: SNESGetSolutionUpdate(snes,&x);
109: if (!viewer) {
110: MPI_Comm comm;
111: PetscObjectGetComm((PetscObject)snes,&comm);
112: viewer = PETSC_VIEWER_DRAW_(comm);
113: }
114: VecView(x,viewer);
116: return(0);
117: }
121: /*@C
122: SNESDefaultMonitor - Monitors progress of the SNES solvers (default).
124: Collective on SNES
126: Input Parameters:
127: + snes - the SNES context
128: . its - iteration number
129: . fgnorm - 2-norm of residual
130: - dummy - unused context
132: Notes:
133: This routine prints the residual norm at each iteration.
135: Level: intermediate
137: .keywords: SNES, nonlinear, default, monitor, norm
139: .seealso: SNESSetMonitor(), SNESVecViewMonitor()
140: @*/
141: PetscErrorCode SNESDefaultMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
142: {
144: PetscViewer viewer = (PetscViewer) dummy;
147: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(snes->comm);
148: PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %14.12e \n",its,fgnorm);
149: return(0);
150: }
152: typedef struct {
153: PetscViewer viewer;
154: PetscReal *history;
155: } SNESRatioMonitorContext;
159: /*@C
160: SNESRatioMonitor - Monitors progress of the SNES solvers by printing the ratio
161: of residual norm at each iteration to the previous.
163: Collective on SNES
165: Input Parameters:
166: + snes - the SNES context
167: . its - iteration number
168: . fgnorm - 2-norm of residual (or gradient)
169: - dummy - context of monitor
171: Level: intermediate
173: .keywords: SNES, nonlinear, monitor, norm
175: .seealso: SNESSetMonitor(), SNESVecViewMonitor()
176: @*/
177: PetscErrorCode SNESRatioMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
178: {
179: PetscErrorCode ierr;
180: PetscInt len;
181: PetscReal *history;
182: SNESRatioMonitorContext *ctx = (SNESRatioMonitorContext*)dummy;
185: SNESGetConvergenceHistory(snes,&history,PETSC_NULL,&len);
186: if (!its || !history || its > len) {
187: PetscViewerASCIIPrintf(ctx->viewer,"%3D SNES Function norm %14.12e \n",its,fgnorm);
188: } else {
189: PetscReal ratio = fgnorm/history[its-1];
190: PetscViewerASCIIPrintf(ctx->viewer,"%3D SNES Function norm %14.12e %G \n",its,fgnorm,ratio);
191: }
192: return(0);
193: }
195: /*
196: If the we set the history monitor space then we must destroy it
197: */
200: PetscErrorCode SNESRatioMonitorDestroy(void *ct)
201: {
202: PetscErrorCode ierr;
203: SNESRatioMonitorContext *ctx = (SNESRatioMonitorContext*)ct;
206: PetscFree(ctx->history);
207: PetscViewerDestroy(ctx->viewer);
208: PetscFree(ctx);
209: return(0);
210: }
214: /*@C
215: SNESSetRatioMonitor - Sets SNES to use a monitor that prints the
216: ratio of the function norm at each iteration.
218: Collective on SNES
220: Input Parameters:
221: + snes - the SNES context
222: - viewer - ASCII viewer to print output
224: Level: intermediate
226: .keywords: SNES, nonlinear, monitor, norm
228: .seealso: SNESSetMonitor(), SNESVecViewMonitor(), SNESDefaultMonitor()
229: @*/
230: PetscErrorCode SNESSetRatioMonitor(SNES snes,PetscViewer viewer)
231: {
232: PetscErrorCode ierr;
233: SNESRatioMonitorContext *ctx;
234: PetscReal *history;
237: if (!viewer) {
238: viewer = PETSC_VIEWER_STDOUT_(snes->comm);
239: PetscObjectReference((PetscObject)viewer);
240: }
241: PetscNew(SNESRatioMonitorContext,&ctx);
242: SNESGetConvergenceHistory(snes,&history,PETSC_NULL,PETSC_NULL);
243: if (!history) {
244: PetscMalloc(100*sizeof(PetscReal),&ctx->history);
245: SNESSetConvergenceHistory(snes,ctx->history,0,100,PETSC_TRUE);
246: }
247: ctx->viewer = viewer;
248: SNESSetMonitor(snes,SNESRatioMonitor,ctx,SNESRatioMonitorDestroy);
249: return(0);
250: }
252: /* ---------------------------------------------------------------- */
255: /*
256: Default (short) SNES Monitor, same as SNESDefaultMonitor() except
257: it prints fewer digits of the residual as the residual gets smaller.
258: This is because the later digits are meaningless and are often
259: different on different machines; by using this routine different
260: machines will usually generate the same output.
261: */
262: PetscErrorCode SNESDefaultSMonitor(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
263: {
265: PetscViewer viewer = (PetscViewer) dummy;
268: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(snes->comm);
269: if (fgnorm > 1.e-9) {
270: PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %G \n",its,fgnorm);
271: } else if (fgnorm > 1.e-11){
272: PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %5.3e \n",its,fgnorm);
273: } else {
274: PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm < 1.e-11\n",its);
275: }
276: return(0);
277: }
278: /* ---------------------------------------------------------------- */
281: /*@C
282: SNESConverged_LS - Monitors the convergence of the solvers for
283: systems of nonlinear equations (default).
285: Collective on SNES
287: Input Parameters:
288: + snes - the SNES context
289: . it - the iteration (0 indicates before any Newton steps)
290: . xnorm - 2-norm of current iterate
291: . pnorm - 2-norm of current step
292: . fnorm - 2-norm of function at current iterate
293: - dummy - unused context
295: Output Parameter:
296: . reason - one of
297: $ SNES_CONVERGED_FNORM_ABS - (fnorm < abstol),
298: $ SNES_CONVERGED_PNORM_RELATIVE - (pnorm < xtol*xnorm),
299: $ SNES_CONVERGED_FNORM_RELATIVE - (fnorm < rtol*fnorm0),
300: $ SNES_DIVERGED_FUNCTION_COUNT - (nfct > maxf),
301: $ SNES_DIVERGED_FNORM_NAN - (fnorm == NaN),
302: $ SNES_CONVERGED_ITERATING - (otherwise),
304: where
305: + maxf - maximum number of function evaluations,
306: set with SNESSetTolerances()
307: . nfct - number of function evaluations,
308: . abstol - absolute function norm tolerance,
309: set with SNESSetTolerances()
310: - rtol - relative function norm tolerance, set with SNESSetTolerances()
312: Level: intermediate
314: .keywords: SNES, nonlinear, default, converged, convergence
316: .seealso: SNESSetConvergenceTest(), SNESEisenstatWalkerConverged()
317: @*/
318: PetscErrorCode SNESConverged_LS(SNES snes,PetscInt it,PetscReal xnorm,PetscReal pnorm,PetscReal fnorm,SNESConvergedReason *reason,void *dummy)
319: {
323: *reason = SNES_CONVERGED_ITERATING;
325: if (!it) {
326: /* set parameter for default relative tolerance convergence test */
327: snes->ttol = fnorm*snes->rtol;
328: }
329: if (fnorm != fnorm) {
330: PetscInfo(snes,"Failed to converged, function norm is NaN\n");
331: *reason = SNES_DIVERGED_FNORM_NAN;
332: } else if (fnorm < snes->abstol) {
333: PetscInfo2(snes,"Converged due to function norm %G < %G\n",fnorm,snes->abstol);
334: *reason = SNES_CONVERGED_FNORM_ABS;
335: } else if (snes->nfuncs >= snes->max_funcs) {
336: PetscInfo2(snes,"Exceeded maximum number of function evaluations: %D > %D\n",snes->nfuncs,snes->max_funcs);
337: *reason = SNES_DIVERGED_FUNCTION_COUNT;
338: }
340: if (it && !*reason) {
341: if (fnorm <= snes->ttol) {
342: PetscInfo2(snes,"Converged due to function norm %G < %G (relative tolerance)\n",fnorm,snes->ttol);
343: *reason = SNES_CONVERGED_FNORM_RELATIVE;
344: } else if (pnorm < snes->xtol*xnorm) {
345: PetscInfo3(snes,"Converged due to small update length: %G < %G * %G\n",pnorm,snes->xtol,xnorm);
346: *reason = SNES_CONVERGED_PNORM_RELATIVE;
347: }
348: }
349: return(0);
350: }
351: /* ------------------------------------------------------------ */
354: /*@
355: SNES_KSP_SetConvergenceTestEW - Sets alternative convergence test
356: for the linear solvers within an inexact Newton method.
358: Collective on SNES
360: Input Parameter:
361: . snes - SNES context
363: Notes:
364: Currently, the default is to use a constant relative tolerance for
365: the inner linear solvers. Alternatively, one can use the
366: Eisenstat-Walker method, where the relative convergence tolerance
367: is reset at each Newton iteration according progress of the nonlinear
368: solver.
370: Level: advanced
372: Reference:
373: S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an
374: inexact Newton method", SISC 17 (1), pp.16-32, 1996.
376: .keywords: SNES, KSP, Eisenstat, Walker, convergence, test, inexact, Newton
377: @*/
378: PetscErrorCode SNES_KSP_SetConvergenceTestEW(SNES snes)
379: {
381: snes->ksp_ewconv = PETSC_TRUE;
382: return(0);
383: }
387: /*@
388: SNES_KSP_SetParametersEW - Sets parameters for Eisenstat-Walker
389: convergence criteria for the linear solvers within an inexact
390: Newton method.
392: Collective on SNES
393:
394: Input Parameters:
395: + snes - SNES context
396: . version - version 1, 2 (default is 2) or 3
397: . rtol_0 - initial relative tolerance (0 <= rtol_0 < 1)
398: . rtol_max - maximum relative tolerance (0 <= rtol_max < 1)
399: . alpha - power for version 2 rtol computation (1 < alpha <= 2)
400: . alpha2 - power for safeguard
401: . gamma2 - multiplicative factor for version 2 rtol computation
402: (0 <= gamma2 <= 1)
403: - threshold - threshold for imposing safeguard (0 < threshold < 1)
405: Note:
406: Version 3 was contributed by Luis Chacon, June 2006.
408: Use PETSC_DEFAULT to retain the default for any of the parameters.
410: Level: advanced
412: Reference:
413: S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an
414: inexact Newton method", Utah State University Math. Stat. Dept. Res.
415: Report 6/94/75, June, 1994, to appear in SIAM J. Sci. Comput.
417: .keywords: SNES, KSP, Eisenstat, Walker, set, parameters
419: .seealso: SNES_KSP_SetConvergenceTestEW()
420: @*/
421: PetscErrorCode SNES_KSP_SetParametersEW(SNES snes,PetscInt version,PetscReal rtol_0,PetscReal rtol_max,PetscReal gamma2,PetscReal alpha,
422: PetscReal alpha2,PetscReal threshold)
423: {
424: SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx*)snes->kspconvctx;
427: if (!kctx) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context existing");
428: if (version != PETSC_DEFAULT) kctx->version = version;
429: if (rtol_0 != PETSC_DEFAULT) kctx->rtol_0 = rtol_0;
430: if (rtol_max != PETSC_DEFAULT) kctx->rtol_max = rtol_max;
431: if (gamma2 != PETSC_DEFAULT) kctx->gamma = gamma2;
432: if (alpha != PETSC_DEFAULT) kctx->alpha = alpha;
433: if (alpha2 != PETSC_DEFAULT) kctx->alpha2 = alpha2;
434: if (threshold != PETSC_DEFAULT) kctx->threshold = threshold;
435: if (kctx->rtol_0 < 0.0 || kctx->rtol_0 >= 1.0) {
436: SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"0.0 <= rtol_0 < 1.0: %G",kctx->rtol_0);
437: }
438: if (kctx->rtol_max < 0.0 || kctx->rtol_max >= 1.0) {
439: SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"0.0 <= rtol_max (%G) < 1.0\n",kctx->rtol_max);
440: }
441: if (kctx->threshold <= 0.0 || kctx->threshold >= 1.0) {
442: SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"0.0 < threshold (%G) < 1.0\n",kctx->threshold);
443: }
444: if (kctx->gamma < 0.0 || kctx->gamma > 1.0) {
445: SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"0.0 <= gamma (%G) <= 1.0\n",kctx->gamma);
446: }
447: if (kctx->alpha <= 1.0 || kctx->alpha > 2.0) {
448: SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"1.0 < alpha (%G) <= 2.0\n",kctx->alpha);
449: }
450: if (kctx->version < 1 || kctx->version > 3) {
451: SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Only versions 1, 2 and 3 are supported: %D",kctx->version);
452: }
453: return(0);
454: }
458: PetscErrorCode SNES_KSP_EW_ComputeRelativeTolerance_Private(SNES snes,KSP ksp)
459: {
460: SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx*)snes->kspconvctx;
461: PetscReal rtol = 0.0,stol;
462: PetscErrorCode ierr;
465: if (!kctx) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context exists");
466: if (!snes->iter) { /* first time in, so use the original user rtol */
467: rtol = kctx->rtol_0;
468: } else {
469: if (kctx->version == 1) {
470: rtol = (snes->norm - kctx->lresid_last)/kctx->norm_last;
471: if (rtol < 0.0) rtol = -rtol;
472: stol = pow(kctx->rtol_last,kctx->alpha2);
473: if (stol > kctx->threshold) rtol = PetscMax(rtol,stol);
474: } else if (kctx->version == 2) {
475: rtol = kctx->gamma * pow(snes->norm/kctx->norm_last,kctx->alpha);
476: stol = kctx->gamma * pow(kctx->rtol_last,kctx->alpha);
477: if (stol > kctx->threshold) rtol = PetscMax(rtol,stol);
478: } else if (kctx->version == 3) {
479: rtol = kctx->gamma * pow(snes->norm/kctx->norm_last,kctx->alpha);
480: /* safeguard: avoid sharp decrease of rtol */
481: rtol = PetscMin(kctx->rtol_0,PetscMax(rtol,kctx->gamma*pow(kctx->rtol_last,kctx->alpha)));
482: /* safeguard: avoid oversolving */
483: rtol = PetscMin(kctx->rtol_0,PetscMax(rtol,kctx->gamma*(snes->ttol)/snes->norm));
485: } else SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Only versions 1, 2 or 3 are supported: %D",kctx->version);
486: }
487: rtol = PetscMin(rtol,kctx->rtol_max);
488: kctx->rtol_last = rtol;
489: PetscInfo3(snes,"iter %D, Eisenstat-Walker (version %D) KSP rtol = %G\n",snes->iter,kctx->version,rtol);
490: KSPSetTolerances(ksp,rtol,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
491: kctx->norm_last = snes->norm;
492: return(0);
493: }
497: PetscErrorCode SNES_KSP_EW_Converged_Private(KSP ksp,PetscInt n,PetscReal rnorm,KSPConvergedReason *reason,void *ctx)
498: {
499: SNES snes = (SNES)ctx;
500: SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx*)snes->kspconvctx;
501: PetscErrorCode ierr;
504: if (!kctx) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context set");
505: if (!n) {SNES_KSP_EW_ComputeRelativeTolerance_Private(snes,ksp);}
506: KSPDefaultConverged(ksp,n,rnorm,reason,ctx);
507: kctx->lresid_last = rnorm;
508: if (*reason) {
509: PetscInfo2(snes,"KSP iterations=%D, rnorm=%G\n",n,rnorm);
510: }
511: return(0);
512: }