Actual source code: snesut.c


  2: #include <petsc/private/snesimpl.h>
  3: #include <petscdm.h>
  4: #include <petscsection.h>
  5: #include <petscblaslapack.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 -  a viewer

 19:    Options Database Keys:
 20: .  -snes_monitor_solution [ascii binary draw][:filename][:viewer format] - plots solution at each iteration

 22:    Level: intermediate

 24: .seealso: SNESMonitorSet(), SNESMonitorDefault(), VecView()
 25: @*/
 26: PetscErrorCode  SNESMonitorSolution(SNES snes,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *vf)
 27: {
 28:   Vec            x;
 29:   PetscViewer    viewer = vf->viewer;

 32:   SNESGetSolution(snes,&x);
 33:   PetscViewerPushFormat(viewer,vf->format);
 34:   VecView(x,viewer);
 35:   PetscViewerPopFormat(viewer);
 36:   return 0;
 37: }

 39: /*@C
 40:    SNESMonitorResidual - Monitors progress of the SNES solvers by calling
 41:    VecView() for the residual at each iteration.

 43:    Collective on SNES

 45:    Input Parameters:
 46: +  snes - the SNES context
 47: .  its - iteration number
 48: .  fgnorm - 2-norm of residual
 49: -  dummy -  a viewer

 51:    Options Database Keys:
 52: .  -snes_monitor_residual [ascii binary draw][:filename][:viewer format] - plots residual (not its norm) at each iteration

 54:    Level: intermediate

 56: .seealso: SNESMonitorSet(), SNESMonitorDefault(), VecView()
 57: @*/
 58: PetscErrorCode  SNESMonitorResidual(SNES snes,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *vf)
 59: {
 60:   Vec            x;
 61:   PetscViewer    viewer = vf->viewer;

 64:   SNESGetFunction(snes,&x,NULL,NULL);
 65:   PetscViewerPushFormat(viewer,vf->format);
 66:   VecView(x,viewer);
 67:   PetscViewerPopFormat(viewer);
 68:   return 0;
 69: }

 71: /*@C
 72:    SNESMonitorSolutionUpdate - Monitors progress of the SNES solvers by calling
 73:    VecView() for the UPDATE to the solution at each iteration.

 75:    Collective on SNES

 77:    Input Parameters:
 78: +  snes - the SNES context
 79: .  its - iteration number
 80: .  fgnorm - 2-norm of residual
 81: -  dummy - a viewer

 83:    Options Database Keys:
 84: .  -snes_monitor_solution_update [ascii binary draw][:filename][:viewer format] - plots update to solution at each iteration

 86:    Level: intermediate

 88: .seealso: SNESMonitorSet(), SNESMonitorDefault(), VecView()
 89: @*/
 90: PetscErrorCode  SNESMonitorSolutionUpdate(SNES snes,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *vf)
 91: {
 92:   Vec            x;
 93:   PetscViewer    viewer = vf->viewer;

 96:   SNESGetSolutionUpdate(snes,&x);
 97:   PetscViewerPushFormat(viewer,vf->format);
 98:   VecView(x,viewer);
 99:   PetscViewerPopFormat(viewer);
100:   return 0;
101: }

103: #include <petscdraw.h>

