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

 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 Key:
 20: .  -snes_monitor_solution [ascii binary draw][:filename][:viewer format] - plots solution at each iteration

 22:    Level: intermediate

 24:    Note:
 25:    This is not called directly by users, rather one calls `SNESMonitorSet()`, with this function as an argument, to cause the monitor
 26:    to be used during the `SNESSolve()`

 28: .seealso: [](chapter_snes), `SNES`, `SNESMonitorSet()`, `SNESMonitorDefault()`, `VecView()`
 29: @*/
 30: PetscErrorCode SNESMonitorSolution(SNES snes, PetscInt its, PetscReal fgnorm, PetscViewerAndFormat *vf)
 31: {
 32:   Vec         x;
 33:   PetscViewer viewer = vf->viewer;

 35:   PetscFunctionBegin;
 37:   PetscCall(SNESGetSolution(snes, &x));
 38:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
 39:   PetscCall(VecView(x, viewer));
 40:   PetscCall(PetscViewerPopFormat(viewer));
 41:   PetscFunctionReturn(PETSC_SUCCESS);
 42: }

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

 48:    Collective

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

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

 59:    Level: intermediate

 61:    Note:
 62:    This is not called directly by users, rather one calls `SNESMonitorSet()`, with this function as an argument, to cause the monitor
 63:    to be used during the `SNES` solve.

 65: .seealso: [](chapter_snes), `SNES`, `SNESMonitorSet()`, `SNESMonitorDefault()`, `VecView()`, `SNESMonitor()`
 66: @*/
 67: PetscErrorCode SNESMonitorResidual(SNES snes, PetscInt its, PetscReal fgnorm, PetscViewerAndFormat *vf)
 68: {
 69:   Vec x;

 71:   PetscFunctionBegin;
 73:   PetscCall(SNESGetFunction(snes, &x, NULL, NULL));
 74:   PetscCall(PetscViewerPushFormat(vf->viewer, vf->format));
 75:   PetscCall(VecView(x, vf->viewer));
 76:   PetscCall(PetscViewerPopFormat(vf->viewer));
 77:   PetscFunctionReturn(PETSC_SUCCESS);
 78: }

 80: /*@C
 81:    SNESMonitorSolutionUpdate - Monitors progress of the `SNES` solvers by calling
 82:    `VecView()` for the UPDATE to the solution at each iteration.

 84:    Collective

 86:    Input Parameters:
 87: +  snes - the `SNES` context
 88: .  its - iteration number
 89: .  fgnorm - 2-norm of residual
 90: -  dummy - a viewer

 92:    Options Database Key:
 93: .  -snes_monitor_solution_update [ascii binary draw][:filename][:viewer format] - plots update to solution at each iteration

 95:    Level: intermediate

 97:    Note:
 98:    This is not called directly by users, rather one calls `SNESMonitorSet()`, with this function as an argument, to cause the monitor
 99:    to be used during the `SNES` solve.

101: .seealso: [](chapter_snes), `SNESMonitorSet()`, `SNESMonitorDefault()`, `VecView()`, `SNESMonitor()`, `SNESMonitor()`
102: @*/
103: PetscErrorCode SNESMonitorSolutionUpdate(SNES snes, PetscInt its, PetscReal fgnorm, PetscViewerAndFormat *vf)
104: {
105:   Vec         x;
106:   PetscViewer viewer = vf->viewer;

108:   PetscFunctionBegin;
110:   PetscCall(SNESGetSolutionUpdate(snes, &x));
111:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
112:   PetscCall(VecView(x, viewer));
113:   PetscCall(PetscViewerPopFormat(viewer));
114:   PetscFunctionReturn(PETSC_SUCCESS);
115: }

117: #include <petscdraw.h>

