Actual source code: snesut.c

petsc-main 2021-04-20
Report Typos and Errors

  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: {
 29:   Vec            x;
 30:   PetscViewer    viewer = vf->viewer;

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

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

 45:    Collective on SNES

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

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

 56:    Level: intermediate

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

 68:   SNESGetFunction(snes,&x,NULL,NULL);
 69:   PetscViewerPushFormat(viewer,vf->format);
 70:   VecView(x,viewer);
 71:   PetscViewerPopFormat(viewer);
 72:   return(0);
 73: }

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

 79:    Collective on SNES

 81:    Input Parameters:
 82: +  snes - the SNES context
 83: .  its - iteration number
 84: .  fgnorm - 2-norm of residual
 85: -  dummy - a viewer

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

 90:    Level: intermediate

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

102:   SNESGetSolutionUpdate(snes,&x);
103:   PetscViewerPushFormat(viewer,vf->format);
104:   VecView(x,viewer);
105:   PetscViewerPopFormat(viewer);
106:   return(0);
107: }

109: #include <petscdraw.h>

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

114:   Collective on ksp

116:   Input Parameters:
117: + ksp   - iterative context
118: . n     - iteration number
119: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
120: - vf    - The viewer context

122:   Options Database Key:
123: . -snes_monitor_ksp - Activates KSPMonitorSNESResidual()

125:   Level: intermediate

127: .seealso: KSPMonitorSet(), KSPMonitorResidual(),KSPMonitorTrueResidualMaxNorm()
128: @*/
129: PetscErrorCode KSPMonitorSNESResidual(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
130: {
131:   PetscViewer       viewer = vf->viewer;
132:   PetscViewerFormat format = vf->format;
133:   SNES              snes   = (SNES) vf->data;
134:   Vec               snes_solution, work1, work2;
135:   PetscReal         snorm;
136:   PetscInt          tablevel;
137:   const char       *prefix;
138:   PetscErrorCode    ierr;

142:   SNESGetSolution(snes, &snes_solution);
143:   VecDuplicate(snes_solution, &work1);
144:   VecDuplicate(snes_solution, &work2);
145:   KSPBuildSolution(ksp, work1, NULL);
146:   VecAYPX(work1, -1.0, snes_solution);
147:   SNESComputeFunction(snes, work1, work2);
148:   VecNorm(work2, NORM_2, &snorm);
149:   VecDestroy(&work1);
150:   VecDestroy(&work2);

152:   PetscObjectGetTabLevel((PetscObject) ksp, &tablevel);
153:   PetscObjectGetOptionsPrefix((PetscObject) ksp, &prefix);
154:   PetscViewerPushFormat(viewer, format);
155:   PetscViewerASCIIAddTab(viewer, tablevel);
156:   if (n == 0 && prefix) {PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);}
157:   PetscViewerASCIIPrintf(viewer, "%3D SNES Residual norm %5.3e KSP Residual norm %5.3e \n", n, (double) snorm, (double) rnorm);
158:   PetscViewerASCIISubtractTab(viewer, tablevel);
159:   PetscViewerPopFormat(viewer);
160:   return(0);
161: }

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

166:   Collective on ksp

168:   Input Parameters:
169: + ksp   - iterative context
170: . n     - iteration number
171: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
172: - vf    - The viewer context

174:   Options Database Key:
175: . -snes_monitor_ksp draw::draw_lg - Activates KSPMonitorSNESResidualDrawLG()

177:   Level: intermediate