105: /*@C
106:   KSPMonitorSNESResidual - Prints the SNES residual norm, as well as the linear residual norm, at each iteration of an iterative solver.

108:   Collective on ksp

110:   Input Parameters:
111: + ksp   - iterative context
112: . n     - iteration number
113: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
114: - vf    - The viewer context

116:   Options Database Key:
117: . -snes_monitor_ksp - Activates KSPMonitorSNESResidual()

119:   Level: intermediate

121: .seealso: KSPMonitorSet(), KSPMonitorResidual(),KSPMonitorTrueResidualMaxNorm()
122: @*/
123: PetscErrorCode KSPMonitorSNESResidual(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
124: {
125:   PetscViewer       viewer = vf->viewer;
126:   PetscViewerFormat format = vf->format;
127:   SNES              snes   = (SNES) vf->data;
128:   Vec               snes_solution, work1, work2;
129:   PetscReal         snorm;
130:   PetscInt          tablevel;
131:   const char       *prefix;

134:   SNESGetSolution(snes, &snes_solution);
135:   VecDuplicate(snes_solution, &work1);
136:   VecDuplicate(snes_solution, &work2);
137:   KSPBuildSolution(ksp, work1, NULL);
138:   VecAYPX(work1, -1.0, snes_solution);
139:   SNESComputeFunction(snes, work1, work2);
140:   VecNorm(work2, NORM_2, &snorm);
141:   VecDestroy(&work1);
142:   VecDestroy(&work2);

144:   PetscObjectGetTabLevel((PetscObject) ksp, &tablevel);
145:   PetscObjectGetOptionsPrefix((PetscObject) ksp, &prefix);
146:   PetscViewerPushFormat(viewer, format);
147:   PetscViewerASCIIAddTab(viewer, tablevel);
148:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
149:   PetscViewerASCIIPrintf(viewer, "%3D SNES Residual norm %5.3e KSP Residual norm %5.3e \n", n, (double) snorm, (double) rnorm);
150:   PetscViewerASCIISubtractTab(viewer, tablevel);
151:   PetscViewerPopFormat(viewer);
152:   return 0;
153: }

155: /*@C
156:   KSPMonitorSNESResidualDrawLG - Plots the linear SNES residual norm at each iteration of an iterative solver.

158:   Collective on ksp

160:   Input Parameters:
161: + ksp   - iterative context
162: . n     - iteration number
163: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
164: - vf    - The viewer context

166:   Options Database Key:
167: . -snes_monitor_ksp draw::draw_lg - Activates KSPMonitorSNESResidualDrawLG()

169:   Level: intermediate

171: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
172: @*/
173: PetscErrorCode KSPMonitorSNESResidualDrawLG(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
174: {
175:   PetscViewer        viewer = vf->viewer;
176:   PetscViewerFormat  format = vf->format;
177:   PetscDrawLG        lg     = vf->lg;
178:   SNES               snes   = (SNES) vf->data;
179:   Vec                snes_solution, work1, work2;
180:   PetscReal          snorm;
181:   KSPConvergedReason reason;
182:   PetscReal          x[2], y[2];

186:   SNESGetSolution(snes, &snes_solution);
187:   VecDuplicate(snes_solution, &work1);
188:   VecDuplicate(snes_solution, &work2);
189:   KSPBuildSolution(ksp, work1, NULL);
190:   VecAYPX(work1, -1.0, snes_solution);
191:   SNESComputeFunction(snes, work1, work2);
192:   VecNorm(work2, NORM_2, &snorm);
193:   VecDestroy(&work1);
194:   VecDestroy(&work2);

196:   PetscViewerPushFormat(viewer, format);
197:   if (!n) PetscDrawLGReset(lg);
198:   x[0] = (PetscReal) n;
199:   if (rnorm > 0.0) y[0] = PetscLog10Real(rnorm);
200:   else y[0] = -15.0;
201:   x[1] = (PetscReal) n;
202:   if (snorm > 0.0) y[1] = PetscLog10Real(snorm);
203:   else y[1] = -15.0;
204:   PetscDrawLGAddPoint(lg, x, y);
205:   KSPGetConvergedReason(ksp, &reason);
206:   if (n <= 20 || !(n % 5) || reason) {
207:     PetscDrawLGDraw(lg);
208:     PetscDrawLGSave(lg);
209:   }
210:   PetscViewerPopFormat(viewer);
211:   return 0;
212: }

214: /*@C
215:   KSPMonitorSNESResidualDrawLGCreate - Creates the plotter for the linear SNES residual.

217:   Collective on ksp

219:   Input Parameters:
220: + viewer - The PetscViewer
221: . format - The viewer format
222: - ctx    - An optional user context

224:   Output Parameter:
225: . vf    - The viewer context

227:   Level: intermediate

229: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
230: @*/
231: PetscErrorCode KSPMonitorSNESResidualDrawLGCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
232: {
233:   const char    *names[] = {"linear", "nonlinear"};

235:   PetscViewerAndFormatCreate(viewer, format, vf);
236:   (*vf)->data = ctx;
237:   KSPMonitorLGCreate(PetscObjectComm((PetscObject) viewer), NULL, NULL, "Log Residual Norm", 2, names, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &(*vf)->lg);
238:   return 0;
239: }

241: PetscErrorCode SNESMonitorDefaultSetUp(SNES snes, PetscViewerAndFormat *vf)
242: {
243:   if (vf->format == PETSC_VIEWER_DRAW_LG) {
244:     KSPMonitorLGCreate(PetscObjectComm((PetscObject) vf->viewer), NULL, NULL, "Log Residual Norm", 1, NULL, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &vf->lg);
245:   }
246:   return 0;
247: }

249: /*@C
250:    SNESMonitorDefault - Monitors progress of the SNES solvers (default).

252:    Collective on SNES

254:    Input Parameters:
255: +  snes - the SNES context
256: .  its - iteration number
257: .  fgnorm - 2-norm of residual
258: -  vf - viewer and format structure

260:    Notes:
261:    This routine prints the residual norm at each iteration.

263:    Level: intermediate

265: .seealso: SNESMonitorSet(), SNESMonitorSolution()
266: @*/
267: PetscErrorCode  SNESMonitorDefault(SNES snes,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *vf)
268: {
269:   PetscViewer       viewer = vf->viewer;
270:   PetscViewerFormat format = vf->format;
271:   PetscBool         isascii, isdraw;

274:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &isascii);
275:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERDRAW,  &isdraw);
276:   PetscViewerPushFormat(viewer,format);
277:   if (isascii) {
278:     PetscViewerASCIIAddTab(viewer,((PetscObject)snes)->tablevel);
279:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %14.12e \n",its,(double)fgnorm);
280:     PetscViewerASCIISubtractTab(viewer,((PetscObject)snes)->tablevel);
281:   } else if (isdraw) {
282:     if (format == PETSC_VIEWER_DRAW_LG) {
283:       PetscDrawLG lg = (PetscDrawLG) vf->lg;
284:       PetscReal   x, y;

287:       if (!its) PetscDrawLGReset(lg);
288:       x = (PetscReal) its;
289:       if (fgnorm > 0.0) y = PetscLog10Real(fgnorm);
290:       else y = -15.0;
291:       PetscDrawLGAddPoint(lg,&x,&y);
292:       if (its <= 20 || !(its % 5) || snes->reason) {
293:         PetscDrawLGDraw(lg);
294:         PetscDrawLGSave(lg);
295:       }
296:     }
297:   }
298:   PetscViewerPopFormat(viewer);
299:   return 0;
300: }