119: /*@C
120:   KSPMonitorSNESResidual - Prints the `SNES` residual norm, as well as the `KSP` residual norm, at each iteration of an iterative solver.

122:   Collective

124:   Input Parameters:
125: + ksp   - iterative context
126: . n     - iteration number
127: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
128: - vf    - The viewer context

130:   Options Database Key:
131: . -snes_monitor_ksp - Activates `KSPMonitorSNESResidual()`

133:   Level: intermediate

135:    Note:
136:    This is not called directly by users, rather one calls `KSPMonitorSet()`, with this function as an argument, to cause the monitor
137:    to be used during the `KSP` solve.

139: .seealso: [](chapter_snes), `KSPMonitorSet()`, `KSPMonitorResidual()`, `KSPMonitorTrueResidualMaxNorm()`, `KSPMonitor()`, `SNESMonitor()`
140: @*/
141: PetscErrorCode KSPMonitorSNESResidual(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
142: {
143:   PetscViewer       viewer = vf->viewer;
144:   PetscViewerFormat format = vf->format;
145:   SNES              snes   = (SNES)vf->data;
146:   Vec               snes_solution, work1, work2;
147:   PetscReal         snorm;
148:   PetscInt          tablevel;
149:   const char       *prefix;

151:   PetscFunctionBegin;
153:   PetscCall(SNESGetSolution(snes, &snes_solution));
154:   PetscCall(VecDuplicate(snes_solution, &work1));
155:   PetscCall(VecDuplicate(snes_solution, &work2));
156:   PetscCall(KSPBuildSolution(ksp, work1, NULL));
157:   PetscCall(VecAYPX(work1, -1.0, snes_solution));
158:   PetscCall(SNESComputeFunction(snes, work1, work2));
159:   PetscCall(VecNorm(work2, NORM_2, &snorm));
160:   PetscCall(VecDestroy(&work1));
161:   PetscCall(VecDestroy(&work2));

163:   PetscCall(PetscObjectGetTabLevel((PetscObject)ksp, &tablevel));
164:   PetscCall(PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix));
165:   PetscCall(PetscViewerPushFormat(viewer, format));
166:   PetscCall(PetscViewerASCIIAddTab(viewer, tablevel));
167:   if (n == 0 && prefix) PetscCall(PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix));
168:   PetscCall(PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " SNES Residual norm %5.3e KSP Residual norm %5.3e \n", n, (double)snorm, (double)rnorm));
169:   PetscCall(PetscViewerASCIISubtractTab(viewer, tablevel));
170:   PetscCall(PetscViewerPopFormat(viewer));
171:   PetscFunctionReturn(PETSC_SUCCESS);
172: }