179: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
180: @*/
181: PetscErrorCode KSPMonitorSNESResidualDrawLG(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
182: {
183:   PetscViewer        viewer = vf->viewer;
184:   PetscViewerFormat  format = vf->format;
185:   PetscDrawLG        lg     = vf->lg;
186:   SNES               snes   = (SNES) vf->data;
187:   Vec                snes_solution, work1, work2;
188:   PetscReal          snorm;
189:   KSPConvergedReason reason;
190:   PetscReal          x[2], y[2];
191:   PetscErrorCode     ierr;

196:   SNESGetSolution(snes, &snes_solution);
197:   VecDuplicate(snes_solution, &work1);
198:   VecDuplicate(snes_solution, &work2);
199:   KSPBuildSolution(ksp, work1, NULL);
200:   VecAYPX(work1, -1.0, snes_solution);
201:   SNESComputeFunction(snes, work1, work2);
202:   VecNorm(work2, NORM_2, &snorm);
203:   VecDestroy(&work1);
204:   VecDestroy(&work2);

206:   PetscViewerPushFormat(viewer, format);
207:   if (!n) {PetscDrawLGReset(lg);}
208:   x[0] = (PetscReal) n;
209:   if (rnorm > 0.0) y[0] = PetscLog10Real(rnorm);
210:   else y[0] = -15.0;
211:   x[1] = (PetscReal) n;
212:   if (snorm > 0.0) y[1] = PetscLog10Real(snorm);
213:   else y[1] = -15.0;
214:   PetscDrawLGAddPoint(lg, x, y);
215:   KSPGetConvergedReason(ksp, &reason);
216:   if (n <= 20 || !(n % 5) || reason) {
217:     PetscDrawLGDraw(lg);
218:     PetscDrawLGSave(lg);
219:   }
220:   PetscViewerPopFormat(viewer);
221:   return(0);
222: }

224: /*@C
225:   KSPMonitorSNESResidualDrawLGCreate - Creates the plotter for the linear SNES residual.

227:   Collective on ksp

229:   Input Parameters:
230: + viewer - The PetscViewer
231: . format - The viewer format
232: - ctx    - An optional user context

234:   Output Parameter:
235: . vf    - The viewer context

237:   Level: intermediate

239: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
240: @*/
241: PetscErrorCode KSPMonitorSNESResidualDrawLGCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
242: {
243:   const char    *names[] = {"linear", "nonlinear"};

247:   PetscViewerAndFormatCreate(viewer, format, vf);
248:   (*vf)->data = ctx;
249:   KSPMonitorLGCreate(PetscObjectComm((PetscObject) viewer), NULL, NULL, "Log Residual Norm", 2, names, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &(*vf)->lg);
250:   return(0);
251: }

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

256:    Collective on SNES

258:    Input Parameters:
259: +  snes - the SNES context
260: .  its - iteration number
261: .  fgnorm - 2-norm of residual
262: -  vf - viewer and format structure

264:    Notes:
265:    This routine prints the residual norm at each iteration.

267:    Level: intermediate

269: .seealso: SNESMonitorSet(), SNESMonitorSolution()
270: @*/
271: PetscErrorCode  SNESMonitorDefault(SNES snes,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *vf)
272: {
273:   PetscViewer       viewer = vf->viewer;
274:   PetscViewerFormat format = vf->format;
275:   PetscBool         isascii, isdraw;
276:   PetscErrorCode    ierr;

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

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

307:   return(0);
308: }

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

313:    Collective on SNES

315:    Input Parameters:
316: +  snes - the SNES context
317: .  its - iteration number
318: .  fgnorm - 2-norm of residual
319: -  vf - viewer and format structure

321:    Notes:
322:    This routine prints the largest value in each row of the Jacobian

324:    Level: intermediate

326: .seealso: SNESMonitorSet(), SNESMonitorSolution()
327: @*/
328: PetscErrorCode  SNESMonitorScaling(SNES snes,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *vf)
329: {
331:   PetscViewer    viewer = vf->viewer;
332:   KSP            ksp;
333:   Mat            J;
334:   Vec            v;

338:   SNESGetKSP(snes,&ksp);
339:   KSPGetOperators(ksp,&J,NULL);
340:   MatCreateVecs(J,&v,NULL);
341:   MatGetRowMaxAbs(J,v,NULL);
342:   PetscViewerPushFormat(viewer,vf->format);
343:   PetscViewerASCIIAddTab(viewer,((PetscObject)snes)->tablevel);
344:   PetscViewerASCIIPrintf(viewer,"%3D SNES Jacobian maximum row entries \n");
345:   VecView(v,viewer);
346:   PetscViewerASCIISubtractTab(viewer,((PetscObject)snes)->tablevel);
347:   PetscViewerPopFormat(viewer);
348:   VecDestroy(&v);
349:   return(0);
350: }

352: PetscErrorCode SNESMonitorJacUpdateSpectrum(SNES snes,PetscInt it,PetscReal fnorm,PetscViewerAndFormat *vf)
353: {
354: #if defined(PETSC_HAVE_ESSL)
355:   SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_SUP,"GEEV - No support for ESSL Lapack Routines");
356: #else
357:   Vec            X;
358:   Mat            J,dJ,dJdense;
360:   PetscErrorCode (*func)(SNES,Vec,Mat,Mat,void*);
361:   PetscInt       n;
362:   PetscBLASInt   nb = 0,lwork;
363:   PetscReal      *eigr,*eigi;
364:   PetscScalar    *work;
365:   PetscScalar    *a;

368:   if (it == 0) return(0);
369:   /* create the difference between the current update and the current jacobian */
370:   SNESGetSolution(snes,&X);
371:   SNESGetJacobian(snes,NULL,&J,&func,NULL);
372:   MatDuplicate(J,MAT_COPY_VALUES,&dJ);
373:   SNESComputeJacobian(snes,X,dJ,dJ);
374:   MatAXPY(dJ,-1.0,J,SAME_NONZERO_PATTERN);

376:   /* compute the spectrum directly */
377:   MatConvert(dJ,MATSEQDENSE,MAT_INITIAL_MATRIX,&dJdense);
378:   MatGetSize(dJ,&n,NULL);
379:   PetscBLASIntCast(n,&nb);
380:   lwork = 3*nb;
381:   PetscMalloc1(n,&eigr);
382:   PetscMalloc1(n,&eigi);
383:   PetscMalloc1(lwork,&work);
384:   MatDenseGetArray(dJdense,&a);
385: #if !defined(PETSC_USE_COMPLEX)
386:   {
387:     PetscBLASInt lierr;
388:     PetscInt     i;
389:     PetscFPTrapPush(PETSC_FP_TRAP_OFF);
390:     PetscStackCallBLAS("LAPACKgeev",LAPACKgeev_("N","N",&nb,a,&nb,eigr,eigi,NULL,&nb,NULL,&nb,work,&lwork,&lierr));
391:     if (lierr) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"geev() error %d",lierr);
392:     PetscFPTrapPop();
393:     PetscPrintf(PetscObjectComm((PetscObject)snes),"Eigenvalues of J_%d - J_%d:\n",it,it-1);
394:     for (i=0;i<n;i++) {
395:       PetscPrintf(PetscObjectComm((PetscObject)snes),"%5d: %20.5g + %20.5gi\n",i,(double)eigr[i],(double)eigi[i]);
396:     }
397:   }
398:   MatDenseRestoreArray(dJdense,&a);
399:   MatDestroy(&dJ);
400:   MatDestroy(&dJdense);
401:   PetscFree(eigr);
402:   PetscFree(eigi);
403:   PetscFree(work);
404:   return(0);
405: #else
406:   SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Not coded for complex");
407: #endif
408: #endif
409: }

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