302: /*@C
303:    SNESMonitorScaling - Monitors the largest value in each row of the Jacobian.

305:    Collective on SNES

307:    Input Parameters:
308: +  snes - the SNES context
309: .  its - iteration number
310: .  fgnorm - 2-norm of residual
311: -  vf - viewer and format structure

313:    Notes:
314:    This routine prints the largest value in each row of the Jacobian

316:    Level: intermediate

318: .seealso: SNESMonitorSet(), SNESMonitorSolution()
319: @*/
320: PetscErrorCode  SNESMonitorScaling(SNES snes,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *vf)
321: {
322:   PetscViewer    viewer = vf->viewer;
323:   KSP            ksp;
324:   Mat            J;
325:   Vec            v;

328:   SNESGetKSP(snes,&ksp);
329:   KSPGetOperators(ksp,&J,NULL);
330:   MatCreateVecs(J,&v,NULL);
331:   MatGetRowMaxAbs(J,v,NULL);
332:   PetscViewerPushFormat(viewer,vf->format);
333:   PetscViewerASCIIAddTab(viewer,((PetscObject)snes)->tablevel);
334:   PetscViewerASCIIPrintf(viewer,"%3D SNES Jacobian maximum row entries \n");
335:   VecView(v,viewer);
336:   PetscViewerASCIISubtractTab(viewer,((PetscObject)snes)->tablevel);
337:   PetscViewerPopFormat(viewer);
338:   VecDestroy(&v);
339:   return 0;
340: }