174: /*@C
175:   KSPMonitorSNESResidualDrawLG - Plots the linear `KSP` residual norm and the `SNES` residual norm at each iteration of an iterative solver.

177:   Collective

179:   Input Parameters:
180: + ksp   - iterative context
181: . n     - iteration number
182: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
183: - vf    - The viewer context, created with `KSPMonitorSNESResidualDrawLGCreate()`

185:   Options Database Key:
186: . -snes_monitor_ksp draw::draw_lg - Activates `KSPMonitorSNESResidualDrawLG()`

188:   Level: intermediate

190:    Note:
191:    This is not called directly by users, rather one calls `SNESMonitorSet()`, with this function as an argument, to cause the monitor
192:    to be used during the `SNESSolve()`

194: .seealso: [](chapter_snes), `KSPMonitorSet()`, `KSPMonitorTrueResidual()`, `SNESMonitor()`, `KSPMonitor()`, `KSPMonitorSNESResidualDrawLGCreate()`
195: @*/
196: PetscErrorCode KSPMonitorSNESResidualDrawLG(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
197: {
198:   PetscViewer        viewer = vf->viewer;
199:   PetscViewerFormat  format = vf->format;
200:   PetscDrawLG        lg     = vf->lg;
201:   SNES               snes   = (SNES)vf->data;
202:   Vec                snes_solution, work1, work2;
203:   PetscReal          snorm;
204:   KSPConvergedReason reason;
205:   PetscReal          x[2], y[2];

207:   PetscFunctionBegin;
210:   PetscCall(SNESGetSolution(snes, &snes_solution));
211:   PetscCall(VecDuplicate(snes_solution, &work1));
212:   PetscCall(VecDuplicate(snes_solution, &work2));
213:   PetscCall(KSPBuildSolution(ksp, work1, NULL));
214:   PetscCall(VecAYPX(work1, -1.0, snes_solution));
215:   PetscCall(SNESComputeFunction(snes, work1, work2));
216:   PetscCall(VecNorm(work2, NORM_2, &snorm));
217:   PetscCall(VecDestroy(&work1));
218:   PetscCall(VecDestroy(&work2));

220:   PetscCall(PetscViewerPushFormat(viewer, format));
221:   if (!n) PetscCall(PetscDrawLGReset(lg));
222:   x[0] = (PetscReal)n;
223:   if (rnorm > 0.0) y[0] = PetscLog10Real(rnorm);
224:   else y[0] = -15.0;
225:   x[1] = (PetscReal)n;
226:   if (snorm > 0.0) y[1] = PetscLog10Real(snorm);
227:   else y[1] = -15.0;
228:   PetscCall(PetscDrawLGAddPoint(lg, x, y));
229:   PetscCall(KSPGetConvergedReason(ksp, &reason));
230:   if (n <= 20 || !(n % 5) || reason) {
231:     PetscCall(PetscDrawLGDraw(lg));
232:     PetscCall(PetscDrawLGSave(lg));
233:   }
234:   PetscCall(PetscViewerPopFormat(viewer));
235:   PetscFunctionReturn(PETSC_SUCCESS);
236: }

238: /*@C
239:   KSPMonitorSNESResidualDrawLGCreate - Creates the `PetscViewer` used by `KSPMonitorSNESResidualDrawLG()`

241:   Collective

243:   Input Parameters:
244: + viewer - The `PetscViewer`
245: . format - The viewer format
246: - ctx    - An optional user context

248:   Output Parameter:
249: . vf    - The viewer context

251:   Level: intermediate

253: .seealso: [](chapter_snes), `KSP`, `SNES`, `KSPMonitorSet()`, `KSPMonitorTrueResidual()`
254: @*/
255: PetscErrorCode KSPMonitorSNESResidualDrawLGCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
256: {
257:   const char *names[] = {"linear", "nonlinear"};

259:   PetscFunctionBegin;
260:   PetscCall(PetscViewerAndFormatCreate(viewer, format, vf));
261:   (*vf)->data = ctx;
262:   PetscCall(KSPMonitorLGCreate(PetscObjectComm((PetscObject)viewer), NULL, NULL, "Log Residual Norm", 2, names, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &(*vf)->lg));
263:   PetscFunctionReturn(PETSC_SUCCESS);
264: }

266: PetscErrorCode SNESMonitorDefaultSetUp(SNES snes, PetscViewerAndFormat *vf)
267: {
268:   PetscFunctionBegin;
269:   if (vf->format == PETSC_VIEWER_DRAW_LG) PetscCall(KSPMonitorLGCreate(PetscObjectComm((PetscObject)vf->viewer), NULL, NULL, "Log Residual Norm", 1, NULL, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &vf->lg));
270:   PetscFunctionReturn(PETSC_SUCCESS);
271: }

273: /*@C
274:    SNESMonitorDefault - Monitors progress of the `SNES` solvers (default).

276:    Collective

278:    Input Parameters:
279: +  snes - the `SNES` context
280: .  its - iteration number
281: .  fgnorm - 2-norm of residual
282: -  vf - viewer and format structure

284:    Options Database Key:
285: .  -snes_monitor - use this function to monitor the convergence of the nonlinear solver

287:    Level: intermediate

289:    Notes:
290:    This routine prints the residual norm at each iteration.

292:    This is not called directly by users, rather one calls `SNESMonitorSet()`, with this function as an argument, to cause the monitor
293:    to be used during the `SNES` solve.

295: .seealso: [](chapter_snes), `SNESMonitorSet()`, `SNESMonitorSolution()`, `SNESMonitorFunction()`, `SNESMonitorSolution()`, `SNESMonitorResidual()`,
296:           `SNESMonitorSolutionUpdate()`, `SNESMonitorDefault()`, `SNESMonitorScaling()`, `SNESMonitorRange()`, `SNESMonitorRatio()`,
297:           `SNESMonitorDefaultField()`
298: @*/
299: PetscErrorCode SNESMonitorDefault(SNES snes, PetscInt its, PetscReal fgnorm, PetscViewerAndFormat *vf)
300: {
301:   PetscViewer       viewer = vf->viewer;
302:   PetscViewerFormat format = vf->format;
303:   PetscBool         isascii, isdraw;

305:   PetscFunctionBegin;
307:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
308:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERDRAW, &isdraw));
309:   PetscCall(PetscViewerPushFormat(viewer, format));
310:   if (isascii) {
311:     PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)snes)->tablevel));
312:     PetscCall(PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " SNES Function norm %14.12e \n", its, (double)fgnorm));
313:     PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)snes)->tablevel));
314:   } else if (isdraw) {
315:     if (format == PETSC_VIEWER_DRAW_LG) {
316:       PetscDrawLG lg = (PetscDrawLG)vf->lg;
317:       PetscReal   x, y;

320:       if (!its) PetscCall(PetscDrawLGReset(lg));
321:       x = (PetscReal)its;
322:       if (fgnorm > 0.0) y = PetscLog10Real(fgnorm);
323:       else y = -15.0;
324:       PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
325:       if (its <= 20 || !(its % 5) || snes->reason) {
326:         PetscCall(PetscDrawLGDraw(lg));
327:         PetscCall(PetscDrawLGSave(lg));
328:       }
329:     }
330:   }
331:   PetscCall(PetscViewerPopFormat(viewer));
332:   PetscFunctionReturn(PETSC_SUCCESS);
333: }

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

338:    Collective

340:    Input Parameters:
341: +  snes - the `SNES` context
342: .  its - iteration number
343: .  fgnorm - 2-norm of residual
344: -  vf - viewer and format structure

346:    Level: intermediate

348:    Notes:
349:    This routine prints the largest value in each row of the Jacobian

351:    This is not called directly by users, rather one calls `SNESMonitorSet()`, with this function as an argument, to cause the monitor
352:    to be used during the `SNES` solve.