413: PetscErrorCode  SNESMonitorRange_Private(SNES snes,PetscInt it,PetscReal *per)
414: {
416:   Vec            resid;
417:   PetscReal      rmax,pwork;
418:   PetscInt       i,n,N;
419:   PetscScalar    *r;

422:   SNESGetFunction(snes,&resid,NULL,NULL);
423:   VecNorm(resid,NORM_INFINITY,&rmax);
424:   VecGetLocalSize(resid,&n);
425:   VecGetSize(resid,&N);
426:   VecGetArray(resid,&r);
427:   pwork = 0.0;
428:   for (i=0; i<n; i++) {
429:     pwork += (PetscAbsScalar(r[i]) > .20*rmax);
430:   }
431:   MPIU_Allreduce(&pwork,per,1,MPIU_REAL,MPIU_SUM,PetscObjectComm((PetscObject)snes));
432:   VecRestoreArray(resid,&r);
433:   *per = *per/N;
434:   return(0);
435: }

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

440:    Collective on SNES

442:    Input Parameters:
443: +  snes   - iterative context
444: .  it    - iteration number
445: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).
446: -  dummy - unused monitor context

448:    Options Database Key:
449: .  -snes_monitor_range - Activates SNESMonitorRange()

451:    Level: intermediate

453: .seealso: SNESMonitorSet(), SNESMonitorDefault(), SNESMonitorLGCreate()
454: @*/
455: PetscErrorCode  SNESMonitorRange(SNES snes,PetscInt it,PetscReal rnorm,PetscViewerAndFormat *vf)
456: {
458:   PetscReal      perc,rel;
459:   PetscViewer    viewer = vf->viewer;
460:   /* should be in a MonitorRangeContext */
461:   static PetscReal prev;

465:   if (!it) prev = rnorm;
466:   SNESMonitorRange_Private(snes,it,&perc);

468:   rel  = (prev - rnorm)/prev;
469:   prev = rnorm;
470:   PetscViewerPushFormat(viewer,vf->format);
471:   PetscViewerASCIIAddTab(viewer,((PetscObject)snes)->tablevel);
472:   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));
473:   PetscViewerASCIISubtractTab(viewer,((PetscObject)snes)->tablevel);
474:   PetscViewerPopFormat(viewer);
475:   return(0);
476: }

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