342: PetscErrorCode SNESMonitorJacUpdateSpectrum(SNES snes,PetscInt it,PetscReal fnorm,PetscViewerAndFormat *vf)
343: {
344:   Vec            X;
345:   Mat            J,dJ,dJdense;
346:   PetscErrorCode (*func)(SNES,Vec,Mat,Mat,void*);
347:   PetscInt       n;
348:   PetscBLASInt   nb = 0,lwork;
349:   PetscReal      *eigr,*eigi;
350:   PetscScalar    *work;
351:   PetscScalar    *a;

353:   if (it == 0) return 0;
354:   /* create the difference between the current update and the current jacobian */
355:   SNESGetSolution(snes,&X);
356:   SNESGetJacobian(snes,NULL,&J,&func,NULL);
357:   MatDuplicate(J,MAT_COPY_VALUES,&dJ);
358:   SNESComputeJacobian(snes,X,dJ,dJ);
359:   MatAXPY(dJ,-1.0,J,SAME_NONZERO_PATTERN);

361:   /* compute the spectrum directly */
362:   MatConvert(dJ,MATSEQDENSE,MAT_INITIAL_MATRIX,&dJdense);
363:   MatGetSize(dJ,&n,NULL);
364:   PetscBLASIntCast(n,&nb);
365:   lwork = 3*nb;
366:   PetscMalloc1(n,&eigr);
367:   PetscMalloc1(n,&eigi);
368:   PetscMalloc1(lwork,&work);
369:   MatDenseGetArray(dJdense,&a);
370: #if !defined(PETSC_USE_COMPLEX)
371:   {
372:     PetscBLASInt lierr;
373:     PetscInt     i;
374:     PetscFPTrapPush(PETSC_FP_TRAP_OFF);
375:     PetscStackCallBLAS("LAPACKgeev",LAPACKgeev_("N","N",&nb,a,&nb,eigr,eigi,NULL,&nb,NULL,&nb,work,&lwork,&lierr));
377:     PetscFPTrapPop();
378:     PetscPrintf(PetscObjectComm((PetscObject)snes),"Eigenvalues of J_%d - J_%d:\n",it,it-1);
379:     for (i=0;i<n;i++) {
380:       PetscPrintf(PetscObjectComm((PetscObject)snes),"%5d: %20.5g + %20.5gi\n",i,(double)eigr[i],(double)eigi[i]);
381:     }
382:   }
383:   MatDenseRestoreArray(dJdense,&a);
384:   MatDestroy(&dJ);
385:   MatDestroy(&dJdense);
386:   PetscFree(eigr);
387:   PetscFree(eigi);
388:   PetscFree(work);
389:   return 0;
390: #else
391:   SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Not coded for complex");
392: #endif
393: }

395: PETSC_INTERN PetscErrorCode  SNESMonitorRange_Private(SNES,PetscInt,PetscReal*);