354: .seealso: [](chapter_snes), `SNESMonitorSet()`, `SNESMonitorSolution()`, `SNESMonitorRange()`, `SNESMonitorJacUpdateSpectrum()`
355: @*/
356: PetscErrorCode SNESMonitorScaling(SNES snes, PetscInt its, PetscReal fgnorm, PetscViewerAndFormat *vf)
357: {
358:   PetscViewer viewer = vf->viewer;
359:   KSP         ksp;
360:   Mat         J;
361:   Vec         v;

363:   PetscFunctionBegin;
365:   PetscCall(SNESGetKSP(snes, &ksp));
366:   PetscCall(KSPGetOperators(ksp, &J, NULL));
367:   PetscCall(MatCreateVecs(J, &v, NULL));
368:   PetscCall(MatGetRowMaxAbs(J, v, NULL));
369:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
370:   PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)snes)->tablevel));
371:   PetscCall(PetscViewerASCIIPrintf(viewer, "SNES Jacobian maximum row entries\n"));
372:   PetscCall(VecView(v, viewer));
373:   PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)snes)->tablevel));
374:   PetscCall(PetscViewerPopFormat(viewer));
375:   PetscCall(VecDestroy(&v));
376:   PetscFunctionReturn(PETSC_SUCCESS);
377: }

379: /*@C
380:    SNESMonitorJacUpdateSpectrum - Monitors the spectrun of the change in the Jacobian from the last Jacobian evaluation

382:    Collective

384:    Input Parameters:
385: +  snes - the `SNES` context
386: .  its - iteration number
387: .  fgnorm - 2-norm of residual
388: -  vf - viewer and format structure

390:    Option Database Key:
391: .  -snes_monitor_jacupdate_spectrum - activates this monitor

393:    Level: intermediate

395:    Notes:
396:    This routine prints the eigenvalues of the difference in the Jacobians

398:    This is not called directly by users, rather one calls `SNESMonitorSet()`, with this function as an argument, to cause the monitor
399:    to be used during the `SNES` solve.

401: .seealso: [](chapter_snes), `SNESMonitorSet()`, `SNESMonitorSolution()`, `SNESMonitorRange()`
402: @*/
403: PetscErrorCode SNESMonitorJacUpdateSpectrum(SNES snes, PetscInt it, PetscReal fnorm, PetscViewerAndFormat *vf)
404: {
405:   Vec X;
406:   Mat J, dJ, dJdense;
407:   PetscErrorCode (*func)(SNES, Vec, Mat, Mat, void *);
408:   PetscInt     n;
409:   PetscBLASInt nb = 0, lwork;
410:   PetscReal   *eigr, *eigi;
411:   PetscScalar *work;
412:   PetscScalar *a;

414:   PetscFunctionBegin;
415:   if (it == 0) PetscFunctionReturn(PETSC_SUCCESS);
416:   /* create the difference between the current update and the current Jacobian */
417:   PetscCall(SNESGetSolution(snes, &X));
418:   PetscCall(SNESGetJacobian(snes, NULL, &J, &func, NULL));
419:   PetscCall(MatDuplicate(J, MAT_COPY_VALUES, &dJ));
420:   PetscCall(SNESComputeJacobian(snes, X, dJ, dJ));
421:   PetscCall(MatAXPY(dJ, -1.0, J, SAME_NONZERO_PATTERN));

423:   /* compute the spectrum directly */
424:   PetscCall(MatConvert(dJ, MATSEQDENSE, MAT_INITIAL_MATRIX, &dJdense));
425:   PetscCall(MatGetSize(dJ, &n, NULL));
426:   PetscCall(PetscBLASIntCast(n, &nb));
427:   lwork = 3 * nb;
428:   PetscCall(PetscMalloc1(n, &eigr));
429:   PetscCall(PetscMalloc1(n, &eigi));
430:   PetscCall(PetscMalloc1(lwork, &work));
431:   PetscCall(MatDenseGetArray(dJdense, &a));
432: #if !defined(PETSC_USE_COMPLEX)
433:   {
434:     PetscBLASInt lierr;
435:     PetscInt     i;
436:     PetscCall(PetscFPTrapPush(PETSC_FP_TRAP_OFF));
437:     PetscCallBLAS("LAPACKgeev", LAPACKgeev_("N", "N", &nb, a, &nb, eigr, eigi, NULL, &nb, NULL, &nb, work, &lwork, &lierr));
438:     PetscCheck(!lierr, PETSC_COMM_SELF, PETSC_ERR_LIB, "geev() error %" PetscBLASInt_FMT, lierr);
439:     PetscCall(PetscFPTrapPop());
440:     PetscCall(PetscPrintf(PetscObjectComm((PetscObject)snes), "Eigenvalues of J_%" PetscInt_FMT " - J_%" PetscInt_FMT ":\n", it, it - 1));
441:     for (i = 0; i < n; i++) PetscCall(PetscPrintf(PetscObjectComm((PetscObject)snes), "%5" PetscInt_FMT ": %20.5g + %20.5gi\n", i, (double)eigr[i], (double)eigi[i]));
442:   }
443:   PetscCall(MatDenseRestoreArray(dJdense, &a));
444:   PetscCall(MatDestroy(&dJ));
445:   PetscCall(MatDestroy(&dJdense));
446:   PetscCall(PetscFree(eigr));
447:   PetscCall(PetscFree(eigi));
448:   PetscCall(PetscFree(work));
449:   PetscFunctionReturn(PETSC_SUCCESS);
450: #else
451:   SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Not coded for complex");
452: #endif
453: }

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