482:    Collective on SNES

484:    Input Parameters:
485: +  snes - the SNES context
486: .  its - iteration number
487: .  fgnorm - 2-norm of residual (or gradient)
488: -  dummy -  context of monitor

490:    Level: intermediate

492:    Notes:
493:     Insure that SNESMonitorRatio() is called when you set this monitor
494: .seealso: SNESMonitorSet(), SNESMonitorSolution(), SNESMonitorRatio()
495: @*/
496: PetscErrorCode  SNESMonitorRatio(SNES snes,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *vf)
497: {
498:   PetscErrorCode          ierr;
499:   PetscInt                len;
500:   PetscReal               *history;
501:   PetscViewer             viewer = vf->viewer;

504:   SNESGetConvergenceHistory(snes,&history,NULL,&len);
505:   PetscViewerPushFormat(viewer,vf->format);
506:   PetscViewerASCIIAddTab(viewer,((PetscObject)snes)->tablevel);
507:   if (!its || !history || its > len) {
508:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %14.12e \n",its,(double)fgnorm);
509:   } else {
510:     PetscReal ratio = fgnorm/history[its-1];
511:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %14.12e %14.12e \n",its,(double)fgnorm,(double)ratio);
512:   }
513:   PetscViewerASCIISubtractTab(viewer,((PetscObject)snes)->tablevel);
514:   PetscViewerPopFormat(viewer);
515:   return(0);
516: }

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

521:    Collective on SNES

523:    Input Parameters:
524: +   snes - the SNES context
525: -   viewer - the PetscViewer object (ignored)

527:    Level: intermediate

529: .seealso: SNESMonitorSet(), SNESMonitorSolution(), SNESMonitorDefault(), SNESMonitorRatio()
530: @*/
531: PetscErrorCode  SNESMonitorRatioSetUp(SNES snes,PetscViewerAndFormat *vf)
532: {
533:   PetscErrorCode          ierr;
534:   PetscReal               *history;

537:   SNESGetConvergenceHistory(snes,&history,NULL,NULL);
538:   if (!history) {
539:     SNESSetConvergenceHistory(snes,NULL,NULL,100,PETSC_TRUE);
540:   }
541:   return(0);
542: }

544: /* ---------------------------------------------------------------- */
545: /*
546:      Default (short) SNES Monitor, same as SNESMonitorDefault() except
547:   it prints fewer digits of the residual as the residual gets smaller.
548:   This is because the later digits are meaningless and are often
549:   different on different machines; by using this routine different
550:   machines will usually generate the same output.

552:   Deprecated: Intentionally has no manual page
553: */
554: PetscErrorCode  SNESMonitorDefaultShort(SNES snes,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *vf)
555: {
557:   PetscViewer    viewer = vf->viewer;

561:   PetscViewerPushFormat(viewer,vf->format);
562:   PetscViewerASCIIAddTab(viewer,((PetscObject)snes)->tablevel);
563:   if (fgnorm > 1.e-9) {
564:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %g \n",its,(double)fgnorm);
565:   } else if (fgnorm > 1.e-11) {
566:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm %5.3e \n",its,(double)fgnorm);
567:   } else {
568:     PetscViewerASCIIPrintf(viewer,"%3D SNES Function norm < 1.e-11\n",its);
569:   }
570:   PetscViewerASCIISubtractTab(viewer,((PetscObject)snes)->tablevel);
571:   PetscViewerPopFormat(viewer);
572:   return(0);
573: }

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

578:   Collective on SNES

580:   Input Parameters:
581: + snes   - the SNES context
582: . its    - iteration number
583: . fgnorm - 2-norm of residual
584: - ctx    - the PetscViewer

586:   Notes:
587:   This routine uses the DM attached to the residual vector

589:   Level: intermediate