397: PetscErrorCode  SNESMonitorRange_Private(SNES snes,PetscInt it,PetscReal *per)
398: {
399:   Vec            resid;
400:   PetscReal      rmax,pwork;
401:   PetscInt       i,n,N;
402:   PetscScalar    *r;

404:   SNESGetFunction(snes,&resid,NULL,NULL);
405:   VecNorm(resid,NORM_INFINITY,&rmax);
406:   VecGetLocalSize(resid,&n);
407:   VecGetSize(resid,&N);
408:   VecGetArray(resid,&r);
409:   pwork = 0.0;
410:   for (i=0; i<n; i++) {
411:     pwork += (PetscAbsScalar(r[i]) > .20*rmax);
412:   }
413:   MPIU_Allreduce(&pwork,per,1,MPIU_REAL,MPIU_SUM,PetscObjectComm((PetscObject)snes));
414:   VecRestoreArray(resid,&r);
415:   *per = *per/N;
416:   return 0;
417: }

419: /*@C
420:    SNESMonitorRange - Prints the percentage of residual elements that are more then 10 percent of the maximum value.

422:    Collective on SNES

424:    Input Parameters:
425: +  snes   - iterative context
426: .  it    - iteration number
427: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).
428: -  dummy - unused monitor context

430:    Options Database Key:
431: .  -snes_monitor_range - Activates SNESMonitorRange()

433:    Level: intermediate

435: .seealso: SNESMonitorSet(), SNESMonitorDefault(), SNESMonitorLGCreate()
436: @*/
437: PetscErrorCode  SNESMonitorRange(SNES snes,PetscInt it,PetscReal rnorm,PetscViewerAndFormat *vf)
438: {
439:   PetscReal      perc,rel;
440:   PetscViewer    viewer = vf->viewer;
441:   /* should be in a MonitorRangeContext */
442:   static PetscReal prev;

445:   if (!it) prev = rnorm;
446:   SNESMonitorRange_Private(snes,it,&perc);

448:   rel  = (prev - rnorm)/prev;
449:   prev = rnorm;
450:   PetscViewerPushFormat(viewer,vf->format);
451:   PetscViewerASCIIAddTab(viewer,((PetscObject)snes)->tablevel);
452:   PetscViewerASCIIPrintf(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,(double)rnorm,(double)(100.0*perc),(double)rel,(double)(rel/perc));
453:   PetscViewerASCIISubtractTab(viewer,((PetscObject)snes)->tablevel);
454:   PetscViewerPopFormat(viewer);
455:   return 0;
456: }

458: /*@C
459:    SNESMonitorRatio - Monitors progress of the SNES solvers by printing the ratio
460:    of residual norm at each iteration to the previous.

462:    Collective on SNES

464:    Input Parameters:
465: +  snes - the SNES context
466: .  its - iteration number
467: .  fgnorm - 2-norm of residual (or gradient)
468: -  dummy -  context of monitor

470:    Level: intermediate

472:    Notes:
473:     Insure that SNESMonitorRatio() is called when you set this monitor
474: .seealso: SNESMonitorSet(), SNESMonitorSolution(), SNESMonitorRatio()
475: @*/
476: PetscErrorCode  SNESMonitorRatio(SNES snes,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *vf)
477: {
478:   PetscInt                len;
479:   PetscReal               *history;
480:   PetscViewer             viewer = vf->viewer;

482:   SNESGetConvergenceHistory(snes,&history,NULL,&len);
483:   PetscViewerPushFormat(viewer,vf->format);
484:   PetscViewerASCIIAddTab(viewer,((PetscObject)snes)->tablevel);
485:   if (!its || !history || its > len) {
486:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %14.12e \n",its,(double)fgnorm);
487:   } else {
488:     PetscReal ratio = fgnorm/history[its-1];
489:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %14.12e %14.12e \n",its,(double)fgnorm,(double)ratio);
490:   }
491:   PetscViewerASCIISubtractTab(viewer,((PetscObject)snes)->tablevel);
492:   PetscViewerPopFormat(viewer);
493:   return 0;
494: }