457: PetscErrorCode SNESMonitorRange_Private(SNES snes, PetscInt it, PetscReal *per)
458: {
459:   Vec          resid;
460:   PetscReal    rmax, pwork;
461:   PetscInt     i, n, N;
462:   PetscScalar *r;

464:   PetscFunctionBegin;
465:   PetscCall(SNESGetFunction(snes, &resid, NULL, NULL));
466:   PetscCall(VecNorm(resid, NORM_INFINITY, &rmax));
467:   PetscCall(VecGetLocalSize(resid, &n));
468:   PetscCall(VecGetSize(resid, &N));
469:   PetscCall(VecGetArray(resid, &r));
470:   pwork = 0.0;
471:   for (i = 0; i < n; i++) pwork += (PetscAbsScalar(r[i]) > .20 * rmax);
472:   PetscCall(MPIU_Allreduce(&pwork, per, 1, MPIU_REAL, MPIU_SUM, PetscObjectComm((PetscObject)snes)));
473:   PetscCall(VecRestoreArray(resid, &r));
474:   *per = *per / N;
475:   PetscFunctionReturn(PETSC_SUCCESS);
476: }

478: /*@C
479:    SNESMonitorRange - Prints the percentage of residual elements that are more then 10 percent of the maximum entry in the residual

481:    Collective

483:    Input Parameters:
484: +  snes   - iterative context
485: .  it    - iteration number
486: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).
487: -  dummy - unused monitor context

489:    Options Database Key:
490: .  -snes_monitor_range - Activates `SNESMonitorRange()`

492:    Level: intermediate

494:    Note:
495:    This is not called directly by users, rather one calls `SNESMonitorSet()`, with this function as an argument, to cause the monitor
496:    to be used during the `SNES` solve.

498: .seealso: [](chapter_snes), `SNESMonitorSet()`, `SNESMonitorDefault()`, `SNESMonitorLGCreate()`, `SNESMonitorScaling()`
499: @*/
500: PetscErrorCode SNESMonitorRange(SNES snes, PetscInt it, PetscReal rnorm, PetscViewerAndFormat *vf)
501: {
502:   PetscReal   perc, rel;
503:   PetscViewer viewer = vf->viewer;
504:   /* should be in a MonitorRangeContext */
505:   static PetscReal prev;

507:   PetscFunctionBegin;
509:   if (!it) prev = rnorm;
510:   PetscCall(SNESMonitorRange_Private(snes, it, &perc));

512:   rel  = (prev - rnorm) / prev;
513:   prev = rnorm;
514:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
515:   PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)snes)->tablevel));
516:   PetscCall(PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " SNES preconditioned resid norm %14.12e Percent values above 20 percent of maximum %5.2g relative decrease %5.2e ratio %5.2e \n", it, (double)rnorm, (double)(100.0 * perc), (double)rel, (double)(rel / perc)));
517:   PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)snes)->tablevel));
518:   PetscCall(PetscViewerPopFormat(viewer));
519:   PetscFunctionReturn(PETSC_SUCCESS);
520: }