591: .seealso: SNESMonitorSet(), SNESMonitorSolution(), SNESMonitorDefault()
592: @*/
593: PetscErrorCode SNESMonitorDefaultField(SNES snes, PetscInt its, PetscReal fgnorm, PetscViewerAndFormat *vf)
594: {
595:   PetscViewer    viewer = vf->viewer;
596:   Vec            r;
597:   DM             dm;
598:   PetscReal      res[256];
599:   PetscInt       tablevel;

604:   SNESGetFunction(snes, &r, NULL, NULL);
605:   VecGetDM(r, &dm);
606:   if (!dm) {SNESMonitorDefault(snes, its, fgnorm, vf);}
607:   else {
608:     PetscSection s, gs;
609:     PetscInt     Nf, f;

611:     DMGetLocalSection(dm, &s);
612:     DMGetGlobalSection(dm, &gs);
613:     if (!s || !gs) {SNESMonitorDefault(snes, its, fgnorm, vf);}
614:     PetscSectionGetNumFields(s, &Nf);
615:     if (Nf > 256) SETERRQ1(PetscObjectComm((PetscObject) snes), PETSC_ERR_SUP, "Do not support %d fields > 256", Nf);
616:     PetscSectionVecNorm(s, gs, r, NORM_2, res);
617:     PetscObjectGetTabLevel((PetscObject) snes, &tablevel);
618:     PetscViewerPushFormat(viewer,vf->format);
619:     PetscViewerASCIIAddTab(viewer, tablevel);
620:     PetscViewerASCIIPrintf(viewer, "%3D SNES Function norm %14.12e [", its, (double) fgnorm);
621:     for (f = 0; f < Nf; ++f) {
622:       if (f) {PetscViewerASCIIPrintf(viewer, ", ");}
623:       PetscViewerASCIIPrintf(viewer, "%14.12e", res[f]);
624:     }
625:     PetscViewerASCIIPrintf(viewer, "] \n");
626:     PetscViewerASCIISubtractTab(viewer, tablevel);
627:     PetscViewerPopFormat(viewer);
628:   }
629:   return(0);
630: }
631: /* ---------------------------------------------------------------- */
632: /*@C
633:    SNESConvergedDefault - Default onvergence test of the solvers for
634:    systems of nonlinear equations.

636:    Collective on SNES

638:    Input Parameters:
639: +  snes - the SNES context
640: .  it - the iteration (0 indicates before any Newton steps)
641: .  xnorm - 2-norm of current iterate
642: .  snorm - 2-norm of current step
643: .  fnorm - 2-norm of function at current iterate
644: -  dummy - unused context

646:    Output Parameter:
647: .   reason  - one of
648: $  SNES_CONVERGED_FNORM_ABS       - (fnorm < abstol),
649: $  SNES_CONVERGED_SNORM_RELATIVE  - (snorm < stol*xnorm),
650: $  SNES_CONVERGED_FNORM_RELATIVE  - (fnorm < rtol*fnorm0),
651: $  SNES_DIVERGED_FUNCTION_COUNT   - (nfct > maxf),
652: $  SNES_DIVERGED_FNORM_NAN        - (fnorm == NaN),
653: $  SNES_CONVERGED_ITERATING       - (otherwise),
654: $  SNES_DIVERGED_DTOL             - (fnorm > divtol*snes->fnorm0)

656:    where
657: +    maxf - maximum number of function evaluations,  set with SNESSetTolerances()
658: .    nfct - number of function evaluations,
659: .    abstol - absolute function norm tolerance, set with SNESSetTolerances()
660: .    rtol - relative function norm tolerance, set with SNESSetTolerances()
661: .    divtol - divergence tolerance, set with SNESSetDivergenceTolerance()
662: -    fnorm0 - 2-norm of the function at the initial solution (initial guess; zeroth iteration)

664:   Options Database Keys:
665: +  -snes_stol - convergence tolerance in terms of the norm  of the change in the solution between steps
666: .  -snes_atol <abstol> - absolute tolerance of residual norm
667: .  -snes_rtol <rtol> - relative decrease in tolerance norm from the initial 2-norm of the solution
668: .  -snes_divergence_tolerance <divtol> - if the residual goes above divtol*rnorm0, exit with divergence
669: .  -snes_max_funcs <max_funcs> - maximum number of function evaluations
670: .  -snes_max_fail <max_fail> - maximum number of line search failures allowed before stopping, default is none
671: -  -snes_max_linear_solve_fail - number of linear solver failures before SNESSolve() stops

673:    Level: intermediate

675: .seealso: SNESSetConvergenceTest(), SNESConvergedSkip(), SNESSetTolerances(), SNESSetDivergenceTolerance()
676: @*/
677: PetscErrorCode  SNESConvergedDefault(SNES snes,PetscInt it,PetscReal xnorm,PetscReal snorm,PetscReal fnorm,SNESConvergedReason *reason,void *dummy)
678: {


685:   *reason = SNES_CONVERGED_ITERATING;

687:   if (!it) {
688:     /* set parameter for default relative tolerance convergence test */
689:     snes->ttol   = fnorm*snes->rtol;
690:     snes->rnorm0 = fnorm;
691:   }
692:   if (PetscIsInfOrNanReal(fnorm)) {
693:     PetscInfo(snes,"Failed to converged, function norm is NaN\n");
694:     *reason = SNES_DIVERGED_FNORM_NAN;
695:   } else if (fnorm < snes->abstol && (it || !snes->forceiteration)) {
696:     PetscInfo2(snes,"Converged due to function norm %14.12e < %14.12e\n",(double)fnorm,(double)snes->abstol);
697:     *reason = SNES_CONVERGED_FNORM_ABS;
698:   } else if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) {
699:     PetscInfo2(snes,"Exceeded maximum number of function evaluations: %D > %D\n",snes->nfuncs,snes->max_funcs);
700:     *reason = SNES_DIVERGED_FUNCTION_COUNT;
701:   }