496: /*@C
497:    SNESMonitorRatioSetUp - Insures the SNES object is saving its history since this monitor needs access to it

499:    Collective on SNES

501:    Input Parameters:
502: +   snes - the SNES context
503: -   viewer - the PetscViewer object (ignored)

505:    Level: intermediate

507: .seealso: SNESMonitorSet(), SNESMonitorSolution(), SNESMonitorDefault(), SNESMonitorRatio()
508: @*/
509: PetscErrorCode  SNESMonitorRatioSetUp(SNES snes,PetscViewerAndFormat *vf)
510: {
511:   PetscReal               *history;

513:   SNESGetConvergenceHistory(snes,&history,NULL,NULL);
514:   if (!history) {
515:     SNESSetConvergenceHistory(snes,NULL,NULL,100,PETSC_TRUE);
516:   }
517:   return 0;
518: }

520: /* ---------------------------------------------------------------- */
521: /*
522:      Default (short) SNES Monitor, same as SNESMonitorDefault() except
523:   it prints fewer digits of the residual as the residual gets smaller.
524:   This is because the later digits are meaningless and are often
525:   different on different machines; by using this routine different
526:   machines will usually generate the same output.

528:   Deprecated: Intentionally has no manual page
529: */
530: PetscErrorCode  SNESMonitorDefaultShort(SNES snes,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *vf)
531: {
532:   PetscViewer    viewer = vf->viewer;

535:   PetscViewerPushFormat(viewer,vf->format);
536:   PetscViewerASCIIAddTab(viewer,((PetscObject)snes)->tablevel);
537:   if (fgnorm > 1.e-9) {
538:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %g \n",its,(double)fgnorm);
539:   } else if (fgnorm > 1.e-11) {
540:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %5.3e \n",its,(double)fgnorm);
541:   } else {
542:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm < 1.e-11\n",its);
543:   }
544:   PetscViewerASCIISubtractTab(viewer,((PetscObject)snes)->tablevel);
545:   PetscViewerPopFormat(viewer);
546:   return 0;
547: }