522: /*@C
523:    SNESMonitorRatio - Monitors progress of the `SNES` solvers by printing the ratio
524:    of residual norm at each iteration to the previous.

526:    Collective

528:    Input Parameters:
529: +  snes - the `SNES` context
530: .  its - iteration number
531: .  fgnorm - 2-norm of residual (or gradient)
532: -  dummy -  context of monitor

534:    Option Database Key:
535: .  -snes_monitor_ratio - activate this monitor

537:    Level: intermediate

539:    Notes:
540:    This is not called directly by users, rather one calls `SNESMonitorSet()`, with this function as an argument, to cause the monitor
541:    to be used during the `SNES` solve.

543:    Be sure to call `SNESMonitorRationSetUp()` before using this monitor.

545: .seealso: [](chapter_snes), `SNESMonitorRationSetUp()`, `SNESMonitorSet()`, `SNESMonitorSolution()`, `SNESMonitorDefault()`
546: @*/
547: PetscErrorCode SNESMonitorRatio(SNES snes, PetscInt its, PetscReal fgnorm, PetscViewerAndFormat *vf)
548: {
549:   PetscInt    len;
550:   PetscReal  *history;
551:   PetscViewer viewer = vf->viewer;

553:   PetscFunctionBegin;
554:   PetscCall(SNESGetConvergenceHistory(snes, &history, NULL, &len));
555:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
556:   PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)snes)->tablevel));
557:   if (!its || !history || its > len) {
558:     PetscCall(PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " SNES Function norm %14.12e \n", its, (double)fgnorm));
559:   } else {
560:     PetscReal ratio = fgnorm / history[its - 1];
561:     PetscCall(PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " SNES Function norm %14.12e %14.12e \n", its, (double)fgnorm, (double)ratio));
562:   }
563:   PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)snes)->tablevel));
564:   PetscCall(PetscViewerPopFormat(viewer));
565:   PetscFunctionReturn(PETSC_SUCCESS);
566: }

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

571:    Collective

573:    Input Parameters:
574: +   snes - the `SNES` context
575: -   viewer - the `PetscViewer` object (ignored)

577:    Level: intermediate

579: .seealso: [](chapter_snes), `SNESMonitorSet()`, `SNESMonitorSolution()`, `SNESMonitorDefault()`, `SNESMonitorRatio()`
580: @*/
581: PetscErrorCode SNESMonitorRatioSetUp(SNES snes, PetscViewerAndFormat *vf)
582: {
583:   PetscReal *history;

585:   PetscFunctionBegin;
586:   PetscCall(SNESGetConvergenceHistory(snes, &history, NULL, NULL));
587:   if (!history) PetscCall(SNESSetConvergenceHistory(snes, NULL, NULL, 100, PETSC_TRUE));
588:   PetscFunctionReturn(PETSC_SUCCESS);
589: }

591: /*
592:      Default (short) SNES Monitor, same as SNESMonitorDefault() except
593:   it prints fewer digits of the residual as the residual gets smaller.
594:   This is because the later digits are meaningless and are often
595:   different on different machines; by using this routine different
596:   machines will usually generate the same output.

598:   Deprecated: Intentionally has no manual page
599: */
600: PetscErrorCode SNESMonitorDefaultShort(SNES snes, PetscInt its, PetscReal fgnorm, PetscViewerAndFormat *vf)
601: {
602:   PetscViewer viewer = vf->viewer;

604:   PetscFunctionBegin;
606:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
607:   PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)snes)->tablevel));
608:   if (fgnorm > 1.e-9) {
609:     PetscCall(PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " SNES Function norm %g \n", its, (double)fgnorm));
610:   } else if (fgnorm > 1.e-11) {
611:     PetscCall(PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " SNES Function norm %5.3e \n", its, (double)fgnorm));
612:   } else {
613:     PetscCall(PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " SNES Function norm < 1.e-11\n", its));
614:   }
615:   PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)snes)->tablevel));
616:   PetscCall(PetscViewerPopFormat(viewer));
617:   PetscFunctionReturn(PETSC_SUCCESS);
618: }

620: /*@C
621:   SNESMonitorDefaultField - Monitors progress of the `SNES` solvers, separated into fields.

623:   Collective

625:   Input Parameters:
626: + snes   - the `SNES` context
627: . its    - iteration number
628: . fgnorm - 2-norm of residual
629: - ctx    - the PetscViewer

631:    Option Database Key:
632: .  -snes_monitor_field - activate this monitor

634:   Level: intermediate

636:   Notes:
637:   This routine uses the `DM` attached to the residual vector to define the fields.

639:   This is not called directly by users, rather one calls `SNESMonitorSet()`, with this function as an argument, to cause the monitor
640:   to be used during the `SNES` solve.

642: .seealso: [](chapter_snes), `SNESMonitorSet()`, `SNESMonitorSolution()`, `SNESMonitorDefault()`
643: @*/
644: PetscErrorCode SNESMonitorDefaultField(SNES snes, PetscInt its, PetscReal fgnorm, PetscViewerAndFormat *vf)
645: {
646:   PetscViewer viewer = vf->viewer;
647:   Vec         r;
648:   DM          dm;
649:   PetscReal   res[256];
650:   PetscInt    tablevel;

652:   PetscFunctionBegin;
654:   PetscCall(SNESGetFunction(snes, &r, NULL, NULL));
655:   PetscCall(VecGetDM(r, &dm));
656:   if (!dm) PetscCall(SNESMonitorDefault(snes, its, fgnorm, vf));
657:   else {
658:     PetscSection s, gs;
659:     PetscInt     Nf, f;

661:     PetscCall(DMGetLocalSection(dm, &s));
662:     PetscCall(DMGetGlobalSection(dm, &gs));
663:     if (!s || !gs) PetscCall(SNESMonitorDefault(snes, its, fgnorm, vf));
664:     PetscCall(PetscSectionGetNumFields(s, &Nf));
665:     PetscCheck(Nf <= 256, PetscObjectComm((PetscObject)snes), PETSC_ERR_SUP, "Do not support %" PetscInt_FMT " fields > 256", Nf);
666:     PetscCall(PetscSectionVecNorm(s, gs, r, NORM_2, res));
667:     PetscCall(PetscObjectGetTabLevel((PetscObject)snes, &tablevel));
668:     PetscCall(PetscViewerPushFormat(viewer, vf->format));
669:     PetscCall(PetscViewerASCIIAddTab(viewer, tablevel));
670:     PetscCall(PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " SNES Function norm %14.12e [", its, (double)fgnorm));
671:     for (f = 0; f < Nf; ++f) {
672:       if (f) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
673:       PetscCall(PetscViewerASCIIPrintf(viewer, "%14.12e", (double)res[f]));
674:     }
675:     PetscCall(PetscViewerASCIIPrintf(viewer, "] \n"));
676:     PetscCall(PetscViewerASCIISubtractTab(viewer, tablevel));
677:     PetscCall(PetscViewerPopFormat(viewer));
678:   }
679:   PetscFunctionReturn(PETSC_SUCCESS);
680: }