703:   if (it && !*reason) {
704:     if (fnorm <= snes->ttol) {
705:       PetscInfo2(snes,"Converged due to function norm %14.12e < %14.12e (relative tolerance)\n",(double)fnorm,(double)snes->ttol);
706:       *reason = SNES_CONVERGED_FNORM_RELATIVE;
707:     } else if (snorm < snes->stol*xnorm) {
708:       PetscInfo3(snes,"Converged due to small update length: %14.12e < %14.12e * %14.12e\n",(double)snorm,(double)snes->stol,(double)xnorm);
709:       *reason = SNES_CONVERGED_SNORM_RELATIVE;
710:     } else if (snes->divtol > 0 && (fnorm > snes->divtol*snes->rnorm0)) {
711:       PetscInfo3(snes,"Diverged due to increase in function norm: %14.12e > %14.12e * %14.12e\n",(double)fnorm,(double)snes->divtol,(double)snes->rnorm0);
712:       *reason = SNES_DIVERGED_DTOL;
713:     }

715:   }
716:   return(0);
717: }

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

723:    Logically Collective on SNES

725:    Input Parameters:
726: +  snes - the SNES context
727: .  it - the iteration (0 indicates before any Newton steps)
728: .  xnorm - 2-norm of current iterate
729: .  snorm - 2-norm of current step
730: .  fnorm - 2-norm of function at current iterate
731: -  dummy - unused context

733:    Output Parameter:
734: .   reason  - SNES_CONVERGED_ITERATING, SNES_CONVERGED_ITS, or SNES_DIVERGED_FNORM_NAN

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

739:    Level: advanced

741: .seealso: SNESConvergedDefault(), SNESSetConvergenceTest()
742: @*/
743: PetscErrorCode  SNESConvergedSkip(SNES snes,PetscInt it,PetscReal xnorm,PetscReal snorm,PetscReal fnorm,SNESConvergedReason *reason,void *dummy)
744: {


751:   *reason = SNES_CONVERGED_ITERATING;

753:   if (fnorm != fnorm) {
754:     PetscInfo(snes,"Failed to converged, function norm is NaN\n");
755:     *reason = SNES_DIVERGED_FNORM_NAN;
756:   } else if (it == snes->max_its) {
757:     *reason = SNES_CONVERGED_ITS;
758:   }
759:   return(0);
760: }

762: /*@C
763:   SNESSetWorkVecs - Gets a number of work vectors.

765:   Input Parameters:
766: + snes  - the SNES context
767: - nw - number of work vectors to allocate

769:   Level: developer

771: @*/
772: PetscErrorCode SNESSetWorkVecs(SNES snes,PetscInt nw)
773: {
774:   DM             dm;
775:   Vec            v;

779:   if (snes->work) {VecDestroyVecs(snes->nwork,&snes->work);}
780:   snes->nwork = nw;

782:   SNESGetDM(snes, &dm);
783:   DMGetGlobalVector(dm, &v);
784:   VecDuplicateVecs(v,snes->nwork,&snes->work);
785:   DMRestoreGlobalVector(dm, &v);
786:   PetscLogObjectParents(snes,nw,snes->work);
787:   return(0);
788: }