549: /*@C
550:   SNESMonitorDefaultField - Monitors progress of the SNES solvers, separated into fields.

552:   Collective on SNES

554:   Input Parameters:
555: + snes   - the SNES context
556: . its    - iteration number
557: . fgnorm - 2-norm of residual
558: - ctx    - the PetscViewer

560:   Notes:
561:   This routine uses the DM attached to the residual vector

563:   Level: intermediate

565: .seealso: SNESMonitorSet(), SNESMonitorSolution(), SNESMonitorDefault()
566: @*/
567: PetscErrorCode SNESMonitorDefaultField(SNES snes, PetscInt its, PetscReal fgnorm, PetscViewerAndFormat *vf)
568: {
569:   PetscViewer    viewer = vf->viewer;
570:   Vec            r;
571:   DM             dm;
572:   PetscReal      res[256];
573:   PetscInt       tablevel;

576:   SNESGetFunction(snes, &r, NULL, NULL);
577:   VecGetDM(r, &dm);
578:   if (!dm) SNESMonitorDefault(snes, its, fgnorm, vf);
579:   else {
580:     PetscSection s, gs;
581:     PetscInt     Nf, f;

583:     DMGetLocalSection(dm, &s);
584:     DMGetGlobalSection(dm, &gs);
585:     if (!s || !gs) SNESMonitorDefault(snes, its, fgnorm, vf);
586:     PetscSectionGetNumFields(s, &Nf);
588:     PetscSectionVecNorm(s, gs, r, NORM_2, res);
589:     PetscObjectGetTabLevel((PetscObject) snes, &tablevel);
590:     PetscViewerPushFormat(viewer,vf->format);
591:     PetscViewerASCIIAddTab(viewer, tablevel);
592:     PetscViewerASCIIPrintf(viewer, "%3D SNES Function norm %14.12e [", its, (double) fgnorm);
593:     for (f = 0; f < Nf; ++f) {
594:       if (f) PetscViewerASCIIPrintf(viewer, ", ");
595:       PetscViewerASCIIPrintf(viewer, "%14.12e", res[f]);
596:     }
597:     PetscViewerASCIIPrintf(viewer, "] \n");
598:     PetscViewerASCIISubtractTab(viewer, tablevel);
599:     PetscViewerPopFormat(viewer);
600:   }
601:   return 0;
602: }
603: /* ---------------------------------------------------------------- */
604: /*@C
605:    SNESConvergedDefault - Default onvergence test of the solvers for
606:    systems of nonlinear equations.

608:    Collective on SNES

610:    Input Parameters:
611: +  snes - the SNES context
612: .  it - the iteration (0 indicates before any Newton steps)
613: .  xnorm - 2-norm of current iterate
614: .  snorm - 2-norm of current step
615: .  fnorm - 2-norm of function at current iterate
616: -  dummy - unused context

618:    Output Parameter:
619: .   reason  - one of
620: $  SNES_CONVERGED_FNORM_ABS       - (fnorm < abstol),
621: $  SNES_CONVERGED_SNORM_RELATIVE  - (snorm < stol*xnorm),
622: $  SNES_CONVERGED_FNORM_RELATIVE  - (fnorm < rtol*fnorm0),
623: $  SNES_DIVERGED_FUNCTION_COUNT   - (nfct > maxf),
624: $  SNES_DIVERGED_FNORM_NAN        - (fnorm == NaN),
625: $  SNES_CONVERGED_ITERATING       - (otherwise),
626: $  SNES_DIVERGED_DTOL             - (fnorm > divtol*snes->fnorm0)

628:    where
629: +    maxf - maximum number of function evaluations,  set with SNESSetTolerances()
630: .    nfct - number of function evaluations,
631: .    abstol - absolute function norm tolerance, set with SNESSetTolerances()
632: .    rtol - relative function norm tolerance, set with SNESSetTolerances()
633: .    divtol - divergence tolerance, set with SNESSetDivergenceTolerance()
634: -    fnorm0 - 2-norm of the function at the initial solution (initial guess; zeroth iteration)

636:   Options Database Keys:
637: +  -snes_convergence_test default - see SNESSetFromOptions()
638: .  -snes_stol - convergence tolerance in terms of the norm  of the change in the solution between steps
639: .  -snes_atol <abstol> - absolute tolerance of residual norm
640: .  -snes_rtol <rtol> - relative decrease in tolerance norm from the initial 2-norm of the solution
641: .  -snes_divergence_tolerance <divtol> - if the residual goes above divtol*rnorm0, exit with divergence
642: .  -snes_max_funcs <max_funcs> - maximum number of function evaluations
643: .  -snes_max_fail <max_fail> - maximum number of line search failures allowed before stopping, default is none
644: -  -snes_max_linear_solve_fail - number of linear solver failures before SNESSolve() stops

646:    Level: intermediate

648: .seealso: SNESSetConvergenceTest(), SNESConvergedSkip(), SNESSetTolerances(), SNESSetDivergenceTolerance()
649: @*/
650: PetscErrorCode  SNESConvergedDefault(SNES snes,PetscInt it,PetscReal xnorm,PetscReal snorm,PetscReal fnorm,SNESConvergedReason *reason,void *dummy)
651: {

655:   *reason = SNES_CONVERGED_ITERATING;

657:   if (!it) {
658:     /* set parameter for default relative tolerance convergence test */
659:     snes->ttol   = fnorm*snes->rtol;
660:     snes->rnorm0 = fnorm;
661:   }
662:   if (PetscIsInfOrNanReal(fnorm)) {
663:     PetscInfo(snes,"Failed to converged, function norm is NaN\n");
664:     *reason = SNES_DIVERGED_FNORM_NAN;
665:   } else if (fnorm < snes->abstol && (it || !snes->forceiteration)) {
666:     PetscInfo(snes,"Converged due to function norm %14.12e < %14.12e\n",(double)fnorm,(double)snes->abstol);
667:     *reason = SNES_CONVERGED_FNORM_ABS;
668:   } else if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) {
669:     PetscInfo(snes,"Exceeded maximum number of function evaluations: %D > %D\n",snes->nfuncs,snes->max_funcs);
670:     *reason = SNES_DIVERGED_FUNCTION_COUNT;
671:   }

