Actual source code: snesut.c
1: #define PETSCSNES_DLL
3: #include private/snesimpl.h
7: /*@C
8: SNESMonitorSolution - 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: SNESMonitorSet(), SNESMonitorDefault(), VecView()
24: @*/
25: PetscErrorCode SNESMonitorSolution(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: SNESMonitorResidual - 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: SNESMonitorSet(), SNESMonitorDefault(), VecView()
62: @*/
63: PetscErrorCode SNESMonitorResidual(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: SNESMonitorSolutionUpdate - 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: SNESMonitorSet(), SNESMonitorDefault(), VecView()
100: @*/
101: PetscErrorCode SNESMonitorSolutionUpdate(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: SNESMonitorDefault - 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: SNESMonitorSet(), SNESMonitorSolution()
140: @*/
141: PetscErrorCode SNESMonitorDefault(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
142: {
143: PetscErrorCode ierr;
144: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;
147: if (!dummy) {
148: PetscViewerASCIIMonitorCreate(((PetscObject)snes)->comm,"stdout",0,&viewer);
149: }
150: PetscViewerASCIIMonitorPrintf(viewer,"%3D SNES Function norm %14.12e \n",its,fgnorm);
151: if (!dummy) {
152: PetscViewerASCIIMonitorDestroy(viewer);
153: }
154: return(0);
155: }
159: PetscErrorCode SNESMonitorRange_Private(SNES snes,PetscInt it,PetscReal *per)
160: {
161: PetscErrorCode ierr;
162: Vec resid;
163: PetscReal rmax,pwork;
164: PetscInt i,n,N;
165: PetscScalar *r;
168: SNESGetFunction(snes,&resid,0,0);
169: VecNorm(resid,NORM_INFINITY,&rmax);
170: VecGetLocalSize(resid,&n);
171: VecGetSize(resid,&N);
172: VecGetArray(resid,&r);
173: pwork = 0.0;
174: for (i=0; i<n; i++) {
175: pwork += (PetscAbsScalar(r[i]) > .20*rmax);
176: }
177: MPI_Allreduce(&pwork,per,1,MPIU_REAL,PetscSum_Op,((PetscObject)snes)->comm);
178: VecRestoreArray(resid,&r);
179: *per = *per/N;
180: return(0);
181: }
185: /*@C
186: SNESMonitorRange - Prints the percentage of residual elements that are more then 10 percent of the maximum value.
188: Collective on SNES
190: Input Parameters:
191: + snes - iterative context
192: . it - iteration number
193: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
194: - dummy - unused monitor context
196: Options Database Key:
197: . -snes_monitor_range - Activates SNESMonitorRange()
200: Notes:
201: When using either ICC or ILU preconditioners in BlockSolve95
202: (via MATMPIROWBS matrix format), then use this monitor will
203: print both the residual norm associated with the original
204: (unscaled) matrix.
206: Level: intermediate
208: .keywords: SNES, default, monitor, residual
210: .seealso: SNESMonitorSet(), SNESMonitorDefault(), SNESMonitorLGCreate()
211: @*/
212: PetscErrorCode SNESMonitorRange(SNES snes,PetscInt it,PetscReal rnorm,void *dummy)
213: {
214: PetscErrorCode ierr;
215: PetscReal perc,rel;
216: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;
217: /* should be in a MonitorRangeContext */
218: static PetscReal prev;
221: if (!it) prev = rnorm;
222: if (!dummy) {PetscViewerASCIIMonitorCreate(((PetscObject)snes)->comm,"stdout",0,&viewer);}
223: SNESMonitorRange_Private(snes,it,&perc);
225: rel = (prev - rnorm)/prev;
226: prev = rnorm;
227: PetscViewerASCIIMonitorPrintf(viewer,"%3D SNES preconditioned resid norm %14.12e Percent values above 20 percent of maximum %5.2f relative decrease %5.2e ratio %5.2e \n",it,rnorm,100.0*perc,rel,rel/perc);
228: if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
229: return(0);
230: }
232: typedef struct {
233: PetscViewerASCIIMonitor viewer;
234: PetscReal *history;
235: } SNESMonitorRatioContext;
239: /*@C
240: SNESMonitorRatio - Monitors progress of the SNES solvers by printing the ratio
241: of residual norm at each iteration to the previous.
243: Collective on SNES
245: Input Parameters:
246: + snes - the SNES context
247: . its - iteration number
248: . fgnorm - 2-norm of residual (or gradient)
249: - dummy - context of monitor
251: Level: intermediate
253: .keywords: SNES, nonlinear, monitor, norm
255: .seealso: SNESMonitorSet(), SNESMonitorSolution()
256: @*/
257: PetscErrorCode SNESMonitorRatio(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
258: {
259: PetscErrorCode ierr;
260: PetscInt len;
261: PetscReal *history;
262: SNESMonitorRatioContext *ctx = (SNESMonitorRatioContext*)dummy;
265: SNESGetConvergenceHistory(snes,&history,PETSC_NULL,&len);
266: if (!its || !history || its > len) {
267: PetscViewerASCIIMonitorPrintf(ctx->viewer,"%3D SNES Function norm %14.12e \n",its,fgnorm);
268: } else {
269: PetscReal ratio = fgnorm/history[its-1];
270: PetscViewerASCIIMonitorPrintf(ctx->viewer,"%3D SNES Function norm %14.12e %G \n",its,fgnorm,ratio);
271: }
272: return(0);
273: }
275: /*
276: If the we set the history monitor space then we must destroy it
277: */
280: PetscErrorCode SNESMonitorRatioDestroy(void *ct)
281: {
282: PetscErrorCode ierr;
283: SNESMonitorRatioContext *ctx = (SNESMonitorRatioContext*)ct;
286: PetscFree(ctx->history);
287: PetscViewerASCIIMonitorDestroy(ctx->viewer);
288: PetscFree(ctx);
289: return(0);
290: }
294: /*@C
295: SNESMonitorSetRatio - Sets SNES to use a monitor that prints the
296: ratio of the function norm at each iteration.
298: Collective on SNES
300: Input Parameters:
301: + snes - the SNES context
302: - viewer - ASCII viewer to print output
304: Level: intermediate
306: .keywords: SNES, nonlinear, monitor, norm
308: .seealso: SNESMonitorSet(), SNESMonitorSolution(), SNESMonitorDefault()
309: @*/
310: PetscErrorCode SNESMonitorSetRatio(SNES snes,PetscViewerASCIIMonitor viewer)
311: {
312: PetscErrorCode ierr;
313: SNESMonitorRatioContext *ctx;
314: PetscReal *history;
317: if (!viewer) {
318: PetscViewerASCIIMonitorCreate(((PetscObject)snes)->comm,"stdout",0,&viewer);
319: PetscObjectReference((PetscObject)viewer);
320: }
321: PetscNewLog(snes,SNESMonitorRatioContext,&ctx);
322: SNESGetConvergenceHistory(snes,&history,PETSC_NULL,PETSC_NULL);
323: if (!history) {
324: PetscMalloc(100*sizeof(PetscReal),&ctx->history);
325: SNESSetConvergenceHistory(snes,ctx->history,0,100,PETSC_TRUE);
326: }
327: ctx->viewer = viewer;
328: SNESMonitorSet(snes,SNESMonitorRatio,ctx,SNESMonitorRatioDestroy);
329: return(0);
330: }
332: /* ---------------------------------------------------------------- */
335: /*
336: Default (short) SNES Monitor, same as SNESMonitorDefault() except
337: it prints fewer digits of the residual as the residual gets smaller.
338: This is because the later digits are meaningless and are often
339: different on different machines; by using this routine different
340: machines will usually generate the same output.
341: */
342: PetscErrorCode SNESMonitorDefaultShort(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
343: {
344: PetscErrorCode ierr;
345: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;
348: if (!dummy) {
349: PetscViewerASCIIMonitorCreate(((PetscObject)snes)->comm,"stdout",0,&viewer);
350: }
351: if (fgnorm > 1.e-9) {
352: PetscViewerASCIIMonitorPrintf(viewer,"%3D SNES Function norm %G \n",its,fgnorm);
353: } else if (fgnorm > 1.e-11){
354: PetscViewerASCIIMonitorPrintf(viewer,"%3D SNES Function norm %5.3e \n",its,fgnorm);
355: } else {
356: PetscViewerASCIIMonitorPrintf(viewer,"%3D SNES Function norm < 1.e-11\n",its);
357: }
358: if (!dummy) {
359: PetscViewerASCIIMonitorDestroy(viewer);
360: }
361: return(0);
362: }
363: /* ---------------------------------------------------------------- */
366: /*@C
367: SNESDefaultConverged - Convergence test of the solvers for
368: systems of nonlinear equations (default).
370: Collective on SNES
372: Input Parameters:
373: + snes - the SNES context
374: . it - the iteration (0 indicates before any Newton steps)
375: . xnorm - 2-norm of current iterate
376: . pnorm - 2-norm of current step
377: . fnorm - 2-norm of function at current iterate
378: - dummy - unused context
380: Output Parameter:
381: . reason - one of
382: $ SNES_CONVERGED_FNORM_ABS - (fnorm < abstol),
383: $ SNES_CONVERGED_PNORM_RELATIVE - (pnorm < xtol*xnorm),
384: $ SNES_CONVERGED_FNORM_RELATIVE - (fnorm < rtol*fnorm0),
385: $ SNES_DIVERGED_FUNCTION_COUNT - (nfct > maxf),
386: $ SNES_DIVERGED_FNORM_NAN - (fnorm == NaN),
387: $ SNES_CONVERGED_ITERATING - (otherwise),
389: where
390: + maxf - maximum number of function evaluations,
391: set with SNESSetTolerances()
392: . nfct - number of function evaluations,
393: . abstol - absolute function norm tolerance,
394: set with SNESSetTolerances()
395: - rtol - relative function norm tolerance, set with SNESSetTolerances()
397: Level: intermediate
399: .keywords: SNES, nonlinear, default, converged, convergence
401: .seealso: SNESSetConvergenceTest()
402: @*/
403: PetscErrorCode SNESDefaultConverged(SNES snes,PetscInt it,PetscReal xnorm,PetscReal pnorm,PetscReal fnorm,SNESConvergedReason *reason,void *dummy)
404: {
410:
411: *reason = SNES_CONVERGED_ITERATING;
413: if (!it) {
414: /* set parameter for default relative tolerance convergence test */
415: snes->ttol = fnorm*snes->rtol;
416: }
417: if (fnorm != fnorm) {
418: PetscInfo(snes,"Failed to converged, function norm is NaN\n");
419: *reason = SNES_DIVERGED_FNORM_NAN;
420: } else if (fnorm < snes->abstol) {
421: PetscInfo2(snes,"Converged due to function norm %G < %G\n",fnorm,snes->abstol);
422: *reason = SNES_CONVERGED_FNORM_ABS;
423: } else if (snes->nfuncs >= snes->max_funcs) {
424: PetscInfo2(snes,"Exceeded maximum number of function evaluations: %D > %D\n",snes->nfuncs,snes->max_funcs);
425: *reason = SNES_DIVERGED_FUNCTION_COUNT;
426: }
428: if (it && !*reason) {
429: if (fnorm <= snes->ttol) {
430: PetscInfo2(snes,"Converged due to function norm %G < %G (relative tolerance)\n",fnorm,snes->ttol);
431: *reason = SNES_CONVERGED_FNORM_RELATIVE;
432: } else if (pnorm < snes->xtol*xnorm) {
433: PetscInfo3(snes,"Converged due to small update length: %G < %G * %G\n",pnorm,snes->xtol,xnorm);
434: *reason = SNES_CONVERGED_PNORM_RELATIVE;
435: }
436: }
437: return(0);
438: }
442: /*@C
443: SNESSkipConverged - Convergence test for SNES that NEVER returns as
444: converged, UNLESS the maximum number of iteration have been reached.
446: Collective on SNES
448: Input Parameters:
449: + snes - the SNES context
450: . it - the iteration (0 indicates before any Newton steps)
451: . xnorm - 2-norm of current iterate
452: . pnorm - 2-norm of current step
453: . fnorm - 2-norm of function at current iterate
454: - dummy - unused context
456: Output Parameter:
457: . reason - SNES_CONVERGED_ITERATING, SNES_CONVERGED_ITS, or SNES_DIVERGED_FNORM_NAN
459: Notes:
460: Convergence is then declared after a fixed number of iterations have been used.
461:
462: Level: advanced
464: .keywords: SNES, nonlinear, skip, converged, convergence
466: .seealso: SNESSetConvergenceTest()
467: @*/
468: PetscErrorCode SNESSkipConverged(SNES snes,PetscInt it,PetscReal xnorm,PetscReal pnorm,PetscReal fnorm,SNESConvergedReason *reason,void *dummy)
469: {
476: *reason = SNES_CONVERGED_ITERATING;
478: if (fnorm != fnorm) {
479: PetscInfo(snes,"Failed to converged, function norm is NaN\n");
480: *reason = SNES_DIVERGED_FNORM_NAN;
481: } else if(it == snes->max_its) {
482: *reason = SNES_CONVERGED_ITS;
483: }
484: return(0);
485: }