682: /*@C
683:    SNESConvergedDefault - Default onvergence test of the solvers for
684:    systems of nonlinear equations.

686:    Collective

688:    Input Parameters:
689: +  snes - the `SNES` context
690: .  it - the iteration (0 indicates before any Newton steps)
691: .  xnorm - 2-norm of current iterate
692: .  snorm - 2-norm of current step
693: .  fnorm - 2-norm of function at current iterate
694: -  dummy - unused context

696:    Output Parameter:
697: .   reason  - one of
698: .vb
699:    SNES_CONVERGED_FNORM_ABS       - (fnorm < abstol),
700:    SNES_CONVERGED_SNORM_RELATIVE  - (snorm < stol*xnorm),
701:    SNES_CONVERGED_FNORM_RELATIVE  - (fnorm < rtol*fnorm0),
702:    SNES_DIVERGED_FUNCTION_COUNT   - (nfct > maxf),
703:    SNES_DIVERGED_FNORM_NAN        - (fnorm == NaN),
704:    SNES_CONVERGED_ITERATING       - (otherwise),
705:    SNES_DIVERGED_DTOL             - (fnorm > divtol*snes->fnorm0)
706: .ve

708:    where
709: +    maxf - maximum number of function evaluations,  set with `SNESSetTolerances()`
710: .    nfct - number of function evaluations,
711: .    abstol - absolute function norm tolerance, set with `SNESSetTolerances()`
712: .    rtol - relative function norm tolerance, set with `SNESSetTolerances()`
713: .    divtol - divergence tolerance, set with `SNESSetDivergenceTolerance()`
714: -    fnorm0 - 2-norm of the function at the initial solution (initial guess; zeroth iteration)

716:   Options Database Keys:
717: +  -snes_convergence_test default - see `SNESSetFromOptions()`
718: .  -snes_stol - convergence tolerance in terms of the norm  of the change in the solution between steps
719: .  -snes_atol <abstol> - absolute tolerance of residual norm
720: .  -snes_rtol <rtol> - relative decrease in tolerance norm from the initial 2-norm of the solution
721: .  -snes_divergence_tolerance <divtol> - if the residual goes above divtol*rnorm0, exit with divergence
722: .  -snes_max_funcs <max_funcs> - maximum number of function evaluations
723: .  -snes_max_fail <max_fail> - maximum number of line search failures allowed before stopping, default is none
724: -  -snes_max_linear_solve_fail - number of linear solver failures before `SNESSolve()` stops

726:    Level: intermediate

728: .seealso: [](chapter_snes), `SNES`, `SNESSolve()`, `SNESSetConvergenceTest()`, `SNESConvergedSkip()`, `SNESSetTolerances()`, `SNESSetDivergenceTolerance()`
729: @*/
730: PetscErrorCode SNESConvergedDefault(SNES snes, PetscInt it, PetscReal xnorm, PetscReal snorm, PetscReal fnorm, SNESConvergedReason *reason, void *dummy)
731: {
732:   PetscFunctionBegin;

736:   *reason = SNES_CONVERGED_ITERATING;
737:   if (!it) {
738:     /* set parameter for default relative tolerance convergence test */
739:     snes->ttol   = fnorm * snes->rtol;
740:     snes->rnorm0 = fnorm;
741:   }
742:   if (PetscIsInfOrNanReal(fnorm)) {
743:     PetscCall(PetscInfo(snes, "Failed to converged, function norm is NaN\n"));
744:     *reason = SNES_DIVERGED_FNORM_NAN;
745:   } else if (fnorm < snes->abstol && (it || !snes->forceiteration)) {
746:     PetscCall(PetscInfo(snes, "Converged due to function norm %14.12e < %14.12e\n", (double)fnorm, (double)snes->abstol));
747:     *reason = SNES_CONVERGED_FNORM_ABS;
748:   } else if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) {
749:     PetscCall(PetscInfo(snes, "Exceeded maximum number of function evaluations: %" PetscInt_FMT " > %" PetscInt_FMT "\n", snes->nfuncs, snes->max_funcs));
750:     *reason = SNES_DIVERGED_FUNCTION_COUNT;
751:   }