673:   if (it && !*reason) {
674:     if (fnorm <= snes->ttol) {
675:       PetscInfo(snes,"Converged due to function norm %14.12e < %14.12e (relative tolerance)\n",(double)fnorm,(double)snes->ttol);
676:       *reason = SNES_CONVERGED_FNORM_RELATIVE;
677:     } else if (snorm < snes->stol*xnorm) {
678:       PetscInfo(snes,"Converged due to small update length: %14.12e < %14.12e * %14.12e\n",(double)snorm,(double)snes->stol,(double)xnorm);
679:       *reason = SNES_CONVERGED_SNORM_RELATIVE;
680:     } else if (snes->divtol > 0 && (fnorm > snes->divtol*snes->rnorm0)) {
681:       PetscInfo(snes,"Diverged due to increase in function norm: %14.12e > %14.12e * %14.12e\n",(double)fnorm,(double)snes->divtol,(double)snes->rnorm0);
682:       *reason = SNES_DIVERGED_DTOL;
683:     }

685:   }
686:   return 0;
687: }

689: /*@C
690:    SNESConvergedSkip - Convergence test for SNES that NEVER returns as
691:    converged, UNLESS the maximum number of iteration have been reached.

693:    Logically Collective on SNES

695:    Input Parameters:
696: +  snes - the SNES context
697: .  it - the iteration (0 indicates before any Newton steps)
698: .  xnorm - 2-norm of current iterate
699: .  snorm - 2-norm of current step
700: .  fnorm - 2-norm of function at current iterate
701: -  dummy - unused context

703:    Output Parameter:
704: .   reason  - SNES_CONVERGED_ITERATING, SNES_CONVERGED_ITS, or SNES_DIVERGED_FNORM_NAN

706:    Notes:
707:    Convergence is then declared after a fixed number of iterations have been used.

709:    Options Database Keys:
710: .  -snes_convergence_test default - see SNESSetFromOptions()

712:    Level: advanced

714: .seealso: SNESConvergedDefault(), SNESSetConvergenceTest()
715: @*/
716: PetscErrorCode  SNESConvergedSkip(SNES snes,PetscInt it,PetscReal xnorm,PetscReal snorm,PetscReal fnorm,SNESConvergedReason *reason,void *dummy)
717: {

721:   *reason = SNES_CONVERGED_ITERATING;

723:   if (fnorm != fnorm) {
724:     PetscInfo(snes,"Failed to converged, function norm is NaN\n");
725:     *reason = SNES_DIVERGED_FNORM_NAN;
726:   } else if (it == snes->max_its) {
727:     *reason = SNES_CONVERGED_ITS;
728:   }
729:   return 0;
730: }

732: /*@C
733:   SNESSetWorkVecs - Gets a number of work vectors.

735:   Input Parameters:
736: + snes  - the SNES context
737: - nw - number of work vectors to allocate

739:   Level: developer

741: @*/
742: PetscErrorCode SNESSetWorkVecs(SNES snes,PetscInt nw)
743: {
744:   DM             dm;
745:   Vec            v;

747:   if (snes->work) VecDestroyVecs(snes->nwork,&snes->work);
748:   snes->nwork = nw;

750:   SNESGetDM(snes, &dm);
751:   DMGetGlobalVector(dm, &v);
752:   VecDuplicateVecs(v,snes->nwork,&snes->work);
753:   DMRestoreGlobalVector(dm, &v);
754:   PetscLogObjectParents(snes,nw,snes->work);
755:   return 0;
756: }