753:   if (it && !*reason) {
754:     if (fnorm <= snes->ttol) {
755:       PetscCall(PetscInfo(snes, "Converged due to function norm %14.12e < %14.12e (relative tolerance)\n", (double)fnorm, (double)snes->ttol));
756:       *reason = SNES_CONVERGED_FNORM_RELATIVE;
757:     } else if (snorm < snes->stol * xnorm) {
758:       PetscCall(PetscInfo(snes, "Converged due to small update length: %14.12e < %14.12e * %14.12e\n", (double)snorm, (double)snes->stol, (double)xnorm));
759:       *reason = SNES_CONVERGED_SNORM_RELATIVE;
760:     } else if (snes->divtol > 0 && (fnorm > snes->divtol * snes->rnorm0)) {
761:       PetscCall(PetscInfo(snes, "Diverged due to increase in function norm: %14.12e > %14.12e * %14.12e\n", (double)fnorm, (double)snes->divtol, (double)snes->rnorm0));
762:       *reason = SNES_DIVERGED_DTOL;
763:     }
764:   }
765:   PetscFunctionReturn(PETSC_SUCCESS);
766: }

768: /*@C
769:    SNESConvergedSkip - Convergence test for `SNES` that NEVER returns as
770:    converged, UNLESS the maximum number of iteration have been reached.

772:    Logically Collective

774:    Input Parameters:
775: +  snes - the `SNES` context
776: .  it - the iteration (0 indicates before any Newton steps)
777: .  xnorm - 2-norm of current iterate
778: .  snorm - 2-norm of current step
779: .  fnorm - 2-norm of function at current iterate
780: -  dummy - unused context

782:    Output Parameter:
783: .   reason  - `SNES_CONVERGED_ITERATING`, `SNES_CONVERGED_ITS`, or `SNES_DIVERGED_FNORM_NAN`

785:    Options Database Key:
786: .  -snes_convergence_test skip - see `SNESSetFromOptions()`

788:    Level: advanced

790: .seealso: [](chapter_snes), `SNES`, `SNESSolve()`, `SNESConvergedDefault()`, `SNESSetConvergenceTest()`
791: @*/
792: PetscErrorCode SNESConvergedSkip(SNES snes, PetscInt it, PetscReal xnorm, PetscReal snorm, PetscReal fnorm, SNESConvergedReason *reason, void *dummy)
793: {
794:   PetscFunctionBegin;

798:   *reason = SNES_CONVERGED_ITERATING;

800:   if (fnorm != fnorm) {
801:     PetscCall(PetscInfo(snes, "Failed to converged, function norm is NaN\n"));
802:     *reason = SNES_DIVERGED_FNORM_NAN;
803:   } else if (it == snes->max_its) {
804:     *reason = SNES_CONVERGED_ITS;
805:   }
806:   PetscFunctionReturn(PETSC_SUCCESS);
807: }

809: /*@C
810:   SNESSetWorkVecs - Gets a number of work vectors to be used internally by `SNES` solvers

812:   Input Parameters:
813: + snes  - the `SNES` context
814: - nw - number of work vectors to allocate

816:   Level: developer

818: .seealso: [](chapter_snes), `SNES`
819: @*/
820: PetscErrorCode SNESSetWorkVecs(SNES snes, PetscInt nw)
821: {
822:   DM  dm;
823:   Vec v;

825:   PetscFunctionBegin;
826:   if (snes->work) PetscCall(VecDestroyVecs(snes->nwork, &snes->work));
827:   snes->nwork = nw;

829:   PetscCall(SNESGetDM(snes, &dm));
830:   PetscCall(DMGetGlobalVector(dm, &v));
831:   PetscCall(VecDuplicateVecs(v, snes->nwork, &snes->work));
832:   PetscCall(DMRestoreGlobalVector(dm, &v));
833:   PetscFunctionReturn(PETSC_SUCCESS);
834: }